AliGeant4/AliTrackingAction.cxx

   1 // $Id$
   2 // Category: event
   3 //
   4 // See the class description in the header file.
   5
   6 #include "AliTrackingAction.h"
   7 #include "AliTrackingActionMessenger.h"
   8 #include "AliTrackInformation.h"
   9 #include "AliRun.h"
  10 #include "AliGlobals.h"
  11 #include "TG4StepManager.h"
  12 #include "TG4PhysicsManager.h"
  13
  14 #include <G4Track.hh>
  15 #include <G4TrackVector.hh>
  16 #include <G4VUserTrackInformation.hh>
  17 #include <G4TrackingManager.hh>
  18 #include <G4SteppingManager.hh>
  19
  20 /*
  21 #include <TTree.h>
  22 #include <TParticle.h>
  23 #include <TClonesArray.h>
  24 */
  25
  26 // static data members
  27 AliTrackingAction* AliTrackingAction::fgInstance = 0;
  28
  29 AliTrackingAction::AliTrackingAction()
  30   : fPrimaryTrackID(0),
  31     fVerboseLevel(2),
  32     fSavePrimaries(true),
  33     fTrackCounter(0)
  34 {
  35 //
  36   if (fgInstance) {
  37     AliGlobals::Exception("AliTrackingAction constructed twice.");
  38   }
  39
  40   fMessenger = new AliTrackingActionMessenger(this);
  41   fgInstance = this;
  42 }
  43
  44 AliTrackingAction::AliTrackingAction(const AliTrackingAction& right) {
  45 //
  46   AliGlobals::Exception("AliTrackingAction is protected from copying.");
  47 }
  48
  49 AliTrackingAction::~AliTrackingAction() {
  50 //
  51   delete fMessenger;
  52 }
  53
  54 // operators
  55
  56 AliTrackingAction&
  57 AliTrackingAction::operator=(const AliTrackingAction &right)
  58 {
  59   // check assignement to self
  60   if (this == &right) return *this;
  61
  62   AliGlobals::Exception("AliTrackingAction is protected from assigning.");
  63
  64   return *this;
  65 }
  66
  67 // private methods
  68
  69 AliTrackInformation* AliTrackingAction::GetTrackInformation(
  70                                            const G4Track* track,
  71                                            const G4String& method) const
  72 {
  73 // Returns user track information.
  74 // ---
  75
  76   G4VUserTrackInformation* trackInfo = track->GetUserInformation();
  77   if (!trackInfo) return 0;
  78
  79   AliTrackInformation* aliTrackInfo
  80     = dynamic_cast<AliTrackInformation*>(trackInfo);
  81   if (!aliTrackInfo) {
  82      G4String text = "AliTrackingAction::" + method + ":\n";
  83      text = text + "   Unknown track information type";
  84      AliGlobals::Exception(text);
  85   }
  86
  87   return aliTrackInfo;
  88 }
  89
  90 // public methods
  91
  92 void AliTrackingAction::PrepareNewEvent()
  93 {
  94 // Called by G4 kernel at the beginning of event.
  95 // ---
  96
  97   fTrackCounter = 0;
  98
  99   // set g4 stepping manager pointer
 100   TG4StepManager* stepManager = TG4StepManager::Instance();
 101   stepManager->SetSteppingManager(fpTrackingManager->GetSteppingManager());
 102 }
 103
 104 void AliTrackingAction::PreTrackingAction(const G4Track* aTrack)
 105 {
 106 // Called by G4 kernel before starting tracking.
 107 // ---
 108
 109   // track index in the particles array
 110   G4int trackID = aTrack->GetTrackID();
 111   G4int parentID = aTrack->GetParentID();
 112   Int_t trackIndex;
 113   if (parentID==0) {
 114     // in AliRoot (from V3.0) track numbering starts from 0
 115     trackIndex = trackID-1;
 116   }
 117   else {
 118     trackIndex = gAlice->GetNtrack();
 119   }
 120
 121   // set track index to track information
 122   AliTrackInformation* trackInfo
 123     = GetTrackInformation(aTrack, "PreTrackingAction");
 124   if (!trackInfo) {
 125     // create track information and set it to G4Track
 126     // if it does not yet exist
 127     trackInfo = new AliTrackInformation(trackIndex);
 128     fpTrackingManager->SetUserTrackInformation(trackInfo);
 129         // the track information is deleted together with its
 130         // G4Track object
 131   }
 132   else
 133     trackInfo->SetTrackParticleID(trackIndex);
 134
 135   // set current track number
 136   gAlice->SetCurrentTrack(trackIndex);
 137
 138   if (parentID == 0) {
 139     // finish previous primary track
 140     FinishPrimaryTrack();
 141     fPrimaryTrackID = aTrack->GetTrackID();
 142   }
 143   else {
 144     // save secondary particles info
 145     SaveTrack(aTrack);
 146   }
 147
 148   // aliroot pre track actions
 149   gAlice->PreTrack();
 150 }
 151
 152 void AliTrackingAction::PostTrackingAction(const G4Track* aTrack)
 153 {
 154 // Called by G4 kernel after finishing tracking.
 155 // ---
 156
 157   fTrackCounter++;
 158
 159   // set parent track particle index to all secondary tracks
 160   G4TrackVector* secondaryTracks
 161     = fpTrackingManager->GetSteppingManager()->GetSecondary();
 162   if (secondaryTracks){
 163     G4int i;
 164     for (i=0; i<secondaryTracks->entries(); i++) {
 165       G4Track* track = (*secondaryTracks)[i];
 166
 167       if (track->GetUserInformation()) {
 168         // this should never happen
 169         G4String text = "AliTrackingAction::PostTrackingAction:\n";
 170         text = text + "    Inconsistent track information.";
 171         AliGlobals::Exception(text);
 172       }
 173
 174       // get parent track index
 175       AliTrackInformation* aliParentInfo
 176         = GetTrackInformation(aTrack, "PostTrackingAction");
 177       G4int parentParticleID
 178         = aliParentInfo->GetTrackParticleID();
 179
 180       // create track information and set it to the G4Track
 181       AliTrackInformation* trackInfo
 182         = new AliTrackInformation(-1, parentParticleID);
 183       track->SetUserInformation(trackInfo);
 184         // the track information is deleted together with its
 185         // G4Track object
 186     }
 187   }
 188
 189   // aliroot post track actions
 190   gAlice->PostTrack();
 191 }
 192
 193 void AliTrackingAction::FinishPrimaryTrack()
 194 {
 195 // Calls AliRun::PurifyKine and fills trees of hits
 196 // after finishing tracking of each primary track.
 197 // !! This method has to be also called from AlEventAction::EndOfEventAction()
 198 // for storing the last primary track of the current event.
 199 // ---
 200
 201   if (fPrimaryTrackID>0) {
 202
 203     // verbose
 204     if (fVerboseLevel == 3) {
 205       G4cout << "$$$ Primary track " << fPrimaryTrackID << G4endl;
 206     }
 207     else if ( fVerboseLevel == 2 &&  fPrimaryTrackID % 10 == 0 ) {
 208       G4cout << "$$$ Primary track " << fPrimaryTrackID  << G4endl;
 209     }
 210     else if ( fVerboseLevel == 1 &&  fPrimaryTrackID % 100 == 0 ) {
 211       G4cout << "$$$ Primary track " << fPrimaryTrackID  << G4endl;
 212     }
 213
 214     // aliroot finish primary track
 215     gAlice->FinishPrimary();
 216   }
 217   fPrimaryTrackID = 0;
 218 }
 219
 220 void AliTrackingAction::SaveTrack(const G4Track* track)
 221 {
 222 // Get all needed parameters from G4track and pass them
 223 // to AliRun::SetTrack() that creates corresponding TParticle
 224 // in the AliRun::fParticles array.
 225 // ----
 226
 227   // parent particle index
 228   G4int parentID = track->GetParentID();
 229   G4int motherIndex;
 230   if (parentID == 0) {
 231     motherIndex = -1;
 232   }
 233   else {
 234     motherIndex
 235       = GetTrackInformation(track,"SaveTrack")->GetParentParticleID();
 236   }
 237
 238   //G4cout << "SaveTrack: TrackID = " << track->GetTrackID()
 239   //       << "  Parent ID = " << track->GetParentID()
 240   //       << "  Index = " << gAlice->CurrentTrack()
 241   //       << "  Parent Index = " << motherIndex
 242   //       << G4endl;
 243
 244   // PDG code
 245   G4int pdg = track->GetDefinition()->GetPDGEncoding();
 246
 247   // track kinematics
 248   G4ThreeVector momentum = track->GetMomentum();
 249
 250   G4double px = momentum.x()/GeV;
 251   G4double py = momentum.y()/GeV;
 252   G4double pz = momentum.z()/GeV;
 253   G4double e = track->GetTotalEnergy()/GeV;
 254
 255   G4ThreeVector position = track->GetPosition();
 256   G4double vx = position.x()/cm;
 257   G4double vy = position.y()/cm;
 258   G4double vz = position.z()/cm;
 259   // time of production - check if ekvivalent with G4
 260   G4double t = track->GetGlobalTime();
 261
 262   G4ThreeVector polarization = track->GetPolarization();
 263   G4double polX = polarization.x();
 264   G4double polY = polarization.y();
 265   G4double polZ = polarization.z();
 266
 267   // production process
 268   AliMCProcess mcProcess;
 269   const G4VProcess* kpProcess = track->GetCreatorProcess();
 270   if (!kpProcess) {
 271     mcProcess = kPPrimary;
 272   }
 273   else {
 274     TG4PhysicsManager* pPhysicsManager = TG4PhysicsManager::Instance();
 275     mcProcess = pPhysicsManager->GetMCProcess(kpProcess);
 276     // distinguish kPDeltaRay from kPEnergyLoss
 277     if (mcProcess == kPEnergyLoss) mcProcess = kPDeltaRay;
 278   }
 279
 280   G4int ntr;
 281   // create particle
 282   gAlice->SetTrack(1, motherIndex, pdg, px, py, pz, e, vx, vy, vz, t,
 283                    polX, polY, polZ, mcProcess, ntr);
 284
 285 }
 286