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 **************************************************************************/
15 //_________________________________________________________________________
16 // An analysis task to check the TOF in simulated data
19 //-Distributions of the matching performance
20 //-TOF Time info and (TOF - expected time) plots
21 //-Summary Plots on TOF Pid
22 //////////////////////////////////////////////////////////////////////////////
33 #include "AliTOFQATask.h"
35 #include "AliESDtrack.h"
37 #include "Riostream.h"
39 //______________________________________________________________________________
40 AliTOFQATask::AliTOFQATask(const char *name) :
41 AliAnalysisTask(name,""),
65 // Input slot #0 works with an Ntuple
66 DefineInput(0, TChain::Class());
67 // Output slot #0 writes into a TH1 container
68 DefineOutput(0, TObjArray::Class()) ;
71 //______________________________________________________________________________
72 AliTOFQATask::AliTOFQATask(const AliTOFQATask &qatask) :
73 AliAnalysisTask("AliTOFQATask",""),
99 fOutputContainer=qatask.fOutputContainer;
100 fhTOFMatch=qatask.fhTOFMatch;
101 fhESDeffPhi=qatask.fhESDeffPhi;
102 fhESDeffTheta=qatask.fhESDeffTheta;
103 fhESDeffMom=qatask.fhESDeffMom;
104 fhTOFeffPhi=qatask.fhTOFeffPhi;
105 fhTOFeffTheta=qatask.fhTOFeffTheta;
106 fhTOFeffMom=qatask.fhTOFeffMom;
107 fhTOFeffPhiMT=qatask.fhTOFeffPhiMT;
108 fhTOFeffThetaMT=qatask.fhTOFeffThetaMT;
109 fhTOFeffMomMT=qatask.fhTOFeffMomMT;
110 fhTOFsector=qatask.fhTOFsector;
111 fhTOFsectorMT=qatask.fhTOFsectorMT;
112 fhTOFTime=qatask.fhTOFTime;
113 fhTOFDeltaTime=qatask.fhTOFDeltaTime;
114 fhTOFDeltaTimeMT=qatask.fhTOFDeltaTimeMT;
115 fhTOFIDSpecies=qatask.fhTOFIDSpecies;
116 fhTOFMassVsMom=qatask.fhTOFMassVsMom;
117 fhTOFMass=qatask.fhTOFMass;
119 //______________________________________________________________________________
120 AliTOFQATask:: ~AliTOFQATask()
122 delete fOutputContainer;
125 delete fhESDeffTheta;
128 delete fhTOFeffTheta;
130 delete fhTOFeffPhiMT;
131 delete fhTOFeffThetaMT;
132 delete fhTOFeffMomMT;
134 delete fhTOFsectorMT;
136 delete fhTOFDeltaTime;
137 delete fhTOFDeltaTimeMT;
138 delete fhTOFIDSpecies;
139 delete fhTOFMassVsMom;
142 //______________________________________________________________________________
143 AliTOFQATask& AliTOFQATask::operator=(const AliTOFQATask &qatask)
145 //assignment operator
146 this->fChain=qatask.fChain;
147 this->fESD=qatask.fESD;
148 this->fOutputContainer=qatask.fOutputContainer;
149 this->fhTOFMatch=qatask.fhTOFMatch;
150 this->fhESDeffPhi=qatask.fhESDeffPhi;
151 this->fhESDeffTheta=qatask.fhESDeffTheta;
152 this->fhESDeffMom=qatask.fhESDeffMom;
153 this->fhTOFeffPhi=qatask.fhTOFeffPhi;
154 this->fhTOFeffTheta=qatask.fhTOFeffTheta;
155 this->fhTOFeffMom=qatask.fhTOFeffMom;
156 this->fhTOFeffPhiMT=qatask.fhTOFeffPhiMT;
157 this->fhTOFeffThetaMT=qatask.fhTOFeffThetaMT;
158 this->fhTOFeffMomMT=qatask.fhTOFeffMomMT;
159 this->fhTOFsector=qatask.fhTOFsector;
160 this->fhTOFsectorMT=qatask.fhTOFsectorMT;
161 this->fhTOFTime=qatask.fhTOFTime;
162 this->fhTOFDeltaTime=qatask.fhTOFDeltaTime;
163 this->fhTOFDeltaTimeMT=qatask.fhTOFDeltaTimeMT;
164 this->fhTOFIDSpecies=qatask.fhTOFIDSpecies;
165 this->fhTOFMassVsMom=qatask.fhTOFMassVsMom;
166 this->fhTOFMass=qatask.fhTOFMass;
169 //______________________________________________________________________________
170 void AliTOFQATask::ConnectInputData(const Option_t*)
172 // Initialisation of branch container and histograms
174 AliInfo(Form("*** Initialization of %s", GetName())) ;
177 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
179 AliError(Form("Input 0 for %s not found\n", GetName()));
183 // One should first check if the branch address was taken by some other task
184 char ** address = (char **)GetBranchAddress(0, "ESD");
186 fESD = (AliESD*)(*address);
189 SetBranchAddress(0, "ESD", &fESD);
193 //________________________________________________________________________
194 void AliTOFQATask::CreateOutputObjects()
200 //______________________________________________________________________________
201 void AliTOFQATask::Exec(Option_t *)
204 //******* The loop over events --------------------------------------------------
209 Int_t pisel=0,kasel=0,prsel=0,elsel=0,musel=0;
210 //Set equal a-priori weights (just charged hadrions)
211 Int_t nCalinSec=8736;
212 Double_t c[5]={0, 0, 1, 1, 1};
213 // Processing of one event
214 Long64_t entry = fChain->GetReadEntry() ;
216 AliError("fESD is not connected to the input!") ;
220 if ( !((entry-1)%100) )
221 AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
223 // ************************ TOF *************************************
226 Int_t ntrk = fESD->GetNumberOfTracks() ;
230 AliESDtrack * t = fESD->GetTrack(ntrk) ;
231 if ( (t->GetStatus() & AliESDtrack::kTIME)==0 )continue;
235 Double_t mom = t->GetP() ;
236 Double_t phi = TMath::ATan2(t->GetX(),t->GetY()) ;
237 Double_t theta = TMath::ACos(t->GetZ()/
238 TMath::Sqrt(t->GetX()*t->GetX()+t->GetY()*t->GetY()+t->GetZ()*t->GetZ()));
239 phi*=180/TMath::Pi();
240 theta*=180/TMath::Pi();
242 fhESDeffPhi->Fill(phi);
243 fhESDeffTheta->Fill(theta);
244 fhESDeffMom->Fill(mom);
247 if(t->GetTOFsignal()<0)continue;
251 Double_t time=t->GetTOFsignal();//TOF time in ps
252 Int_t detid=t->GetTOFCalChannel();//which pad was hit
253 Int_t sector = detid/nCalinSec;
254 fhTOFTime->Fill(time*1.E-3);
255 fhTOFeffPhi->Fill(phi);
256 fhTOFeffTheta->Fill(theta);
257 fhTOFeffMom->Fill(mom);
258 fhTOFsector->Fill(sector);
261 Int_t label=TMath::Abs(t->GetLabel());
262 Int_t clab[3]; t->GetTOFLabel(clab);
263 if(label==clab[0] || label==clab[1] || label==clab[2]) {
264 fhTOFeffPhiMT->Fill(phi);
265 fhTOFeffThetaMT->Fill(theta);
266 fhTOFeffMomMT->Fill(mom);
267 fhTOFsectorMT->Fill(sector);
272 UInt_t status=AliESDtrack::kESDpid;status|=AliESDtrack::kTOFpid;
273 if (!((t->GetStatus()&status) == status))continue;
276 t->GetIntegratedTimes(times);//ps
277 Double_t l =t->GetIntegratedLength()/100.; // (m)
279 Double_t invBetaGamma= (0.299*time*1.E-3/l)*(0.299*time*1.E-3/l) -1.;
280 if(invBetaGamma<0){mass = -mom*TMath::Sqrt(-invBetaGamma);}
281 else{mass = mom*TMath::Sqrt(invBetaGamma);}
283 //The Mass/ vs Momentum Plot:
284 fhTOFMassVsMom->Fill(mass,mom);
285 fhTOFMass->Fill(mass);
288 Double_t r[10]; t->GetTOFpid(r);
291 for (i=0; i<AliPID::kSPECIES; i++) rcc+=(c[i]*r[i]);
292 if (rcc==0.) continue;
294 for (i=0; i<AliPID::kSPECIES; i++) w[i]=c[i]*r[i]/rcc;
296 fhTOFDeltaTime->Fill((time-times[2])*1.E-3);
298 if(label==clab[0] || label==clab[1] || label==clab[2]) {
299 fhTOFDeltaTimeMT->Fill((time-times[2])*1.E-3);
302 if (w[4]>w[0] && w[4]>w[1] && w[4]>w[2] && w[4]>w[3]){
304 fhTOFIDSpecies->Fill(4);
306 if (w[3]>w[0] && w[3]>w[1] && w[3]>w[2] && w[3]>w[4]){
308 fhTOFIDSpecies->Fill(3);
310 if (w[2]>w[0] && w[2]>w[1] && w[2]>w[3] && w[2]>w[4]){
312 fhTOFIDSpecies->Fill(2);
314 if (w[1]>w[0] && w[1]>w[2] && w[1]>w[3] && w[1]>w[4]){
316 fhTOFIDSpecies->Fill(1);
318 if (w[0]>w[1] && w[0]>w[2] && w[0]>w[3] && w[0]>w[4]){
320 fhTOFIDSpecies->Fill(0);
325 if(nselESD>10)fracM=((Float_t) nmatchTOF)/((Float_t) nselESD);
326 fhTOFMatch->Fill(fracM);
328 PostData(0, fOutputContainer);
329 AliInfo("Finishing event processing...") ;
332 //______________________________________________________________________________
333 void AliTOFQATask::BookHistos()
335 // Construct histograms:
338 new TH1F("hTOFMatch","Fraction of Matched TOF tracks",101,-0.005,1.005);
340 new TH1F("hESDeffPhi","ESD tracks Phi(vtx)", 180, -180., 180.) ;
342 new TH1F("hESDeffTheta","ESD tracks Theta (vtx)",90, 45., 135.) ;
344 new TH1F("hESDeffMom","ESD tracks Momentum (vtx)",40, 0., 6.) ;
346 new TH1F("hTOFeffPhi","TOF, Matched vs Phi(vtx)", 180,-180, 180.);
348 new TH1F("hTOFeffPhiMT","TOF, Well Matched vs Phi(vtx)",180,-180,180.);
350 new TH1F("hTOFeffTheta","TOF, Matched vs Theta(vtx)",90,45.,135.);
352 new TH1F("hTOFeffThetaMT","TOF, Well Matched vs Theta(vtx)",90,45.,135.);
354 new TH1F("hTOFeffMom","TOF, Matched vs Momentum ", 40, 0.,6.);
356 new TH1F("hTOFeffMomMT","TOF, Well Matched vs Momentum ", 40, 0.,6.);
358 new TH1F("hTOFsector","TOF, Matched vs Sector ", 18,0.,18.);
360 new TH1F("hTOFsectorMT","TOF, Well Matched vs Sector", 18, 0.,18.);
362 fhESDeffMom->Sumw2(); fhTOFeffMom->Sumw2(); fhTOFeffMomMT->Sumw2();
363 fhESDeffTheta->Sumw2(); fhTOFeffTheta->Sumw2(); fhTOFeffThetaMT->Sumw2();
364 fhESDeffPhi->Sumw2(); fhTOFeffPhi->Sumw2(); fhTOFeffPhiMT->Sumw2();
365 fhTOFsector->Sumw2(); fhTOFsectorMT->Sumw2();
368 new TH1F("hTOFTime","TOF, t(TOF)in ns ",1000,0,100.);
370 new TH1F("hTOFDeltaTime","TOF,t(TOF)-t(exp,pion), ns ",1000,-4.4,20.);
372 new TH1F("hTOFDeltaTimeMT","TOF, t(TOF)-t(exp,pion) for Well Matched, ns ",1000,-4.4,20.);
374 new TH1F("hTOFIDSpecies","TOF, ID Sample Composition ",5,-0.5,4.5);
376 new TH2F("hTOFMassVsMom","TOF, Mass Vs Momentum ",280,-0.2,1.2,600,0.,6.);
378 new TH1F("hTOFMass","TOF, reconstructed mass ",280,-0.2,1.2);
381 // create the output container
383 fOutputContainer = new TObjArray(18) ;
384 fOutputContainer->SetName(GetName()) ;
386 fOutputContainer->AddAt(fhTOFMatch, 0) ;
387 fOutputContainer->AddAt(fhESDeffPhi, 1) ;
388 fOutputContainer->AddAt(fhESDeffTheta, 2) ;
389 fOutputContainer->AddAt(fhESDeffMom, 3) ;
390 fOutputContainer->AddAt(fhTOFeffPhi, 4) ;
391 fOutputContainer->AddAt(fhTOFeffPhiMT, 5) ;
392 fOutputContainer->AddAt(fhTOFeffTheta, 6) ;
393 fOutputContainer->AddAt(fhTOFeffThetaMT, 7) ;
394 fOutputContainer->AddAt(fhTOFeffMom, 8) ;
395 fOutputContainer->AddAt(fhTOFeffMomMT, 9) ;
396 fOutputContainer->AddAt(fhTOFsector, 10) ;
397 fOutputContainer->AddAt(fhTOFsectorMT, 11) ;
398 fOutputContainer->AddAt(fhTOFTime, 12) ;
399 fOutputContainer->AddAt(fhTOFDeltaTime , 13) ;
400 fOutputContainer->AddAt(fhTOFDeltaTimeMT, 14) ;
401 fOutputContainer->AddAt(fhTOFIDSpecies, 15) ;
402 fOutputContainer->AddAt(fhTOFMassVsMom, 16) ;
403 fOutputContainer->AddAt(fhTOFMass, 17) ;
406 //______________________________________________________________________________
407 void AliTOFQATask::GetEfficiency()
409 // calculates the efficiency
412 fhTOFeffPhiMT->Divide(fhTOFeffPhiMT,fhTOFeffPhi,1,1,"B");
413 fhTOFeffThetaMT->Divide(fhTOFeffThetaMT,fhTOFeffTheta,1,1,"B");
414 fhTOFeffMomMT->Divide(fhTOFeffMomMT,fhTOFeffMom,1,1,"B");
416 fhTOFeffPhi->Divide(fhTOFeffPhi,fhESDeffPhi,1,1,"B");
417 fhTOFeffTheta->Divide(fhTOFeffTheta,fhESDeffTheta,1,1,"B");
418 fhTOFeffMom->Divide(fhTOFeffMom,fhESDeffMom,1,1,"B");
422 //______________________________________________________________________________
423 void AliTOFQATask::DrawHistos()
427 AliInfo("Plotting....") ;
429 gROOT->SetStyle("Plain");
430 gStyle->SetPalette(1);
431 gStyle->SetOptStat(110010);
432 gStyle->SetPalette(1);
433 gStyle->SetCanvasColor(0);
434 gStyle->SetFrameFillColor(0);
436 TGaxis::SetMaxDigits(3);
437 gStyle->SetLabelFont(52, "XYZ");
438 gStyle->SetTitleFont(62, "XYZ");
439 gStyle->SetPadRightMargin(0.02);
442 TCanvas * cTOFeff = new TCanvas("cTOFeff", "TOF ESD General", 400, 30, 550, 630) ;
444 cTOFeff->Divide(1,2) ;
446 fhTOFMatch->GetXaxis()->SetTitle("Fraction of matched ESD tracks/event");
447 fhTOFMatch->GetYaxis()->SetTitle("Events");
448 fhTOFMatch->SetFillColor(4);
451 TPad *b = (TPad*)gPad;
454 fhTOFsector->GetXaxis()->SetTitle("Sector index of matched TOF Cluster");
455 fhTOFsector->GetYaxis()->SetTitle("Entries");
456 fhTOFsector->SetFillColor(4);
457 fhTOFsector->SetMinimum(0.);
458 fhTOFsector->GetYaxis()->SetNdivisions(205);
459 fhTOFsector->GetYaxis()->SetTitleOffset(1.2);
460 fhTOFsector->Draw("histo");
462 fhTOFeffMom->SetMaximum(1.2);
463 fhTOFeffMom->GetXaxis()->SetTitle("Track momentum (GeV/c)");
464 fhTOFeffMom->GetYaxis()->SetNdivisions(205);
465 fhTOFeffMom->GetYaxis()->SetTitleOffset(1.2);
466 fhTOFeffMom->GetYaxis()->SetTitle("Fraction of matched ESD tracks");
467 fhTOFeffMom->SetFillColor(4);
469 fhTOFeffMom->Draw("histo,same");
471 cTOFeff->Print("TOF_eff.gif");
472 cTOFeff->Print("TOF_eff.eps");
476 TCanvas * cTOFtime = new TCanvas("cTOFtime", "TOF measured Times ", 400, 30, 550, 630) ;
477 cTOFtime->Divide(1,2) ;
479 cTOFtime->GetPad(1)->SetLogy(1);
480 fhTOFTime->GetXaxis()->SetTitle("TOF time (ns)");
481 fhTOFTime->GetYaxis()->SetTitle("Entries");
482 fhTOFTime->SetFillColor(4);
485 cTOFtime->GetPad(2)->SetLogy(1);
486 fhTOFDeltaTime->GetXaxis()->SetTitle("t^{TOF}-t^{exp}_{#pi} (ns)");
487 fhTOFDeltaTime->GetYaxis()->SetTitle("Entries");
488 fhTOFDeltaTime->SetFillColor(4);
489 fhTOFDeltaTime->Draw();
491 cTOFtime->Print("TOF_time.gif");
492 cTOFtime->Print("TOF_time.eps");
495 TCanvas * cTOFpid = new TCanvas("cTOFpid", "TOF PID ", 400, 30, 550, 630) ;
496 cTOFpid->Divide(1,3) ;
499 fhTOFMass->GetXaxis()->SetTitle("Reconstructed Mass (GeV/c^{2})");
500 fhTOFMass->GetYaxis()->SetTitle("Entries");
501 fhTOFMass->SetFillColor(4);
505 fhTOFMassVsMom->GetYaxis()->SetRange(0,400);
506 fhTOFMassVsMom->GetXaxis()->SetTitle("Reconstructed Mass (GeV/c^{2})");
507 fhTOFMassVsMom->GetYaxis()->SetTitle("Track Momentum (GeV/c)");
508 fhTOFMassVsMom->GetXaxis()->SetTitleSize(0.05);
509 fhTOFMassVsMom->GetYaxis()->SetTitleSize(0.05);
510 fhTOFMassVsMom->SetMarkerStyle(20);
511 fhTOFMassVsMom->SetMarkerSize(0.05);
512 fhTOFMassVsMom->SetMarkerColor(2);
513 fhTOFMassVsMom->Draw();
516 TLatex * tex = new TLatex(1., 1.25, "Bayesian PID: a-priori concentrations: [0,0,1,1,1]");
517 tex->SetTextColor(1);
518 tex->SetTextSize(0.045);
519 tex->SetLineWidth(2);
521 Float_t norm=1./fhTOFIDSpecies->GetEntries();
522 fhTOFIDSpecies->Scale(norm);
523 fhTOFIDSpecies->SetMaximum(1.2);
524 fhTOFIDSpecies->GetXaxis()->SetTitle("Particle Type");
525 fhTOFIDSpecies->GetYaxis()->SetTitle("ID Fractions");
526 fhTOFIDSpecies->GetXaxis()->SetTitleSize(0.05);
527 fhTOFIDSpecies->GetYaxis()->SetTitleSize(0.05);
528 fhTOFIDSpecies->SetFillColor(4);
529 fhTOFIDSpecies->Draw();
534 Float_t pifrac=fhTOFIDSpecies->GetBinContent(3);
535 Float_t kafrac=fhTOFIDSpecies->GetBinContent(4);
536 Float_t prfrac=fhTOFIDSpecies->GetBinContent(5);
538 sprintf(ch,"pion fraction = %5.3f",pifrac);
539 TLatex * texpi = new TLatex(-0.3, 0.9, ch);
540 texpi->SetTextColor(1);
541 texpi->SetTextSize(0.05);
542 texpi->SetLineWidth(2);
544 sprintf(ch,"kaon fraction = %5.3f",kafrac);
545 TLatex * texka = new TLatex(-0.3, 0.8, ch);
546 texka->SetTextColor(1);
547 texka->SetTextSize(0.05);
548 texka->SetLineWidth(2);
550 sprintf(ch,"proton fraction = %5.3f",prfrac);
551 TLatex * texpr = new TLatex(-0.3, 0.7, ch);
552 texpr->SetTextColor(1);
553 texpr->SetTextSize(0.05);
554 texpr->SetLineWidth(2);
558 cTOFpid->Print("TOF_pid.gif");
559 cTOFpid->Print("TOF_pid.eps");
568 //______________________________________________________________________________
569 void AliTOFQATask::Terminate(Option_t *)
571 // Processing when the event loop is ended
575 AliInfo("TOF QA Task: End of events loop");
576 fOutputContainer = (TObjArray*)GetOutputData(0);
577 fhTOFMatch = (TH1F*)fOutputContainer->At(0);
578 fhESDeffPhi = (TH1F*)fOutputContainer->At(1);
579 fhESDeffTheta = (TH1F*)fOutputContainer->At(2);
580 fhESDeffMom = (TH1F*)fOutputContainer->At(3);
581 fhTOFeffPhi = (TH1F*)fOutputContainer->At(4);
582 fhTOFeffPhiMT = (TH1F*)fOutputContainer->At(5);
583 fhTOFeffTheta = (TH1F*)fOutputContainer->At(6);
584 fhTOFeffThetaMT = (TH1F*)fOutputContainer->At(7);
585 fhTOFeffMom = (TH1F*)fOutputContainer->At(8);
586 fhTOFeffMomMT = (TH1F*)fOutputContainer->At(9);
587 fhTOFsector = (TH1F*)fOutputContainer->At(10);
588 fhTOFsectorMT = (TH1F*)fOutputContainer->At(11);
589 fhTOFTime = (TH1F*)fOutputContainer->At(12);
590 fhTOFDeltaTime = (TH1F*)fOutputContainer->At(13);
591 fhTOFDeltaTimeMT = (TH1F*)fOutputContainer->At(14);
592 fhTOFIDSpecies = (TH1F*)fOutputContainer->At(15);
593 fhTOFMassVsMom = (TH2F*)fOutputContainer->At(16);
594 fhTOFMass = (TH1F*)fOutputContainer->At(17);
596 // PostData(0, fOutputContainer);