Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4RichTrajectory.cc
이 파일의 문서화 페이지로 가기
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 110262 2018-05-17 14:25:55Z 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  return _instance;
65 }
66 
68  fpRichPointsContainer(0),
69  fpCreatorProcess(0),
70  fCreatorModelID(0),
71  fpEndingProcess(0),
72  fFinalKineticEnergy(0.)
73 {
74 }
75 
77  G4Trajectory(aTrack) // Note: this initialises the base class data
78  // members and, unfortunately but never mind,
79  // creates a G4TrajectoryPoint in
80  // TrajectoryPointContainer that we cannot
81  // access because it's private. We store the
82  // same information (plus more) in a
83  // G4RichTrajectoryPoint in the
84  // RichTrajectoryPointsContainer
85 {
90  // On construction, set final values to initial values.
91  // Final values are updated at the addition of every step - see AppendStep.
96  // Insert the first rich trajectory point (see note above)...
98  fpRichPointsContainer->push_back(new G4RichTrajectoryPoint(aTrack));
99 }
100 
102  G4Trajectory(right)
103 {
113  for(size_t i=0;i<right.fpRichPointsContainer->size();i++)
114  {
115  G4RichTrajectoryPoint* rightPoint =
117  fpRichPointsContainer->push_back(new G4RichTrajectoryPoint(*rightPoint));
118  }
119 }
120 
122 {
123  if (fpRichPointsContainer) {
124  // fpRichPointsContainer->clearAndDestroy();
125  size_t i;
126  for(i=0;i<fpRichPointsContainer->size();i++){
127  delete (*fpRichPointsContainer)[i];
128  }
129  fpRichPointsContainer->clear();
130  delete fpRichPointsContainer;
131  }
132 }
133 
135 {
136  fpRichPointsContainer->push_back(new G4RichTrajectoryPoint(aStep));
137  // Except for first step, which is a sort of virtual step to start
138  // the track, compute the final values...
139  const G4Track* track = aStep->GetTrack();
140  const G4StepPoint* postStepPoint = aStep->GetPostStepPoint();
141  if (track->GetCurrentStepNumber() > 0) {
144  fpEndingProcess = postStepPoint->GetProcessDefinedStep();
146  aStep->GetPreStepPoint()->GetKineticEnergy() -
147  aStep->GetTotalEnergyDeposit();
148  }
149 }
150 
152 {
153  if(!secondTrajectory) return;
154 
155  G4RichTrajectory* seco = (G4RichTrajectory*)secondTrajectory;
156  G4int ent = seco->GetPointEntries();
157  for(G4int i=1;i<ent;i++) {
158  // initial point of the second trajectory should not be merged
159  fpRichPointsContainer->push_back((*(seco->fpRichPointsContainer))[i]);
160  // fpRichPointsContainer->push_back(seco->fpRichPointsContainer->removeAt(1));
161  }
162  delete (*seco->fpRichPointsContainer)[0];
163  seco->fpRichPointsContainer->clear();
164 }
165 
166 void G4RichTrajectory::ShowTrajectory(std::ostream& os) const
167 {
168  // Invoke the default implementation in G4VTrajectory...
170  // ... or override with your own code here.
171 }
172 
174 {
175  // Invoke the default implementation in G4VTrajectory...
177  // ... or override with your own code here.
178 }
179 
180 const std::map<G4String,G4AttDef>* G4RichTrajectory::GetAttDefs() const
181 {
182  G4bool isNew;
183  std::map<G4String,G4AttDef>* store
184  = G4AttDefStore::GetInstance("G4RichTrajectory",isNew);
185  if (isNew) {
186 
187  // Get att defs from base class...
188  *store = *(G4Trajectory::GetAttDefs());
189 
190  G4String ID;
191 
192  ID = "IVPath";
193  (*store)[ID] = G4AttDef(ID,"Initial Volume Path",
194  "Physics","","G4String");
195 
196  ID = "INVPath";
197  (*store)[ID] = G4AttDef(ID,"Initial Next Volume Path",
198  "Physics","","G4String");
199 
200  ID = "CPN";
201  (*store)[ID] = G4AttDef(ID,"Creator Process Name",
202  "Physics","","G4String");
203 
204  ID = "CPTN";
205  (*store)[ID] = G4AttDef(ID,"Creator Process Type Name",
206  "Physics","","G4String");
207 
208  ID = "CMID";
209  (*store)[ID] = G4AttDef(ID,"Creator Model ID",
210  "Physics","","G4int");
211 
212  ID = "CMN";
213  (*store)[ID] = G4AttDef(ID,"Creator Model Name",
214  "Physics","","G4String");
215 
216  ID = "FVPath";
217  (*store)[ID] = G4AttDef(ID,"Final Volume Path",
218  "Physics","","G4String");
219 
220  ID = "FNVPath";
221  (*store)[ID] = G4AttDef(ID,"Final Next Volume Path",
222  "Physics","","G4String");
223 
224  ID = "EPN";
225  (*store)[ID] = G4AttDef(ID,"Ending Process Name",
226  "Physics","","G4String");
227 
228  ID = "EPTN";
229  (*store)[ID] = G4AttDef(ID,"Ending Process Type Name",
230  "Physics","","G4String");
231 
232  ID = "FKE";
233  (*store)[ID] = G4AttDef(ID,"Final kinetic energy",
234  "Physics","G4BestUnit","G4double");
235 
236  }
237 
238  return store;
239 }
240 
241 static G4String Path(const G4TouchableHandle& th)
242 {
243  std::ostringstream oss;
244  G4int depth = th->GetHistoryDepth();
245  for (G4int i = depth; i >= 0; --i) {
246  oss << th->GetVolume(i)->GetName()
247  << ':' << th->GetCopyNumber(i);
248  if (i != 0) oss << '/';
249  }
250  return oss.str();
251 }
252 
253 std::vector<G4AttValue>* G4RichTrajectory::CreateAttValues() const
254 {
255  // Create base class att values...
256  std::vector<G4AttValue>* values = G4Trajectory::CreateAttValues();
257 
259  values->push_back(G4AttValue("IVPath",Path(fpInitialVolume),""));
260  } else {
261  values->push_back(G4AttValue("IVPath","None",""));
262  }
263 
265  values->push_back(G4AttValue("INVPath",Path(fpInitialNextVolume),""));
266  } else {
267  values->push_back(G4AttValue("INVPath","None",""));
268  }
269 
270  if (fpCreatorProcess) {
271  values->push_back
274  values->push_back
275  (G4AttValue("CPTN",G4VProcess::GetProcessTypeName(type),""));
276  values->push_back
278  const G4String& creatorModelName =
280  values->push_back(G4AttValue("CMN",creatorModelName,""));
281  } else {
282  values->push_back(G4AttValue("CPN","None",""));
283  values->push_back(G4AttValue("CPTN","None",""));
284  values->push_back(G4AttValue("CMID","None",""));
285  values->push_back(G4AttValue("CMN","None",""));
286  }
287 
289  values->push_back(G4AttValue("FVPath",Path(fpFinalVolume),""));
290  } else {
291  values->push_back(G4AttValue("FVPath","None",""));
292  }
293 
295  values->push_back(G4AttValue("FNVPath",Path(fpFinalNextVolume),""));
296  } else {
297  values->push_back(G4AttValue("FNVPath","None",""));
298  }
299 
300  if (fpEndingProcess) {
301  values->push_back(G4AttValue("EPN",fpEndingProcess->GetProcessName(),""));
303  values->push_back(G4AttValue("EPTN",G4VProcess::GetProcessTypeName(type),""));
304  } else {
305  values->push_back(G4AttValue("EPN","None",""));
306  values->push_back(G4AttValue("EPTN","None",""));
307  }
308 
309  values->push_back
310  (G4AttValue("FKE",G4BestUnit(fFinalKineticEnergy,"Energy"),""));
311 
312 #ifdef G4ATTDEBUG
313  G4cout << G4AttCheck(values,GetAttDefs());
314 #endif
315 
316  return values;
317 }
G4double GetKineticEnergy() const
const G4VProcess * fpEndingProcess
void AppendStep(const G4Step *aStep)
int GetPointEntries() const
G4double fFinalKineticEnergy
G4StepPoint * GetPreStepPoint() const
G4TouchableHandle fpFinalNextVolume
const G4VProcess * fpCreatorProcess
const G4TouchableHandle & GetTouchableHandle() const
virtual void ShowTrajectory(std::ostream &os=G4cout) const
virtual G4VPhysicalVolume * GetVolume(G4int depth=0) const
Definition: G4VTouchable.cc:44
G4int GetCurrentStepNumber() const
#define G4ThreadLocalStatic
Definition: tls.hh:68
void DrawTrajectory() const
static const G4String & GetProcessTypeName(G4ProcessType)
Definition: G4VProcess.cc:142
void MergeTrajectory(G4VTrajectory *secondTrajectory)
std::vector< G4VTrajectoryPoint * > RichTrajectoryPointsContainer
bool G4bool
Definition: G4Types.hh:79
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
static const G4String & GetModelName(G4int)
const G4TouchableHandle & GetNextTouchableHandle() const
G4TRACKING_DLL G4Allocator< G4RichTrajectory > *& aRichTrajectoryAllocator()
static G4String Path(const G4TouchableHandle &th)
G4ProcessType
G4Track * GetTrack() const
const G4String & GetProcessName() const
Definition: G4VProcess.hh:411
G4double GetTotalEnergyDeposit() const
G4ProcessType GetProcessType() const
Definition: G4VProcess.hh:417
Definition: G4Step.hh:76
G4StepPoint * GetPostStepPoint() const
G4double GetKineticEnergy() const
G4int GetCopyNumber(G4int depth=0) const
virtual std::vector< G4AttValue > * CreateAttValues() const
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
static G4String ConvertToString(G4bool boolVal)
Definition: G4UIcommand.cc:374
void ShowTrajectory(std::ostream &os=G4cout) const
std::map< G4String, G4AttDef > * GetInstance(const G4String &storeKey, G4bool &isNew)
G4TouchableHandle fpInitialNextVolume
int G4int
Definition: G4Types.hh:78
G4int GetCreatorModelID() const
RichTrajectoryPointsContainer * fpRichPointsContainer
virtual void DrawTrajectory() const
const G4VProcess * GetCreatorProcess() const
virtual ~G4RichTrajectory()
G4GLOB_DLL std::ostream G4cout
virtual std::vector< G4AttValue > * CreateAttValues() const
virtual G4int GetHistoryDepth() const
Definition: G4VTouchable.cc:79
const G4VProcess * GetProcessDefinedStep() const
G4TouchableHandle fpFinalVolume
const G4String & GetName() const
G4TouchableHandle fpInitialVolume