#include <TVirtualMC.h>
#include <TParticle.h>
#include <TROOT.h>
+#include <TFile.h>
#include "AliLog.h"
#include "AliDetector.h"
fTransParName("\0"),
fMCQA(0),
fHitLists(0),
- fTrackReferences(0)
+ fTmpTreeTR(0),
+ fTmpFileTR(0),
+ fTrackReferences(0),
+ fTmpTrackReferences(0)
{
//default constructor
fTransParName("\0"),
fMCQA(0),
fHitLists(new TList()),
- fTrackReferences(new TClonesArray("AliTrackReference", 100))
+ fTmpTreeTR(0),
+ fTmpFileTR(0),
+ fTrackReferences(new TClonesArray("AliTrackReference", 100)),
+ fTmpTrackReferences(new TClonesArray("AliTrackReference", 100))
{
//constructor
// Set transport parameters
fTransParName("\0"),
fMCQA(0),
fHitLists(0),
- fTrackReferences(0)
+ fTmpTreeTR(0),
+ fTmpFileTR(0),
+ fTrackReferences(0),
+ fTmpTrackReferences(0)
{
//
// Copy constructor for AliMC
fTrackReferences = 0;
}
+if (fTmpTrackReferences) {
+ fTmpTrackReferences->Delete();
+ delete fTmpTrackReferences;
+ fTmpTrackReferences = 0;
+ }
+
}
//_______________________________________________________________________
// Reset Hits info
ResetHits();
ResetTrackReferences();
-
}
//_______________________________________________________________________
AddEnergyDeposit(gMC->CurrentVolID(copy),gMC->Edep());
//
// write tracke reference for track which is dissapearing - MI
- if (gMC->IsTrackDisappeared()) {
- if (gMC->Etot()>0.05) AddTrackReference(GetCurrentTrackNumber(),
+
+ if (gMC->IsTrackDisappeared() && !(gMC->IsTrackAlive())) {
+ if (gMC->Etot() > 0.05) AddTrackReference(GetCurrentTrackNumber(),
AliTrackReference::kDisappeared);
+
+
}
-
+
//Call the appropriate stepping routine;
AliModule *det = dynamic_cast<AliModule*>(gAlice->Modules()->At(id));
if(det && det->StepManagerIsEnabled()) {
AliDebug(1, ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
AliDebug(1, ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>");
- AliRunLoader *runloader=gAlice->GetRunLoader();
+ AliRunLoader *runloader=gAlice->GetRunLoader();
/*******************************/
/* Clean after eventual */
// Clean detector information
if (runloader->Stack())
- runloader->Stack()->Reset();//clean stack -> tree is unloaded
+ runloader->Stack()->Reset();//clean stack -> tree is unloaded
else
- runloader->MakeStack();//or make a new one
+ runloader->MakeStack();//or make a new one
-// runloader->Stack()->BeginEvent();
-if(gAlice->Lego() == 0x0)
- {
- AliDebug(1, "fRunLoader->MakeTree(K)");
- runloader->MakeTree("K");
- }
-
+ if(gAlice->Lego() == 0x0)
+ {
+ AliDebug(1, "fRunLoader->MakeTree(K)");
+ runloader->MakeTree("K");
+ }
+
AliDebug(1, "gMC->SetStack(fRunLoader->Stack())");
gMC->SetStack(gAlice->GetRunLoader()->Stack());//Was in InitMC - but was moved here
//because we don't have guarantee that
//
runloader->GetHeader()->Reset(gAlice->GetRunNumber(),gAlice->GetEvNumber(),
- gAlice->GetEventNrInRun());
+ gAlice->GetEventNrInRun());
// fRunLoader->WriteKinematics("OVERWRITE"); is there any reason to rewrite here since MakeTree does so
if(gAlice->Lego())
- {
- gAlice->Lego()->BeginEvent();
- return;
- }
+ {
+ gAlice->Lego()->BeginEvent();
+ return;
+ }
AliDebug(1, "ResetHits()");
AliDebug(1, "fRunLoader->MakeTree(H)");
runloader->MakeTree("H");
- AliDebug(1, "fRunLoader->MakeTrackRefsContainer()");
- runloader->MakeTrackRefsContainer();//for insurance
+ MakeTmpTrackRefsTree();
//create new branches and SetAdresses
TIter next(gAlice->Modules());
AliModule *detector;
while((detector = (AliModule*)next()))
{
- AliDebug(2, Form("%s->MakeBranch(H)",detector->GetName()));
- detector->MakeBranch("H");
+ AliDebug(2, Form("%s->MakeBranch(H)",detector->GetName()));
+ detector->MakeBranch("H");
}
- // make branch for AliRun track References
- TTree * treeTR = runloader->TreeTR();
- if (treeTR){
- // make branch for central track references
- if (!fTrackReferences) fTrackReferences = new TClonesArray("AliTrackReference",0);
- TBranch *branch;
- branch = treeTR->Branch("TrackReferences",&fTrackReferences);
- branch->SetAddress(&fTrackReferences);
- }
- //
}
//_______________________________________________________________________
// This primary is finished, purify stack
#if ROOT_VERSION_CODE > 262152
if (!(gMC->SecondariesAreOrdered())) {
- runloader->Stack()->ReorderKine();
- RemapHits();
+ if (runloader->Stack()->ReorderKine()) RemapHits();
}
#endif
- runloader->Stack()->PurifyKine();
- RemapHits();
+ if (runloader->Stack()->PurifyKine()) RemapHits();
TIter next(gAlice->Modules());
AliModule *detector;
}
// Write out track references if any
- if (runloader->TreeTR()) runloader->TreeTR()->Fill();
+ if (fTmpTreeTR) fTmpTreeTR->Fill();
}
void AliMC::RemapHits()
// This for detectors which have a special mapping mechanism
// for hits, such as TPC and TRD
//
+
TObjArray* modules = gAlice->Modules();
TIter nextmod(modules);
- AliDetector *detector;
- while((detector = dynamic_cast<AliDetector*>(nextmod()))) {
- detector->RemapTrackHitIDs(stack->TrackLabelMap());
+ AliModule *module;
+ while((module = (AliModule*) nextmod())) {
+ AliDetector* det = dynamic_cast<AliDetector*> (module);
+ if (det) det->RemapTrackHitIDs(stack->TrackLabelMap());
}
//
RemapTrackReferencesIDs(stack->TrackLabelMap());
header->SetNprimary(stack->GetNprimary());
header->SetNtrack(stack->GetNtrack());
-
// Write out the kinematics
if (!gAlice->Lego()) stack->FinishEvent();
+
+ // Synchronize the TreeTR with TreeK
+ if (fTmpTreeTR) ReorderAndExpandTreeTR();
// Write out the event Header information
TTree* treeE = runloader->TreeE();
Int_t primary = GetPrimary(label);
Particle(primary)->SetBit(kKeepBit);
- Int_t nref = fTrackReferences->GetEntriesFast();
- TClonesArray &lref = *fTrackReferences;
+ Int_t nref = fTmpTrackReferences->GetEntriesFast();
+ TClonesArray &lref = *fTmpTrackReferences;
return new(lref[nref]) AliTrackReference(label, id);
}
//
// Reset all references
//
- if (fTrackReferences) fTrackReferences->Clear();
+ if (fTmpTrackReferences) fTmpTrackReferences->Clear();
}
void AliMC::RemapTrackReferencesIDs(Int_t *map)
// Remapping track reference
// Called at finish primary
//
- if (!fTrackReferences) return;
- Int_t nEntries = fTrackReferences->GetEntries();
+ if (!fTmpTrackReferences) return;
+ Int_t nEntries = fTmpTrackReferences->GetEntries();
for (Int_t i=0; i < nEntries; i++){
- AliTrackReference * ref = dynamic_cast<AliTrackReference*>(fTrackReferences->UncheckedAt(i));
+ AliTrackReference * ref = dynamic_cast<AliTrackReference*>(fTmpTrackReferences->UncheckedAt(i));
if (ref) {
Int_t newID = map[ref->GetTrack()];
if (newID>=0) ref->SetTrack(newID);
else {
ref->SetBit(kNotDeleted,kFALSE);
- fTrackReferences->RemoveAt(i);
+ fTmpTrackReferences->RemoveAt(i);
}
} // if ref
}
- fTrackReferences->Compress();
+ fTmpTrackReferences->Compress();
}
void AliMC::FixParticleDecaytime()
//
gMC->ForceDecayTime(t / 2.99792458e10);
}
+
+void AliMC::MakeTmpTrackRefsTree()
+{
+ // Make the temporary track reference tree
+ fTmpFileTR = new TFile("TrackRefsTmp.root", "recreate");
+ fTmpTreeTR = new TTree("TreeTR", "Track References");
+ if (!fTmpTrackReferences) fTmpTrackReferences = new TClonesArray("AliTrackReference", 100);
+ fTmpTreeTR->Branch("TrackReferences", "TClonesArray", &fTmpTrackReferences, 4000);
+}
+
+void AliMC::ReorderAndExpandTreeTR()
+{
+//
+// Reorder and expand the temporary track reference tree in order to match the kinematics tree
+//
+
+ AliRunLoader *rl = gAlice->GetRunLoader();
+//
+// TreeTR
+ AliDebug(1, "fRunLoader->MakeTrackRefsContainer()");
+ rl->MakeTrackRefsContainer();
+ TTree * treeTR = rl->TreeTR();
+ if (treeTR){
+ // make branch for central track references
+ if (!fTrackReferences) fTrackReferences = new TClonesArray("AliTrackReference",0);
+ TBranch *branch;
+ branch = treeTR->Branch("TrackReferences",&fTrackReferences);
+ branch->SetAddress(&fTrackReferences);
+ }
+
+ AliStack* stack = rl->Stack();
+ Int_t np = stack->GetNprimary();
+ Int_t nt = fTmpTreeTR->GetEntries();
+ //
+ // Loop over tracks and find the secondaries with the help of the kine tree
+ Int_t ifills = 0;
+ Int_t it = 0;
+ for (Int_t ip = np - 1; ip > -1; ip--) {
+ TParticle *part = stack->Particle(ip);
+// printf("Particle %5d %5d %5d %5d %5d \n", ip, part->GetPdgCode(), part->GetFirstMother(), part->GetFirstDaughter(), part->GetLastDaughter());
+
+ // Skip primaries that have not been transported
+ Int_t dau1 = part->GetFirstDaughter();
+ Int_t dau2 = -1;
+ // if ((dau1 > -1 && dau1 < np) || part->GetStatusCode() > 1) continue;
+ if (!part->TestBit(kTransportBit)) continue;
+ //
+ fTmpTreeTR->GetEntry(it++);
+ Int_t nh = fTmpTrackReferences->GetEntries();
+ // Determine range of secondaries produced by this primary
+ if (dau1 > -1) {
+ Int_t inext = ip - 1;
+ while (dau2 < 0) {
+ if (inext >= 0) {
+ part = stack->Particle(inext);
+ dau2 = part->GetFirstDaughter();
+ if (dau2 == -1 || dau2 < np) {
+ dau2 = -1;
+ } else {
+ dau2--;
+ }
+ } else {
+ dau2 = stack->GetNtrack() - 1;
+ }
+ inext--;
+ } // find upper bound
+ } // dau2 < 0
+// printf("Check (1) %5d %5d %5d %5d %5d \n", ip, np, it, dau1, dau2);
+
+ //
+ // Loop over reference hits and find secondary label
+ for (Int_t id = dau1; (id <= dau2) && (dau1 > -1); id++) {
+ for (Int_t ih = 0; ih < nh; ih++) {
+ AliTrackReference* tr = (AliTrackReference*) fTmpTrackReferences->At(ih);
+ Int_t label = tr->Label();
+ // Skip primaries
+ if (label == ip) continue;
+ if (label > dau2 || label < dau1)
+ AliWarning(Form("Track Reference Label out of range !: %5d %5d %5d \n", label, dau1, dau2));
+ if (label == id) {
+ // secondary found
+ Int_t nref = fTrackReferences->GetEntriesFast();
+ TClonesArray &lref = *fTrackReferences;
+ new(lref[nref]) AliTrackReference(*tr);
+ }
+ } // hits
+ treeTR->Fill();
+ fTrackReferences->Clear();
+ ifills++;
+ } // daughters
+ } // tracks
+ //
+ // Now loop again and write the primaries
+ it = nt - 1;
+ for (Int_t ip = 0; ip < np; ip++) {
+ TParticle* part = stack->Particle(ip);
+// if ((part->GetFirstDaughter() == -1 && part->GetStatusCode() <= 1) || part->GetFirstDaughter() >= np)
+ if (part->TestBit(kTransportBit))
+ {
+ // Skip particles that have not been transported
+ fTmpTreeTR->GetEntry(it--);
+ Int_t nh = fTmpTrackReferences->GetEntries();
+ //
+ // Loop over reference hits and find primary labels
+ for (Int_t ih = 0; ih < nh; ih++) {
+ AliTrackReference* tr = (AliTrackReference*) fTmpTrackReferences->At(ih);
+ Int_t label = tr->Label();
+ if (label == ip) {
+ Int_t nref = fTrackReferences->GetEntriesFast();
+ TClonesArray &lref = *fTrackReferences;
+ new(lref[nref]) AliTrackReference(*tr);
+ }
+ }
+ }
+ treeTR->Fill();
+ fTrackReferences->Clear();
+ ifills++;
+ } // tracks
+ // Check
+ if (ifills != stack->GetNtrack())
+ AliWarning(Form("Number of entries in TreeTR (%5d) unequal to TreeK (%5d) \n", ifills, stack->GetNtrack()));
+//
+// Clean-up
+ delete fTmpTreeTR;
+ fTmpFileTR->Close();
+ delete fTmpFileTR;
+ delete fTmpTrackReferences;
+ fTmpTrackReferences = 0;
+ gSystem->Exec("rm -rf TrackRefsTmp.root");
+}
+