Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
RE05StackingAction Class Reference

#include <RE05StackingAction.hh>

Inheritance diagram for RE05StackingAction:
Collaboration diagram for RE05StackingAction:

Public Member Functions

 RE05StackingAction ()
 
virtual ~RE05StackingAction ()
 
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 Member Functions inherited from G4UserStackingAction
 G4UserStackingAction ()
 
virtual ~G4UserStackingAction ()
 
void SetStackManager (G4StackManager *value)
 

Additional Inherited Members

- Protected Attributes inherited from G4UserStackingAction
G4StackManagerstackManager
 

Detailed Description

Definition at line 45 of file RE05StackingAction.hh.

Constructor & Destructor Documentation

RE05StackingAction::RE05StackingAction ( )

Definition at line 47 of file RE05StackingAction.cc.

49  fTrkHits(0), fMuonHits(0), fMessenger(0),
50  fStage(0),
51  fReqMuon(2),
52  fReqIso(10),
53  fAngRoI(30.0*deg)
54 {
55  fMessenger = new RE05StackingActionMessenger(this);
56 }
static constexpr double deg
Definition: G4SIunits.hh:152
RE05StackingAction::~RE05StackingAction ( )
virtual

Definition at line 60 of file RE05StackingAction.cc.

61 { delete fMessenger; }

Member Function Documentation

G4ClassificationOfNewTrack RE05StackingAction::ClassifyNewTrack ( const G4Track aTrack)
virtual

Reimplemented from G4UserStackingAction.

Definition at line 66 of file RE05StackingAction.cc.

67 {
68  G4ClassificationOfNewTrack classification = fWaiting;
69  switch(fStage)
70  {
71  case 0: // Stage 0 : Primary muons only
72  if(aTrack->GetParentID()==0)
73  {
74  G4ParticleDefinition * particleType = aTrack->GetDefinition();
75  if((particleType==G4MuonPlus::MuonPlusDefinition())
76  ||(particleType==G4MuonMinus::MuonMinusDefinition()))
77  { classification = fUrgent; }
78  }
79  break;
80 
81  case 1: // Stage 1 : Charged primaries only
82  // Suspended tracks will be sent to the waiting stack
83  if(aTrack->GetParentID()!=0) { break; }
84  if(aTrack->GetTrackStatus()==fSuspend) { break; }
85  if(aTrack->GetDefinition()->GetPDGCharge()==0.) { break; }
86  classification = fUrgent;
87  break;
88 
89  default: // Stage 2 : Accept all primaries
90  // Accept all secondaries in RoI
91  // Kill secondaries outside RoI
92  if(aTrack->GetParentID()==0)
93  {
94  classification = fUrgent;
95  break;
96  }
97  if((fAngRoI<0.)||InsideRoI(aTrack,fAngRoI))
98  {
99  classification = fUrgent;
100  break;
101  }
102  classification = fKill;
103  }
104  return classification;
105 }
G4ParticleDefinition * GetDefinition() const
G4int GetParentID() const
static G4MuonPlus * MuonPlusDefinition()
Definition: G4MuonPlus.cc:94
G4TrackStatus GetTrackStatus() const
static G4MuonMinus * MuonMinusDefinition()
Definition: G4MuonMinus.cc:95
G4double GetPDGCharge() const

Here is the call graph for this function:

G4int RE05StackingAction::GetNIsolation ( ) const
inline

Definition at line 76 of file RE05StackingAction.hh.

76 { return fReqIso; }

Here is the caller graph for this function:

G4int RE05StackingAction::GetNRequestIsoMuon ( ) const
inline

Definition at line 74 of file RE05StackingAction.hh.

74 { return fReqIsoMuon; }

Here is the caller graph for this function:

G4int RE05StackingAction::GetNRequestMuon ( ) const
inline

Definition at line 72 of file RE05StackingAction.hh.

72 { return fReqMuon; }

Here is the caller graph for this function:

G4double RE05StackingAction::GetRoIAngle ( ) const
inline

Definition at line 78 of file RE05StackingAction.hh.

78 { return fAngRoI; }

Here is the caller graph for this function:

void RE05StackingAction::NewStage ( )
virtual

Reimplemented from G4UserStackingAction.

Definition at line 148 of file RE05StackingAction.cc.

149 {
150  fStage++;
151  G4int nhits;
152  if(fStage==1)
153  {
154  // Stage 0->1 : check if at least "fReqMuon" hits on muon chamber
155  // otherwise abort current event
156  if(!fMuonHits)
157  { fMuonHits = (RE05MuonHitsCollection*)GetCollection("muonCollection"); }
158  if(!fMuonHits)
159  { G4cerr << "muonCollection NOT FOUND" << G4endl;
160  return; }
161  nhits = fMuonHits->entries();
164  if(nhits<fReqMuon)
165  {
166  stackManager->clear();
168  return;
169  }
171  return;
172  }
173 
174  else if(fStage==2)
175  {
176  // Stage 1->2 : check the isolation of muon tracks
177  // at least "fReqIsoMuon" isolated muons
178  // otherwise abort current event.
179  // Isolation requires "fReqIso" or less hits
180  // (including own hits) in the RoI region
181  // in the tracker layers.
182  nhits = fMuonHits->entries();
183  if(!fTrkHits)
184  { fTrkHits = (RE05TrackerHitsCollection*)GetCollection("trackerCollection"); }
185  if(!fTrkHits)
186  { G4cerr << "trackerCollection NOT FOUND" << G4endl;
187  return; }
188  G4int nTrkhits = fTrkHits->entries();
189  G4int isoMuon = 0;
190  for(G4int j=0;j<nhits;j++)
191  {
192  G4ThreeVector hitPos = (*fMuonHits)[j]->GetPos();
193  G4int nhitIn = 0;
194  for(G4int jj=0;(jj<nTrkhits)&&(nhitIn<=fReqIso);jj++)
195  {
196  G4ThreeVector trkhitPos = (*fTrkHits)[jj]->GetPos();
197  if(trkhitPos.angle(hitPos)<fAngRoI) nhitIn++;
198  }
199  if(nhitIn<=fReqIso) isoMuon++;
200  }
202  if(isoMuon<fReqIsoMuon)
203  {
204  stackManager->clear();
206  return;
207  }
209  return;
210  }
211 
212  else
213  {
214  // Other fStage change : just re-classify
216  }
217 }
double angle(const Hep3Vector &) const
G4int entries() const
int G4int
Definition: G4Types.hh:78
#define G4endl
Definition: G4ios.hh:61
G4StackManager * stackManager
G4GLOB_DLL std::ostream G4cerr

Here is the call graph for this function:

void RE05StackingAction::PrepareNewEvent ( )
virtual

Reimplemented from G4UserStackingAction.

Definition at line 221 of file RE05StackingAction.cc.

222 {
223  fStage = 0;
224  fTrkHits = 0;
225  fMuonHits = 0;
226 }
void RE05StackingAction::SetNIsolation ( G4int  val)
inline

Definition at line 75 of file RE05StackingAction.hh.

75 { fReqIso = val; }

Here is the caller graph for this function:

void RE05StackingAction::SetNRequestIsoMuon ( G4int  val)
inline

Definition at line 73 of file RE05StackingAction.hh.

73 { fReqIsoMuon = val; }

Here is the caller graph for this function:

void RE05StackingAction::SetNRequestMuon ( G4int  val)
inline

Definition at line 71 of file RE05StackingAction.hh.

71 { fReqMuon = val; }

Here is the caller graph for this function:

void RE05StackingAction::SetRoIAngle ( G4double  val)
inline

Definition at line 77 of file RE05StackingAction.hh.

77 { fAngRoI = val; }

Here is the caller graph for this function:


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