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__)
41 //#include "TString.h"
42 #include "TStopwatch.h"
43 #include "TParticle.h"
44 //#include "TSystem.h"
45 //#include "TCanvas.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"
60 //#include "AliMagF.h"
61 //#include "AliHelix.h"
62 #include "AliTrackPointArray.h"
65 #include "AliMCInfo.h"
66 #include "AliGenV0Info.h"
67 #include "AliGenKinkInfo.h"
68 #include "AliGenInfoMaker.h"
72 ClassImp(AliGenInfoMaker)
80 ////////////////////////////////////////////////////////////////////////
81 AliGenInfoMaker::AliGenInfoMaker():
82 fDebug(0), //! debug flag
83 fEventNr(0), //! current event number
84 fLabel(0), //! track label
85 fNEvents(0), //! number of events to process
86 fFirstEventNr(0), //! first event to process
87 fNParticles(0), //! number of particles in TreeK
88 fTreeGenTracks(0), //! output tree with generated tracks
89 fTreeKinks(0), //! output tree with Kinks
90 fTreeV0(0), //! output tree with V0
91 fTreeHitLines(0), //! tree with hit lines
92 fFileGenTracks(0), //! output file with stored fTreeGenTracks
93 fLoader(0), //! pointer to the run loader
94 fTreeD(0), //! current tree with digits
95 fTreeTR(0), //! current tree with TR
96 fStack(0), //! current stack
97 fGenInfo(0), //! array with pointers to gen info
98 fNInfos(0), //! number of tracks with infos
99 fParamTPC(0), //! AliTPCParam
107 ////////////////////////////////////////////////////////////////////////
108 AliGenInfoMaker::AliGenInfoMaker(const char * fnGalice, const char* fnRes,
109 Int_t nEvents, Int_t firstEvent):
110 fDebug(0), //! debug flag
111 fEventNr(0), //! current event number
112 fLabel(0), //! track label
113 fNEvents(0), //! number of events to process
114 fFirstEventNr(0), //! first event to process
115 fNParticles(0), //! number of particles in TreeK
116 fTreeGenTracks(0), //! output tree with generated tracks
117 fTreeKinks(0), //! output tree with Kinks
118 fTreeV0(0), //! output tree with V0
119 fTreeHitLines(0), //! tree with hit lines
120 fFileGenTracks(0), //! output file with stored fTreeGenTracks
121 fLoader(0), //! pointer to the run loader
122 fTreeD(0), //! current tree with digits
123 fTreeTR(0), //! current tree with TR
124 fStack(0), //! current stack
125 fGenInfo(0), //! array with pointers to gen info
126 fNInfos(0), //! number of tracks with infos
127 fParamTPC(0), //! AliTPCParam
137 fFirstEventNr = firstEvent;
138 fEventNr = firstEvent;
140 sprintf(fFnRes,"%s",fnRes);
142 fLoader = AliRunLoader::Open(fnGalice);
144 delete gAlice->GetRunLoader();
148 if (fLoader->LoadgAlice()){
149 cerr<<"Error occured while l"<<endl;
151 Int_t nall = fLoader->GetNumberOfEvents();
159 cerr<<"no events available"<<endl;
163 if (firstEvent+nEvents>nall) {
164 fEventNr = nall-firstEvent;
165 cerr<<"restricted number of events availaible"<<endl;
167 AliMagF * magf = gAlice->Field();
168 AliTracker::SetFieldMap(magf,0);
172 AliMCInfo * AliGenInfoMaker::MakeInfo(UInt_t i)
175 // Make info structure for given particle index
178 if (fGenInfo[i]) return fGenInfo[i];
179 fGenInfo[i] = new AliMCInfo;
187 ////////////////////////////////////////////////////////////////////////
188 AliGenInfoMaker::~AliGenInfoMaker()
195 fLoader->UnloadgAlice();
201 Int_t AliGenInfoMaker::SetIO()
204 // Set IO for given event
206 CreateTreeGenTracks();
207 if (!fTreeGenTracks) return 1;
208 // AliTracker::SetFieldFactor();
210 fParamTPC = GetTPCParam();
215 ////////////////////////////////////////////////////////////////////////
216 Int_t AliGenInfoMaker::SetIO(Int_t eventNr)
221 fLoader->SetEventNumber(eventNr);
223 fLoader->LoadHeader();
224 fLoader->LoadKinematics();
225 fStack = fLoader->Stack();
227 fLoader->LoadTrackRefs();
229 fTreeTR = fLoader->TreeTR();
231 AliTPCLoader * tpcl = (AliTPCLoader*)fLoader->GetLoader("TPCLoader");
233 fTreeD = tpcl->TreeD();
237 Int_t AliGenInfoMaker::CloseIOEvent()
240 // Close IO for current event
242 fLoader->UnloadHeader();
243 fLoader->UnloadKinematics();
244 fLoader->UnloadTrackRefs();
245 AliTPCLoader * tpcl = (AliTPCLoader*)fLoader->GetLoader("TPCLoader");
246 tpcl->UnloadDigits();
250 Int_t AliGenInfoMaker::CloseIO()
252 fLoader->UnloadgAlice();
258 ////////////////////////////////////////////////////////////////////////
259 Int_t AliGenInfoMaker::Exec(Int_t nEvents, Int_t firstEventNr)
262 // Execute action for
264 // firstEventNr - first event number
267 fFirstEventNr = firstEventNr;
271 ////////////////////////////////////////////////////////////////////////
272 Int_t AliGenInfoMaker::Exec()
275 // Make a comparision MC tree
276 // Connect MC information -TPArticle - AliTrackRefernces ...
280 Int_t status =SetIO();
281 if (status>0) return status;
284 for (fEventNr = fFirstEventNr; fEventNr < fFirstEventNr+fNEvents;
287 fNParticles = fStack->GetNtrack();
289 fGenInfo = new AliMCInfo*[fNParticles];
290 for (UInt_t i = 0; i<fNParticles; i++) {
294 cout<<"Start to process event "<<fEventNr<<endl;
295 cout<<"\tfNParticles = "<<fNParticles<<endl;
296 if (fDebug>2) cout<<"\n\n\n\tStart loop over TreeTR"<<endl;
297 if (TreeTRLoop()>0) return 1;
299 if (fDebug>2) cout<<"\n\n\n\tStart loop over TreeD"<<endl;
300 if (TreeDLoop()>0) return 1;
302 if (fDebug>2) cout<<"\n\n\n\tStart loop over TreeK"<<endl;
303 if (TreeKLoop()>0) return 1;
304 if (fDebug>2) cout<<"\tEnd loop over TreeK"<<endl;
306 if (BuildKinkInfo()>0) return 1;
307 if (BuildV0Info()>0) return 1;
308 //if (BuildHitLines()>0) return 1;
309 if (fDebug>2) cout<<"\tEnd loop over TreeK"<<endl;
311 for (UInt_t i = 0; i<fNParticles; i++) {
312 if (fGenInfo[i]) delete fGenInfo[i];
321 cerr<<"Exec finished"<<endl;
328 ////////////////////////////////////////////////////////////////////////
329 void AliGenInfoMaker::CreateTreeGenTracks()
334 fFileGenTracks = TFile::Open(fFnRes,"RECREATE");
335 if (!fFileGenTracks) {
336 cerr<<"Error in CreateTreeGenTracks: cannot open file "<<fFnRes<<endl;
339 fTreeGenTracks = new TTree("genTracksTree","genTracksTree");
340 AliMCInfo * info = new AliMCInfo;
341 fTreeGenTracks->Branch("MC","AliMCInfo",&info);
344 AliGenKinkInfo *kinkinfo = new AliGenKinkInfo;
345 fTreeKinks = new TTree("genKinksTree","genKinksTree");
346 fTreeKinks->Branch("MC","AliGenKinkInfo",&kinkinfo);
349 AliGenV0Info *v0info = new AliGenV0Info;
350 fTreeV0 = new TTree("genV0Tree","genV0Tree");
351 fTreeV0->Branch("MC","AliGenV0Info",&v0info);
355 AliTrackPointArray * points0 = new AliTrackPointArray;
356 AliTrackPointArray * points1 = new AliTrackPointArray;
357 AliTrackPointArray * points2 = new AliTrackPointArray;
358 fTreeHitLines = new TTree("HitLines","HitLines");
359 fTreeHitLines->Branch("TPC.","AliTrackPointArray",&points0,32000,10);
360 fTreeHitLines->Branch("TRD.","AliTrackPointArray",&points1,32000,10);
361 fTreeHitLines->Branch("ITS.","AliTrackPointArray",&points2,32000,10);
363 fTreeHitLines->Branch("Label",&label,"label/I");
365 fTreeGenTracks->AutoSave();
366 fTreeKinks->AutoSave();
368 fTreeHitLines->AutoSave();
370 ////////////////////////////////////////////////////////////////////////
371 void AliGenInfoMaker::CloseOutputFile()
374 // Close Output files
376 if (!fFileGenTracks) {
377 cerr<<"File "<<fFnRes<<" not found as an open file."<<endl;
380 fFileGenTracks->cd();
381 fTreeGenTracks->Write();
382 delete fTreeGenTracks;
388 fFileGenTracks->Close();
389 delete fFileGenTracks;
393 ////////////////////////////////////////////////////////////////////////
394 Int_t AliGenInfoMaker::TreeKLoop()
397 // open the file with treeK
398 // loop over all entries there and save information about some tracks
401 AliStack * stack = fStack;
402 if (!stack) {cerr<<"Stack was not found!\n"; return 1;}
405 cout<<"There are "<<fNParticles<<" primary and secondary particles in event "
409 TParticlePDG *ppdg = 0;
410 // not all generators give primary vertex position. Take the vertex
411 // of the particle 0 as primary vertex.
412 TDatabasePDG pdg; //get pdg table
413 //thank you very much root for this
414 TBranch * br = fTreeGenTracks->GetBranch("MC");
415 TParticle *particle = stack->ParticleFromTreeK(0);
416 fVPrim[0] = particle->Vx();
417 fVPrim[1] = particle->Vy();
418 fVPrim[2] = particle->Vz();
419 for (UInt_t iParticle = 0; iParticle < fNParticles; iParticle++) {
420 // load only particles with TR
421 AliMCInfo * info = GetInfo(iParticle);
423 //////////////////////////////////////////////////////////////////////
424 info->fLabel = iParticle;
426 info->fParticle = *(stack->Particle(iParticle));
427 info->fVDist[0] = info->fParticle.Vx()-fVPrim[0];
428 info->fVDist[1] = info->fParticle.Vy()-fVPrim[1];
429 info->fVDist[2] = info->fParticle.Vz()-fVPrim[2];
430 info->fVDist[3] = TMath::Sqrt(info->fVDist[0]*info->fVDist[0]+
431 info->fVDist[1]*info->fVDist[1]+info->fVDist[2]*info->fVDist[2]);
434 ipdg = info->fParticle.GetPdgCode();
436 ppdg = pdg.GetParticle(ipdg);
437 info->fEventNr = fEventNr;
439 //////////////////////////////////////////////////////////////////////
440 br->SetAddress(&info);
441 fTreeGenTracks->Fill();
443 fTreeGenTracks->AutoSave();
444 if (fDebug > 2) cerr<<"end of TreeKLoop"<<endl;
449 ////////////////////////////////////////////////////////////////////////////////////
450 Int_t AliGenInfoMaker::BuildKinkInfo()
453 // Build tree with Kink Information
455 AliStack * stack = fStack;
456 if (!stack) {cerr<<"Stack was not found!\n"; return 1;}
459 cout<<"There are "<<fNParticles<<" primary and secondary particles in event "
463 //TParticlePDG *ppdg = 0;
464 // not all generators give primary vertex position. Take the vertex
465 // of the particle 0 as primary vertex.
466 TDatabasePDG pdg; //get pdg table
467 //thank you very much root for this
468 TBranch * br = fTreeKinks->GetBranch("MC");
470 AliGenKinkInfo * kinkinfo = new AliGenKinkInfo;
472 for (UInt_t iParticle = 0; iParticle < fNParticles; iParticle++) {
473 // load only particles with TR
474 AliMCInfo * info = GetInfo(iParticle);
476 if (info->fCharge==0) continue;
477 if (info->fTRdecay.P()<0.13) continue; //momenta cut
478 if (info->fTRdecay.R()>500) continue; //R cut - decay outside barrel
479 TParticle & particle = info->fParticle;
480 Int_t first = particle.GetDaughter(0);
481 if (first ==0) continue;
483 Int_t last = particle.GetDaughter(1);
484 if (last ==0) last=first;
485 AliMCInfo * info2 = 0;
486 AliMCInfo * dinfo = 0;
488 for (Int_t id2=first;id2<=last;id2++){
489 info2 = GetInfo(id2);
490 if (!info2) continue;
491 TParticle &p2 = info2->fParticle;
492 Double_t r = TMath::Sqrt(p2.Vx()*p2.Vx()+p2.Vy()*p2.Vy());
493 if ( TMath::Abs(info->fTRdecay.R()-r)>0.01) continue;
494 if (!(p2.GetPDG())) continue;
497 kinkinfo->SetInfoMother(*info);
498 kinkinfo->SetInfoDaughter(*dinfo);
499 if (kinkinfo->GetMinus().GetParticle().GetPDG()==0 || kinkinfo->GetPlus().GetParticle().GetPDG()==0) continue;
501 br->SetAddress(&kinkinfo);
506 kinkinfo->fMCm = (*info);
507 kinkinfo->GetPlus() = (*dinfo);
509 br->SetAddress(&kinkinfo);
514 fTreeGenTracks->AutoSave();
515 if (fDebug > 2) cerr<<"end of Kink Loop"<<endl;
520 ////////////////////////////////////////////////////////////////////////////////////
521 Int_t AliGenInfoMaker::BuildV0Info()
524 // Build tree with V0 Information
526 AliStack * stack = fStack;
527 if (!stack) {cerr<<"Stack was not found!\n"; return 1;}
530 cout<<"There are "<<fNParticles<<" primary and secondary particles in event "
534 //TParticlePDG *ppdg = 0;
535 // not all generators give primary vertex position. Take the vertex
536 // of the particle 0 as primary vertex.
537 TDatabasePDG pdg; //get pdg table
538 //thank you very much root for this
539 TBranch * br = fTreeV0->GetBranch("MC");
541 AliGenV0Info * v0info = new AliGenV0Info;
544 for (UInt_t iParticle = 0; iParticle < fNParticles; iParticle++) {
545 // load only particles with TR
546 AliMCInfo * info = GetInfo(iParticle);
548 if (info->fCharge==0) continue;
551 TParticle & particle = info->fParticle;
552 Int_t mother = particle.GetMother(0);
553 if (mother <=0) continue;
555 TParticle * motherparticle = stack->Particle(mother);
556 if (!motherparticle) continue;
558 Int_t last = motherparticle->GetDaughter(1);
559 if (last==(int)iParticle) continue;
560 AliMCInfo * info2 = info;
561 AliMCInfo * dinfo = GetInfo(last);
562 if (!dinfo) continue;
563 if (!dinfo->fParticle.GetPDG()) continue;
564 if (!info2->fParticle.GetPDG()) continue;
566 v0info->SetInfoP(*info);
567 v0info->SetInfoM(*dinfo);
568 v0info->SetMother(*motherparticle);
569 v0info->Update(fVPrim);
570 br->SetAddress(&v0info);
575 if (fDebug > 2) cerr<<"end of V0 Loop"<<endl;
582 ////////////////////////////////////////////////////////////////////////
583 Int_t AliGenInfoMaker::BuildHitLines()
587 // open the file with treeK
588 // loop over all entries there and save information about some tracks
591 AliStack * stack = fStack;
592 if (!stack) {cerr<<"Stack was not found!\n"; return 1;}
595 cout<<"There are "<<fNParticles<<" primary and secondary particles in event "
599 // // TParticlePDG *ppdg = 0;
600 // // not all generators give primary vertex position. Take the vertex
601 // // of the particle 0 as primary vertex.
602 // TDatabasePDG pdg; //get pdg table
603 // //thank you very much root for this
604 // AliTrackPointArray *tpcp = new AliTrackPointArray;
605 // AliTrackPointArray *trdp = new AliTrackPointArray;
606 // AliTrackPointArray *itsp = new AliTrackPointArray;
609 // TBranch * brtpc = fTreeHitLines->GetBranch("TPC.");
610 // TBranch * brtrd = fTreeHitLines->GetBranch("TRD.");
611 // TBranch * brits = fTreeHitLines->GetBranch("ITS.");
612 // TBranch * brlabel = fTreeHitLines->GetBranch("Label");
613 // brlabel->SetAddress(&label);
614 // brtpc->SetAddress(&tpcp);
615 // brtrd->SetAddress(&trdp);
616 // brits->SetAddress(&itsp);
618 // AliDetector *dtpc = gAlice->GetDetector("TPC");
619 // AliDetector *dtrd = gAlice->GetDetector("TRD");
620 // AliDetector *dits = gAlice->GetDetector("ITS");
622 // for (UInt_t iParticle = 0; iParticle < fNParticles; iParticle++) {
623 // // load only particles with TR
624 // AliMCInfo * info = GetInfo(iParticle);
625 // if (!info) continue;
626 // Int_t primpart = info->fPrimPart;
627 // ipdg = info->fParticle.GetPdgCode();
628 // label = iParticle;
630 // gAlice->ResetHits();
634 // tpcp->fLabel1 = ipdg;
635 // trdp->fLabel1 = ipdg;
636 // itsp->fLabel1 = ipdg;
637 // if (dtpc->TreeH()->GetEvent(primpart)){
638 // dtpc->LoadPoints(primpart);
639 // tpcp->Reset(dtpc,iParticle);
641 // if (dtrd->TreeH()->GetEvent(primpart)){
642 // dtrd->LoadPoints(primpart);
643 // trdp->Reset(dtrd,iParticle);
645 // if (dits->TreeH()->GetEvent(primpart)){
646 // dits->LoadPoints(primpart);
647 // itsp->Reset(dits,iParticle);
650 // fTreeHitLines->Fill();
651 // dtpc->ResetPoints();
652 // dtrd->ResetPoints();
653 // dits->ResetPoints();
658 // fTreeHitLines->AutoSave();
659 // if (fDebug > 2) cerr<<"end of TreeKLoop"<<endl;
664 ////////////////////////////////////////////////////////////////////////
665 Int_t AliGenInfoMaker::TreeDLoop()
668 // open the file with treeD
669 // loop over all entries there and save information about some tracks
672 Int_t nInnerSector = fParamTPC->GetNInnerSector();
674 Int_t zero=fParamTPC->GetZeroSup()+6;
677 AliSimDigits digitsAddress, *digits=&digitsAddress;
678 fTreeD->GetBranch("Segment")->SetAddress(&digits);
680 Int_t sectorsByRows=(Int_t)fTreeD->GetEntries();
681 if (fDebug > 1) cout<<"\tsectorsByRows = "<<sectorsByRows<<endl;
682 for (Int_t i=0; i<sectorsByRows; i++) {
683 if (!fTreeD->GetEvent(i)) continue;
685 fParamTPC->AdjustSectorRow(digits->GetID(),sec,row);
686 if (fDebug > 1) cout<<sec<<' '<<row<<" \r";
687 // here I expect that upper sectors follow lower sectors
688 if (sec > nInnerSector) rowShift = fParamTPC->GetNRowLow();
690 digits->ExpandTrackBuffer();
693 Int_t iRow=digits->CurrentRow();
694 Int_t iColumn=digits->CurrentColumn();
695 Short_t digitValue = digits->CurrentDigit();
696 if (digitValue >= zero) {
698 for (Int_t j = 0; j<3; j++) {
699 // label = digits->GetTrackID(iRow,iColumn,j);
700 label = digits->GetTrackIDFast(iRow,iColumn,j)-2;
701 if (label >= (Int_t)fNParticles) { //don't label from bakground event
704 if (label >= 0 && label <= (Int_t)fNParticles) {
706 cout<<"Inner loop: sector, iRow, iColumn, label, value, row "
708 <<iRow<<" "<<iColumn<<" "<<label<<" "<<digitValue
711 AliMCInfo * info = GetInfo(label);
713 info->fTPCRow.SetRow(row+rowShift);
718 } while (digits->Next());
721 if (fDebug > 2) cerr<<"end of TreeDLoop"<<endl;
726 ////////////////////////////////////////////////////////////////////////
727 Int_t AliGenInfoMaker::TreeTRLoop()
730 // loop over TrackReferences and store the first one for each track
732 TTree * treeTR = fTreeTR;
733 Int_t nPrimaries = (Int_t) treeTR->GetEntries();
734 if (fDebug > 1) cout<<"There are "<<nPrimaries<<" entries in TreeTR"<<endl;
737 //track references for TPC
738 TClonesArray* tpcArrayTR = new TClonesArray("AliTrackReference");
739 TClonesArray* itsArrayTR = new TClonesArray("AliTrackReference");
740 TClonesArray* trdArrayTR = new TClonesArray("AliTrackReference");
741 TClonesArray* tofArrayTR = new TClonesArray("AliTrackReference");
742 TClonesArray* runArrayTR = new TClonesArray("AliTrackReference");
744 if (treeTR->GetBranch("TPC")) treeTR->GetBranch("TPC")->SetAddress(&tpcArrayTR);
745 if (treeTR->GetBranch("ITS")) treeTR->GetBranch("ITS")->SetAddress(&itsArrayTR);
746 if (treeTR->GetBranch("TRD")) treeTR->GetBranch("TRD")->SetAddress(&trdArrayTR);
747 if (treeTR->GetBranch("TOF")) treeTR->GetBranch("TOF")->SetAddress(&tofArrayTR);
748 if (treeTR->GetBranch("AliRun")) treeTR->GetBranch("AliRun")->SetAddress(&runArrayTR);
752 for (Int_t iPrimPart = 0; iPrimPart<nPrimaries; iPrimPart++) {
753 treeTR->GetEntry(iPrimPart);
755 // Loop over TPC references
757 for (Int_t iTrackRef = 0; iTrackRef < tpcArrayTR->GetEntriesFast(); iTrackRef++) {
758 AliTrackReference *trackRef = (AliTrackReference*)tpcArrayTR->At(iTrackRef);
760 if (trackRef->TestBit(BIT(2))){
762 if (trackRef->P()<fTPCPtCut) continue;
763 Int_t label = trackRef->GetTrack();
764 AliMCInfo * info = GetInfo(label);
765 if (!info) info = MakeInfo(label);
766 info->fTRdecay = *trackRef;
769 if (trackRef->P()<fTPCPtCut) continue;
770 Int_t label = trackRef->GetTrack();
771 AliMCInfo * info = GetInfo(label);
772 if (!info) info = MakeInfo(label);
774 info->fPrimPart = iPrimPart;
775 TClonesArray & arr = *(info->fTPCReferences);
776 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
779 // Loop over ITS references
781 for (Int_t iTrackRef = 0; iTrackRef < itsArrayTR->GetEntriesFast(); iTrackRef++) {
782 AliTrackReference *trackRef = (AliTrackReference*)itsArrayTR->At(iTrackRef);
785 if (trackRef->P()<fTPCPtCut) continue;
786 Int_t label = trackRef->GetTrack();
787 AliMCInfo * info = GetInfo(label);
788 if ( (!info) && trackRef->Pt()<fITSPtCut) continue;
789 if (!info) info = MakeInfo(label);
791 info->fPrimPart = iPrimPart;
792 TClonesArray & arr = *(info->fITSReferences);
793 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
796 // Loop over TRD references
798 for (Int_t iTrackRef = 0; iTrackRef < trdArrayTR->GetEntriesFast(); iTrackRef++) {
799 AliTrackReference *trackRef = (AliTrackReference*)trdArrayTR->At(iTrackRef);
801 if (trackRef->P()<fTRDPtCut) continue;
802 Int_t label = trackRef->GetTrack();
803 AliMCInfo * info = GetInfo(label);
804 if ( (!info) && trackRef->Pt()<fTRDPtCut) continue;
805 if (!info) info = MakeInfo(label);
807 info->fPrimPart = iPrimPart;
808 TClonesArray & arr = *(info->fTRDReferences);
809 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
812 // Loop over TOF references
814 for (Int_t iTrackRef = 0; iTrackRef < tofArrayTR->GetEntriesFast(); iTrackRef++) {
815 AliTrackReference *trackRef = (AliTrackReference*)tofArrayTR->At(iTrackRef);
816 Int_t label = trackRef->GetTrack();
817 AliMCInfo * info = GetInfo(label);
819 if (trackRef->Pt()<fTPCPtCut) continue;
820 if ( (!info) && trackRef->Pt()<fTOFPtCut) continue;
822 if (!info) info = MakeInfo(label);
824 info->fPrimPart = iPrimPart;
825 TClonesArray & arr = *(info->fTOFReferences);
826 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
830 // get dacay position
832 for (Int_t iTrackRef = 0; iTrackRef < runArrayTR->GetEntriesFast(); iTrackRef++) {
833 AliTrackReference *trackRef = (AliTrackReference*)runArrayTR->At(iTrackRef);
835 Int_t label = trackRef->GetTrack();
836 AliMCInfo * info = GetInfo(label);
838 if (!trackRef->TestBit(BIT(2))) continue; //if not decay
839 // if (TMath::Abs(trackRef.X());
840 info->fTRdecay = *trackRef;
844 tpcArrayTR->Delete();
846 trdArrayTR->Delete();
848 tofArrayTR->Delete();
850 itsArrayTR->Delete();
852 runArrayTR->Delete();
855 return TreeTRLoopNew();
862 ////////////////////////////////////////////////////////////////////////
863 Int_t AliGenInfoMaker::TreeTRLoopNew()
866 // loop over TrackReferences and store the first one for each track
868 TTree * treeTR = fTreeTR;
869 Int_t nPrimaries = (Int_t) treeTR->GetEntries();
870 if (fDebug > 1) cout<<"There are "<<nPrimaries<<" entries in TreeTR"<<endl;
874 TClonesArray* runArrayTR = new TClonesArray("AliTrackReference");
876 if (treeTR->GetBranch("TrackReferences")) treeTR->GetBranch("TrackReferences")->SetAddress(&runArrayTR);
880 for (Int_t ipart = 0; ipart<fStack->GetNtrack(); ipart++) {
881 TParticle * part = fStack->Particle(ipart);
883 if (part->Pt()<fITSPtCut) continue;
884 if (part->R()>250.) continue;
885 if (TMath::Abs(part->Vz())>250.) continue;
886 if (TMath::Abs(part->Pz()/part->Pt())>2.5) continue;
892 for (Int_t iPrimPart = 0; iPrimPart<nPrimaries; iPrimPart++) {
893 treeTR->GetEntry(iPrimPart);
895 // gettrack references
897 for (Int_t iTrackRef = 0; iTrackRef < runArrayTR->GetEntriesFast(); iTrackRef++) {
898 AliTrackReference *trackRef = (AliTrackReference*)runArrayTR->At(iTrackRef);
900 Int_t label = trackRef->GetTrack();
904 if (trackRef->DetectorId()== AliTrackReference::kTPC){
906 Int_t label = trackRef->GetTrack();
907 AliMCInfo * info = GetInfo(label);
908 if (!info && trackRef->Pt()<fTPCPtCut) continue;
909 if (!info) info = MakeInfo(label);
911 info->fPrimPart = iPrimPart;
912 TClonesArray & arr = *(info->fTPCReferences);
913 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
918 if (trackRef->DetectorId()== AliTrackReference::kITS){
920 Int_t label = trackRef->GetTrack();
921 AliMCInfo * info = GetInfo(label);
922 if (!info && trackRef->Pt()<fITSPtCut) continue;
923 if (!info) info = MakeInfo(label);
925 info->fPrimPart = iPrimPart;
926 TClonesArray & arr = *(info->fITSReferences);
927 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
932 if (trackRef->DetectorId()== AliTrackReference::kTRD){
934 Int_t label = trackRef->GetTrack();
935 AliMCInfo * info = GetInfo(label);
936 if (!info&&trackRef->P()<fTRDPtCut) continue;
937 if (!info) info = MakeInfo(label);
939 info->fPrimPart = iPrimPart;
940 TClonesArray & arr = *(info->fTRDReferences);
941 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
946 if (trackRef->DetectorId()== AliTrackReference::kTOF){
948 Int_t label = trackRef->GetTrack();
949 AliMCInfo * info = GetInfo(label);
950 if (!info&&trackRef->P()<fTOFPtCut) continue;
951 if (!info) info = MakeInfo(label);
953 info->fPrimPart = iPrimPart;
954 TClonesArray & arr = *(info->fTOFReferences);
955 new (arr[arr.GetEntriesFast()]) AliTrackReference(*trackRef);
960 if (trackRef->DetectorId()== AliTrackReference::kDisappeared){
962 AliMCInfo * info = GetInfo(label);
964 //if (!trackRef->TestBit(BIT(2))) continue; //if not decay
965 info->fTRdecay = *trackRef;
975 ////////////////////////////////////////////////////////////////////////
976 Float_t AliGenInfoMaker::TR2LocalX(AliTrackReference *trackRef,
977 AliTPCParam *paramTPC) const {
980 // Transformation from Gloabal to local tacking system
982 Float_t x[3] = { trackRef->X(),trackRef->Y(),trackRef->Z()};
984 paramTPC->Transform0to1(x,index);
985 paramTPC->Transform1to2(x,index);
988 ////////////////////////////////////////////////////////////////////////
992 AliTPCParam * AliGenInfoMaker::GetTPCParam(){
994 // create default TPC parameters
996 AliTPCParamSR * par = new AliTPCParamSR;