Structure modules as built
[u/mrichter/AliRoot.git] / test / genkine / sim / test.C
CommitLineData
3928b038 1// Usage in compiled mode
2// gSystem->SetIncludePath("-I$ROOTSYS/include -I$ALICE_ROOT/include");
3// gROOT->LoadMacro("test.C+");
4// test()
5
6#if !defined(__CINT__) || defined(__MAKECINT__)
7
8// Root include files
9#include <Riostream.h>
10#include <TFile.h>
11#include <TTree.h>
12#include <TBranch.h>
13#include <TStopwatch.h>
14#include <TObject.h>
15#include <TParticle.h>
16
17// AliRoot include files
18#include "AliESD.h"
19#include "AliRunLoader.h"
20#include "AliRun.h"
21#include "AliStack.h"
22
23#endif
24
25void test() {
26
27 TStopwatch timer;
28 timer.Start();
29
30 TString name;
31
32 // Signal file, tree, and branch
33 name = "AliESDs.root";
34 TFile * fSig = TFile::Open(name.Data());
35 TTree * tSig = (TTree*)fSig->Get("esdTree");
36 TBranch * bSig = tSig->GetBranch("ESD");
37
38 AliESD * esdSig = 0; // The signal ESD object is put here
39 bSig->SetAddress(&esdSig);
40
41 // Run loader
42 name = "galice.root";
43 AliRunLoader* rlSig = AliRunLoader::Open(name.Data());
44
45 // gAlice
46 rlSig->LoadgAlice();
47 gAlice = rlSig->GetAliRun();
48
49 // Now load kinematics and event header
50 rlSig->LoadKinematics();
51 rlSig->LoadHeader();
52
53 // Loop on events
54 Long64_t nevSig = rlSig->GetNumberOfEvents();
55
56 cout << nevSig << " events" << endl;
57
58 for (Int_t iev=0; iev<nevSig; iev++) {
59 cout << "Signal event " << iev << endl;
60
61 // Get ESD
62 bSig->GetEntry(iev);
63
64 // Get MC
65 rlSig->GetEvent(iev);
66
67 // Particle stack
68 AliStack * stackSig = rlSig->Stack();
69
70 Int_t nrec = esdSig->GetNumberOfTracks();
71
72 //-------------------------------------------------------------------------------------------------
73 Int_t nstack = stackSig->GetNtrack();
74 for(Int_t istack=0; istack < nstack; istack++){
75
76 TParticle * part = stackSig->Particle(istack);
77
78 // Loop on particles: check if the D* decay products are reconstructed
79
80 if(!part) continue;
81
82 if(TMath::Abs(part->GetPdgCode())== 413){
83 cout<<"particle "<< istack << " is D*"<<endl;
84
85
86 Int_t iDaughter1 = part->GetFirstDaughter(); //id of the Daughter = D^0
87 if( iDaughter1<0) continue;
88 TParticle* daughter1 = stackSig->Particle(iDaughter1);
89 if(!daughter1) continue;
90 cout<<"first daughter: "<<daughter1->GetPdgCode()<<endl;
91
92 Int_t iDaughter2 = part->GetLastDaughter(); //id of the Daughter = pi+
93 if( iDaughter2<0) continue;
94 TParticle* daughter2 = stackSig->Particle(iDaughter2);
95 if(!daughter2) continue;
96 cout<<"last daughter: "<<daughter2->GetPdgCode()<<endl;
97
98 Int_t iD0 = -1;
99 Int_t iPi = -1;
100
101 if(TMath::Abs(daughter1->GetPdgCode())== 421){
102 iD0=iDaughter1;
103 iPi=iDaughter2;
104 }
105 else if(TMath::Abs(daughter2->GetPdgCode())== 421){
106 iD0=iDaughter2;
107 iPi=iDaughter1;
108 }
109
110 if (iD0<0) continue;
111
112 TParticle* secondmother = stackSig->Particle(iD0);
113
114 Int_t iD0Daughter1 = secondmother->GetFirstDaughter();
115 TParticle* D0Daughter1 = stackSig->Particle(iD0Daughter1);
116 Int_t iD0Daughter2 = secondmother->GetLastDaughter();
117 TParticle* D0Daughter2 = stackSig->Particle(iD0Daughter2);
118
119 for(Int_t irec=0; irec<nrec; irec++) {//loop on the ESDTree;
120 AliESDtrack * track = esdSig->GetTrack(irec);
121 UInt_t label = TMath::Abs(track->GetLabel());
122 if(label<10000000) {
123 if(label == iPi) cout<<label<< " We found the Pi from the D* decay"<<endl;
124 if(label == iD0Daughter1) cout<<label<<" We found the K from the D0 decay"<<endl;
125 if(label == iD0Daughter2) cout<<label<<" We found the Pi from the D0 decay"<<endl;
126 }
127 }
128 }
129 }
130 }
131
132 // end loop on kine tree
133
134 fSig->Close();
135
136 timer.Stop();
137 timer.Print();
138}