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

#include <RE01StackingAction.hh>

Inheritance diagram for RE01StackingAction:
Collaboration diagram for RE01StackingAction:

Public Member Functions

 RE01StackingAction ()
 
virtual ~RE01StackingAction ()
 
virtual G4ClassificationOfNewTrack ClassifyNewTrack (const G4Track *aTrack)
 
virtual void NewStage ()
 
virtual void PrepareNewEvent ()
 
- 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 41 of file RE01StackingAction.hh.

Constructor & Destructor Documentation

RE01StackingAction::RE01StackingAction ( )

Definition at line 48 of file RE01StackingAction.cc.

50  fStage(0),fCalorimeterHitsColID(-1)
51 {;}
RE01StackingAction::~RE01StackingAction ( )
virtual

Definition at line 54 of file RE01StackingAction.cc.

55 {;}

Member Function Documentation

G4ClassificationOfNewTrack RE01StackingAction::ClassifyNewTrack ( const G4Track aTrack)
virtual

Reimplemented from G4UserStackingAction.

Definition at line 59 of file RE01StackingAction.cc.

60 {
61  G4ClassificationOfNewTrack classification = fUrgent;
62 
63  if(fStage==0)
64  {
65  RE01TrackInformation* trackInfo;
66  if(aTrack->GetTrackStatus()==fSuspend) // Track reached to calorimeter
67  {
68  trackInfo = (RE01TrackInformation*)(aTrack->GetUserInformation());
69  trackInfo->SetTrackingStatus(0);
70  trackInfo->SetSourceTrackInformation(aTrack);
71 // G4cout << "Track " << aTrack->GetTrackID()
72 // << " (parentID " << aTrack->GetParentID()
73 // << ") has reached to calorimeter and has been suspended at "
74 // << aTrack->GetPosition() << G4endl;
75  classification = fWaiting;
76  }
77  else if(aTrack->GetParentID()==0) // Primary particle
78  {
79  trackInfo = new RE01TrackInformation(aTrack);
80  trackInfo->SetTrackingStatus(1);
81  G4Track* theTrack = (G4Track*)aTrack;
82  theTrack->SetUserInformation(trackInfo);
83  }
84  }
85  return classification;
86 }
G4int GetParentID() const
G4TrackStatus GetTrackStatus() const
void SetSourceTrackInformation(const G4Track *aTrack)
G4VUserTrackInformation * GetUserInformation() const
void SetUserInformation(G4VUserTrackInformation *aValue) const

Here is the call graph for this function:

void RE01StackingAction::NewStage ( )
virtual

Reimplemented from G4UserStackingAction.

Definition at line 105 of file RE01StackingAction.cc.

106 {
107  // G4cout << "+++++++++++ Stage " << fStage << G4endl;
108  if(fStage==0)
109  {
110  // display trajetory information in the tracking region
111  G4cout << G4endl;
112  G4cout << "Tracks in tracking region have been processed. -- Stage 0 over."
113  << G4endl;
114  G4cout << G4endl;
115  }
116  else
117  {
118  // display calorimeter information caused by a "source" track
119  // in the tracker region
120  // G4cout << G4endl;
121  // G4cout << "Processing one shower originated by a source track "
122  // << "in tracker region is over." << G4endl;
123  // G4cout << G4endl;
125  (RE01CalorimeterHitsCollection*)GetCalCollection();
126  if(CHC)
127  {
128  int n_hit = CHC->entries();
129  G4double totE = 0;
130  G4int n_hitByATrack = 0;
131  for(int i=0;i<n_hit;i++)
132  {
133  G4double edepByATrack = (*CHC)[i]->GetEdepByATrack();
134  if(edepByATrack>0.)
135  {
136  totE += edepByATrack;
137  if(n_hitByATrack==0)
138  { (*CHC)[i]->GetTrackInformation()->Print(); }
139  n_hitByATrack++;
140  G4cout << "Cell[" << (*CHC)[i]->GetZ() << ","
141  << (*CHC)[i]->GetPhi() << "] "
142  << edepByATrack/GeV << " [GeV]" << G4endl;
143  (*CHC)[i]->ClearEdepByATrack();
144  }
145  }
146  if(n_hitByATrack>0)
147  {
148  G4cout <<
149  "### Total energy deposition in calorimeter by a source track in "
150  << n_hitByATrack << " cells : " << totE / GeV << " (GeV)"
151  << G4endl;
152  G4cout << G4endl;
153  }
154  }
155  }
156 
158  {
159  // Transfer all tracks in Urgent stack to Waiting stack, since all tracks
160  // in Waiting stack have already been transfered to Urgent stack before
161  // invokation of this method.
163 
164  // Then, transfer only one track to Urgent stack.
166 
167  fStage++;
168  }
169 }
G4int entries() const
int G4int
Definition: G4Types.hh:78
G4int GetNUrgentTrack() const
G4GLOB_DLL std::ostream G4cout
void TransferOneStackedTrack(G4ClassificationOfNewTrack origin, G4ClassificationOfNewTrack destination)
void TransferStackedTracks(G4ClassificationOfNewTrack origin, G4ClassificationOfNewTrack destination)
static constexpr double GeV
Definition: G4SIunits.hh:217
#define G4endl
Definition: G4ios.hh:61
G4StackManager * stackManager
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

void RE01StackingAction::PrepareNewEvent ( )
virtual

Reimplemented from G4UserStackingAction.

Definition at line 172 of file RE01StackingAction.cc.

173 {
174  fStage = 0;
175 }

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