Geant4  10.01.p02
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 69003 2013-04-15 09:25:23Z 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 "G4UnitsTable.hh"
51 #include "G4VProcess.hh"
52 
53 //#define G4ATTDEBUG
54 #ifdef G4ATTDEBUG
55 #include "G4AttCheck.hh"
56 #endif
57 
58 #include <sstream>
59 
61 
63  fpRichPointsContainer(0),
64  fpCreatorProcess(0),
65  fpEndingProcess(0),
66  fFinalKineticEnergy(0.)
67 {
68 }
69 
71  G4Trajectory(aTrack) // Note: this initialises the base class data
72  // members and, unfortunately but never mind,
73  // creates a G4TrajectoryPoint in
74  // TrajectoryPointContainer that we cannot
75  // access because it's private. We store the
76  // same information (plus more) in a
77  // G4RichTrajectoryPoint in the
78  // RichTrajectoryPointsContainer
79 {
83  // On construction, set final values to initial values.
84  // Final values are updated at the addition of every step - see AppendStep.
89  // Insert the first rich trajectory point (see note above)...
91  fpRichPointsContainer->push_back(new G4RichTrajectoryPoint(aTrack));
92 }
93 
95  G4Trajectory(right)
96 {
105  for(size_t i=0;i<right.fpRichPointsContainer->size();i++)
106  {
107  G4RichTrajectoryPoint* rightPoint =
109  fpRichPointsContainer->push_back(new G4RichTrajectoryPoint(*rightPoint));
110  }
111 }
112 
114 {
115  if (fpRichPointsContainer) {
116  // fpRichPointsContainer->clearAndDestroy();
117  size_t i;
118  for(i=0;i<fpRichPointsContainer->size();i++){
119  delete (*fpRichPointsContainer)[i];
120  }
121  fpRichPointsContainer->clear();
122  delete fpRichPointsContainer;
123  }
124 }
125 
127 {
128  fpRichPointsContainer->push_back(new G4RichTrajectoryPoint(aStep));
129  // Except for first step, which is a sort of virtual step to start
130  // the track, compute the final values...
131  const G4Track* track = aStep->GetTrack();
132  const G4StepPoint* postStepPoint = aStep->GetPostStepPoint();
133  if (track->GetCurrentStepNumber() > 0) {
136  fpEndingProcess = postStepPoint->GetProcessDefinedStep();
138  aStep->GetPreStepPoint()->GetKineticEnergy() -
139  aStep->GetTotalEnergyDeposit();
140  }
141 }
142 
144 {
145  if(!secondTrajectory) return;
146 
147  G4RichTrajectory* seco = (G4RichTrajectory*)secondTrajectory;
148  G4int ent = seco->GetPointEntries();
149  for(G4int i=1;i<ent;i++) {
150  // initial point of the second trajectory should not be merged
151  fpRichPointsContainer->push_back((*(seco->fpRichPointsContainer))[i]);
152  // fpRichPointsContainer->push_back(seco->fpRichPointsContainer->removeAt(1));
153  }
154  delete (*seco->fpRichPointsContainer)[0];
155  seco->fpRichPointsContainer->clear();
156 }
157 
158 void G4RichTrajectory::ShowTrajectory(std::ostream& os) const
159 {
160  // Invoke the default implementation in G4VTrajectory...
162  // ... or override with your own code here.
163 }
164 
166 {
167  // Invoke the default implementation in G4VTrajectory...
169  // ... or override with your own code here.
170 }
171 
172 const std::map<G4String,G4AttDef>* G4RichTrajectory::GetAttDefs() const
173 {
174  G4bool isNew;
175  std::map<G4String,G4AttDef>* store
176  = G4AttDefStore::GetInstance("G4RichTrajectory",isNew);
177  if (isNew) {
178 
179  // Get att defs from base class...
180  *store = *(G4Trajectory::GetAttDefs());
181 
182  G4String ID;
183 
184  ID = "IVPath";
185  (*store)[ID] = G4AttDef(ID,"Initial Volume Path",
186  "Physics","","G4String");
187 
188  ID = "INVPath";
189  (*store)[ID] = G4AttDef(ID,"Initial Next Volume Path",
190  "Physics","","G4String");
191 
192  ID = "CPN";
193  (*store)[ID] = G4AttDef(ID,"Creator Process Name",
194  "Physics","","G4String");
195 
196  ID = "CPTN";
197  (*store)[ID] = G4AttDef(ID,"Creator Process Type Name",
198  "Physics","","G4String");
199 
200  ID = "FVPath";
201  (*store)[ID] = G4AttDef(ID,"Final Volume Path",
202  "Physics","","G4String");
203 
204  ID = "FNVPath";
205  (*store)[ID] = G4AttDef(ID,"Final Next Volume Path",
206  "Physics","","G4String");
207 
208  ID = "EPN";
209  (*store)[ID] = G4AttDef(ID,"Ending Process Name",
210  "Physics","","G4String");
211 
212  ID = "EPTN";
213  (*store)[ID] = G4AttDef(ID,"Ending Process Type Name",
214  "Physics","","G4String");
215 
216  ID = "FKE";
217  (*store)[ID] = G4AttDef(ID,"Final kinetic energy",
218  "Physics","G4BestUnit","G4double");
219 
220  }
221 
222  return store;
223 }
224 
225 static G4String Path(const G4TouchableHandle& th)
226 {
227  std::ostringstream oss;
228  G4int depth = th->GetHistoryDepth();
229  for (G4int i = depth; i >= 0; --i) {
230  oss << th->GetVolume(i)->GetName()
231  << ':' << th->GetCopyNumber(i);
232  if (i != 0) oss << '/';
233  }
234  return oss.str();
235 }
236 
237 std::vector<G4AttValue>* G4RichTrajectory::CreateAttValues() const
238 {
239  // Create base class att values...
240  std::vector<G4AttValue>* values = G4Trajectory::CreateAttValues();
241 
243  values->push_back(G4AttValue("IVPath",Path(fpInitialVolume),""));
244  } else {
245  values->push_back(G4AttValue("IVPath","None",""));
246  }
247 
249  values->push_back(G4AttValue("INVPath",Path(fpInitialNextVolume),""));
250  } else {
251  values->push_back(G4AttValue("INVPath","None",""));
252  }
253 
254  if (fpCreatorProcess) {
255  values->push_back(G4AttValue("CPN",fpCreatorProcess->GetProcessName(),""));
257  values->push_back(G4AttValue("CPTN",G4VProcess::GetProcessTypeName(type),""));
258  } else {
259  values->push_back(G4AttValue("CPN","None",""));
260  values->push_back(G4AttValue("CPTN","None",""));
261  }
262 
264  values->push_back(G4AttValue("FVPath",Path(fpFinalVolume),""));
265  } else {
266  values->push_back(G4AttValue("FVPath","None",""));
267  }
268 
270  values->push_back(G4AttValue("FNVPath",Path(fpFinalNextVolume),""));
271  } else {
272  values->push_back(G4AttValue("FNVPath","None",""));
273  }
274 
275  if (fpEndingProcess) {
276  values->push_back(G4AttValue("EPN",fpEndingProcess->GetProcessName(),""));
278  values->push_back(G4AttValue("EPTN",G4VProcess::GetProcessTypeName(type),""));
279  } else {
280  values->push_back(G4AttValue("EPN","None",""));
281  values->push_back(G4AttValue("EPTN","None",""));
282  }
283 
284  values->push_back
285  (G4AttValue("FKE",G4BestUnit(fFinalKineticEnergy,"Energy"),""));
286 
287 #ifdef G4ATTDEBUG
288  G4cout << G4AttCheck(values,GetAttDefs());
289 #endif
290 
291  return values;
292 }
G4ThreadLocal G4Allocator< G4RichTrajectory > * aRichTrajectoryAllocator
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
#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)
G4TouchableHandle fpInitialNextVolume
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
G4double fFinalKineticEnergy
G4int GetCurrentStepNumber() const
const G4String & GetName() const
G4TouchableHandle fpInitialVolume
bool G4bool
Definition: G4Types.hh:79
virtual std::vector< G4AttValue > * CreateAttValues() const
const G4TouchableHandle & GetNextTouchableHandle() const
Definition: G4Step.hh:76
RichTrajectoryPointsContainer * fpRichPointsContainer
virtual G4int GetHistoryDepth() const
Definition: G4VTouchable.cc:79
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
G4TouchableHandle fpFinalVolume
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)
G4TouchableHandle fpFinalNextVolume
G4StepPoint * GetPostStepPoint() const
void DrawTrajectory() const
virtual ~G4RichTrajectory()
std::map< G4String, G4AttDef > * GetInstance(const G4String &storeKey, G4bool &isNew)
const G4VProcess * fpEndingProcess
G4double GetKineticEnergy() const
G4Track * GetTrack() const
virtual void ShowTrajectory(std::ostream &os=G4cout) const
const G4VProcess * fpCreatorProcess
G4ProcessType