TGeant4/TG4StepManager.cxx

   1 // $Id$
   2 // Category: event
   3 //
   4 // See the class description in the header file.
   5
   6 #include "TG4StepManager.h"
   7 #include "TG4GeometryManager.h"
   8 #include "TG4PhysicsManager.h"
   9 #include "TG4VSensitiveDetector.h"
  10 #include "TG4Globals.h"
  11 #include "TG3Units.h"
  12
  13 #include <G4SteppingManager.hh>
  14 #include <G4UserLimits.hh>
  15 #include <G4UImanager.hh>
  16 #include <G4AffineTransform.hh>
  17 #include <G4TransportationManager.hh>
  18 #include <G4Navigator.hh>
  19
  20 #include <Randomize.hh>
  21 #include <TLorentzVector.h>
  22
  23 TG4StepManager* TG4StepManager::fgInstance = 0;
  24
  25 TG4StepManager::TG4StepManager()
  26   : fTrack(0),
  27     fStep(0),
  28     fStepStatus(kNormalStep),
  29     fSteppingManager(0)
  30 {
  31 //
  32   if (fgInstance) {
  33     TG4Globals::Exception(
  34       "TG4StepManager: attempt to create two instances of singleton.");
  35   }
  36
  37   fgInstance = this;
  38 }
  39
  40 TG4StepManager::TG4StepManager(const TG4StepManager& right) {
  41 //
  42   TG4Globals::Exception(
  43     "Attempt to copy TG4StepManager singleton.");
  44 }
  45
  46 TG4StepManager::~TG4StepManager() {
  47 //
  48 }
  49
  50 // operators
  51
  52 TG4StepManager& TG4StepManager::operator=(const TG4StepManager& right)
  53 {
  54   // check assignement to self
  55   if (this == &right) return *this;
  56
  57   TG4Globals::Exception(
  58     "Attempt to assign TG4StepManager singleton.");
  59
  60   return *this;
  61 }
  62
  63 // private methods
  64
  65 void TG4StepManager::CheckTrack() const
  66 {
  67 // Gives exception in case the track is not defined.
  68 // ---
  69
  70   if (!fTrack)
  71     TG4Globals::Exception("TG4StepManager: Track is not defined.");
  72 }
  73
  74
  75 void TG4StepManager::CheckStep(const G4String& method) const
  76 {
  77 // Gives exception in case the step is not defined.
  78 // ---
  79
  80   if (!fStep) {
  81     G4String text = "TG4StepManager::";
  82     text = text + method + ": Step is not defined.";
  83     TG4Globals::Exception(text);
  84   }
  85 }
  86
  87
  88 void TG4StepManager::CheckSteppingManager() const
  89 {
  90 // Gives exception in case the step is not defined.
  91 // ---
  92
  93   if (!fSteppingManager)
  94     TG4Globals::Exception("TG4StepManager: Stepping manager is not defined.");
  95 }
  96
  97
  98 void TG4StepManager::SetTLorentzVector(G4ThreeVector xyz, G4double t,
  99                                        TLorentzVector& lv) const
 100 {
 101 // Fills TLorentzVector with G4ThreeVector and G4double.
 102 // ---
 103
 104    lv[0] = xyz.x();
 105    lv[1] = xyz.y();
 106    lv[2] = xyz.z();
 107    lv[3] = t;
 108 }
 109
 110 G4VPhysicalVolume* TG4StepManager::GetCurrentOffPhysicalVolume(G4int off) const
 111 {
 112 // Returns the physical volume of the off-th mother's
 113 // of the current volume.
 114 // ---
 115
 116   G4VPhysicalVolume* physVolume = GetCurrentPhysicalVolume();
 117
 118   G4VPhysicalVolume* mother = physVolume;
 119
 120   Int_t level = off;
 121   while (level > 0) {
 122     if (mother) mother = mother->GetMother();
 123     level--;
 124   }
 125
 126   if (!mother) {
 127     G4String text = "TG4StepManager::CurrentVolOff: \n";
 128     text = text + "    Volume ";
 129     text = text + physVolume->GetName();
 130     text = text + " has not defined mother in the required level.";
 131     TG4Globals::Warning(text);
 132   }
 133
 134   return mother;
 135 }
 136
 137 G4int TG4StepManager::GetVolumeID(G4VPhysicalVolume* physVolume) const
 138 {
 139 // Returns the sensitive detector ID of the specified
 140 // physical volume.
 141 // ---
 142
 143   // sensitive detector ID
 144   G4VSensitiveDetector* sd
 145     = physVolume->GetLogicalVolume()->GetSensitiveDetector();
 146   if (sd) {
 147     TG4VSensitiveDetector* tsd = dynamic_cast<TG4VSensitiveDetector*>(sd);
 148     if (tsd)
 149       return tsd->GetID();
 150     else {
 151       TG4Globals::Exception(
 152         "TG4StepManager::GetVolumeID: Unknown sensitive detector type");
 153       return 0;
 154     }
 155   }
 156   else {
 157     G4String text = "TG4StepManager::GetVolumeID: \n";
 158     text = text + "    Volume " + physVolume->GetName();
 159     text = text + " has not a sensitive detector.";
 160     //TG4Globals::Exception(text);
 161     TG4Globals::Warning(text);
 162     return 0;
 163   }
 164 }
 165
 166
 167 // public methods
 168
 169 void TG4StepManager::StopTrack()
 170 {
 171 // Stops the current track and skips to the next.
 172 // ?? do we want to invoke rest processes?
 173 // ?? do we want to stop secondaries too?
 174 //   possible "stop" track status from G4:
 175 //   fStopButAlive,      // Invoke active rest physics processes and
 176 //                       // and kill the current track afterward
 177 //   fStopAndKill,       // Kill the current track
 178 //   fKillTrackAndSecondaries,
 179 //                       // Kill the current track and also associated
 180 //                       // secondaries.
 181 // ---
 182
 183 #ifdef TGEANT4_DEBUG
 184   CheckTrack();
 185 #endif
 186
 187   fTrack->SetTrackStatus(fStopAndKill);
 188   // fTrack->SetTrackStatus(fStopButAlive);
 189   // fTrack->SetTrackStatus(fKillTrackAndSecondaries);
 190 }
 191
 192 void TG4StepManager::StopEvent()
 193 {
 194 // Aborts the current event processing.
 195 // ---
 196
 197 #ifdef TGEANT4_DEBUG
 198   CheckTrack();
 199 #endif
 200
 201   fTrack->SetTrackStatus(fKillTrackAndSecondaries);
 202           //StopTrack();   // cannot be used as it keeps secondaries
 203   G4UImanager::GetUIpointer()->ApplyCommand("/event/abort");
 204   G4UImanager::GetUIpointer()->ApplyCommand("/alStacking/clearStack");
 205 }
 206
 207 void TG4StepManager::Rndm(Float_t* array, const Int_t size) const
 208 {
 209 // Random numbers array of the specified size.
 210 // ---
 211
 212   G4double* const kpDoubleArray = new G4double[size];
 213   RandFlat::shootArray(size, kpDoubleArray);
 214   for (G4int i=0; i<size; i++) {
 215     array[i] = kpDoubleArray[i];
 216   }
 217   delete [] kpDoubleArray;
 218 }
 219
 220 void TG4StepManager::SetMaxStep(Float_t step)
 221 {
 222 // Maximum step allowed in the current logical volume.
 223 // The maximum value is kept for following tracks - is it ok ??
 224 // ---
 225
 226   // check this
 227   G4LogicalVolume* curLogVolume
 228     = GetCurrentPhysicalVolume()->GetLogicalVolume();
 229   G4UserLimits* userLimits
 230     = curLogVolume->GetUserLimits();
 231   if (userLimits == 0) {
 232     userLimits = new G4UserLimits(step);
 233     curLogVolume->SetUserLimits(userLimits);
 234   }
 235   else
 236     userLimits->SetMaxAllowedStep(step);
 237 }
 238
 239 void TG4StepManager::SetMaxNStep(Int_t maxNofSteps)
 240 {
 241 // Not yet implemented.
 242 // ---
 243
 244   TG4Globals::Warning(
 245     "TG4StepManager::SetMaxNStep(..) is not yet implemented.");
 246 }
 247
 248 void TG4StepManager::SetUserDecay(Int_t pdg)
 249 {
 250 // Not yet implemented.
 251 // ---
 252
 253   TG4Globals::Exception(
 254     "TG4StepManager::SetUserDecay(..) is not yet implemented.");
 255 }
 256
 257 G4VPhysicalVolume* TG4StepManager::GetCurrentPhysicalVolume() const
 258 {
 259 // Returns the current physical volume.
 260 // According to fStepStatus the volume from track vertex,
 261 // pre step point or post step point is returned.
 262 // ---
 263
 264   G4VPhysicalVolume* physVolume;
 265   if (fStepStatus == kNormalStep) {
 266
 267 #ifdef TGEANT4_DEBUG
 268     CheckStep("GetCurrentPhysicalVolume");
 269 #endif
 270
 271     physVolume = fStep->GetPreStepPoint()->GetPhysicalVolume();
 272   }
 273   else if (fStepStatus == kBoundary) {
 274
 275 #ifdef TGEANT4_DEBUG
 276     CheckStep("GetCurrentPhysicalVolume");
 277 #endif
 278
 279     physVolume = fStep->GetPostStepPoint()->GetPhysicalVolume();
 280   }
 281   else {
 282
 283 #ifdef TGEANT4_DEBUG
 284     CheckTrack();
 285 #endif
 286
 287     G4ThreeVector position = fTrack->GetPosition();
 288     G4Navigator* navigator =
 289       G4TransportationManager::GetTransportationManager()->
 290         GetNavigatorForTracking();
 291     physVolume
 292      = navigator->LocateGlobalPointAndSetup(position);
 293   }
 294
 295   return physVolume;
 296 }
 297
 298 Int_t TG4StepManager::CurrentVolID(Int_t& copyNo) const
 299 {
 300 // Returns the current sensitive detector ID
 301 // and the copy number of the current physical volume.
 302 // ---
 303
 304   G4VPhysicalVolume* physVolume = GetCurrentPhysicalVolume();
 305   copyNo = physVolume->GetCopyNo() + 1;
 306
 307   // sensitive detector ID
 308   return GetVolumeID(physVolume);
 309 }
 310
 311 Int_t TG4StepManager::CurrentVolOffID(Int_t off, Int_t&  copyNo) const
 312 {
 313 // Returns the off-th mother's of the current volume
 314 // the sensitive detector ID and the copy number.
 315 // ---
 316
 317   if (off == 0) return CurrentVolID(copyNo);
 318
 319   G4VPhysicalVolume* mother = GetCurrentOffPhysicalVolume(off);
 320
 321   if (mother) {
 322     copyNo = mother->GetCopyNo() + 1;
 323
 324     // sensitive detector ID
 325     return GetVolumeID(mother);
 326   }
 327   else {
 328     copyNo = 0;
 329     return 0;
 330   }
 331 }
 332
 333 const char* TG4StepManager::CurrentVolName() const
 334 {
 335 // Returns the current physical volume name.
 336 // ---
 337
 338   return GetCurrentPhysicalVolume()->GetName();
 339 }
 340
 341 const char* TG4StepManager::CurrentVolOffName(Int_t off) const
 342 {
 343 // Returns the off-th mother's physical volume name.
 344 // ---
 345
 346   if (off == 0) return CurrentVolName();
 347
 348   G4VPhysicalVolume* mother = GetCurrentOffPhysicalVolume(off);
 349
 350   if (mother)
 351     return mother->GetName();
 352   else
 353     return 0;
 354 }
 355
 356 Int_t TG4StepManager::CurrentMaterial(Float_t &a, Float_t &z, Float_t &dens,
 357                           Float_t &radl, Float_t &absl) const
 358 {
 359 // Returns the parameters of the current material during transport
 360 // the return value is the number of elements in the mixture.
 361 // ---
 362
 363   G4VPhysicalVolume* physVolume = GetCurrentPhysicalVolume();
 364
 365   G4Material* material
 366     = physVolume->GetLogicalVolume()->GetMaterial();
 367
 368   if (material) {
 369     G4int nofElements = material->GetNumberOfElements();
 370     TG4GeometryManager* pGeometryManager = TG4GeometryManager::Instance();
 371     a = pGeometryManager->GetEffA(material);
 372     z = pGeometryManager->GetEffZ(material);
 373
 374     // density
 375     dens = material->GetDensity();
 376     dens /= TG3Units::MassDensity();
 377
 378     // radiation length
 379     radl = material->GetRadlen();
 380     radl /= TG3Units::Length();
 381
 382     absl = 0.;  // this parameter is not defined in Geant4
 383     return nofElements;
 384   }
 385   else {
 386     TG4Globals::Exception(
 387      "TG4StepManager::CurrentMaterial(..): material is not defined.");
 388     return 0;
 389   }
 390 }
 391
 392 void TG4StepManager::Gmtod(Float_t* xm, Float_t* xd, Int_t iflag)
 393 {
 394 // Transforms a position from the world reference frame
 395 // to the current volume reference frame.
 396 //
 397 //  Geant3 desription:
 398 //  ==================
 399 //       Computes coordinates XD (in DRS)
 400 //       from known coordinates XM in MRS
 401 //       The local reference system can be initialized by
 402 //         - the tracking routines and GMTOD used in GUSTEP
 403 //         - a call to GMEDIA(XM,NUMED)
 404 //         - a call to GLVOLU(NLEVEL,NAMES,NUMBER,IER)
 405 //             (inverse routine is GDTOM)
 406 //
 407 //        If IFLAG=1  convert coordinates
 408 //           IFLAG=2  convert direction cosinus
 409 //
 410 // ---
 411
 412   G4AffineTransform affineTransform;
 413
 414   if (fStepStatus == kVertex) {
 415     G4Navigator* navigator =
 416       G4TransportationManager::GetTransportationManager()->
 417         GetNavigatorForTracking();
 418
 419     affineTransform = navigator->GetGlobalToLocalTransform();
 420   }
 421   else {
 422
 423 #ifdef TGEANT4_DEBUG
 424     CheckStep("Gmtod");
 425 #endif
 426
 427     affineTransform
 428       = fStep->GetPreStepPoint()->GetTouchable()->GetHistory()
 429         ->GetTopTransform();
 430   }
 431
 432   G4ThreeVector theGlobalPoint(xm[0],xm[1],xm[2]);
 433   G4ThreeVector theLocalPoint;
 434   if(iflag == 1)
 435        theLocalPoint = affineTransform.TransformPoint(theGlobalPoint);
 436   else if ( iflag == 2)
 437        theLocalPoint = affineTransform.TransformAxis(theGlobalPoint);
 438   else
 439     TG4Globals::Exception(
 440       "TG4StepManager::Gmtod(..,iflag): iflag is not in 1..2");
 441
 442   xd[0] = theLocalPoint.x();
 443   xd[1] = theLocalPoint.y();
 444   xd[2] = theLocalPoint.z();
 445 }
 446
 447 void TG4StepManager::Gdtom(Float_t* xd, Float_t* xm, Int_t iflag)
 448 {
 449 // Transforms a position from the current volume reference frame
 450 // to the world reference frame.
 451 //
 452 //  Geant3 desription:
 453 //  ==================
 454 //  Computes coordinates XM (Master Reference System
 455 //  knowing the coordinates XD (Detector Ref System)
 456 //  The local reference system can be initialized by
 457 //    - the tracking routines and GDTOM used in GUSTEP
 458 //    - a call to GSCMED(NLEVEL,NAMES,NUMBER)
 459 //        (inverse routine is GMTOD)
 460 //
 461 //   If IFLAG=1  convert coordinates
 462 //      IFLAG=2  convert direction cosinus
 463 //
 464 // ---
 465
 466   G4AffineTransform affineTransform;
 467
 468   if (fStepStatus == kVertex) {
 469     G4Navigator* navigator =
 470       G4TransportationManager::GetTransportationManager()->
 471         GetNavigatorForTracking();
 472
 473     affineTransform = navigator->GetLocalToGlobalTransform();
 474   }
 475   else {
 476
 477 #ifdef TGEANT4_DEBUG
 478     CheckStep("Gdtom");
 479 #endif
 480
 481     // check this
 482
 483     affineTransform
 484       = fStep->GetPreStepPoint()->GetTouchable()->GetHistory()
 485         ->GetTopTransform().Inverse();
 486   }
 487
 488   G4ThreeVector theLocalPoint(xd[0],xd[1],xd[2]);
 489   G4ThreeVector theGlobalPoint;
 490   if(iflag == 1)
 491        theGlobalPoint = affineTransform.TransformPoint(theLocalPoint);
 492   else if( iflag == 2)
 493        theGlobalPoint = affineTransform.TransformAxis(theLocalPoint);
 494   else
 495     TG4Globals::Warning(
 496       "TG4StepManager::Gdtom(...,iflag): iflag is not in 1..2");
 497
 498   xm[0] = theGlobalPoint.x();
 499   xm[1] = theGlobalPoint.y();
 500   xm[2] = theGlobalPoint.z();
 501 }
 502
 503 Float_t TG4StepManager::MaxStep() const
 504 {
 505 // Returns maximum step allowed in the current logical volume
 506 // by User Limits.
 507 // ---
 508
 509   G4LogicalVolume* curLogVolume
 510     = GetCurrentPhysicalVolume()->GetLogicalVolume();
 511
 512   // check this
 513   G4UserLimits* userLimits
 514     = curLogVolume->GetUserLimits();
 515
 516   G4double maxStep;
 517   if (userLimits == 0) {
 518     G4String text = "User Limits are not defined for log volume ";
 519     text = text + curLogVolume->GetName();
 520     TG4Globals::Warning(text);
 521     return FLT_MAX;
 522   }
 523   else {
 524     const G4Track& trackRef = *(fTrack);
 525     maxStep = userLimits->GetMaxAllowedStep(trackRef);
 526     maxStep /= TG3Units::Length();
 527     return maxStep;
 528   }
 529 }
 530
 531 Int_t TG4StepManager::GetMaxNStep() const
 532 {
 533 // Not yet implemented.
 534 // ---
 535
 536   TG4Globals::Warning(
 537     "Method GetMaxNStep is not yet implemented in TG4StepManager.");
 538   return 0;
 539 }
 540
 541 void TG4StepManager::TrackPosition(TLorentzVector& position) const
 542 {
 543 // Current particle position (in the world reference frame)
 544 // and the local time since the current track is created
 545 // (position of the PostStepPoint).
 546 // ---
 547
 548 #ifdef TGEANT4_DEBUG
 549   CheckTrack();
 550 #endif
 551
 552   // get position
 553   // check if this is == to PostStepPoint position !!
 554   G4ThreeVector positionVector = fTrack->GetPosition();
 555   positionVector *= 1./(TG3Units::Length());
 556
 557   // local time
 558   G4double time = fTrack->GetLocalTime();
 559   time /= TG3Units::Time();
 560
 561   SetTLorentzVector(positionVector, time, position);
 562 }
 563
 564 Int_t TG4StepManager::GetMedium() const
 565 {
 566 // Returns the second index of the current material (corresponding to
 567 // G3 tracking medium index).
 568 // ---
 569
 570   // current material
 571   G4Material* curMaterial
 572     = GetCurrentPhysicalVolume()->GetLogicalVolume()->GetMaterial();
 573
 574   // medium index
 575   TG4GeometryManager* pGeometryManager = TG4GeometryManager::Instance();
 576   return pGeometryManager->GetMediumId(curMaterial);
 577 }
 578
 579 void TG4StepManager::TrackMomentum(TLorentzVector& momentum) const
 580 {
 581 // Current particle "momentum" (px, py, pz, Etot).
 582 // ---
 583
 584 #ifdef TGEANT4_DEBUG
 585   CheckTrack();
 586 #endif
 587
 588   G4ThreeVector momentumVector = fTrack->GetMomentum();
 589   momentumVector *= 1./(TG3Units::Energy());
 590
 591   G4double energy = fTrack->GetDynamicParticle()->GetTotalEnergy();
 592   energy /= TG3Units::Energy();
 593
 594   SetTLorentzVector(momentumVector, energy, momentum);
 595 }
 596
 597 void TG4StepManager::TrackVertexPosition(TLorentzVector& position) const
 598 {
 599 // The vertex particle position (in the world reference frame)
 600 // and the local time since the current track is created.
 601 // ---
 602
 603 #ifdef TGEANT4_DEBUG
 604   CheckTrack();
 605 #endif
 606
 607   // position
 608   G4ThreeVector positionVector = fTrack->GetVertexPosition();
 609   positionVector *= 1./(TG3Units::Length());
 610
 611   // local time
 612   // to be checked
 613   G4double time = fTrack->GetLocalTime();
 614   time /= TG3Units::Time();
 615
 616   SetTLorentzVector(positionVector, time, position);
 617 }
 618
 619 void TG4StepManager::TrackVertexMomentum(TLorentzVector& momentum) const
 620 {
 621 // The vertex particle "momentum" (px, py, pz, Ekin)
 622 // to do: change Ekin -> Etot
 623 // ---
 624
 625 #ifdef TGEANT4_DEBUG
 626   CheckTrack();
 627 #endif
 628
 629   G4ThreeVector momentumVector = fTrack->GetVertexMomentumDirection();
 630   momentumVector *= 1./(TG3Units::Energy());
 631
 632   G4double energy = fTrack->GetVertexKineticEnergy();
 633   energy /= TG3Units::Energy();
 634
 635   SetTLorentzVector(momentumVector, energy, momentum);
 636 }
 637
 638 Float_t TG4StepManager::TrackStep() const
 639 {
 640 // Returns the current step length.
 641 // ---
 642
 643   G4double length;
 644   if (fStepStatus == kNormalStep) {
 645
 646 #ifdef TGEANT4_DEBUG
 647     CheckStep("TrackStep");
 648 #endif
 649
 650     length = fStep->GetStepLength();
 651     length /= TG3Units::Length();
 652   }
 653   else
 654     length = 0;
 655
 656   return length;
 657 }
 658
 659 Float_t TG4StepManager::TrackLength() const
 660 {
 661 // Returns the length of the current track from its origin.
 662 // ---
 663
 664 #ifdef TGEANT4_DEBUG
 665   CheckTrack();
 666 #endif
 667
 668   G4double length = fTrack->GetTrackLength();
 669   length /= TG3Units::Length();
 670   return length;
 671 }
 672
 673 Float_t TG4StepManager::TrackTime() const
 674 {
 675 // Returns the local time since the current track is created.
 676 // Comment:
 677 // in Geant4: there is also defined proper time as
 678 // the proper time of the dynamical particle of the current track.
 679 // ---
 680
 681 #ifdef TGEANT4_DEBUG
 682   CheckTrack();
 683 #endif
 684
 685   G4double time = fTrack->GetLocalTime();
 686   time /= TG3Units::Time();
 687   return time;
 688 }
 689
 690 Float_t TG4StepManager::Edep() const
 691 {
 692 // Returns total energy deposit in this step.
 693 // ---
 694
 695   G4double energyDeposit;
 696   if (fStepStatus == kNormalStep) {
 697
 698 #ifdef TGEANT4_DEBUG
 699     CheckStep("Edep");
 700 #endif
 701
 702     energyDeposit = fStep->GetTotalEnergyDeposit();
 703     energyDeposit /= TG3Units::Energy();
 704   }
 705   else
 706     energyDeposit = 0;
 707
 708   return energyDeposit;
 709 }
 710
 711 Int_t TG4StepManager::TrackPid() const
 712 {
 713 // Returns the current particle PDG encoding.
 714 // ---
 715
 716 #ifdef TGEANT4_DEBUG
 717   CheckTrack();
 718 #endif
 719
 720   G4ParticleDefinition* particle
 721     = fTrack->GetDynamicParticle()->GetDefinition();
 722
 723   // ask TG4PhysicsManager to get PDG encoding
 724   // (in order to get PDG from extended TDatabasePDG
 725   // in case the standard PDG code is not defined)
 726   TG4PhysicsManager* pPhysicsManager = TG4PhysicsManager::Instance();
 727   G4int pdgEncoding = pPhysicsManager->GetPDGEncodingFast(particle);
 728
 729   return pdgEncoding;
 730 }
 731
 732 Float_t TG4StepManager::TrackCharge() const
 733 {
 734 // Returns the current particle charge.
 735 // ---
 736
 737 #ifdef TGEANT4_DEBUG
 738   CheckTrack();
 739 #endif
 740
 741   G4double charge
 742     = fTrack->GetDynamicParticle()->GetDefinition()
 743       ->GetPDGCharge();
 744   charge /= TG3Units::Charge();
 745   return charge;
 746 }
 747
 748 Float_t TG4StepManager::TrackMass() const
 749 {
 750 // Returns current particle rest mass.
 751 // ---
 752
 753 #ifdef TGEANT4_DEBUG
 754   CheckTrack();
 755 #endif
 756
 757   G4double mass
 758     = fTrack->GetDynamicParticle()->GetDefinition()
 759       ->GetPDGMass();
 760   mass /= TG3Units::Mass();
 761   return mass;
 762 }
 763
 764 Float_t TG4StepManager::Etot() const
 765 {
 766 // Returns total energy of the current particle.
 767 // ---
 768
 769 #ifdef TGEANT4_DEBUG
 770   CheckTrack();
 771 #endif
 772
 773   G4double energy
 774     = fTrack->GetDynamicParticle()->GetTotalEnergy();
 775   energy /= TG3Units::Energy();
 776   return energy;
 777 }
 778
 779 Bool_t TG4StepManager::IsTrackInside() const
 780 {
 781 // Returns true if particle does not cross geometrical boundary
 782 // and is not in vertex.
 783 // ---
 784
 785   if (fStepStatus == kNormalStep  && !(IsTrackExiting()) ) {
 786     // track is always inside during a normal step
 787     return true;
 788   }
 789
 790   return false;
 791 }
 792
 793 Bool_t TG4StepManager::IsTrackEntering() const
 794 {
 795 // Returns true if particle cross a geometrical boundary
 796 // or is in vertex.
 797 // ---
 798
 799   if (fStepStatus != kNormalStep) {
 800     // track is entering during a vertex or boundary step
 801     return true;
 802   }
 803
 804   return false;
 805 }
 806
 807 Bool_t TG4StepManager::IsTrackExiting() const
 808 {
 809 // Returns true if particle cross a geometrical boundary.
 810 // ---
 811
 812   if (fStepStatus == kNormalStep) {
 813
 814 #ifdef TGEANT4_DEBUG
 815     CheckStep("IsTrackExiting");
 816 #endif
 817
 818     if (fStep->GetPostStepPoint()->GetStepStatus() == fGeomBoundary)
 819        return true;
 820   }
 821
 822   return false;
 823 }
 824
 825 Bool_t TG4StepManager::IsTrackOut() const
 826 {
 827 // Returns true if particle cross the world boundary
 828 // at post-step point.
 829 // ---
 830
 831   if (fStepStatus == kVertex) return false;
 832
 833 #ifdef TGEANT4_DEBUG
 834   CheckStep("IsTrackCut");
 835 #endif
 836
 837   // check
 838   G4StepStatus status
 839     = fStep->GetPostStepPoint()->GetStepStatus();
 840   if (status == fWorldBoundary)
 841     return true;
 842   else
 843     return false;
 844 }
 845
 846 Bool_t TG4StepManager::IsTrackStop() const
 847 {
 848 // Returns true if particle has stopped
 849 // or has been killed, suspended or postponed to next event.
 850 //
 851 // Possible track status from G4:
 852 //   fAlive,             // Continue the tracking
 853 //   fStopButAlive,      // Invoke active rest physics processes and
 854 //                            // and kill the current track afterward
 855 //   fStopAndKill,       // Kill the current track
 856 //   fKillTrackAndSecondaries,
 857 //                       // Kill the current track and also associated
 858 //                       // secondaries.
 859 //   fSuspend,           // Suspend the current track
 860 //   fPostponeToNextEvent
 861 //                       // Postpones the tracking of thecurrent track
 862 //                       // to the next event.
 863 // ---
 864
 865 #ifdef TGEANT4_DEBUG
 866   CheckTrack();
 867 #endif
 868
 869   // check
 870   G4TrackStatus status
 871      = fTrack->GetTrackStatus();
 872   if ((status == fStopAndKill) ||
 873       (status == fKillTrackAndSecondaries) ||
 874       (status == fSuspend) ||
 875       (status == fPostponeToNextEvent)) {
 876     return true;
 877   }
 878   else
 879     return false;
 880 }
 881
 882 Bool_t TG4StepManager::IsTrackDisappeared() const
 883 {
 884 // Returns true if particle has disappeared
 885 // (due to any physical process)
 886 // or has been killed, suspended or postponed to next event.
 887 // ---
 888
 889 #ifdef TGEANT4_DEBUG
 890   CheckTrack();
 891 #endif
 892
 893   // check
 894   G4TrackStatus status
 895      = fTrack->GetTrackStatus();
 896   if ((status == fStopButAlive) ||
 897       (status == fKillTrackAndSecondaries) ||
 898       (status == fSuspend) ||
 899       (status == fPostponeToNextEvent)) {
 900     return true;
 901   }
 902   else
 903     return false;
 904 }
 905
 906 Bool_t TG4StepManager::IsTrackAlive() const
 907 {
 908 // Returns true if particle continues tracking.
 909 // ---
 910
 911 #ifdef TGEANT4_DEBUG
 912   CheckTrack();
 913 #endif
 914
 915   G4TrackStatus status
 916      = fTrack->GetTrackStatus();
 917   if (status == fAlive)
 918     return true;
 919   else
 920     return false;
 921 }
 922
 923 Bool_t TG4StepManager::IsNewTrack() const
 924 {
 925 // Returns true when track performs the first step.
 926 // ---
 927
 928   if (fStepStatus == kVertex)
 929     return true;
 930   else
 931     return false;
 932 }
 933
 934 Int_t TG4StepManager::NSecondaries() const
 935 {
 936 // Returns the number of secondary particles generated
 937 // in the current step.
 938 // ---
 939
 940 #ifdef TGEANT4_DEBUG
 941   CheckSteppingManager();
 942 #endif
 943
 944   G4int nofSecondaries = 0;
 945   nofSecondaries += fSteppingManager->GetfN2ndariesAtRestDoIt();
 946   nofSecondaries += fSteppingManager->GetfN2ndariesAlongStepDoIt();
 947   nofSecondaries += fSteppingManager->GetfN2ndariesPostStepDoIt();
 948
 949   return nofSecondaries;
 950 }
 951
 952 void TG4StepManager::GetSecondary(Int_t index, Int_t& particleId,
 953                           TLorentzVector& position, TLorentzVector& momentum)
 954 {
 955 // Fills the parameters (particle encoding, position, momentum)
 956 // of the generated secondary particle which is specified by index.
 957 // !! Check if indexing of secondaries is same !!
 958 // ---
 959
 960 #ifdef TGEANT4_DEBUG
 961   CheckSteppingManager();
 962 #endif
 963
 964   G4int nofSecondaries = NSecondaries();
 965   G4TrackVector* secondaryTracks = fSteppingManager->GetSecondary();
 966
 967   if (secondaryTracks){
 968     if (index < nofSecondaries) {
 969
 970       // the index of the first secondary of this step
 971       G4int startIndex
 972         = secondaryTracks->entries() - nofSecondaries;
 973              // (the secondaryTracks vector contains secondaries
 974              // produced by the track at previous steps, too)
 975       G4Track* track
 976         = (*secondaryTracks)[startIndex + index];
 977
 978       // particle encoding
 979       particleId
 980         = track->GetDynamicParticle()->GetDefinition()->GetPDGEncoding();
 981
 982       // position & time
 983       G4ThreeVector positionVector = track->GetPosition();
 984       positionVector *= 1./(TG3Units::Length());
 985       G4double time = track->GetLocalTime();
 986       time /= TG3Units::Time();
 987       SetTLorentzVector(positionVector, time, position);
 988
 989       // momentum & energy
 990       G4ThreeVector momentumVector = track->GetMomentum();
 991       G4double energy = track->GetDynamicParticle()->GetTotalEnergy();
 992       energy /= TG3Units::Energy();
 993       SetTLorentzVector(momentumVector, energy, momentum);
 994     }
 995     else {
 996       TG4Globals::Exception(
 997         "TG4StepManager::GetSecondary(): wrong secondary track index.");
 998     }
 999   }
1000   else {
1001     TG4Globals::Exception(
1002       "TG4StepManager::GetSecondary(): secondary tracks vector is empty");
1003   }
1004 }
1005
1006 const char* TG4StepManager::ProdProcess() const
1007 {
1008 // Returns the name of the process that defined current step
1009 // (and may produce the secondary particles).
1010 // ---
1011
1012   if (fStepStatus == kVertex) return "NONE";
1013
1014   G4int nofSecondaries = NSecondaries();
1015   if (nofSecondaries == 0) return "NONE";
1016
1017 #ifdef TGEANT4_DEBUG
1018   CheckStep("ProdProcess");
1019 #endif
1020
1021   const G4VProcess* curProcess
1022     = fStep->GetPostStepPoint()->GetProcessDefinedStep();
1023
1024   G4String g4Name = curProcess->GetProcessName();
1025   return g4Name;
1026 }