1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
17 ///////////////////////////////////////////////////////////////////////////
20 Origin: marian.ivanov@cern.ch
22 Generate complex MC information - used for Comparison later on
25 gSystem->Load("libPWG1.so")
26 AliGenInfoMaker *t = new AliGenInfoMaker("galice.root","genTracks.root",0,0)
31 #if !defined(__CINT__) || defined(__MAKECINT__)
42 #include "TStopwatch.h"
43 #include "TParticle.h"
46 #include "TGeometry.h"
47 #include "TPolyLine3D.h"
52 #include "AliSimDigits.h"
53 #include "AliTPCParam.h"
55 #include "AliTPCLoader.h"
56 #include "AliDetector.h"
57 #include "AliTrackReference.h"
58 #include "AliTPCParamSR.h"
59 #include "AliTracker.h"
62 #include "AliTrackPointArray.h"
65 #include "AliGenInfo.h"
66 #include "AliGenInfoMaker.h"
70 ClassImp(AliGenInfoMaker)
78 ////////////////////////////////////////////////////////////////////////
79 AliGenInfoMaker::AliGenInfoMaker():
80 fDebug(0), //! debug flag
81 fEventNr(0), //! current event number
82 fLabel(0), //! track label
83 fNEvents(0), //! number of events to process
84 fFirstEventNr(0), //! first event to process
85 fNParticles(0), //! number of particles in TreeK
86 fTreeGenTracks(0), //! output tree with generated tracks
87 fTreeKinks(0), //! output tree with Kinks
88 fTreeV0(0), //! output tree with V0
89 fTreeHitLines(0), //! tree with hit lines
90 fFileGenTracks(0), //! output file with stored fTreeGenTracks
91 fLoader(0), //! pointer to the run loader
92 fTreeD(0), //! current tree with digits
93 fTreeTR(0), //! current tree with TR
94 fStack(0), //! current stack
95 fGenInfo(0), //! array with pointers to gen info
96 fNInfos(0), //! number of tracks with infos
97 fParamTPC(0), //! AliTPCParam
105 ////////////////////////////////////////////////////////////////////////
106 AliGenInfoMaker::AliGenInfoMaker(const char * fnGalice, const char* fnRes,
107 Int_t nEvents, Int_t firstEvent):
108 fDebug(0), //! debug flag
109 fEventNr(0), //! current event number
110 fLabel(0), //! track label
111 fNEvents(0), //! number of events to process
112 fFirstEventNr(0), //! first event to process
113 fNParticles(0), //! number of particles in TreeK
114 fTreeGenTracks(0), //! output tree with generated tracks
115 fTreeKinks(0), //! output tree with Kinks
116 fTreeV0(0), //! output tree with V0
117 fTreeHitLines(0), //! tree with hit lines
118 fFileGenTracks(0), //! output file with stored fTreeGenTracks
119 fLoader(0), //! pointer to the run loader
120 fTreeD(0), //! current tree with digits
121 fTreeTR(0), //! current tree with TR
122 fStack(0), //! current stack
123 fGenInfo(0), //! array with pointers to gen info
124 fNInfos(0), //! number of tracks with infos
125 fParamTPC(0), //! AliTPCParam
135 fFirstEventNr = firstEvent;
136 fEventNr = firstEvent;
138 sprintf(fFnRes,"%s",fnRes);
140 fLoader = AliRunLoader::Open(fnGalice);
142 delete gAlice->GetRunLoader();
146 if (fLoader->LoadgAlice()){
147 cerr<<"Error occured while l"<<endl;
149 Int_t nall = fLoader->GetNumberOfEvents();
157 cerr<<"no events available"<<endl;
161 if (firstEvent+nEvents>nall) {
162 fEventNr = nall-firstEvent;
163 cerr<<"restricted number of events availaible"<<endl;
165 AliMagF * magf = gAlice->Field();
166 AliTracker::SetFieldMap(magf,0);
170 AliMCInfo * AliGenInfoMaker::MakeInfo(UInt_t i)
174 if (fGenInfo[i]) return fGenInfo[i];
175 fGenInfo[i] = new AliMCInfo;
183 ////////////////////////////////////////////////////////////////////////
184 AliGenInfoMaker::~AliGenInfoMaker()
188 fLoader->UnloadgAlice();
194 Int_t AliGenInfoMaker::SetIO()
198 CreateTreeGenTracks();
199 if (!fTreeGenTracks) return 1;
200 // AliTracker::SetFieldFactor();
202 fParamTPC = GetTPCParam();
207 ////////////////////////////////////////////////////////////////////////
208 Int_t AliGenInfoMaker::SetIO(Int_t eventNr)
213 fLoader->SetEventNumber(eventNr);
215 fLoader->LoadHeader();
216 fLoader->LoadKinematics();
217 fStack = fLoader->Stack();
219 fLoader->LoadTrackRefs();
221 fTreeTR = fLoader->TreeTR();
223 AliTPCLoader * tpcl = (AliTPCLoader*)fLoader->GetLoader("TPCLoader");
225 fTreeD = tpcl->TreeD();
229 Int_t AliGenInfoMaker::CloseIOEvent()
231 fLoader->UnloadHeader();
232 fLoader->UnloadKinematics();
233 fLoader->UnloadTrackRefs();
234 AliTPCLoader * tpcl = (AliTPCLoader*)fLoader->GetLoader("TPCLoader");
235 tpcl->UnloadDigits();
239 Int_t AliGenInfoMaker::CloseIO()
241 fLoader->UnloadgAlice();
247 ////////////////////////////////////////////////////////////////////////
248 Int_t AliGenInfoMaker::Exec(Int_t nEvents, Int_t firstEventNr)
251 fFirstEventNr = firstEventNr;
255 ////////////////////////////////////////////////////////////////////////
256 Int_t AliGenInfoMaker::Exec()
260 Int_t status =SetIO();
261 if (status>0) return status;
264 for (fEventNr = fFirstEventNr; fEventNr < fFirstEventNr+fNEvents;
267 fNParticles = fStack->GetNtrack();
269 fGenInfo = new AliMCInfo*[fNParticles];
270 for (UInt_t i = 0; i<fNParticles; i++) {
274 cout<<"Start to process event "<<fEventNr<<endl;
275 cout<<"\tfNParticles = "<<fNParticles<<endl;
276 if (fDebug>2) cout<<"\n\n\n\tStart loop over TreeTR"<<endl;
277 if (TreeTRLoop()>0) return 1;
279 if (fDebug>2) cout<<"\n\n\n\tStart loop over TreeD"<<endl;
280 if (TreeDLoop()>0) return 1;
282 if (fDebug>2) cout<<"\n\n\n\tStart loop over TreeK"<<endl;
283 if (TreeKLoop()>0) return 1;
284 if (fDebug>2) cout<<"\tEnd loop over TreeK"<<endl;
286 if (BuildKinkInfo()>0) return 1;
287 if (BuildV0Info()>0) return 1;
288 //if (BuildHitLines()>0) return 1;
289 if (fDebug>2) cout<<"\tEnd loop over TreeK"<<endl;
291 for (UInt_t i = 0; i<fNParticles; i++) {
292 if (fGenInfo[i]) delete fGenInfo[i];
301 cerr<<"Exec finished"<<endl;
308 ////////////////////////////////////////////////////////////////////////
309 void AliGenInfoMaker::CreateTreeGenTracks()
311 fFileGenTracks = TFile::Open(fFnRes,"RECREATE");
312 if (!fFileGenTracks) {
313 cerr<<"Error in CreateTreeGenTracks: cannot open file "<<fFnRes<<endl;
316 fTreeGenTracks = new TTree("genTracksTree","genTracksTree");
317 AliMCInfo * info = new AliMCInfo;
318 fTreeGenTracks->Branch("MC","AliMCInfo",&info);
321 AliGenKinkInfo *kinkinfo = new AliGenKinkInfo;
322 fTreeKinks = new TTree("genKinksTree","genKinksTree");
323 fTreeKinks->Branch("MC","AliGenKinkInfo",&kinkinfo);
326 AliGenV0Info *v0info = new AliGenV0Info;
327 fTreeV0 = new TTree("genV0Tree","genV0Tree");
328 fTreeV0->Branch("MC","AliGenV0Info",&v0info);
332 AliTrackPointArray * points0 = new AliTrackPointArray;
333 AliTrackPointArray * points1 = new AliTrackPointArray;
334 AliTrackPointArray * points2 = new AliTrackPointArray;
335 fTreeHitLines = new TTree("HitLines","HitLines");
336 fTreeHitLines->Branch("TPC.","AliTrackPointArray",&points0,32000,10);
337 fTreeHitLines->Branch("TRD.","AliTrackPointArray",&points1,32000,10);
338 fTreeHitLines->Branch("ITS.","AliTrackPointArray",&points2,32000,10);
340 fTreeHitLines->Branch("Label",&label,"label/I");
342 fTreeGenTracks->AutoSave();
343 fTreeKinks->AutoSave();
345 fTreeHitLines->AutoSave();
347 ////////////////////////////////////////////////////////////////////////
348 void AliGenInfoMaker::CloseOutputFile()
350 if (!fFileGenTracks) {
351 cerr<<"File "<<fFnRes<<" not found as an open file."<<endl;
354 fFileGenTracks->cd();
355 fTreeGenTracks->Write();
356 delete fTreeGenTracks;
362 fFileGenTracks->Close();
363 delete fFileGenTracks;
367 ////////////////////////////////////////////////////////////////////////
368 Int_t AliGenInfoMaker::TreeKLoop()
371 // open the file with treeK
372 // loop over all entries there and save information about some tracks
375 AliStack * stack = fStack;
376 if (!stack) {cerr<<"Stack was not found!\n"; return 1;}
379 cout<<"There are "<<fNParticles<<" primary and secondary particles in event "
383 TParticlePDG *ppdg = 0;
384 // not all generators give primary vertex position. Take the vertex
385 // of the particle 0 as primary vertex.
386 TDatabasePDG pdg; //get pdg table
387 //thank you very much root for this
388 TBranch * br = fTreeGenTracks->GetBranch("MC");
389 TParticle *particle = stack->ParticleFromTreeK(0);
390 fVPrim[0] = particle->Vx();
391 fVPrim[1] = particle->Vy();
392 fVPrim[2] = particle->Vz();
393 for (UInt_t iParticle = 0; iParticle < fNParticles; iParticle++) {
394 // load only particles with TR
395 AliMCInfo * info = GetInfo(iParticle);
397 //////////////////////////////////////////////////////////////////////
398 info->fLabel = iParticle;
400 info->fParticle = *(stack->Particle(iParticle));
401 info->fVDist[0] = info->fParticle.Vx()-fVPrim[0];
402 info->fVDist[1] = info->fParticle.Vy()-fVPrim[1];
403 info->fVDist[2] = info->fParticle.Vz()-fVPrim[2];
404 info->fVDist[3] = TMath::Sqrt(info->fVDist[0]*info->fVDist[0]+
405 info->fVDist[1]*info->fVDist[1]+info->fVDist[2]*info->fVDist[2]);
408 ipdg = info->fParticle.GetPdgCode();
410 ppdg = pdg.GetParticle(ipdg);
411 info->fEventNr = fEventNr;
413 //////////////////////////////////////////////////////////////////////
414 br->SetAddress(&info);
415 fTreeGenTracks->Fill();
417 fTreeGenTracks->AutoSave();
418 if (fDebug > 2) cerr<<"end of TreeKLoop"<<endl;
423 ////////////////////////////////////////////////////////////////////////////////////
424 Int_t AliGenInfoMaker::BuildKinkInfo()
427 // Build tree with Kink Information
429 AliStack * stack = fStack;
430 if (!stack) {cerr<<"Stack was not found!\n"; return 1;}
433 cout<<"There are "<<fNParticles<<" primary and secondary particles in event "
437 //TParticlePDG *ppdg = 0;
438 // not all generators give primary vertex position. Take the vertex
439 // of the particle 0 as primary vertex.
440 TDatabasePDG pdg; //get pdg table
441 //thank you very much root for this
442 TBranch * br = fTreeKinks->GetBranch("MC");
444 AliGenKinkInfo * kinkinfo = new AliGenKinkInfo;
446 for (UInt_t iParticle = 0; iParticle < fNParticles; iParticle++) {
447 // load only particles with TR
448 AliMCInfo * info = GetInfo(iParticle);
450 if (info->fCharge==0) continue;
451 if (info->fTRdecay.P()<0.13) continue; //momenta cut
452 if (info->fTRdecay.R()>500) continue; //R cut - decay outside barrel
453 TParticle & particle = info->fParticle;
454 Int_t first = particle.GetDaughter(0);
455 if (first ==0) continue;
457 Int_t last = particle.GetDaughter(1);
458 if (last ==0) last=first;
459 AliMCInfo * info2 = 0;
460 AliMCInfo * dinfo = 0;
462 for (Int_t id2=first;id2<=last;id2++){
463 info2 = GetInfo(id2);
464 if (!info2) continue;
465 TParticle &p2 = info2->fParticle;
466 Double_t r = TMath::Sqrt(p2.Vx()*p2.Vx()+p2.Vy()*p2.Vy());
467 if ( TMath::Abs(info->fTRdecay.R()-r)>0.01) continue;
468 if (!(p2.GetPDG())) continue;
471 kinkinfo->SetInfoMother(*info);
472 kinkinfo->SetInfoDaughter(*dinfo);
473 if (kinkinfo->GetMinus().GetParticle().GetPDG()==0 || kinkinfo->GetPlus().GetParticle().GetPDG()==0) continue;
475 br->SetAddress(&kinkinfo);
480 kinkinfo->fMCm = (*info);
481 kinkinfo->GetPlus() = (*dinfo);
483 br->SetAddress(&kinkinfo);
488 fTreeGenTracks->AutoSave();
489 if (fDebug > 2) cerr<<"end of Kink Loop"<<endl;
494 ////////////////////////////////////////////////////////////////////////////////////
495 Int_t AliGenInfoMaker::BuildV0Info()
498 // Build tree with V0 Information
500 AliStack * stack = fStack;
501 if (!stack) {cerr<<"Stack was not found!\n"; return 1;}
504 cout<<"There are "<<fNParticles<<" primary and secondary particles in event "
508 //TParticlePDG *ppdg = 0;
509 // not all generators give primary vertex position. Take the vertex
510 // of the particle 0 as primary vertex.
511 TDatabasePDG pdg; //get pdg table
512 //thank you very much root for this
513 TBranch * br = fTreeV0->GetBranch("MC");
515 AliGenV0Info * v0info = new AliGenV0Info;
518 for (UInt_t iParticle = 0; iParticle < fNParticles; iParticle++) {
519 // load only particles with TR
520 AliMCInfo * info = GetInfo(iParticle);
522 if (info->fCharge==0) continue;
525 TParticle & particle = info->fParticle;
526 Int_t mother = particle.GetMother(0);
527 if (mother <=0) continue;
529 TParticle * motherparticle = stack->Particle(mother);
530 if (!motherparticle) continue;
532 Int_t last = motherparticle->GetDaughter(1);
533 if (last==(int)iParticle) continue;
534 AliMCInfo * info2 = info;
535 AliMCInfo * dinfo = GetInfo(last);
536 if (!dinfo) continue;
537 if (!dinfo->fParticle.GetPDG()) continue;
538 if (!info2->fParticle.GetPDG()) continue;
540 v0info->SetInfoP(*info);
541 v0info->SetInfoM(*dinfo);
542 v0info->SetMother(*motherparticle);
543 v0info->Update(fVPrim);
544 br->SetAddress(&v0info);
549 if (fDebug > 2) cerr<<"end of V0 Loop"<<endl;
556 ////////////////////////////////////////////////////////////////////////
557 Int_t AliGenInfoMaker::BuildHitLines()
561 // open the file with treeK
562 // loop over all entries there and save information about some tracks
565 AliStack * stack = fStack;
566 if (!stack) {cerr<<"Stack was not found!\n"; return 1;}
569 cout<<"There are "<<fNParticles<<" primary and secondary particles in event "
573 // // TParticlePDG *ppdg = 0;
574 // // not all generators give primary vertex position. Take the vertex
575 // // of the particle 0 as primary vertex.
576 // TDatabasePDG pdg; //get pdg table
577 // //thank you very much root for this
578 // AliTrackPointArray *tpcp = new AliTrackPointArray;
579 // AliTrackPointArray *trdp = new AliTrackPointArray;
580 // AliTrackPointArray *itsp = new AliTrackPointArray;
583 // TBranch * brtpc = fTreeHitLines->GetBranch("TPC.");
584 // TBranch * brtrd = fTreeHitLines->GetBranch("TRD.");
585 // TBranch * brits = fTreeHitLines->GetBranch("ITS.");
586 // TBranch * brlabel = fTreeHitLines->GetBranch("Label");
587 // brlabel->SetAddress(&label);
588 // brtpc->SetAddress(&tpcp);
589 // brtrd->SetAddress(&trdp);
590 // brits->SetAddress(&itsp);
592 // AliDetector *dtpc = gAlice->GetDetector("TPC");
593 // AliDetector *dtrd = gAlice->GetDetector("TRD");
594 // AliDetector *dits = gAlice->GetDetector("ITS");
596 // for (UInt_t iParticle = 0; iParticle < fNParticles; iParticle++) {
597 // // load only particles with TR
598 // AliMCInfo * info = GetInfo(iParticle);
599 // if (!info) continue;
600 // Int_t primpart = info->fPrimPart;
601 // ipdg = info->fParticle.GetPdgCode();
602 // label = iParticle;
604 // gAlice->ResetHits();
608 // tpcp->fLabel1 = ipdg;
609 // trdp->fLabel1 = ipdg;
610 // itsp->fLabel1 = ipdg;
611 // if (dtpc->TreeH()->GetEvent(primpart)){
612 // dtpc->LoadPoints(primpart);
613 // tpcp->Reset(dtpc,iParticle);
615 // if (dtrd->TreeH()->GetEvent(primpart)){
616 // dtrd->LoadPoints(primpart);
617 // trdp->Reset(dtrd,iParticle);
619 // if (dits->TreeH()->GetEvent(primpart)){
620 // dits->LoadPoints(primpart);
621 // itsp->Reset(dits,iParticle);
624 // fTreeHitLines->Fill();
625 // dtpc->ResetPoints();
626 // dtrd->ResetPoints();
627 // dits->ResetPoints();
632 // fTreeHitLines->AutoSave();
633 // if (fDebug > 2) cerr<<"end of TreeKLoop"<<endl;
638 ////////////////////////////////////////////////////////////////////////
639 Int_t AliGenInfoMaker::TreeDLoop()
642 // open the file with treeD
643 // loop over all entries there and save information about some tracks
646 Int_t nInnerSector = fParamTPC->GetNInnerSector();
648 Int_t zero=fParamTPC->GetZeroSup()+6;
651 AliSimDigits digitsAddress, *digits=&digitsAddress;
652 fTreeD->GetBranch("Segment")->SetAddress(&digits);
654 Int_t sectorsByRows=(Int_t)fTreeD->GetEntries();
655 if (fDebug > 1) cout<<"\tsectorsByRows = "<<sectorsByRows<<endl;
656 for (Int_t i=0; i<sectorsByRows; i++) {
657 if (!fTreeD->GetEvent(i)) continue;
659 fParamTPC->AdjustSectorRow(digits->GetID(),sec,row);
660 if (fDebug > 1) cout<<sec<<' '<<row<<" \r";
661 // here I expect that upper sectors follow lower sectors
662 if (sec > nInnerSector) rowShift = fParamTPC->GetNRowLow();
664 digits->ExpandTrackBuffer();
667 Int_t iRow=digits->CurrentRow();
668 Int_t iColumn=digits->CurrentColumn();
669 Short_t digitValue = digits->CurrentDigit();
670 if (digitValue >= zero) {
672 for (Int_t j = 0; j<3; j++) {
673 // label = digits->GetTrackID(iRow,iColumn,j);
674 label = digits->GetTrackIDFast(iRow,iColumn,j)-2;
675 if (label >= (Int_t)fNParticles) { //don't label from bakground event
678 if (label >= 0 && label <= (Int_t)fNParticles) {
680 cout<<"Inner loop: sector, iRow, iColumn, label, value, row "
682 <<iRow<<" "<<iColumn<<" "<<label<<" "<<digitValue
685 AliMCInfo * info = GetInfo(label);
687 info->fTPCRow.SetRow(row+rowShift);
692 } while (digits->Next());
695 if (fDebug > 2) cerr<<"end of TreeDLoop"<<endl;
700 ////////////////////////////////////////////////////////////////////////
701 Int_t AliGenInfoMaker::TreeTRLoop()
704 // loop over TrackReferences and store the first one for each track
706 TTree * treeTR = fTreeTR;
707 Int_t nPrimaries = (Int_t) treeTR->GetEntries();
708 if (fDebug > 1) cout<<"There are "<<nPrimaries<<" entries in TreeTR"<<endl;
711 //track references for TPC
712 TClonesArray* tpcArrayTR = new TClonesArray("AliTrackReference");
713 TClonesArray* itsArrayTR = new TClonesArray("AliTrackReference");
714 TClonesArray* trdArrayTR = new TClonesArray("AliTrackReference");
715 TClonesArray* tofArrayTR = new TClonesArray("AliTrackReference");
716 TClonesArray* runArrayTR = new TClonesArray("AliTrackReference");
718 if (treeTR->GetBranch("TPC")) treeTR->GetBranch("TPC")->SetAddress(&tpcArrayTR);
719 if (treeTR->GetBranch("ITS")) treeTR->GetBranch("ITS")->SetAddress(&itsArrayTR);
720 if (treeTR->GetBranch("TRD")) treeTR->GetBranch("TRD")->SetAddress(&trdArrayTR);
721 if (treeTR->GetBranch("TOF")) treeTR->GetBranch("TOF")->SetAddress(&tofArrayTR);
722 if (treeTR->GetBranch("AliRun")) treeTR->GetBranch("AliRun")->SetAddress(&runArrayTR);
726 for (Int_t iPrimPart = 0; iPrimPart<nPrimaries; iPrimPart++) {
727 treeTR->GetEntry(iPrimPart);
729 // Loop over TPC references
731 for (Int_t iTrackRef = 0; iTrackRef < tpcArrayTR->GetEntriesFast(); iTrackRef++) {
732 AliTrackReference *trackRef = (AliTrackReference*)tpcArrayTR->At(iTrackRef);
734 if (trackRef->TestBit(BIT(2))){
736 if (trackRef->P()<fTPCPtCut) continue;
737 Int_t label = trackRef->GetTrack();
738 AliMCInfo * info = GetInfo(label);
739 if (!info) info = MakeInfo(label);
740 info->fTRdecay = *trackRef;
743 if (trackRef->P()<fTPCPtCut) continue;
744 Int_t label = trackRef->GetTrack();
745 AliMCInfo * info = GetInfo(label);
746 if (!info) info = MakeInfo(label);
748 info->fPrimPart = iPrimPart;
749 TClonesArray & arr = *(info->fTPCReferences);
750 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
753 // Loop over ITS references
755 for (Int_t iTrackRef = 0; iTrackRef < itsArrayTR->GetEntriesFast(); iTrackRef++) {
756 AliTrackReference *trackRef = (AliTrackReference*)itsArrayTR->At(iTrackRef);
759 if (trackRef->P()<fTPCPtCut) continue;
760 Int_t label = trackRef->GetTrack();
761 AliMCInfo * info = GetInfo(label);
762 if ( (!info) && trackRef->Pt()<fITSPtCut) continue;
763 if (!info) info = MakeInfo(label);
765 info->fPrimPart = iPrimPart;
766 TClonesArray & arr = *(info->fITSReferences);
767 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
770 // Loop over TRD references
772 for (Int_t iTrackRef = 0; iTrackRef < trdArrayTR->GetEntriesFast(); iTrackRef++) {
773 AliTrackReference *trackRef = (AliTrackReference*)trdArrayTR->At(iTrackRef);
775 if (trackRef->P()<fTPCPtCut) continue;
776 Int_t label = trackRef->GetTrack();
777 AliMCInfo * info = GetInfo(label);
778 if ( (!info) && trackRef->Pt()<fTRDPtCut) continue;
779 if (!info) info = MakeInfo(label);
781 info->fPrimPart = iPrimPart;
782 TClonesArray & arr = *(info->fTRDReferences);
783 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
786 // Loop over TOF references
788 for (Int_t iTrackRef = 0; iTrackRef < tofArrayTR->GetEntriesFast(); iTrackRef++) {
789 AliTrackReference *trackRef = (AliTrackReference*)tofArrayTR->At(iTrackRef);
790 Int_t label = trackRef->GetTrack();
791 AliMCInfo * info = GetInfo(label);
793 if (trackRef->Pt()<fTPCPtCut) continue;
794 if ( (!info) && trackRef->Pt()<fTOFPtCut) continue;
796 if (!info) info = MakeInfo(label);
798 info->fPrimPart = iPrimPart;
799 TClonesArray & arr = *(info->fTOFReferences);
800 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
804 // get dacay position
806 for (Int_t iTrackRef = 0; iTrackRef < runArrayTR->GetEntriesFast(); iTrackRef++) {
807 AliTrackReference *trackRef = (AliTrackReference*)runArrayTR->At(iTrackRef);
809 Int_t label = trackRef->GetTrack();
810 AliMCInfo * info = GetInfo(label);
812 if (!trackRef->TestBit(BIT(2))) continue; //if not decay
813 // if (TMath::Abs(trackRef.X());
814 info->fTRdecay = *trackRef;
818 tpcArrayTR->Delete();
820 trdArrayTR->Delete();
822 tofArrayTR->Delete();
824 itsArrayTR->Delete();
826 runArrayTR->Delete();
829 return TreeTRLoopNew();
832 ////////////////////////////////////////////////////////////////////////
833 Int_t AliGenInfoMaker::TreeTRLoopNew()
836 // loop over TrackReferences and store the first one for each track
838 TTree * treeTR = fTreeTR;
839 Int_t nPrimaries = (Int_t) treeTR->GetEntries();
840 if (fDebug > 1) cout<<"There are "<<nPrimaries<<" entries in TreeTR"<<endl;
844 TClonesArray* runArrayTR = new TClonesArray("AliTrackReference");
846 if (treeTR->GetBranch("TrackReferences")) treeTR->GetBranch("TrackReferences")->SetAddress(&runArrayTR);
850 for (Int_t iPrimPart = 0; iPrimPart<nPrimaries; iPrimPart++) {
851 treeTR->GetEntry(iPrimPart);
853 // gettrack references
855 for (Int_t iTrackRef = 0; iTrackRef < runArrayTR->GetEntriesFast(); iTrackRef++) {
856 AliTrackReference *trackRef = (AliTrackReference*)runArrayTR->At(iTrackRef);
858 Int_t label = trackRef->GetTrack();
862 if (trackRef->DetectorId()== AliTrackReference::kTPC){
864 if (trackRef->P()<fTPCPtCut) continue;
865 Int_t label = trackRef->GetTrack();
866 AliMCInfo * info = GetInfo(label);
867 if (!info) info = MakeInfo(label);
869 info->fPrimPart = iPrimPart;
870 TClonesArray & arr = *(info->fTPCReferences);
871 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
876 if (trackRef->DetectorId()== AliTrackReference::kITS){
878 if (trackRef->P()<fITSPtCut) continue;
879 Int_t label = trackRef->GetTrack();
880 AliMCInfo * info = GetInfo(label);
881 if (!info) info = MakeInfo(label);
883 info->fPrimPart = iPrimPart;
884 TClonesArray & arr = *(info->fITSReferences);
885 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
890 if (trackRef->DetectorId()== AliTrackReference::kTRD){
892 if (trackRef->P()<fTRDPtCut) continue;
893 Int_t label = trackRef->GetTrack();
894 AliMCInfo * info = GetInfo(label);
895 if (!info) info = MakeInfo(label);
897 info->fPrimPart = iPrimPart;
898 TClonesArray & arr = *(info->fTRDReferences);
899 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
904 if (trackRef->DetectorId()== AliTrackReference::kTOF){
906 if (trackRef->P()<fTOFPtCut) continue;
907 Int_t label = trackRef->GetTrack();
908 AliMCInfo * info = GetInfo(label);
909 if (!info) info = MakeInfo(label);
911 info->fPrimPart = iPrimPart;
912 TClonesArray & arr = *(info->fTOFReferences);
913 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
918 if (trackRef->DetectorId()== AliTrackReference::kDisappeared){
920 if (!trackRef->TestBit(BIT(2))) continue; //if not decay
921 AliMCInfo * info = GetInfo(label);
923 info->fTRdecay = *trackRef;
933 ////////////////////////////////////////////////////////////////////////
934 Float_t AliGenInfoMaker::TR2LocalX(AliTrackReference *trackRef,
935 AliTPCParam *paramTPC) const {
937 Float_t x[3] = { trackRef->X(),trackRef->Y(),trackRef->Z()};
939 paramTPC->Transform0to1(x,index);
940 paramTPC->Transform1to2(x,index);
943 ////////////////////////////////////////////////////////////////////////
947 AliTPCParam * AliGenInfoMaker::GetTPCParam(){
948 AliTPCParamSR * par = new AliTPCParamSR;