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