Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4RichTrajectoryPoint Class Reference

#include <G4RichTrajectoryPoint.hh>

+ Inheritance diagram for G4RichTrajectoryPoint:

Public Member Functions

 G4RichTrajectoryPoint ()
 
 G4RichTrajectoryPoint (const G4Track *)
 
 G4RichTrajectoryPoint (const G4Step *)
 
 G4RichTrajectoryPoint (const G4RichTrajectoryPoint &right)
 
virtual ~G4RichTrajectoryPoint ()
 
const std::vector< G4ThreeVector > * GetAuxiliaryPoints () const
 
void * operator new (size_t)
 
void operator delete (void *aRichTrajectoryPoint)
 
int operator== (const G4RichTrajectoryPoint &right) const
 
virtual const std::map< G4String, G4AttDef > * GetAttDefs () const
 
virtual std::vector< G4AttValue > * CreateAttValues () const
 
- Public Member Functions inherited from G4TrajectoryPoint
 G4TrajectoryPoint ()
 
 G4TrajectoryPoint (G4ThreeVector pos)
 
 G4TrajectoryPoint (const G4TrajectoryPoint &right)
 
virtual ~G4TrajectoryPoint ()
 
void * operator new (size_t)
 
void operator delete (void *aTrajectoryPoint)
 
int operator== (const G4TrajectoryPoint &right) const
 
const G4ThreeVector GetPosition () const
 
virtual const std::map< G4String, G4AttDef > * GetAttDefs () const
 
virtual std::vector< G4AttValue > * CreateAttValues () const
 
- Public Member Functions inherited from G4VTrajectoryPoint
 G4VTrajectoryPoint ()
 
virtual ~G4VTrajectoryPoint ()
 
G4bool operator== (const G4VTrajectoryPoint &right) const
 
virtual const G4ThreeVector GetPosition () const =0
 
virtual const std::vector< G4ThreeVector > * GetAuxiliaryPoints () const
 
virtual const std::map< G4String, G4AttDef > * GetAttDefs () const
 
virtual std::vector< G4AttValue > * CreateAttValues () const
 

Detailed Description

Definition at line 71 of file G4RichTrajectoryPoint.hh.

Constructor & Destructor Documentation

◆ G4RichTrajectoryPoint() [1/4]

G4RichTrajectoryPoint::G4RichTrajectoryPoint ( )

Definition at line 66 of file G4RichTrajectoryPoint.cc.

66 :
67 fpAuxiliaryPointVector(0),
68 fTotEDep(0.),
69 fRemainingEnergy(0.),
70 fpProcess(0),
71 fPreStepPointStatus(fUndefined),
72 fPostStepPointStatus(fUndefined),
73 fPreStepPointGlobalTime(0),
74 fPostStepPointGlobalTime(0),
75 fPreStepPointWeight(1.),
76 fPostStepPointWeight(1.)
77{}
@ fUndefined
Definition: G4StepStatus.hh:66

◆ G4RichTrajectoryPoint() [2/4]

G4RichTrajectoryPoint::G4RichTrajectoryPoint ( const G4Track aTrack)

Definition at line 79 of file G4RichTrajectoryPoint.cc.

79 :
81 fpAuxiliaryPointVector(0),
82 fTotEDep(0.),
83 fRemainingEnergy(aTrack->GetKineticEnergy()),
84 fpProcess(0),
85 fPreStepPointStatus(fUndefined),
86 fPostStepPointStatus(fUndefined),
87 fPreStepPointGlobalTime(aTrack->GetGlobalTime()),
88 fPostStepPointGlobalTime(aTrack->GetGlobalTime()),
89 fpPreStepPointVolume(aTrack->GetTouchableHandle()),
90 fpPostStepPointVolume(aTrack->GetNextTouchableHandle()),
91 fPreStepPointWeight(aTrack->GetWeight()),
92 fPostStepPointWeight(aTrack->GetWeight())
93{}
const G4TouchableHandle & GetNextTouchableHandle() const
G4double GetWeight() const
const G4ThreeVector & GetPosition() const
G4double GetGlobalTime() const
const G4TouchableHandle & GetTouchableHandle() const
G4double GetKineticEnergy() const

◆ G4RichTrajectoryPoint() [3/4]

G4RichTrajectoryPoint::G4RichTrajectoryPoint ( const G4Step aStep)

Definition at line 95 of file G4RichTrajectoryPoint.cc.

95 :
97 fpAuxiliaryPointVector(aStep->GetPointerToVectorOfAuxiliaryPoints()),
98 fTotEDep(aStep->GetTotalEnergyDeposit())
99{
100 G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
101 G4StepPoint* postStepPoint = aStep->GetPostStepPoint();
102 if (aStep->GetTrack()->GetCurrentStepNumber() <= 0) { // First step
103 fRemainingEnergy = aStep->GetTrack()->GetKineticEnergy();
104 } else {
105 fRemainingEnergy = preStepPoint->GetKineticEnergy() - fTotEDep;
106 }
107 fpProcess = postStepPoint->GetProcessDefinedStep();
108 fPreStepPointStatus = preStepPoint->GetStepStatus();
109 fPostStepPointStatus = postStepPoint->GetStepStatus();
110 fPreStepPointGlobalTime = preStepPoint->GetGlobalTime();
111 fPostStepPointGlobalTime = postStepPoint->GetGlobalTime();
112 fpPreStepPointVolume = preStepPoint->GetTouchableHandle();
113 fpPostStepPointVolume = postStepPoint->GetTouchableHandle();
114 fPreStepPointWeight = preStepPoint->GetWeight();
115 fPostStepPointWeight = postStepPoint->GetWeight();
116
117 /*
118 G4cout << "fpAuxiliaryPointVector "
119 << (void*) fpAuxiliaryPointVector;
120 G4cout << ": ";
121 if (fpAuxiliaryPointVector) {
122 G4cout << "size: " << fpAuxiliaryPointVector->size();
123 for (size_t i = 0; i < fpAuxiliaryPointVector->size(); ++i)
124 G4cout << "\n " << (*fpAuxiliaryPointVector)[i];
125 } else {
126 G4cout << "non-existent";
127 }
128 G4cout << G4endl;
129
130 static const G4Step* lastStep = 0;
131 if (aStep && aStep == lastStep) {
132 G4cout << "********* aStep is same as last" << G4endl;
133 }
134 lastStep = aStep;
135
136 static std::vector<G4ThreeVector>* lastAuxiliaryPointVector = 0;
137 if (fpAuxiliaryPointVector &&
138 fpAuxiliaryPointVector == lastAuxiliaryPointVector) {
139 G4cout << "********* fpAuxiliaryPointVector is same as last" << G4endl;
140 }
141 lastAuxiliaryPointVector = fpAuxiliaryPointVector;
142 */
143}
G4StepStatus GetStepStatus() const
G4double GetGlobalTime() const
const G4VProcess * GetProcessDefinedStep() const
const G4ThreeVector & GetPosition() const
const G4TouchableHandle & GetTouchableHandle() const
G4double GetKineticEnergy() const
G4double GetWeight() const
std::vector< G4ThreeVector > * GetPointerToVectorOfAuxiliaryPoints() const
Definition: G4Step.hh:240
G4Track * GetTrack() const
G4StepPoint * GetPreStepPoint() const
G4double GetTotalEnergyDeposit() const
G4StepPoint * GetPostStepPoint() const
G4int GetCurrentStepNumber() const

◆ G4RichTrajectoryPoint() [4/4]

G4RichTrajectoryPoint::G4RichTrajectoryPoint ( const G4RichTrajectoryPoint right)

Definition at line 145 of file G4RichTrajectoryPoint.cc.

146 :
147 G4TrajectoryPoint(right),
148 fpAuxiliaryPointVector(right.fpAuxiliaryPointVector),
149 fTotEDep(right.fTotEDep),
150 fRemainingEnergy(right.fRemainingEnergy),
151 fpProcess(right.fpProcess),
152 fPreStepPointStatus(right.fPreStepPointStatus),
153 fPostStepPointStatus(right.fPostStepPointStatus),
154 fPreStepPointGlobalTime(right.fPreStepPointGlobalTime),
155 fPostStepPointGlobalTime(right.fPostStepPointGlobalTime),
156 fpPreStepPointVolume(right.fpPreStepPointVolume),
157 fpPostStepPointVolume(right.fpPostStepPointVolume),
158 fPreStepPointWeight(right.fPreStepPointWeight),
159 fPostStepPointWeight(right.fPostStepPointWeight)
160{}

◆ ~G4RichTrajectoryPoint()

G4RichTrajectoryPoint::~G4RichTrajectoryPoint ( )
virtual

Definition at line 162 of file G4RichTrajectoryPoint.cc.

163{
164 if(fpAuxiliaryPointVector) {
165 /*
166 G4cout << "Deleting fpAuxiliaryPointVector at "
167 << (void*) fpAuxiliaryPointVector
168 << G4endl;
169 */
170 delete fpAuxiliaryPointVector;
171 }
172}

Member Function Documentation

◆ CreateAttValues()

std::vector< G4AttValue > * G4RichTrajectoryPoint::CreateAttValues ( ) const
virtual

Reimplemented from G4TrajectoryPoint.

Definition at line 259 of file G4RichTrajectoryPoint.cc.

260{
261 // Create base class att values...
262 std::vector<G4AttValue>* values = G4TrajectoryPoint::CreateAttValues();
263
264 if (fpAuxiliaryPointVector) {
265 std::vector<G4ThreeVector>::iterator iAux;
266 for (iAux = fpAuxiliaryPointVector->begin();
267 iAux != fpAuxiliaryPointVector->end(); ++iAux) {
268 values->push_back(G4AttValue("Aux",G4BestUnit(*iAux,"Length"),""));
269 }
270 }
271
272 values->push_back(G4AttValue("TED",G4BestUnit(fTotEDep,"Energy"),""));
273
274 values->push_back(G4AttValue("RE",G4BestUnit(fRemainingEnergy,"Energy"),""));
275
276 if (fpProcess) {
277 values->push_back
278 (G4AttValue("PDS",fpProcess->GetProcessName(),""));
279 values->push_back
282 ""));
283 } else {
284 values->push_back(G4AttValue("PDS","None",""));
285 values->push_back(G4AttValue("PTDS","None",""));
286 }
287
288 values->push_back
289 (G4AttValue("PreStatus",Status(fPreStepPointStatus),""));
290
291 values->push_back
292 (G4AttValue("PostStatus",Status(fPostStepPointStatus),""));
293
294 values->push_back
295 (G4AttValue("PreT",G4BestUnit(fPreStepPointGlobalTime,"Time"),""));
296
297 values->push_back
298 (G4AttValue("PostT",G4BestUnit(fPostStepPointGlobalTime,"Time"),""));
299
300 if (fpPreStepPointVolume && fpPreStepPointVolume->GetVolume()) {
301 values->push_back(G4AttValue("PreVPath",Path(fpPreStepPointVolume),""));
302 } else {
303 values->push_back(G4AttValue("PreVPath","None",""));
304 }
305
306 if (fpPostStepPointVolume && fpPostStepPointVolume->GetVolume()) {
307 values->push_back(G4AttValue("PostVPath",Path(fpPostStepPointVolume),""));
308 } else {
309 values->push_back(G4AttValue("PostVPath","None",""));
310 }
311
312 {
313 std::ostringstream oss;
314 oss << fPreStepPointWeight;
315 values->push_back
316 (G4AttValue("PreW",oss.str(),""));
317 }
318
319 {
320 std::ostringstream oss;
321 oss << fPostStepPointWeight;
322 values->push_back
323 (G4AttValue("PostW",oss.str(),""));
324 }
325
326#ifdef G4ATTDEBUG
327 G4cout << G4AttCheck(values,GetAttDefs());
328#endif
329
330 return values;
331}
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
G4DLLIMPORT std::ostream G4cout
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
virtual std::vector< G4AttValue > * CreateAttValues() const
static const G4String & GetProcessTypeName(G4ProcessType)
Definition: G4VProcess.cc:150
G4ProcessType GetProcessType() const
Definition: G4VProcess.hh:385
const G4String & GetProcessName() const
Definition: G4VProcess.hh:379
virtual G4VPhysicalVolume * GetVolume(G4int depth=0) const
Definition: G4VTouchable.cc:44

◆ GetAttDefs()

const std::map< G4String, G4AttDef > * G4RichTrajectoryPoint::GetAttDefs ( ) const
virtual

Reimplemented from G4TrajectoryPoint.

Definition at line 175 of file G4RichTrajectoryPoint.cc.

176{
177 G4bool isNew;
178 std::map<G4String,G4AttDef>* store
179 = G4AttDefStore::GetInstance("G4RichTrajectoryPoint",isNew);
180 if (isNew) {
181
182 // Copy base class att defs...
183 *store = *(G4TrajectoryPoint::GetAttDefs());
184
185 G4String ID;
186
187 ID = "Aux";
188 (*store)[ID] = G4AttDef(ID, "Auxiliary Point Position",
189 "Physics","G4BestUnit","G4ThreeVector");
190 ID = "TED";
191 (*store)[ID] = G4AttDef(ID,"Total Energy Deposit",
192 "Physics","G4BestUnit","G4double");
193 ID = "RE";
194 (*store)[ID] = G4AttDef(ID,"Remaining Energy",
195 "Physics","G4BestUnit","G4double");
196 ID = "PDS";
197 (*store)[ID] = G4AttDef(ID,"Process Defined Step",
198 "Physics","","G4String");
199 ID = "PTDS";
200 (*store)[ID] = G4AttDef(ID,"Process Type Defined Step",
201 "Physics","","G4String");
202 ID = "PreStatus";
203 (*store)[ID] = G4AttDef(ID,"Pre-step-point status",
204 "Physics","","G4String");
205 ID = "PostStatus";
206 (*store)[ID] = G4AttDef(ID,"Post-step-point status",
207 "Physics","","G4String");
208 ID = "PreT";
209 (*store)[ID] = G4AttDef(ID,"Pre-step-point global time",
210 "Physics","G4BestUnit","G4double");
211 ID = "PostT";
212 (*store)[ID] = G4AttDef(ID,"Post-step-point global time",
213 "Physics","G4BestUnit","G4double");
214 ID = "PreVPath";
215 (*store)[ID] = G4AttDef(ID,"Pre-step Volume Path",
216 "Physics","","G4String");
217 ID = "PostVPath";
218 (*store)[ID] = G4AttDef(ID,"Post-step Volume Path",
219 "Physics","","G4String");
220 ID = "PreW";
221 (*store)[ID] = G4AttDef(ID,"Pre-step-point weight",
222 "Physics","","G4double");
223 ID = "PostW";
224 (*store)[ID] = G4AttDef(ID,"Post-step-point weight",
225 "Physics","","G4double");
226 }
227 return store;
228}
bool G4bool
Definition: G4Types.hh:67
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
std::map< G4String, G4AttDef > * GetInstance(G4String storeKey, G4bool &isNew)

Referenced by CreateAttValues(), G4VisCommandList::SetNewValue(), and G4VisCommandSceneAddTrajectories::SetNewValue().

◆ GetAuxiliaryPoints()

const std::vector< G4ThreeVector > * G4RichTrajectoryPoint::GetAuxiliaryPoints ( ) const
inlinevirtual

Reimplemented from G4VTrajectoryPoint.

Definition at line 89 of file G4RichTrajectoryPoint.hh.

90 { return fpAuxiliaryPointVector; }

◆ operator delete()

void G4RichTrajectoryPoint::operator delete ( void *  aRichTrajectoryPoint)
inline

Definition at line 135 of file G4RichTrajectoryPoint.hh.

137{
139 ((G4RichTrajectoryPoint *) aRichTrajectoryPoint);
140}
G4DLLIMPORT G4Allocator< G4RichTrajectoryPoint > aRichTrajectoryPointAllocator

◆ operator new()

void * G4RichTrajectoryPoint::operator new ( size_t  )
inline

Definition at line 127 of file G4RichTrajectoryPoint.hh.

128{
129 void *aRichTrajectoryPoint;
130 aRichTrajectoryPoint =
131 (void *) aRichTrajectoryPointAllocator.MallocSingle();
132 return aRichTrajectoryPoint;
133}

◆ operator==()

int G4RichTrajectoryPoint::operator== ( const G4RichTrajectoryPoint right) const
inline

Definition at line 95 of file G4RichTrajectoryPoint.hh.

96 { return (this==&right); }

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