Geant4 Cross Reference |
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 // G4Track inline methods implementation 27 // 28 // Author: Katsuya Amako, KEK - 1996 29 // Revisions: Hisaya Kurashige, 1998-2011 30 // -------------------------------------------------------------------- 31 32 extern G4TRACK_DLL G4Allocator<G4Track>*& aTrackAllocator(); 33 34 //------------------------------------------------------------- 35 // To implement bi-directional association between G4Step and 36 // and G4Track, a combined usage of 'forward declaration' and 37 // 'include' is necessary. 38 //------------------------------------------------------------- 39 #include "G4Step.hh" 40 41 inline void* G4Track::operator new(std::size_t) 42 { 43 if(aTrackAllocator() == nullptr) 44 { 45 aTrackAllocator() = new G4Allocator<G4Track>; 46 } 47 return (void*) aTrackAllocator()->MallocSingle(); 48 } 49 50 inline void G4Track::operator delete(void* aTrack) 51 { 52 aTrackAllocator()->FreeSingle((G4Track*) aTrack); 53 } 54 55 inline G4bool G4Track::operator==(const G4Track& trk) 56 { 57 return (this == &trk); 58 } 59 60 inline G4bool G4Track::operator!=(const G4Track& trk) 61 { 62 return (this != &trk); 63 } 64 65 // dynamic particle 66 inline const G4DynamicParticle* G4Track::GetDynamicParticle() const 67 { 68 return fpDynamicParticle; 69 } 70 71 // particle definition 72 inline G4ParticleDefinition* G4Track::GetDefinition() const 73 { 74 return fpDynamicParticle->GetDefinition(); 75 } 76 77 // particle definition 78 inline const G4ParticleDefinition* G4Track::GetParticleDefinition() const 79 { 80 return fpDynamicParticle->GetParticleDefinition(); 81 } 82 83 // parent track ID 84 inline G4int G4Track::GetParentID() const 85 { 86 return fParentID; 87 } 88 89 inline void G4Track::SetParentID(const G4int aValue) 90 { 91 fParentID = aValue; 92 } 93 94 // current track ID 95 inline G4int G4Track::GetTrackID() const 96 { 97 return fTrackID; 98 } 99 100 inline void G4Track::SetTrackID(const G4int aValue) 101 { 102 fTrackID = aValue; 103 } 104 105 // position 106 inline const G4ThreeVector& G4Track::GetPosition() const 107 { 108 return fPosition; 109 } 110 111 inline void G4Track::SetPosition(const G4ThreeVector& aValue) 112 { 113 fPosition = aValue; 114 } 115 116 // global time 117 inline G4double G4Track::GetGlobalTime() const 118 { 119 return fGlobalTime; 120 } 121 122 inline void G4Track::SetGlobalTime(const G4double aValue) 123 { 124 fGlobalTime = aValue; 125 } 126 127 // local time 128 inline G4double G4Track::GetLocalTime() const 129 { 130 return fLocalTime; 131 } 132 133 inline void G4Track::SetLocalTime(const G4double aValue) 134 { 135 fLocalTime = aValue; 136 } 137 138 // proper time 139 inline G4double G4Track::GetProperTime() const 140 { 141 return fpDynamicParticle->GetProperTime(); 142 } 143 144 inline void G4Track::SetProperTime(const G4double aValue) 145 { 146 fpDynamicParticle->SetProperTime(aValue); 147 } 148 149 // velocity 150 inline G4double G4Track::GetVelocity() const 151 { 152 return (useGivenVelocity) ? fVelocity 153 : ((!is_OpticalPhoton) 154 ? CLHEP::c_light * fpDynamicParticle->GetBeta() 155 : CalculateVelocityForOpticalPhoton()); 156 } 157 158 inline G4double G4Track::CalculateVelocity() const 159 { 160 return GetVelocity(); 161 } 162 163 inline void G4Track::SetVelocity(G4double val) 164 { 165 fVelocity = val; 166 } 167 168 inline G4bool G4Track::UseGivenVelocity() const 169 { 170 return useGivenVelocity; 171 } 172 173 inline void G4Track::UseGivenVelocity(G4bool val) 174 { 175 useGivenVelocity = val; 176 } 177 178 // volume 179 inline G4VPhysicalVolume* G4Track::GetVolume() const 180 { 181 return (fpTouchable) ? fpTouchable->GetVolume() : nullptr; 182 } 183 184 inline G4VPhysicalVolume* G4Track::GetNextVolume() const 185 { 186 return (fpNextTouchable) ? fpNextTouchable->GetVolume() : nullptr; 187 } 188 189 // material - assuming that the pointer to G4Step is defined 190 inline const G4MaterialCutsCouple* G4Track::GetMaterialCutsCouple() const 191 { 192 return fpStep->GetPreStepPoint()->GetMaterialCutsCouple(); 193 } 194 195 inline const G4MaterialCutsCouple* G4Track::GetNextMaterialCutsCouple() const 196 { 197 return fpStep->GetPostStepPoint()->GetMaterialCutsCouple(); 198 } 199 200 inline G4Material* G4Track::GetMaterial() const 201 { 202 return fpStep->GetPreStepPoint()->GetMaterial(); 203 } 204 205 inline G4Material* G4Track::GetNextMaterial() const 206 { 207 return fpStep->GetPostStepPoint()->GetMaterial(); 208 } 209 210 // touchable 211 inline const G4VTouchable* G4Track::GetTouchable() const 212 { 213 return fpTouchable(); 214 } 215 216 inline const G4TouchableHandle& G4Track::GetTouchableHandle() const 217 { 218 return fpTouchable; 219 } 220 221 inline void G4Track::SetTouchableHandle(const G4TouchableHandle& apValue) 222 { 223 fpTouchable = apValue; 224 } 225 226 inline const G4VTouchable* G4Track::GetNextTouchable() const 227 { 228 return fpNextTouchable(); 229 } 230 231 inline const G4TouchableHandle& G4Track::GetNextTouchableHandle() const 232 { 233 return fpNextTouchable; 234 } 235 236 inline void G4Track::SetNextTouchableHandle(const G4TouchableHandle& apValue) 237 { 238 fpNextTouchable = apValue; 239 } 240 241 inline const G4VTouchable* G4Track::GetOriginTouchable() const 242 { 243 return fpOriginTouchable(); 244 } 245 246 inline const G4TouchableHandle& G4Track::GetOriginTouchableHandle() const 247 { 248 return fpOriginTouchable; 249 } 250 251 inline void G4Track::SetOriginTouchableHandle(const G4TouchableHandle& apValue) 252 { 253 fpOriginTouchable = apValue; 254 } 255 256 // kinetic energy 257 inline G4double G4Track::GetKineticEnergy() const 258 { 259 return fpDynamicParticle->GetKineticEnergy(); 260 } 261 262 inline void G4Track::SetKineticEnergy(const G4double aValue) 263 { 264 fpDynamicParticle->SetKineticEnergy(aValue); 265 } 266 267 // total energy 268 inline G4double G4Track::GetTotalEnergy() const 269 { 270 return fpDynamicParticle->GetTotalEnergy(); 271 } 272 273 // momentum 274 inline G4ThreeVector G4Track::GetMomentum() const 275 { 276 return fpDynamicParticle->GetMomentum(); 277 } 278 279 // momentum (direction) 280 inline const G4ThreeVector& G4Track::GetMomentumDirection() const 281 { 282 return fpDynamicParticle->GetMomentumDirection(); 283 } 284 285 inline void G4Track::SetMomentumDirection(const G4ThreeVector& aValue) 286 { 287 fpDynamicParticle->SetMomentumDirection(aValue); 288 } 289 290 // polarization 291 inline const G4ThreeVector& G4Track::GetPolarization() const 292 { 293 return fpDynamicParticle->GetPolarization(); 294 } 295 296 inline void G4Track::SetPolarization(const G4ThreeVector& aValue) 297 { 298 fpDynamicParticle->SetPolarization(aValue); 299 } 300 301 // track status 302 inline G4TrackStatus G4Track::GetTrackStatus() const 303 { 304 return fTrackStatus; 305 } 306 307 inline void G4Track::SetTrackStatus(const G4TrackStatus aTrackStatus) 308 { 309 fTrackStatus = aTrackStatus; 310 } 311 312 // track length 313 inline G4double G4Track::GetTrackLength() const 314 { 315 return fTrackLength; 316 } 317 318 inline void G4Track::AddTrackLength(const G4double aValue) 319 { 320 fTrackLength += aValue; 321 } 322 // Accumulated track length 323 324 // step number 325 inline G4int G4Track::GetCurrentStepNumber() const 326 { 327 return fCurrentStepNumber; 328 } 329 330 inline void G4Track::IncrementCurrentStepNumber() 331 { 332 ++fCurrentStepNumber; 333 } 334 335 // step length 336 inline G4double G4Track::GetStepLength() const 337 { 338 return fStepLength; 339 } 340 341 inline void G4Track::SetStepLength(G4double value) 342 { 343 fStepLength = value; 344 } 345 346 // vertex (where this track was created) information 347 inline const G4ThreeVector& G4Track::GetVertexPosition() const 348 { 349 return fVtxPosition; 350 } 351 352 inline void G4Track::SetVertexPosition(const G4ThreeVector& aValue) 353 { 354 fVtxPosition = aValue; 355 } 356 357 inline const G4ThreeVector& G4Track::GetVertexMomentumDirection() const 358 { 359 return fVtxMomentumDirection; 360 } 361 362 inline void G4Track::SetVertexMomentumDirection(const G4ThreeVector& aValue) 363 { 364 fVtxMomentumDirection = aValue; 365 } 366 367 inline G4double G4Track::GetVertexKineticEnergy() const 368 { 369 return fVtxKineticEnergy; 370 } 371 372 inline void G4Track::SetVertexKineticEnergy(const G4double aValue) 373 { 374 fVtxKineticEnergy = aValue; 375 } 376 377 inline const G4LogicalVolume* G4Track::GetLogicalVolumeAtVertex() const 378 { 379 return fpLVAtVertex; 380 } 381 382 inline void G4Track::SetLogicalVolumeAtVertex(const G4LogicalVolume* aValue) 383 { 384 fpLVAtVertex = aValue; 385 } 386 387 inline const G4VProcess* G4Track::GetCreatorProcess() const 388 { 389 // If the pointer is null, this means the track is created 390 // by the event generator, i.e. the primary track. If it is not 391 // null, it points to the process which created this track. 392 393 return fpCreatorProcess; 394 } 395 396 inline void G4Track::SetCreatorProcess(const G4VProcess* aValue) 397 { 398 fpCreatorProcess = aValue; 399 } 400 401 inline void G4Track::SetCreatorModelID(const G4int id) 402 { 403 fCreatorModelID = id; 404 } 405 406 inline G4int G4Track::GetCreatorModelID() const 407 { 408 return fCreatorModelID; 409 } 410 411 inline G4int G4Track::GetCreatorModelIndex() const 412 { 413 return G4PhysicsModelCatalog::GetModelIndex(fCreatorModelID); 414 } 415 416 inline const G4String G4Track::GetCreatorModelName() const 417 { 418 return G4PhysicsModelCatalog::GetModelNameFromID(fCreatorModelID); 419 } 420 421 inline const G4ParticleDefinition* G4Track::GetParentResonanceDef() const 422 { 423 return fParentResonanceDef; 424 } 425 426 inline void G4Track::SetParentResonanceDef(const G4ParticleDefinition* parentDef) 427 { 428 fParentResonanceDef = parentDef; 429 } 430 431 inline G4int G4Track::GetParentResonanceID() const 432 { 433 return fParentResonanceID; 434 } 435 436 inline void G4Track::SetParentResonanceID(const G4int parentID) 437 { 438 fParentResonanceID = parentID; 439 } 440 441 inline G4bool G4Track::HasParentResonance() const { 442 return ( fParentResonanceDef != nullptr ); 443 } 444 445 inline G4int G4Track::GetParentResonancePDGEncoding() const { 446 return ( fParentResonanceDef == nullptr ? 0 : fParentResonanceDef->GetPDGEncoding() ); 447 } 448 449 inline G4String G4Track::GetParentResonanceName() const { 450 return ( fParentResonanceDef == nullptr ? G4String() : fParentResonanceDef->GetParticleName() ); 451 } 452 453 inline G4double G4Track::GetParentResonanceMass() const { 454 return fParentResonanceID * CLHEP::keV; // the ID is the mass of the resonance in eV 455 } 456 457 // flag for "Below Threshold" 458 inline G4bool G4Track::IsBelowThreshold() const 459 { 460 return fBelowThreshold; 461 } 462 463 inline void G4Track::SetBelowThresholdFlag(G4bool value) 464 { 465 fBelowThreshold = value; 466 } 467 468 // flag for " Good for Tracking" 469 inline G4bool G4Track::IsGoodForTracking() const 470 { 471 return fGoodForTracking; 472 } 473 474 inline void G4Track::SetGoodForTrackingFlag(G4bool value) 475 { 476 fGoodForTracking = value; 477 } 478 479 // track weight 480 inline void G4Track::SetWeight(G4double aValue) 481 { 482 fWeight = aValue; 483 } 484 485 inline G4double G4Track::GetWeight() const 486 { 487 return fWeight; 488 } 489 490 // user information 491 inline G4VUserTrackInformation* G4Track::GetUserInformation() const 492 { 493 return fpUserInformation; 494 } 495 496 inline void G4Track::SetUserInformation(G4VUserTrackInformation* aValue) const 497 { 498 fpUserInformation = aValue; 499 } 500 501 inline const G4Step* G4Track::GetStep() const 502 { 503 return fpStep; 504 } 505 506 inline void G4Track::SetStep(const G4Step* aValue) 507 { 508 fpStep = aValue; 509 } 510 511 inline std::map<G4int, G4VAuxiliaryTrackInformation*>* 512 G4Track::GetAuxiliaryTrackInformationMap() const 513 { 514 return fpAuxiliaryTrackInformationMap; 515 } 516