+
+void AliMC::MakeTmpTrackRefsTree()
+{
+ // Make the temporary track reference tree
+ fTmpFileTR = new TFile("TrackRefsTmp.root", "recreate");
+ fTmpTreeTR = new TTree("TreeTR", "Track References");
+ TClonesArray* pRef = &fTmpTrackReferences;
+ fTmpTreeTR->Branch("TrackReferences", "TClonesArray", &pRef, 4000);
+}
+
+void AliMC::ReorderAndExpandTreeTR()
+{
+//
+// Reorder and expand the temporary track reference tree in order to match the kinematics tree
+//
+
+ AliRunLoader *rl = AliRunLoader::Instance();
+//
+// TreeTR
+ AliDebug(1, "fRunLoader->MakeTrackRefsContainer()");
+ rl->MakeTrackRefsContainer();
+ TTree * treeTR = rl->TreeTR();
+ if (treeTR){
+ // make branch for central track references
+ TBranch *branch;
+ TClonesArray* pRef = &fTrackReferences;
+ branch = treeTR->Branch("TrackReferences", &pRef);
+ branch->SetAddress(&pRef);
+ }
+
+ 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 (!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 (!(part->TestBit(kTransportBit)) || dau2 == -1 || dau2 < np) {
+// 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();
+ new(fTrackReferences[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();
+ new(fTrackReferences[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;
+ fTmpTrackReferences.Clear();
+ gSystem->Exec("rm -rf TrackRefsTmp.root");
+}
+