LArSoft  v10_04_05
Liquid Argon Software toolkit - https://larsoft.org/
LArStackingAction Class Reference

#include "LArStackingAction.h"

Inheritance diagram for LArStackingAction:

Public Member Functions

 LArStackingAction (int)
 
virtual G4ClassificationOfNewTrack ClassifyNewTrack (const G4Track *aTrack)
 
virtual void NewStage ()
 
virtual void PrepareNewEvent ()
 
void SetNRequestMuon (G4int val)
 
G4int GetNRequestMuon () const
 
void SetNRequestIsoMuon (G4int val)
 
G4int GetNRequestIsoMuon () const
 
void SetNIsolation (G4int val)
 
G4int GetNIsolation () const
 
void SetRoIAngle (G4double val)
 
G4double GetRoIAngle () const
 

Public Attributes

G4int fStack
 

Private Member Functions

std::string InsideTPC (const G4Track *aTrack)
 

Private Attributes

G4int fstage
 
G4int freqMuon
 
G4int freqIsoMuon
 
G4int freqIso
 
G4double fangRoI
 

Detailed Description

Definition at line 39 of file LArStackingAction.h.

Constructor & Destructor Documentation

LArStackingAction::LArStackingAction ( int  )

Definition at line 50 of file LArStackingAction.cxx.

References fStack.

51  : fstage(0), freqMuon(2), freqIsoMuon(0), freqIso(10), fangRoI(30. * CLHEP::deg)
52 {
53  fStack = dum;
54  // Positive values effect action in this routine. Negative values effect action in
55  // G4BadIdeaAction.
56 }

Member Function Documentation

G4ClassificationOfNewTrack LArStackingAction::ClassifyNewTrack ( const G4Track *  aTrack)
virtual

Definition at line 58 of file LArStackingAction.cxx.

References fStack, fstage, InsideTPC(), geo::GeometryCore::TPC(), and x4.

59 {
60  G4ClassificationOfNewTrack classification = fWaiting;
62  TString volName(InsideTPC(aTrack));
63  Double_t buffer = 500; // Keep muNucl neutrals within 5m (for now) of larTPC.
64 
65  // These 3 for now for investigation in gdb.
66  int ppdg = aTrack->GetParentID();
67  TString process("NA");
68  if (ppdg) process = (TString)aTrack->GetCreatorProcess()->GetProcessName();
69 
70  switch (fstage) {
71  case 0: // Fstage 0 : Primary muons only
72  if (aTrack->GetParentID() == 0) {
73  G4ParticleDefinition* particleType = aTrack->GetDefinition();
74  if (((particleType == G4MuonPlus::MuonPlusDefinition()) ||
75  (particleType == G4MuonMinus::MuonMinusDefinition())) &&
76  !volName.Contains("unknown")) {
77  classification = fUrgent;
78  }
79  }
80  if (volName.Contains("unknown")) classification = fKill;
81  break;
82 
83  case 1: // Stage 1 : K0,Lambda,n's made urgent here.
84  // Suspended tracks will be sent to the waiting stack
85  if (aTrack->GetTrackStatus() == fSuspend) { break; }
86 
87  if ((aTrack->GetDefinition()->GetPDGEncoding() == 2112 ||
88  aTrack->GetDefinition()->GetPDGEncoding() == 130 ||
89  aTrack->GetDefinition()->GetPDGEncoding() == 310 ||
90  aTrack->GetDefinition()->GetPDGEncoding() == 311 ||
91  aTrack->GetDefinition()->GetPDGEncoding() == 3122) &&
92  (aTrack->GetParentID() == 1) && !volName.Contains("unknown")) {
93 
94  const G4ThreeVector tr4Pos = aTrack->GetPosition();
95  // G4 returns positions in mm, have to convert to cm for LArSoft coordinate systems
96  const TVector3 trPos(tr4Pos.x() / CLHEP::cm, tr4Pos.y() / CLHEP::cm, tr4Pos.z() / CLHEP::cm);
97  classification = fUrgent;
98 
99  auto const& tpc = geom->TPC({0, 0});
100  if (trPos.X() < (tpc.HalfWidth() * 2.0 + buffer) && trPos.X() > (-buffer) &&
101  trPos.Y() < (tpc.HalfHeight() * 2.0 + buffer) &&
102  trPos.Y() > (-tpc.HalfHeight() * 2.0 - buffer) && trPos.Z() < (tpc.Length() + buffer) &&
103  trPos.Z() > (-buffer)) {
104  classification = fUrgent;
105  break;
106  }
107  // These tracks need to be "scored" cuz every now and then they might get to the
108  // LAr.
109  else {
110  classification = fKill;
111  break;
112  }
113  }
114 
115  if (volName.Contains("unknown")) classification = fKill;
116  break;
117 
118  default:
119  // Track all other Primaries. Accept all secondaries in TPC. Kill muon ionization
120  // electrons outside TPC ignore primaries since they have no creator process
121  if (aTrack->GetParentID() == 0 && !volName.Contains("unknown")) {
122  classification = fUrgent;
123  break;
124  }
125 
126  auto const& tpc = geom->TPC({0, 0});
127  if (volName.Contains(tpc.ActiveVolume()->GetName()) && aTrack->GetParentID() != 0) {
128  classification = fUrgent;
129  if (fStack & 0x4 && aTrack->GetCreatorProcess()->GetProcessName().contains("muIoni")) {
130  classification = fKill;
131  }
132  break;
133  }
134  else if (volName.Contains("unknown")) {
135  classification = fKill;
136  break;
137  }
138  // Leave this here, even though I claim we've Killed these in stage 2.
139  if (aTrack->GetDefinition()->GetPDGEncoding() == 11 &&
140  aTrack->GetCreatorProcess()->GetProcessName().contains("muIoni")) {
141  classification = fKill;
142  break;
143  }
144  // For now, kill every other thing, no matter where it is.
145  classification = fKill;
146 
147  } // end switch
148 
149  return classification;
150 }
std::string InsideTPC(const G4Track *aTrack)
Float_t x4[n_points_geant4]
TPCGeo const & TPC(TPCID const &tpcid=details::tpc_zero) const
Returns the specified TPC.
Definition: GeometryCore.h:448
G4int LArStackingAction::GetNIsolation ( ) const
inline

Definition at line 71 of file LArStackingAction.h.

References freqIso.

71 { return freqIso; }
G4int LArStackingAction::GetNRequestIsoMuon ( ) const
inline

Definition at line 69 of file LArStackingAction.h.

References freqIsoMuon.

69 { return freqIsoMuon; }
G4int LArStackingAction::GetNRequestMuon ( ) const
inline

Definition at line 67 of file LArStackingAction.h.

References freqMuon.

67 { return freqMuon; }
G4double LArStackingAction::GetRoIAngle ( ) const
inline

Definition at line 73 of file LArStackingAction.h.

References fangRoI.

73 { return fangRoI; }
std::string LArStackingAction::InsideTPC ( const G4Track *  aTrack)
private

Definition at line 152 of file LArStackingAction.cxx.

References geo::GeometryCore::VolumeName().

Referenced by ClassifyNewTrack().

153 {
155  const G4ThreeVector tr4Pos = aTrack->GetPosition();
156 
157  // G4 returns positions in mm, have to convert to cm for LArSoft coordinate systems
158  const geo::Point_t trPos(tr4Pos.x() / CLHEP::cm, tr4Pos.y() / CLHEP::cm, tr4Pos.z() / CLHEP::cm);
159 
160  return geom->VolumeName(trPos);
161 }
std::string VolumeName(Point_t const &point) const
Returns the name of the deepest volume containing specified point.
ROOT::Math::PositionVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Point_t
Type for representation of position in physical 3D space.
Definition: geo_vectors.h:180
void LArStackingAction::NewStage ( )
virtual

Definition at line 163 of file LArStackingAction.cxx.

References fstage.

164 {
165  // Here when Urgent stack is empty. Waiting stack about to be made Urgent, upon saying
166  // ReClassify().
167  fstage++;
168 
169  // I yanked the ExN04's use here of stackManager->clear(), which clears stack and
170  // prepares to end the event. Think I may wanna do something like this if muon has been
171  // tracked, its doca is large, there are no hit voxels in the TPC, and I'm in further
172  // need of optimization.
173 
174  if (fstage == 1) {
175  // Stage 0->1 : check if at least "reqMuon" hits on muon chamber otherwise abort
176  // current event
177 
178  stackManager->ReClassify();
179  return;
180  }
181 
182  else if (fstage == 2) {
183  // Stage 1->2 : check the isolation of muon tracks at least "reqIsoMuon" isolated
184  // muons otherwise abort current event. Isolation requires "reqIso" or
185  // less hits (including own hits) in the RoI region in the tracker
186  // layers.
187  stackManager->ReClassify();
188  return;
189  }
190 
191  else {
192  // Other stage change : just re-classify
193  stackManager->ReClassify();
194  }
195 }
void LArStackingAction::PrepareNewEvent ( )
virtual

Definition at line 197 of file LArStackingAction.cxx.

References fstage.

198 {
199  fstage = 0;
200 }
void LArStackingAction::SetNIsolation ( G4int  val)
inline

Definition at line 70 of file LArStackingAction.h.

70 { freqIso = val; }
void LArStackingAction::SetNRequestIsoMuon ( G4int  val)
inline

Definition at line 68 of file LArStackingAction.h.

68 { freqIsoMuon = val; }
void LArStackingAction::SetNRequestMuon ( G4int  val)
inline

Definition at line 66 of file LArStackingAction.h.

66 { freqMuon = val; }
void LArStackingAction::SetRoIAngle ( G4double  val)
inline

Definition at line 72 of file LArStackingAction.h.

72 { fangRoI = val; }

Member Data Documentation

G4double LArStackingAction::fangRoI
private

Definition at line 63 of file LArStackingAction.h.

Referenced by GetRoIAngle().

G4int LArStackingAction::freqIso
private

Definition at line 62 of file LArStackingAction.h.

Referenced by GetNIsolation().

G4int LArStackingAction::freqIsoMuon
private

Definition at line 61 of file LArStackingAction.h.

Referenced by GetNRequestIsoMuon().

G4int LArStackingAction::freqMuon
private

Definition at line 60 of file LArStackingAction.h.

Referenced by GetNRequestMuon().

G4int LArStackingAction::fStack

Definition at line 48 of file LArStackingAction.h.

Referenced by ClassifyNewTrack(), and LArStackingAction().

G4int LArStackingAction::fstage
private

Definition at line 59 of file LArStackingAction.h.

Referenced by ClassifyNewTrack(), NewStage(), and PrepareNewEvent().


The documentation for this class was generated from the following files: