Additional protection
[u/mrichter/AliRoot.git] / TOF / AliTOFComparison.C
CommitLineData
3f83f224 1/****************************************************************************
2 * This macro estimates efficiency of matching with the TOF. *
3 * TOF "Good" tracks are those originating from the primary vertex, *
4 * being "good" in the ITS and having at least one digit in the TOF. *
5 * (To get the list of "good" tracks one should first run *
6 * AliTPCComparison.C and AliITSComparisonV2.C macros) *
7 * Origin: I.Belikov, CERN, Jouri.Belikov@cern.ch *
8 ****************************************************************************/
9
10#if !defined(__CINT__) || defined(__MAKECINT__)
4fa1e629 11 #include <TMath.h>
12 #include <TError.h>
3f83f224 13 #include <Riostream.h>
4fa1e629 14 #include <TH1F.h>
15 #include <TTree.h>
16 #include <TParticle.h>
17 #include <TCanvas.h>
18 #include <TLine.h>
19 #include <TText.h>
20 #include <TBenchmark.h>
21 #include <TStyle.h>
8b462fd8 22 #include <TFile.h>
4fa1e629 23 #include <TROOT.h>
3f83f224 24
3f83f224 25 #include "AliStack.h"
4fa1e629 26 #include "AliHeader.h"
27 #include "AliTrackReference.h"
3f83f224 28 #include "AliRunLoader.h"
4fa1e629 29 #include "AliRun.h"
4514157e 30 #include "AliESDEvent.h"
0c73e091 31 #include "AliESDtrack.h"
4fa1e629 32
f4af68da 33 #include "AliTOFcluster.h"
4fa1e629 34 #include "AliLoader.h"
3f83f224 35
3f83f224 36 #include "TClonesArray.h"
3f83f224 37#endif
38
4fa1e629 39Int_t GoodTracksTOF(const Char_t *dir=".");
3f83f224 40
41extern AliRun *gAlice;
4fa1e629 42extern TBenchmark *gBenchmark;
43extern TROOT *gROOT;
44
45static Int_t allgood=0;
46static Int_t allmatched=0;
47static Int_t allmismatched=0;
48
49Int_t AliTOFComparison(const Char_t *dir=".") {
50 gBenchmark->Start("AliTOFComparison");
51
52 ::Info("AliTOFComparison.C","Doing comparison...");
53
54
55 Double_t pmin=0.2;
56 Double_t pmax=3.0;
57
58 TH1F *hgood=(TH1F*)gROOT->FindObject("hgood");
59 if (!hgood) hgood=new TH1F("hgood","Good tracks",30,pmin,pmax);
60
61 TH1F *hfound=(TH1F*)gROOT->FindObject("hfound");
62 if (!hfound) hfound=new TH1F("hfound","Matched tracks",30,pmin,pmax);
63
64 TH1F *hfake=(TH1F*)gROOT->FindObject("hfake");
65 if (!hfake) hfake=new TH1F("hfake","Mismatched tracks",30,pmin,pmax);
66
67 TH1F *hgp=(TH1F*)gROOT->FindObject("hgp");
68 if (!hgp) hgp=new TH1F("hgp","",30,pmin,pmax);
69 hgp->SetLineColor(4); hgp->SetLineWidth(2);
70
71 TH1F *hfp=(TH1F*)gROOT->FindObject("hfp");
72 if (!hfp) hfp=new TH1F("hfp","Probability of mismatching",30,pmin,pmax);
73 hfp->SetFillColor(1); hfp->SetFillStyle(3013); hfp->SetLineWidth(2);
74
75 TH1F *hgoo=(TH1F*)gROOT->FindObject("hgoo");
76 if (!hgoo) hgoo=new TH1F("hgoo","Good tracks",30,-1,1);
77
78 TH1F *hfoun=(TH1F*)gROOT->FindObject("hfoun");
79 if (!hfoun) hfoun=new TH1F("hfoun","Matched tracks",30,-1,1);
80
81 TH1F *hfak=(TH1F*)gROOT->FindObject("hfak");
82 if (!hfak) hfak=new TH1F("hfak","Mismatched tracks",30,-1,1);
83
84 TH1F *hgl=(TH1F*)gROOT->FindObject("hgl");
85 if (!hgl) hgl=new TH1F("hgl","",30,-1,1);
86 hgl->SetLineColor(4); hgl->SetLineWidth(2);
87
88 TH1F *hfl=(TH1F*)gROOT->FindObject("hfl");
89 if (!hfl) hfl=new TH1F("hfl","Probability of mismatching",30,-1,1);
90 hfl->SetFillColor(1); hfl->SetFillStyle(3013); hfl->SetLineWidth(2);
91
3f83f224 92
3f83f224 93
4fa1e629 94 Char_t fname[100];
95 sprintf(fname,"%s/GoodTracksTOF.root",dir);
96
97 TFile *refFile=TFile::Open(fname,"old");
98 if (!refFile || !refFile->IsOpen()) {
99 ::Info("AliTOFComparison.C","Marking good tracks (will take a while)...");
100 if (GoodTracksTOF(dir)) {
101 ::Error("AliTOFComparison.C","Can't generate the reference file !");
102 return 1;
103 }
104 }
105 refFile=TFile::Open(fname,"old");
106 if (!refFile || !refFile->IsOpen()) {
107 ::Error("AliTOFComparison.C","Can't open the reference file !");
108 return 1;
109 }
110
111 TTree *tofTree=(TTree*)refFile->Get("tofTree");
112 if (!tofTree) {
113 ::Error("AliTOFComparison.C","Can't get the reference tree !");
114 return 2;
115 }
116 TBranch *branch=tofTree->GetBranch("TOF");
117 if (!branch) {
118 ::Error("AliTOFComparison.C","Can't get the TOF branch !");
119 return 3;
120 }
121 TClonesArray dummy("AliTrackReference",1000), *refs=&dummy;
122 branch->SetAddress(&refs);
123
3f83f224 124
125 if (gAlice) {
880d6abe 126 delete AliRunLoader::GetRunLoader();
3f83f224 127 delete gAlice;//if everything was OK here it is already NULL
128 gAlice = 0x0;
129 }
4fa1e629 130 sprintf(fname,"%s/galice.root",dir);
131 AliRunLoader *rl = AliRunLoader::Open(fname,"COMPARISON");
3f83f224 132 if (rl == 0x0) {
133 cerr<<"Can not open session"<<endl;
134 return 1;
135 }
136 AliLoader* tofl = rl->GetLoader("TOFLoader");
137 if (tofl == 0x0) {
138 cerr<<"Can not get the TOF loader"<<endl;
139 return 2;
140 }
f4af68da 141 tofl->LoadRecPoints("read");
3f83f224 142
3f83f224 143
4fa1e629 144 sprintf(fname,"%s/AliESDs.root",dir);
145 TFile *ef=TFile::Open(fname);
146 if ((!ef)||(!ef->IsOpen())) {
147 ::Error("AliTOFComparison.C","Can't open AliESDs.root !");
148 delete rl;
3f83f224 149 return 4;
150 }
4514157e 151 AliESDEvent* event = new AliESDEvent();
8b462fd8 152 TTree* esdTree = (TTree*) ef->Get("esdTree");
153 if (!esdTree) {
154 ::Error("AliTOFComparison.C", "no ESD tree found");
155 return 5;
156 }
4514157e 157 event->ReadFromTree(esdTree);
3f83f224 158
3f83f224 159
3f83f224 160
4fa1e629 161 //******* Loop over events *********
162 Int_t e=0;
8b462fd8 163 while (esdTree->GetEvent(e)) {
4fa1e629 164 cout<<endl<<endl<<"********* Processing event number: "<<e<<"*******\n";
165
8b462fd8 166 rl->GetEvent(e);
4fa1e629 167
f4af68da 168 TTree *clsTree=tofl->TreeR();
169 if (!clsTree) {
4fa1e629 170 cerr<<"Can't get the TOF cluster tree !\n";
171 return 3;
172 }
f4af68da 173 TBranch *branch=clsTree->GetBranch("TOF");
4fa1e629 174 if (!branch) {
175 cerr<<"Can't get the branch with the TOF digits !\n";
176 return 4;
177 }
f4af68da 178 TClonesArray dummy("AliTOFcluster",10000), *clusters=&dummy;
179 branch->SetAddress(&clusters);
3f83f224 180
f4af68da 181 clsTree->GetEvent(0);
182 Int_t nd=clusters->GetEntriesFast();
183 cerr<<"Number of the TOF clusters: "<<nd<<endl;
1c3e95d8 184
3f83f224 185
3f83f224 186
3f83f224 187 Int_t ntrk=event->GetNumberOfTracks();
188 cerr<<"Number of ESD tracks : "<<ntrk<<endl;
189
4fa1e629 190
191 if (tofTree->GetEvent(e++)==0) {
192 cerr<<"No reconstructable tracks !\n";
193 continue;
194 }
195
196 Int_t ngood=refs->GetEntriesFast();
197
3f83f224 198 Int_t matched=0;
199 Int_t mismatched=0;
200 for (Int_t i=0; i<ngood; i++) {
4fa1e629 201 AliTrackReference *ref=(AliTrackReference*)refs->UncheckedAt(i);
202 Int_t lab=ref->Label();
203 Float_t ptg=TMath::Sqrt(ref->Px()*ref->Px() + ref->Py()*ref->Py());
3f83f224 204
4fa1e629 205 Double_t tgl=ref->Pz()/ptg; //tan(lambda)
3f83f224 206
207 if (ptg>pmin) { hgood->Fill(ptg); hgoo->Fill(tgl); }
208
209 Int_t j;
210 AliESDtrack *t=0;
211 for (j=0; j<ntrk; j++) {
212 AliESDtrack *tt=event->GetTrack(j);
213 if (lab!=TMath::Abs(tt->GetLabel())) continue;
214 t=tt;
215 //if ((tt->GetStatus()&AliESDtrack::kTOFpid) == 0) continue;
216 if (tt->GetTOFsignal() < 0) continue;
217 UInt_t idx=tt->GetTOFcluster();
218 if ((Int_t)idx>=nd) {
f4af68da 219 cerr<<"Wrong cluster index ! "<<idx<<endl;
3f83f224 220 return 5;
221 }
bad86b82 222 AliTOFcluster *cls=(AliTOFcluster*)clusters->At(idx);
223 if (cls) {
224 if (cls->GetLabel(0)!=lab)
225 if (cls->GetLabel(1)!=lab)
226 if (cls->GetLabel(2)!=lab) {
227 mismatched++;
228 if (ptg>pmin) { hfake->Fill(ptg); hfak->Fill(tgl); }
229 break;
230 }
231 if (ptg>pmin) { hfound->Fill(ptg); hfoun->Fill(tgl); }
232 matched++;
233 break;
234 }
3f83f224 235 }
236 if (j==ntrk) {
237 cerr<<"Not matched: "<<lab<<" ";
238 if (t) {
239 cerr<<(t->GetStatus()&AliESDtrack::kITSout)<<' '
240 <<(t->GetStatus()&AliESDtrack::kTPCout)<<' '
241 <<(t->GetStatus()&AliESDtrack::kTRDout)<<' '
242 <<(t->GetStatus()&AliESDtrack::kTIME);
243 } else cerr<<"No ESD track !";
244 cerr<<endl;
245 }
246 }
4fa1e629 247 cout<<"Number of good tracks: "<<ngood<<endl;
248 cout<<"Number of matched tracks: "<<matched<<endl;
249 cout<<"Number of mismatched tracks: "<<mismatched<<endl;
3f83f224 250
4fa1e629 251 allgood+=ngood; allmatched+=matched; allmismatched+=mismatched;
252
253 refs->Clear();
4fa1e629 254 } //***** End of the loop over events
3f83f224 255
8b462fd8 256 delete event;
4514157e 257 delete esdTree;
4fa1e629 258 ef->Close();
259
260 delete tofTree;
261 refFile->Close();
262
263 if (allgood!=0) cerr<<"\n\nEfficiency: "<<Float_t(allmatched)/allgood<<endl;
264 cout<<"Total number of good tracks: "<<allgood<<endl;
265 cout<<"Total number of matched tracks: "<<allmatched<<endl;
266 cout<<"Total number of mismatched tracks: "<<allmismatched<<endl;
267
268 TCanvas *c1=(TCanvas*)gROOT->FindObject("c1");
269 if (!c1) {
270 c1=new TCanvas("c1","",0,0,600,900);
271 c1->Divide(1,2);
272 }
273 hfound->Sumw2(); hgood->Sumw2(); hfake->Sumw2();
274 hgp->Divide(hfound,hgood,1,1.,"b");
275 hfp->Divide(hfake,hgood,1,1.,"b");
276 hgp->SetMaximum(1.4);
277 hgp->SetYTitle("Matching efficiency");
278 hgp->SetXTitle("Pt (GeV/c)");
279
280 hfoun->Sumw2(); hgoo->Sumw2(); hfak->Sumw2();
281 hgl->Divide(hfoun,hgoo,1,1.,"b");
282 hfl->Divide(hfak,hgoo,1,1.,"b");
283 hgl->SetMaximum(1.4);
284 hgl->SetYTitle("Matching efficiency");
285 hgl->SetXTitle("Tan(lambda)");
286
287 c1->cd(1);
288
289 hgp->Draw();
290 hfp->Draw("histsame");
291 TLine *line1 = new TLine(pmin,1.0,pmax,1.0); line1->SetLineStyle(4);
292 line1->Draw("same");
293 TLine *line2 = new TLine(pmin,0.9,pmax,0.9); line2->SetLineStyle(4);
294 line2->Draw("same");
295
296 c1->cd(2);
297
298 hgl->Draw();
299 hfl->Draw("histsame");
300 TLine *line3 = new TLine(-1,1.0,1,1.0); line3->SetLineStyle(4);
301 line3->Draw("same");
302 TLine *line4 = new TLine(-1,0.9,1,0.9); line4->SetLineStyle(4);
303 line4->Draw("same");
304
305 c1->Update();
306
3f7a17bb 307 TFile fc("AliTOFComparison.root","RECREATE");
308 c1->Write();
309 fc.Close();
310
4fa1e629 311 gBenchmark->Stop("AliTOFComparison");
312 gBenchmark->Show("AliTOFComparison");
313
314 delete rl;
3f83f224 315 return 0;
316}
317
4fa1e629 318Int_t GoodTracksTOF(const Char_t *dir) {
319 if (gAlice) {
880d6abe 320 delete AliRunLoader::GetRunLoader();
4fa1e629 321 delete gAlice;//if everything was OK here it is already NULL
322 gAlice = 0x0;
3f83f224 323 }
324
4fa1e629 325 Char_t fname[100];
326 sprintf(fname,"%s/galice.root",dir);
3f83f224 327
4fa1e629 328 AliRunLoader *rl = AliRunLoader::Open(fname,"COMPARISON");
329 if (!rl) {
330 ::Error("GoodTracksTOF","Can't start session !");
331 return 1;
332 }
3f83f224 333
334 rl->LoadgAlice();
335 rl->LoadHeader();
4fa1e629 336 rl->LoadKinematics();
3f83f224 337
3f83f224 338
339 AliLoader* tofl = rl->GetLoader("TOFLoader");
340 if (tofl == 0x0) {
4fa1e629 341 ::Error("GoodTracksTOF","Can not get the TOF loader !");
342 delete rl;
343 return 2;
3f83f224 344 }
f4af68da 345 tofl->LoadRecPoints("read");
3f83f224 346
4fa1e629 347 Int_t nev=rl->GetNumberOfEvents();
348 ::Info("GoodTracksTOF","Number of events : %d\n",nev);
3f83f224 349
4fa1e629 350 sprintf(fname,"%s/GoodTracksITS.root",dir);
351 TFile *itsFile=TFile::Open(fname);
352 if ((!itsFile)||(!itsFile->IsOpen())) {
353 ::Error("GoodTracksTOF","Can't open the GoodTracksITS.root !");
354 delete rl;
355 return 5;
3f83f224 356 }
4fa1e629 357 TClonesArray dum("AliTrackReference",1000), *itsRefs=&dum;
358 TTree *itsTree=(TTree*)itsFile->Get("itsTree");
359 if (!itsTree) {
360 ::Error("GoodTracksTOF","Can't get the ITS reference tree !");
361 delete rl;
362 return 6;
363 }
364 TBranch *itsBranch=itsTree->GetBranch("ITS");
365 if (!itsBranch) {
366 ::Error("GoodTracksTOF","Can't get the ITS reference branch !");
367 delete rl;
368 return 7;
369 }
370 itsBranch->SetAddress(&itsRefs);
371
372
373 sprintf(fname,"%s/GoodTracksTOF.root",dir);
374 TFile *tofFile=TFile::Open(fname,"recreate");
375 TClonesArray dummy("AliTrackReference",1000), *tofRefs=&dummy;
376 TTree tofTree("tofTree","Tree with info about the reconstructable TOF tracks");
377 tofTree.Branch("TOF",&tofRefs);
378
379
380 //******** Loop over generated events
381 for (Int_t e=0; e<nev; e++) {
3f83f224 382
4fa1e629 383 rl->GetEvent(e); tofFile->cd();
384
385 Int_t np = rl->GetHeader()->GetNtrack();
386 cout<<"Event "<<e<<" Number of particles: "<<np<<endl;
387
388 //******** Fill the "good" masks
389 Int_t *good=new Int_t[np]; Int_t k; for (k=0; k<np; k++) good[k]=0;
390
f4af68da 391 TTree *cTree=tofl->TreeR();
392 if (!cTree) {
4fa1e629 393 ::Error("GoodTracksTOF","Can't get the TOF cluster tree !");
394 delete rl;
395 return 8;
396 }
f4af68da 397 TBranch *branch=cTree->GetBranch("TOF");
4fa1e629 398 if (!branch) {
399 ::Error("GoodTracksTOF","Can't get the branch with the TOF digits !");
400 return 9;
401 }
f4af68da 402 TClonesArray dummy("AliTOFcluster",10000), *clusters=&dummy;
403 branch->SetAddress(&clusters);
3f83f224 404
f4af68da 405 cTree->GetEvent(0);
406 Int_t nd=clusters->GetEntriesFast();
4fa1e629 407
408 for (Int_t i=0; i<nd; i++) {
f4af68da 409 AliTOFcluster *c=(AliTOFcluster*)clusters->UncheckedAt(i);
410 Int_t l0=c->GetLabel(0);
4fa1e629 411 if (l0>=np) {cerr<<"Wrong label: "<<l0<<endl; continue;}
f4af68da 412 Int_t l1=c->GetLabel(1);
4fa1e629 413 if (l1>=np) {cerr<<"Wrong label: "<<l1<<endl; continue;}
f4af68da 414 Int_t l2=c->GetLabel(2);
4fa1e629 415 if (l2>=np) {cerr<<"Wrong label: "<<l2<<endl; continue;}
416 if (l0>=0) good[l0]++;
417 if (l1>=0) good[l1]++;
418 if (l2>=0) good[l2]++;
419 }
f4af68da 420 clusters->Clear();
4fa1e629 421
422
423 //****** select tracks which are "good" enough
424 AliStack* stack = rl->Stack();
425
426 itsTree->GetEvent(e);
427 Int_t nk=itsRefs->GetEntriesFast();
428
429 Int_t nt=0;
430 for (k=0; k<nk; k++) {
431 AliTrackReference *itsRef=(AliTrackReference *)itsRefs->UncheckedAt(k);
432 Int_t lab=itsRef->Label();
433 if (good[lab] == 0) continue;
434 TParticle *p = (TParticle*)stack->Particle(lab);
435 if (p == 0x0) {
436 cerr<<"Can not get particle "<<lab<<endl;
437 continue;
438 }
439
440 if (TMath::Abs(p->Vx())>0.1) continue;
441 if (TMath::Abs(p->Vy())>0.1) continue;
442 //if (TMath::Abs(p->Vz())>0.1) continue;
443
444 new((*tofRefs)[nt]) AliTrackReference(*itsRef);
445 nt++;
446 }
447 itsRefs->Clear();
448
449 tofTree.Fill();
450 tofRefs->Clear();
451
452 delete[] good;
453
454 } //*** end of the loop over generated events
455
456 tofTree.Write();
457 tofFile->Close();
458
459 delete itsTree;
460 itsFile->Close();
461
462 delete rl;
463 return 0;
3f83f224 464}