Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4RichTrajectoryPoint.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: G4RichTrajectoryPoint.cc 110262 2018-05-17 14:25:55Z gcosmo $
28 //
29 //
30 // ---------------------------------------------------------------
31 //
32 // G4RichTrajectoryPoint.cc
33 //
34 // Contact:
35 // Questions and comments on G4TrajectoryPoint, on which this is based,
36 // should be sent to
37 // Katsuya Amako (e-mail: Katsuya.Amako@kek.jp)
38 // Makoto Asai (e-mail: asai@kekvax.kek.jp)
39 // Takashi Sasaki (e-mail: Takashi.Sasaki@kek.jp)
40 // and on the extended code to:
41 // John Allison (e-mail: John.Allison@manchester.ac.uk)
42 // Joseph Perl (e-mail: perl@slac.stanford.edu)
43 //
44 // ---------------------------------------------------------------
45 
46 #include "G4RichTrajectoryPoint.hh"
47 
48 #include "G4Track.hh"
49 #include "G4Step.hh"
50 #include "G4VProcess.hh"
51 
52 #include "G4AttDefStore.hh"
53 #include "G4AttDef.hh"
54 #include "G4AttValue.hh"
55 #include "G4UnitsTable.hh"
56 
57 //#define G4ATTDEBUG
58 #ifdef G4ATTDEBUG
59 #include "G4AttCheck.hh"
60 #endif
61 
62 #include <sstream>
63 
65 {
67  return _instance;
68 }
69 
71  fpAuxiliaryPointVector(0),
72  fTotEDep(0.),
73  fRemainingEnergy(0.),
74  fpProcess(0),
75  fPreStepPointStatus(fUndefined),
76  fPostStepPointStatus(fUndefined),
77  fPreStepPointGlobalTime(0),
78  fPostStepPointGlobalTime(0),
79  fPreStepPointWeight(1.),
80  fPostStepPointWeight(1.)
81 {}
82 
84  G4TrajectoryPoint(aTrack->GetPosition()),
85  fpAuxiliaryPointVector(0),
86  fTotEDep(0.),
87  fRemainingEnergy(aTrack->GetKineticEnergy()),
88  fpProcess(0),
89  fPreStepPointStatus(fUndefined),
90  fPostStepPointStatus(fUndefined),
91  fPreStepPointGlobalTime(aTrack->GetGlobalTime()),
92  fPostStepPointGlobalTime(aTrack->GetGlobalTime()),
93  fpPreStepPointVolume(aTrack->GetTouchableHandle()),
94  fpPostStepPointVolume(aTrack->GetNextTouchableHandle()),
95  fPreStepPointWeight(aTrack->GetWeight()),
96  fPostStepPointWeight(aTrack->GetWeight())
97 {}
98 
100  G4TrajectoryPoint(aStep->GetPostStepPoint()->GetPosition()),
101  fpAuxiliaryPointVector(aStep->GetPointerToVectorOfAuxiliaryPoints()),
102  fTotEDep(aStep->GetTotalEnergyDeposit())
103 {
104  G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
105  G4StepPoint* postStepPoint = aStep->GetPostStepPoint();
106  if (aStep->GetTrack()->GetCurrentStepNumber() <= 0) { // First step
108  } else {
109  fRemainingEnergy = preStepPoint->GetKineticEnergy() - fTotEDep;
110  }
111  fpProcess = postStepPoint->GetProcessDefinedStep();
112  fPreStepPointStatus = preStepPoint->GetStepStatus();
113  fPostStepPointStatus = postStepPoint->GetStepStatus();
114  fPreStepPointGlobalTime = preStepPoint->GetGlobalTime();
115  fPostStepPointGlobalTime = postStepPoint->GetGlobalTime();
116  fpPreStepPointVolume = preStepPoint->GetTouchableHandle();
117  fpPostStepPointVolume = postStepPoint->GetTouchableHandle();
118  fPreStepPointWeight = preStepPoint->GetWeight();
119  fPostStepPointWeight = postStepPoint->GetWeight();
120 }
121 
124  G4TrajectoryPoint(right),
125  fpAuxiliaryPointVector(right.fpAuxiliaryPointVector),
126  fTotEDep(right.fTotEDep),
127  fRemainingEnergy(right.fRemainingEnergy),
128  fpProcess(right.fpProcess),
129  fPreStepPointStatus(right.fPreStepPointStatus),
130  fPostStepPointStatus(right.fPostStepPointStatus),
131  fPreStepPointGlobalTime(right.fPreStepPointGlobalTime),
132  fPostStepPointGlobalTime(right.fPostStepPointGlobalTime),
133  fpPreStepPointVolume(right.fpPreStepPointVolume),
134  fpPostStepPointVolume(right.fpPostStepPointVolume),
135  fPreStepPointWeight(right.fPreStepPointWeight),
136  fPostStepPointWeight(right.fPostStepPointWeight)
137 {}
138 
140 {
142  /*
143  G4cout << "Deleting fpAuxiliaryPointVector at "
144  << (void*) fpAuxiliaryPointVector
145  << G4endl;
146  */
147  delete fpAuxiliaryPointVector;
148  }
149 }
150 
151 const std::map<G4String,G4AttDef>*
153 {
154  G4bool isNew;
155  std::map<G4String,G4AttDef>* store
156  = G4AttDefStore::GetInstance("G4RichTrajectoryPoint",isNew);
157  if (isNew) {
158 
159  // Copy base class att defs...
160  *store = *(G4TrajectoryPoint::GetAttDefs());
161 
162  G4String ID;
163 
164  ID = "Aux";
165  (*store)[ID] = G4AttDef(ID, "Auxiliary Point Position",
166  "Physics","G4BestUnit","G4ThreeVector");
167  ID = "TED";
168  (*store)[ID] = G4AttDef(ID,"Total Energy Deposit",
169  "Physics","G4BestUnit","G4double");
170  ID = "RE";
171  (*store)[ID] = G4AttDef(ID,"Remaining Energy",
172  "Physics","G4BestUnit","G4double");
173  ID = "PDS";
174  (*store)[ID] = G4AttDef(ID,"Process Defined Step",
175  "Physics","","G4String");
176  ID = "PTDS";
177  (*store)[ID] = G4AttDef(ID,"Process Type Defined Step",
178  "Physics","","G4String");
179  ID = "PreStatus";
180  (*store)[ID] = G4AttDef(ID,"Pre-step-point status",
181  "Physics","","G4String");
182  ID = "PostStatus";
183  (*store)[ID] = G4AttDef(ID,"Post-step-point status",
184  "Physics","","G4String");
185  ID = "PreT";
186  (*store)[ID] = G4AttDef(ID,"Pre-step-point global time",
187  "Physics","G4BestUnit","G4double");
188  ID = "PostT";
189  (*store)[ID] = G4AttDef(ID,"Post-step-point global time",
190  "Physics","G4BestUnit","G4double");
191  ID = "PreVPath";
192  (*store)[ID] = G4AttDef(ID,"Pre-step Volume Path",
193  "Physics","","G4String");
194  ID = "PostVPath";
195  (*store)[ID] = G4AttDef(ID,"Post-step Volume Path",
196  "Physics","","G4String");
197  ID = "PreW";
198  (*store)[ID] = G4AttDef(ID,"Pre-step-point weight",
199  "Physics","","G4double");
200  ID = "PostW";
201  (*store)[ID] = G4AttDef(ID,"Post-step-point weight",
202  "Physics","","G4double");
203  }
204  return store;
205 }
206 
208 {
209  G4String status;
210  switch (stps) {
211  case fWorldBoundary: status = "fWorldBoundary"; break;
212  case fGeomBoundary: status = "fGeomBoundary"; break;
213  case fAtRestDoItProc: status = "fAtRestDoItProc"; break;
214  case fAlongStepDoItProc: status = "fAlongStepDoItProc"; break;
215  case fPostStepDoItProc: status = "fPostStepDoItProc"; break;
216  case fUserDefinedLimit: status = "fUserDefinedLimit"; break;
217  case fExclusivelyForcedProc: status = "fExclusivelyForcedProc"; break;
218  case fUndefined: status = "fUndefined"; break;
219  default: status = "Not recognised"; break;
220  }
221  return status;
222 }
223 
224 static G4String Path(const G4TouchableHandle& th)
225 {
226  std::ostringstream oss;
227  G4int depth = th->GetHistoryDepth();
228  for (G4int i = depth; i >= 0; --i) {
229  oss << th->GetVolume(i)->GetName()
230  << ':' << th->GetCopyNumber(i);
231  if (i != 0) oss << '/';
232  }
233  return oss.str();
234 }
235 
236 std::vector<G4AttValue>* G4RichTrajectoryPoint::CreateAttValues() const
237 {
238  // Create base class att values...
239  std::vector<G4AttValue>* values = G4TrajectoryPoint::CreateAttValues();
240 
242  std::vector<G4ThreeVector>::iterator iAux;
243  for (iAux = fpAuxiliaryPointVector->begin();
244  iAux != fpAuxiliaryPointVector->end(); ++iAux) {
245  values->push_back(G4AttValue("Aux",G4BestUnit(*iAux,"Length"),""));
246  }
247  }
248 
249  values->push_back(G4AttValue("TED",G4BestUnit(fTotEDep,"Energy"),""));
250 
251  values->push_back(G4AttValue("RE",G4BestUnit(fRemainingEnergy,"Energy"),""));
252 
253  if (fpProcess) {
254  values->push_back
255  (G4AttValue("PDS",fpProcess->GetProcessName(),""));
256  values->push_back
257  (G4AttValue
259  ""));
260  } else {
261  values->push_back(G4AttValue("PDS","None",""));
262  values->push_back(G4AttValue("PTDS","None",""));
263  }
264 
265  values->push_back
266  (G4AttValue("PreStatus",Status(fPreStepPointStatus),""));
267 
268  values->push_back
269  (G4AttValue("PostStatus",Status(fPostStepPointStatus),""));
270 
271  values->push_back
272  (G4AttValue("PreT",G4BestUnit(fPreStepPointGlobalTime,"Time"),""));
273 
274  values->push_back
275  (G4AttValue("PostT",G4BestUnit(fPostStepPointGlobalTime,"Time"),""));
276 
278  values->push_back(G4AttValue("PreVPath",Path(fpPreStepPointVolume),""));
279  } else {
280  values->push_back(G4AttValue("PreVPath","None",""));
281  }
282 
284  values->push_back(G4AttValue("PostVPath",Path(fpPostStepPointVolume),""));
285  } else {
286  values->push_back(G4AttValue("PostVPath","None",""));
287  }
288 
289  {
290  std::ostringstream oss;
291  oss << fPreStepPointWeight;
292  values->push_back
293  (G4AttValue("PreW",oss.str(),""));
294  }
295 
296  {
297  std::ostringstream oss;
298  oss << fPostStepPointWeight;
299  values->push_back
300  (G4AttValue("PostW",oss.str(),""));
301  }
302 
303 #ifdef G4ATTDEBUG
304  G4cout << G4AttCheck(values,GetAttDefs());
305 #endif
306 
307  return values;
308 }
G4double GetKineticEnergy() const
G4StepPoint * GetPreStepPoint() const
G4TRACKING_DLL G4Allocator< G4RichTrajectoryPoint > *& aRichTrajectoryPointAllocator()
virtual G4VPhysicalVolume * GetVolume(G4int depth=0) const
Definition: G4VTouchable.cc:44
std::vector< G4ThreeVector > * fpAuxiliaryPointVector
G4int GetCurrentStepNumber() const
#define G4ThreadLocalStatic
Definition: tls.hh:68
const G4TouchableHandle & GetTouchableHandle() const
static const G4String & GetProcessTypeName(G4ProcessType)
Definition: G4VProcess.cc:142
bool G4bool
Definition: G4Types.hh:79
G4StepStatus GetStepStatus() const
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
G4Track * GetTrack() const
const G4String & GetProcessName() const
Definition: G4VProcess.hh:411
G4StepStatus
Definition: G4StepStatus.hh:49
G4double GetGlobalTime() 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
std::map< G4String, G4AttDef > * GetInstance(const G4String &storeKey, G4bool &isNew)
int G4int
Definition: G4Types.hh:78
static G4String Path(const G4TouchableHandle &th)
static G4String Status(G4StepStatus stps)
G4GLOB_DLL std::ostream G4cout
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
virtual std::vector< G4AttValue > * CreateAttValues() const
G4double GetWeight() const
virtual G4int GetHistoryDepth() const
Definition: G4VTouchable.cc:79
const G4VProcess * GetProcessDefinedStep() const
const G4VProcess * fpProcess
const G4String & GetName() const
G4TouchableHandle fpPreStepPointVolume
G4TouchableHandle fpPostStepPointVolume