Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4RichTrajectory.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // $Id: G4RichTrajectory.cc 91269 2015-06-29 07:05:59Z gcosmo $
28 //
29 // ---------------------------------------------------------------
30 //
31 // G4RichTrajectory.cc
32 //
33 // Contact:
34 // Questions and comments on G4Trajectory, on which this is based,
35 // should be sent to
36 // Katsuya Amako (e-mail: Katsuya.Amako@kek.jp)
37 // Makoto Asai (e-mail: asai@kekvax.kek.jp)
38 // Takashi Sasaki (e-mail: Takashi.Sasaki@kek.jp)
39 // and on the extended code to:
40 // John Allison (e-mail: John.Allison@manchester.ac.uk)
41 // Joseph Perl (e-mail: perl@slac.stanford.edu)
42 //
43 // ---------------------------------------------------------------
44 
45 #include "G4RichTrajectory.hh"
46 #include "G4RichTrajectoryPoint.hh"
47 #include "G4AttDefStore.hh"
48 #include "G4AttDef.hh"
49 #include "G4AttValue.hh"
50 #include "G4UIcommand.hh"
51 #include "G4UnitsTable.hh"
52 #include "G4VProcess.hh"
53 
54 //#define G4ATTDEBUG
55 #ifdef G4ATTDEBUG
56 #include "G4AttCheck.hh"
57 #endif
58 
59 #include <sstream>
60 
62 
64  fpRichPointsContainer(0),
65  fpCreatorProcess(0),
66  fCreatorModelID(0),
67  fpEndingProcess(0),
68  fFinalKineticEnergy(0.)
69 {
70 }
71 
73  G4Trajectory(aTrack) // Note: this initialises the base class data
74  // members and, unfortunately but never mind,
75  // creates a G4TrajectoryPoint in
76  // TrajectoryPointContainer that we cannot
77  // access because it's private. We store the
78  // same information (plus more) in a
79  // G4RichTrajectoryPoint in the
80  // RichTrajectoryPointsContainer
81 {
82  fpInitialVolume = aTrack->GetTouchableHandle();
83  fpInitialNextVolume = aTrack->GetNextTouchableHandle();
84  fpCreatorProcess = aTrack->GetCreatorProcess();
85  fCreatorModelID = aTrack->GetCreatorModelID();
86  // On construction, set final values to initial values.
87  // Final values are updated at the addition of every step - see AppendStep.
88  fpFinalVolume = aTrack->GetTouchableHandle();
89  fpFinalNextVolume = aTrack->GetNextTouchableHandle();
90  fpEndingProcess = aTrack->GetCreatorProcess();
91  fFinalKineticEnergy = aTrack->GetKineticEnergy();
92  // Insert the first rich trajectory point (see note above)...
93  fpRichPointsContainer = new RichTrajectoryPointsContainer;
94  fpRichPointsContainer->push_back(new G4RichTrajectoryPoint(aTrack));
95 }
96 
98  G4Trajectory(right)
99 {
100  fpInitialVolume = right.fpInitialVolume;
101  fpInitialNextVolume = right.fpInitialNextVolume;
102  fpCreatorProcess = right.fpCreatorProcess;
103  fCreatorModelID = right.fCreatorModelID;
104  fpFinalVolume = right.fpFinalVolume;
105  fpFinalNextVolume = right.fpFinalNextVolume;
106  fpEndingProcess = right.fpEndingProcess;
107  fFinalKineticEnergy = right.fFinalKineticEnergy;
108  fpRichPointsContainer = new RichTrajectoryPointsContainer;
109  for(size_t i=0;i<right.fpRichPointsContainer->size();i++)
110  {
111  G4RichTrajectoryPoint* rightPoint =
112  (G4RichTrajectoryPoint*)((*(right.fpRichPointsContainer))[i]);
113  fpRichPointsContainer->push_back(new G4RichTrajectoryPoint(*rightPoint));
114  }
115 }
116 
118 {
119  if (fpRichPointsContainer) {
120  // fpRichPointsContainer->clearAndDestroy();
121  size_t i;
122  for(i=0;i<fpRichPointsContainer->size();i++){
123  delete (*fpRichPointsContainer)[i];
124  }
125  fpRichPointsContainer->clear();
126  delete fpRichPointsContainer;
127  }
128 }
129 
131 {
132  fpRichPointsContainer->push_back(new G4RichTrajectoryPoint(aStep));
133  // Except for first step, which is a sort of virtual step to start
134  // the track, compute the final values...
135  const G4Track* track = aStep->GetTrack();
136  const G4StepPoint* postStepPoint = aStep->GetPostStepPoint();
137  if (track->GetCurrentStepNumber() > 0) {
138  fpFinalVolume = track->GetTouchableHandle();
139  fpFinalNextVolume = track->GetNextTouchableHandle();
140  fpEndingProcess = postStepPoint->GetProcessDefinedStep();
141  fFinalKineticEnergy =
142  aStep->GetPreStepPoint()->GetKineticEnergy() -
143  aStep->GetTotalEnergyDeposit();
144  }
145 }
146 
148 {
149  if(!secondTrajectory) return;
150 
151  G4RichTrajectory* seco = (G4RichTrajectory*)secondTrajectory;
152  G4int ent = seco->GetPointEntries();
153  for(G4int i=1;i<ent;i++) {
154  // initial point of the second trajectory should not be merged
155  fpRichPointsContainer->push_back((*(seco->fpRichPointsContainer))[i]);
156  // fpRichPointsContainer->push_back(seco->fpRichPointsContainer->removeAt(1));
157  }
158  delete (*seco->fpRichPointsContainer)[0];
159  seco->fpRichPointsContainer->clear();
160 }
161 
162 void G4RichTrajectory::ShowTrajectory(std::ostream& os) const
163 {
164  // Invoke the default implementation in G4VTrajectory...
166  // ... or override with your own code here.
167 }
168 
170 {
171  // Invoke the default implementation in G4VTrajectory...
173  // ... or override with your own code here.
174 }
175 
176 const std::map<G4String,G4AttDef>* G4RichTrajectory::GetAttDefs() const
177 {
178  G4bool isNew;
179  std::map<G4String,G4AttDef>* store
180  = G4AttDefStore::GetInstance("G4RichTrajectory",isNew);
181  if (isNew) {
182 
183  // Get att defs from base class...
184  *store = *(G4Trajectory::GetAttDefs());
185 
186  G4String ID;
187 
188  ID = "IVPath";
189  (*store)[ID] = G4AttDef(ID,"Initial Volume Path",
190  "Physics","","G4String");
191 
192  ID = "INVPath";
193  (*store)[ID] = G4AttDef(ID,"Initial Next Volume Path",
194  "Physics","","G4String");
195 
196  ID = "CPN";
197  (*store)[ID] = G4AttDef(ID,"Creator Process Name",
198  "Physics","","G4String");
199 
200  ID = "CPTN";
201  (*store)[ID] = G4AttDef(ID,"Creator Process Type Name",
202  "Physics","","G4String");
203 
204  ID = "CMID";
205  (*store)[ID] = G4AttDef(ID,"Creator Model ID",
206  "Physics","","G4int");
207 
208  ID = "CMN";
209  (*store)[ID] = G4AttDef(ID,"Creator Model Name",
210  "Physics","","G4String");
211 
212  ID = "FVPath";
213  (*store)[ID] = G4AttDef(ID,"Final Volume Path",
214  "Physics","","G4String");
215 
216  ID = "FNVPath";
217  (*store)[ID] = G4AttDef(ID,"Final Next Volume Path",
218  "Physics","","G4String");
219 
220  ID = "EPN";
221  (*store)[ID] = G4AttDef(ID,"Ending Process Name",
222  "Physics","","G4String");
223 
224  ID = "EPTN";
225  (*store)[ID] = G4AttDef(ID,"Ending Process Type Name",
226  "Physics","","G4String");
227 
228  ID = "FKE";
229  (*store)[ID] = G4AttDef(ID,"Final kinetic energy",
230  "Physics","G4BestUnit","G4double");
231 
232  }
233 
234  return store;
235 }
236 
237 static G4String Path(const G4TouchableHandle& th)
238 {
239  std::ostringstream oss;
240  G4int depth = th->GetHistoryDepth();
241  for (G4int i = depth; i >= 0; --i) {
242  oss << th->GetVolume(i)->GetName()
243  << ':' << th->GetCopyNumber(i);
244  if (i != 0) oss << '/';
245  }
246  return oss.str();
247 }
248 
249 std::vector<G4AttValue>* G4RichTrajectory::CreateAttValues() const
250 {
251  // Create base class att values...
252  std::vector<G4AttValue>* values = G4Trajectory::CreateAttValues();
253 
254  if (fpInitialVolume && fpInitialVolume->GetVolume()) {
255  values->push_back(G4AttValue("IVPath",Path(fpInitialVolume),""));
256  } else {
257  values->push_back(G4AttValue("IVPath","None",""));
258  }
259 
260  if (fpInitialNextVolume && fpInitialNextVolume->GetVolume()) {
261  values->push_back(G4AttValue("INVPath",Path(fpInitialNextVolume),""));
262  } else {
263  values->push_back(G4AttValue("INVPath","None",""));
264  }
265 
266  if (fpCreatorProcess) {
267  values->push_back
268  (G4AttValue("CPN",fpCreatorProcess->GetProcessName(),""));
269  G4ProcessType type = fpCreatorProcess->GetProcessType();
270  values->push_back
271  (G4AttValue("CPTN",G4VProcess::GetProcessTypeName(type),""));
272  values->push_back
273  (G4AttValue("CMID",G4UIcommand::ConvertToString(fCreatorModelID),""));
274  const G4String& creatorModelName =
275  G4PhysicsModelCatalog::GetModelName(fCreatorModelID);
276  values->push_back(G4AttValue("CMN",creatorModelName,""));
277  } else {
278  values->push_back(G4AttValue("CPN","None",""));
279  values->push_back(G4AttValue("CPTN","None",""));
280  values->push_back(G4AttValue("CMID","None",""));
281  values->push_back(G4AttValue("CMN","None",""));
282  }
283 
284  if (fpFinalVolume && fpFinalVolume->GetVolume()) {
285  values->push_back(G4AttValue("FVPath",Path(fpFinalVolume),""));
286  } else {
287  values->push_back(G4AttValue("FVPath","None",""));
288  }
289 
290  if (fpFinalNextVolume && fpFinalNextVolume->GetVolume()) {
291  values->push_back(G4AttValue("FNVPath",Path(fpFinalNextVolume),""));
292  } else {
293  values->push_back(G4AttValue("FNVPath","None",""));
294  }
295 
296  if (fpEndingProcess) {
297  values->push_back(G4AttValue("EPN",fpEndingProcess->GetProcessName(),""));
298  G4ProcessType type = fpEndingProcess->GetProcessType();
299  values->push_back(G4AttValue("EPTN",G4VProcess::GetProcessTypeName(type),""));
300  } else {
301  values->push_back(G4AttValue("EPN","None",""));
302  values->push_back(G4AttValue("EPTN","None",""));
303  }
304 
305  values->push_back
306  (G4AttValue("FKE",G4BestUnit(fFinalKineticEnergy,"Energy"),""));
307 
308 #ifdef G4ATTDEBUG
309  G4cout << G4AttCheck(values,GetAttDefs());
310 #endif
311 
312  return values;
313 }
G4int GetCreatorModelID() const
static const G4String & GetProcessTypeName(G4ProcessType)
Definition: G4VProcess.cc:141
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
virtual std::vector< G4AttValue > * CreateAttValues() const
void AppendStep(const G4Step *aStep)
void ShowTrajectory(std::ostream &os=G4cout) const
static G4String ConvertToString(G4bool boolVal)
Definition: G4UIcommand.cc:372
static const G4String & GetModelName(G4int)
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
int GetPointEntries() const
G4ProcessType GetProcessType() const
Definition: G4VProcess.hh:414
#define G4ThreadLocal
Definition: tls.hh:89
G4int GetCopyNumber(G4int depth=0) const
static G4String Path(const G4TouchableHandle &th)
int G4int
Definition: G4Types.hh:78
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
virtual void DrawTrajectory() const
G4StepPoint * GetPreStepPoint() const
const G4VProcess * GetCreatorProcess() const
G4double GetKineticEnergy() const
std::vector< G4VTrajectoryPoint * > RichTrajectoryPointsContainer
G4GLOB_DLL std::ostream G4cout
G4int GetCurrentStepNumber() const
const G4String & GetName() const
bool G4bool
Definition: G4Types.hh:79
virtual std::vector< G4AttValue > * CreateAttValues() const
const G4TouchableHandle & GetNextTouchableHandle() const
Definition: G4Step.hh:76
virtual G4int GetHistoryDepth() const
Definition: G4VTouchable.cc:79
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
const G4TouchableHandle & GetTouchableHandle() const
G4double GetTotalEnergyDeposit() const
const G4VProcess * GetProcessDefinedStep() const
virtual G4VPhysicalVolume * GetVolume(G4int depth=0) const
Definition: G4VTouchable.cc:44
void MergeTrajectory(G4VTrajectory *secondTrajectory)
G4StepPoint * GetPostStepPoint() const
void DrawTrajectory() const
virtual ~G4RichTrajectory()
std::map< G4String, G4AttDef > * GetInstance(const G4String &storeKey, G4bool &isNew)
G4TRACKING_DLL G4ThreadLocal G4Allocator< G4RichTrajectory > * aRichTrajectoryAllocator
G4double GetKineticEnergy() const
G4Track * GetTrack() const
virtual void ShowTrajectory(std::ostream &os=G4cout) const
G4ProcessType