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 **************************************************************************/
16 //==============================================================================
17 // AliHMPIDAnalysysTask - Class representing a basic analysis tool of HMPID data
18 // A set of histograms is created.
19 //==============================================================================
21 // By means of AliHMPIDAnalysisTask.C macro it is possible to use this class
22 // to perform the analysis on local data, on data on alien using local machine
25 #ifndef AliHMPIDAnalysisTASK_CXX
26 #define AliHMPIDAnalysisTASK_CXX
33 #include "TGraphErrors.h"
34 #include "AliAnalysisManager.h"
35 #include "AliESDInputHandler.h"
36 #include "AliMCEventHandler.h"
37 #include "AliMCEvent.h"
38 #include "AliESDtrack.h"
41 #include "AliHMPIDAnalysisTask.h"
43 ClassImp(AliHMPIDAnalysisTask)
45 //__________________________________________________________________________
46 AliHMPIDAnalysisTask::AliHMPIDAnalysisTask() :
47 fESD(0x0),fMC(0x0),fUseMC(kTRUE),
49 fHmpPesdPhmp(0x0),fHmpCkovPesd(0x0),fHmpCkovPhmp(0x0),
50 fHmpMipTrkDist(0x0),fHmpMipTrkDistX(0x0),fHmpMipTrkDistY(0x0),
51 fHmpMipCharge3cm(0x0),fHmpMipCharge1cm(0x0),
52 fHmpNumPhots(0x0),fHmpTrkFlags(0x0),
80 for (Int_t i=0; i<28; i++) fVar[i]=0;
83 //___________________________________________________________________________
84 AliHMPIDAnalysisTask::AliHMPIDAnalysisTask(const Char_t* name) :
85 AliAnalysisTaskSE(name),
86 fESD(0x0), fMC(0x0), fUseMC(kTRUE),
88 fHmpPesdPhmp(0x0), fHmpCkovPesd(0x0), fHmpCkovPhmp(0x0),
89 fHmpMipTrkDist(0x0), fHmpMipTrkDistX(0x0), fHmpMipTrkDistY(0x0),
90 fHmpMipCharge3cm(0x0), fHmpMipCharge1cm(0x0),
91 fHmpNumPhots(0x0), fHmpTrkFlags(0x0),
106 fThetavsPiFromK(0x0),
117 // Constructor. Initialization of Inputs and Outputs
119 for (Int_t i=0; i<28; i++) fVar[i]=0;
121 DefineOutput(1,TList::Class());
122 DefineOutput(2,TTree::Class());
125 //___________________________________________________________________________
126 AliHMPIDAnalysisTask& AliHMPIDAnalysisTask::operator=(const AliHMPIDAnalysisTask& c)
129 // Assignment operator
132 AliAnalysisTaskSE::operator=(c);
136 fHmpHistList = c.fHmpHistList;
137 fHmpPesdPhmp = c.fHmpPesdPhmp;
138 fHmpCkovPesd = c.fHmpCkovPesd;
139 fHmpCkovPhmp = c.fHmpCkovPhmp;
140 fHmpMipTrkDist = c.fHmpMipTrkDist;
141 fHmpMipTrkDistX = c.fHmpMipTrkDistX;
142 fHmpMipTrkDistY = c.fHmpMipTrkDistY;
143 fHmpMipCharge3cm = c.fHmpMipCharge3cm;
144 fHmpMipCharge1cm = c.fHmpMipCharge1cm;
145 fHmpNumPhots = c.fHmpNumPhots;
146 fHmpTrkFlags = c.fHmpTrkFlags;
149 fPionEff = c.fPionEff;
150 fKaonEff = c.fKaonEff;
151 fProtEff = c.fProtEff;
152 fPionTot = c.fPionTot;
153 fKaonTot = c.fKaonTot;
154 fProtTot = c.fProtTot;
155 fPionNot = c.fPionNot;
156 fKaonNot = c.fKaonNot;
157 fProtNot = c.fProtNot;
158 fPionCon = c.fPionCon;
159 fKaonCon = c.fKaonCon;
160 fProtCon = c.fProtCon;
161 fThetavsPiFromK = c.fThetavsPiFromK;
162 fThetapivsPesd = c.fThetapivsPesd;
163 fThetaKvsPesd = c.fThetaKvsPesd;
164 fThetaPvsPesd = c.fThetaPvsPesd;
165 fProtGen = c.fProtGen;
166 fPbarGen = c.fPbarGen;
167 fProtHmp = c.fProtHmp;
168 fPbarHmp = c.fPbarHmp;
170 for (Int_t i=0; i<28; i++) fVar[i]=c.fVar[i];
175 //___________________________________________________________________________
176 AliHMPIDAnalysisTask::AliHMPIDAnalysisTask(const AliHMPIDAnalysisTask& c) :
177 AliAnalysisTaskSE(c),
178 fESD(c.fESD),fMC(c.fMC),fUseMC(c.fUseMC),
179 fHmpHistList(c.fHmpHistList),
180 fHmpPesdPhmp(c.fHmpPesdPhmp),fHmpCkovPesd(c.fHmpCkovPesd),fHmpCkovPhmp(c.fHmpCkovPhmp),
181 fHmpMipTrkDist(c.fHmpMipTrkDist),fHmpMipTrkDistX(c.fHmpMipTrkDistX),fHmpMipTrkDistY(c.fHmpMipTrkDistY),
182 fHmpMipCharge3cm(c.fHmpMipCharge3cm), fHmpMipCharge1cm(c.fHmpMipCharge1cm),
183 fHmpNumPhots(c.fHmpNumPhots), fHmpTrkFlags(c.fHmpTrkFlags),
186 fPionEff(c.fPionEff),
187 fKaonEff(c.fKaonEff),
188 fProtEff(c.fProtEff),
189 fPionTot(c.fPionTot),
190 fKaonTot(c.fKaonTot),
191 fProtTot(c.fProtTot),
192 fPionNot(c.fPionNot),
193 fKaonNot(c.fKaonNot),
194 fProtNot(c.fProtNot),
195 fPionCon(c.fPionCon),
196 fKaonCon(c.fKaonCon),
197 fProtCon(c.fProtCon),
198 fThetavsPiFromK(c.fThetavsPiFromK),
199 fThetapivsPesd(c.fThetapivsPesd),
200 fThetaKvsPesd(c.fThetaKvsPesd),
201 fThetaPvsPesd(c.fThetaPvsPesd),
202 fProtGen(c.fProtGen),
203 fPbarGen(c.fPbarGen),
204 fProtHmp(c.fProtHmp),
205 fPbarHmp(c.fPbarHmp),
211 for (Int_t i=0; i<28; i++) fVar[i]=c.fVar[i];
214 //___________________________________________________________________________
215 AliHMPIDAnalysisTask::~AliHMPIDAnalysisTask() {
219 Info("~AliHMPIDAnalysisTask","Calling Destructor");
220 if (fHmpHistList && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fHmpHistList;
223 //___________________________________________________________________________
224 void AliHMPIDAnalysisTask::ConnectInputData(Option_t *)
228 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
230 AliDebug(2,Form("ERROR: Could not get ESDInputHandler"));
232 fESD = esdH->GetEvent();
236 AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
238 AliDebug(2,Form("ERROR: Could not get MCEventHandler"));
241 fMC = mcH->MCEvent();
242 if (!fMC) AliDebug(2,Form("ERROR: Could not get MCEvent"));
246 //___________________________________________________________________________
247 void AliHMPIDAnalysisTask::UserExec(Option_t *)
249 Double_t priors[5]={1.,1.,1.,1.,1.}; //{0.01,0.01,0.83,0.10,0.5};
251 AliPID *pPid = new AliPID();
252 pPid->SetPriors(priors);
253 AliESDtrack *track = 0;
254 TParticle *pPart = 0;
255 AliStack* pStack = 0;
258 pStack = fMC->Stack();
261 Double_t ktol = 0.001;
262 Double_t rin[3], rout[3];
263 Float_t distCut = 0.7;
264 Float_t thetaCut = 0.175;
267 // Main loop function, executed on Event basis
269 for (Int_t iTrack = 0; iTrack<fESD->GetNumberOfTracks(); iTrack++) {
271 track = fESD->GetTrack(iTrack);
273 track->GetInnerXYZ(rin);
274 track->GetOuterXYZ(rout);
276 if(Equal(track->GetHMPIDsignal(),-20.,ktol)) fHmpTrkFlags->Fill(0);
277 else if(Equal(track->GetHMPIDsignal(),-9.,ktol)) fHmpTrkFlags->Fill(1);
278 else if(Equal(track->GetHMPIDsignal(),-5.,ktol)) fHmpTrkFlags->Fill(2);
279 else if(Equal(track->GetHMPIDsignal(),-11.,ktol)) fHmpTrkFlags->Fill(3);
280 else if(Equal(track->GetHMPIDsignal(),-22.,ktol)) fHmpTrkFlags->Fill(4);
281 else fHmpTrkFlags->Fill(5);
283 if(Equal(track->GetHMPIDsignal(),-20.,ktol)) continue;
284 if(track->GetHMPIDcluIdx() < 0) continue;
288 Float_t xpc, ypc, th, ph;
289 track->GetHMPIDmip(x,y,q,nph);
290 track->GetHMPIDtrk(xpc,ypc,th,ph);
292 if(Equal(x,0.,ktol) && Equal(y,0.,ktol) && Equal(xpc,0.,ktol) && Equal(ypc,0.,ktol)) continue;
294 Float_t dist = TMath::Sqrt((xpc-x)*(xpc-x)+(ypc-y)*(ypc-y));
295 fHmpMipTrkDist->Fill(dist);
296 fHmpMipTrkDistX->Fill(xpc - x);
297 fHmpMipTrkDistY->Fill(ypc - y);
298 Double_t pHmp[3] = {0}, pHmp3 = 0;
299 if (track->GetOuterHmpPxPyPz(pHmp)) pHmp3 = TMath::Sqrt(pHmp[0]*pHmp[0]+pHmp[1]*pHmp[1]+pHmp[2]*pHmp[2]);
300 if (dist <= 3.0) fHmpMipCharge3cm->Fill(q);
304 track->GetImpactParameters(b,bCov);
306 track->GetHMPIDpid(probs);
307 pPid->SetProbabilities(probs);
309 if ((label = track->GetLabel()) < 0) continue;
310 pPart = pStack->Particle(label);
313 if(track->GetHMPIDsignal() > 0 ){
314 Double_t thetaCkov = track->GetHMPIDsignal() - (Int_t)track->GetHMPIDsignal();
315 if (pHmp3) fHmpPesdPhmp->Fill(track->P(),pHmp3);
316 if (dist<=1.0) fHmpMipCharge1cm->Fill(q);
317 fHmpNumPhots->Fill(nph);
318 fHmpCkovPesd->Fill(track->P(),thetaCkov);
319 if (pHmp3) fHmpCkovPhmp->Fill(pHmp3,thetaCkov);
321 if (fUseMC && dist<distCut && TMath::Abs(th)<thetaCut){
322 if (!pStack->IsPhysicalPrimary(label)) continue;
323 Int_t pdgCode = TMath::Abs(pPart->GetPdgCode());
325 fThetapivsPesd->Fill(track->P(),thetaCkov);
326 Int_t mot=pPart->GetFirstMother();
328 TParticle *pMot=pStack->Particle(mot);
329 TString str=pMot->GetName();
330 if (str.Contains("K")) fThetavsPiFromK->Fill(pHmp3,thetaCkov);
333 if (pdgCode==321) fThetaKvsPesd->Fill(track->P(),thetaCkov);
334 if (pdgCode==2212) fThetaPvsPesd->Fill(track->P(),thetaCkov);
336 if (track->P()<1. || track->P()>5.) continue;
337 Int_t pBin=(Int_t) (2*(track->P()-1));
338 if (pdgCode!=2212) fProtCon->Fill(pBin);
340 fProtTot->Fill(pBin);
341 fProtEff->Fill(pBin,pPid->GetProbability(AliPID::kProton));
342 fPionNot->Fill(pBin,pPid->GetProbability(AliPID::kPion));
343 fKaonNot->Fill(pBin,pPid->GetProbability(AliPID::kKaon));
345 if (pdgCode!=211) fPionCon->Fill(pBin);
346 if (pdgCode!=321) fKaonCon->Fill(pBin);
349 Float_t weight=pPid->GetProbability(AliPID::kElectron)+
350 pPid->GetProbability(AliPID::kMuon)+
351 pPid->GetProbability(AliPID::kPion);
352 fPionTot->Fill(pBin);
353 fPionEff->Fill(pBin,weight);
354 fKaonNot->Fill(pBin,pPid->GetProbability(AliPID::kKaon));
356 fProtNot->Fill(pBin,pPid->GetProbability(AliPID::kProton));
360 fKaonTot->Fill(pBin);
361 fKaonEff->Fill(pBin,pPid->GetProbability(AliPID::kKaon));
362 fPionNot->Fill(pBin,(pPid->GetProbability(AliPID::kPion)));
364 fProtNot->Fill(pBin,(pPid->GetProbability(AliPID::kProton)));
368 fVar[0] = track->GetHMPIDcluIdx()/1000000;
370 fVar[2] = (Float_t)track->P();
375 fVar[7] = (Float_t)thetaCkov;
379 fVar[11] = (Float_t)track->GetSign();
380 fVar[12] = (Float_t)nph;
381 fVar[13] = (Float_t)track->GetNcls(1);
382 fVar[14] = (Float_t)probs[0];
383 fVar[15] = (Float_t)probs[1];
384 fVar[16] = (Float_t)probs[2];
385 fVar[17] = (Float_t)probs[3];
386 fVar[18] = (Float_t)probs[4];
387 fVar[19] = (Float_t)track->GetTOFsignal();
388 fVar[20] = (Float_t)track->GetKinkIndex(0);
389 fVar[21] = (Float_t)track->Xv();
390 fVar[22] = (Float_t)track->Yv();
391 fVar[23] = (Float_t)track->Zv();
392 fVar[24] = (Float_t)track->GetTPCchi2();
395 fVar[27] = track->GetHMPIDcluIdx()%1000000/1000;
401 Float_t phiMin=(5./180.)*TMath::Pi() , phiMax=(55./180.)*TMath::Pi();
402 for (Int_t iPart=0; iPart<pStack->GetNprimary(); iPart++){
403 TString namepart=pStack->Particle(iPart)->GetName();
404 if (!namepart.Contains("proton")) continue;
405 pPart=pStack->Particle(iPart);
406 if (!pStack->IsPhysicalPrimary(iPart) || pPart->P()<1. || pPart->P()>5.) continue;
407 if (pPart->Phi()< phiMin || pPart->Phi()>phiMax) continue;
408 if (TMath::Abs(pPart->Eta()) > 0.5) continue;
409 Int_t pBin=(Int_t) (2*(pPart->P()-1));
410 if (namepart == "proton") fProtGen->Fill(pBin);
411 if (namepart == "antiproton") fPbarGen->Fill(pBin);
412 for (Int_t iTrack=0; iTrack<fESD->GetNumberOfTracks(); iTrack++){
413 if (fESD->GetTrack(iTrack)->GetLabel() != iPart) continue;
414 if (fESD->GetTrack(iTrack)->GetHMPIDsignal() < 0 ) continue;
415 if (namepart == "proton") fProtHmp->Fill(pBin);
416 if (namepart == "antiproton") fPbarHmp->Fill(pBin);
421 /* PostData(0) is taken care of by AliAnalysisTaskSE */
422 PostData(1,fHmpHistList);
426 //___________________________________________________________________________
427 void AliHMPIDAnalysisTask::Terminate(Option_t*)
429 // The Terminate() function is the last function to be called during
430 // a query. It always runs on the client, it can be used to present
431 // the results graphically or save the results to file.
433 Info("Terminate"," ");
437 fHmpHistList = dynamic_cast<TList*> (GetOutputData(1));
440 AliError("Histogram List is not available");
444 fPionEff = dynamic_cast<TH1F*> (fHmpHistList->FindObject("PionEff"));
445 fPionTot = dynamic_cast<TH1I*> (fHmpHistList->FindObject("PionTot"));
446 fPionNot = dynamic_cast<TH1F*> (fHmpHistList->FindObject("PionNot"));
447 fPionCon = dynamic_cast<TH1I*> (fHmpHistList->FindObject("PionCon"));
448 fKaonEff = dynamic_cast<TH1F*> (fHmpHistList->FindObject("KaonEff"));
449 fKaonTot = dynamic_cast<TH1I*> (fHmpHistList->FindObject("KaonTot"));
450 fKaonNot = dynamic_cast<TH1F*> (fHmpHistList->FindObject("KaonNot"));
451 fKaonCon = dynamic_cast<TH1I*> (fHmpHistList->FindObject("KaonCon"));
452 fProtEff = dynamic_cast<TH1F*> (fHmpHistList->FindObject("ProtEff"));
453 fProtTot = dynamic_cast<TH1I*> (fHmpHistList->FindObject("ProtTot"));
454 fProtNot = dynamic_cast<TH1F*> (fHmpHistList->FindObject("ProtNot"));
455 fProtCon = dynamic_cast<TH1I*> (fHmpHistList->FindObject("ProtCon"));
457 Float_t *pionEff = fPionEff->GetArray();
458 Int_t *pionTot = fPionTot->GetArray();
459 Float_t *pionNot = fPionNot->GetArray();
460 Int_t *pionCon = fPionCon->GetArray();
461 Float_t *kaonEff = fKaonEff->GetArray();
462 Int_t *kaonTot = fKaonTot->GetArray();
463 Float_t *kaonNot = fKaonNot->GetArray();
464 Int_t *kaonCon = fKaonCon->GetArray();
465 Float_t *protEff = fProtEff->GetArray();
466 Int_t *protTot = fProtTot->GetArray();
467 Float_t *protNot = fProtNot->GetArray();
468 Int_t *protCon = fProtCon->GetArray();
470 TGraphErrors *effPi = new TGraphErrors(fN1);
471 TGraphErrors *effKa = new TGraphErrors(fN1);
472 TGraphErrors *effPr = new TGraphErrors(fN2);
473 TGraphErrors *conPi = new TGraphErrors(fN1);
474 TGraphErrors *conKa = new TGraphErrors(fN1);
475 TGraphErrors *conPr = new TGraphErrors(fN2);
477 Float_t pt[8]={1.25,1.75,2.25,2.75,3.25,3.75,4.25,4.75};
478 Float_t eff=0, effErr=0, con=0, conErr=0;
479 for (Int_t i=0; i< fN2; i++){
480 eff = protEff[i+1]/TMath::Max(protTot[i+1],1);
481 effErr = TMath::Sqrt(eff*(1.-eff)/TMath::Max(protTot[i+1],1));
482 con = protNot[i+1]/TMath::Max(protCon[i+1],1);
483 conErr = TMath::Sqrt(con*(1.-con)/protCon[i+1]);
484 effPr->SetPoint(i,pt[i],eff);
485 effPr->SetPointError(i,0,effErr);
486 conPr->SetPoint(i,pt[i],con);
487 conPr->SetPointError(i,0,conErr);
489 if (i>=fN1) continue;
490 eff = pionEff[i+1]/pionTot[i+1];
491 effErr = TMath::Sqrt(eff*(1.-eff)/pionTot[i+1]);
492 con=pionNot[i+1]/pionCon[i+1];
493 conErr = TMath::Sqrt(con*(1.-con)/pionCon[i+1]);
494 effPi->SetPoint(i,pt[i],(Float_t)pionEff[i+1]/(Float_t)pionTot[i+1]);
495 effPi->SetPointError(i,0,effErr);
496 conPi->SetPoint(i,pt[i],(Float_t)pionNot[i+1]/(Float_t)pionCon[i+1]);
497 conPi->SetPointError(i,0,conErr);
499 eff = kaonEff[i+1]/TMath::Max(kaonTot[i+1],1);
500 effErr = TMath::Sqrt(eff*(1.-eff)/kaonTot[i+1]);
501 con = kaonNot[i+1]/TMath::Max(kaonCon[i+1],1);
502 conErr = TMath::Sqrt(con*(1.-con)/kaonCon[i+1]);
503 effKa->SetPoint(i,pt[i],(Float_t)kaonEff[i+1]/TMath::Max(kaonTot[i+1],1));
504 effKa->SetPointError(i,0,effErr);
505 conKa->SetPoint(i,pt[i],(Float_t)kaonNot[i+1]/TMath::Max(kaonCon[i+1],1));
506 conKa->SetPointError(i,0,conErr);
509 fProtGen = dynamic_cast<TH1I*> (fHmpHistList->FindObject("ProtGen"));
510 fPbarGen = dynamic_cast<TH1I*> (fHmpHistList->FindObject("PbarGen"));
511 fProtHmp = dynamic_cast<TH1I*> (fHmpHistList->FindObject("ProtHmp"));
512 fPbarHmp = dynamic_cast<TH1I*> (fHmpHistList->FindObject("PbarHmp"));
514 Int_t *protGen = fProtGen->GetArray();
515 Int_t *pbarGen = fPbarGen->GetArray();
516 Int_t *protHmp = fProtHmp->GetArray();
517 Int_t *pbarHmp = fPbarHmp->GetArray();
519 TGraphErrors *ratPr = new TGraphErrors(fN2);
520 TGraphErrors *ratPb = new TGraphErrors(fN2);
522 Float_t rat1=0, rat1Err=0, rat2=0, rat2Err=0;
523 for (Int_t i=0; i< fN2; i++){
524 rat1 = (Float_t)protHmp[i+1]/TMath::Max(protGen[i+1],1);
525 rat1Err = TMath::Sqrt(rat1*(1.-rat1)/TMath::Max(protGen[i+1],1));
526 rat2 = (Float_t)pbarHmp[i+1]/TMath::Max(pbarGen[i+1],1);
527 rat2Err = TMath::Sqrt(rat2*(1.-rat2)/TMath::Max(pbarGen[i+1],1));
528 ratPr->SetPoint(i,pt[i],rat1);
529 ratPr->SetPointError(i,0,rat1Err);
530 ratPb->SetPoint(i,pt[i],rat2);
531 ratPb->SetPointError(i,0,rat2Err);
534 TCanvas *pCan=new TCanvas("Hmp","Efficiency and contamination",500,900);
538 effPi->SetTitle("Pions");
539 effPi->GetXaxis()->SetTitle("p_{T} (GeV/c)");
540 effPi->GetYaxis()->SetRangeUser(0.,1.);
541 effPi->SetMarkerStyle(20);
543 conPi->SetMarkerStyle(21);
544 conPi->Draw("sameLP");
547 effKa->SetTitle("Kaons");
548 effKa->GetXaxis()->SetTitle("p_{T} (GeV/c)");
549 effKa->GetYaxis()->SetRangeUser(0.,1.);
550 effKa->SetMarkerStyle(20);
552 conKa->SetMarkerStyle(21);
553 conKa->Draw("sameLP");
556 effPr->SetTitle("Protons");
557 effPr->GetXaxis()->SetTitle("p_{T} (GeV/c)");
558 effPr->GetYaxis()->SetRangeUser(0.,1.);
559 effPr->SetMarkerStyle(20);
561 conPr->SetMarkerStyle(21);
562 conPr->Draw("sameLP");
565 ratPr->SetTitle("Ratios");
566 ratPr->GetXaxis()->SetTitle("p_{T} (GeV/c)");
567 ratPr->GetYaxis()->SetRangeUser(0.,1.);
568 ratPr->SetMarkerStyle(20);
570 ratPb->SetMarkerStyle(21);
571 ratPb->Draw("sameLP");
573 TFile *outFile = new TFile("HmpidGraphs.root","recreate");
577 AliAnalysisTaskSE::Terminate();
581 //___________________________________________________________________________
582 void AliHMPIDAnalysisTask::UserCreateOutputObjects() {
584 //HERE ONE CAN CREATE OUTPUT OBJECTS
585 //TO BE SET BEFORE THE EXECUTION OF THE TASK
590 fHmpHistList = new TList();
591 fHmpHistList->SetOwner();
593 fHmpPesdPhmp = new TH2F("fHmpPesdPhmp","HMPID: ESD p - running p;HMP p (GeV/c);ESD p (Gev/c)",100,0,10,100,0,10);
594 fHmpHistList->Add(fHmpPesdPhmp);
596 fHmpCkovPesd = new TH2F("fHmpCkovPesd","HMPID: ThetaCherenkov vs P;p_esd (GeV/c);#Theta_C;Entries",100,0,10,110,0,1.1);
597 fHmpHistList->Add(fHmpCkovPesd);
599 fHmpCkovPhmp = new TH2F("fHmpCkovPhmp","HMPID: ThetaCherenkov vs P;p_hmp (GeV/c);#Theta_C;Entries",100,0,10,110,0,1.1);
600 fHmpHistList->Add(fHmpCkovPhmp);
602 fHmpMipTrkDist = new TH1F("fHmpMipTrkDist","HMPID MIP-Track distance;distance (cm);Entries",400,0,20);
603 fHmpHistList->Add(fHmpMipTrkDist);
605 fHmpMipTrkDistX = new TH1F("fHmpMipTrkDistX","HMPID MIP-Track distance in local X;distance (cm);Entries",800,-20,20);
606 fHmpHistList->Add(fHmpMipTrkDistX);
608 fHmpMipTrkDistY = new TH1F("fHmpMipTrkDistY","HMPID MIP-Track distance in local Y;distance (cm);Entries",800,-20,20);
609 fHmpHistList->Add(fHmpMipTrkDistY);
611 fHmpMipCharge3cm = new TH1F("fHmpMipCharge3cm","HMPID MIP Charge;MIP Charge (ADC);Entries",5001,-0.5,5000.5);
612 fHmpHistList->Add(fHmpMipCharge3cm);
614 fHmpMipCharge1cm = new TH1F("fHmpMipCharge1cm","HMPID MIP Charge;MIP Charge (ADC);Entries",5001,-0.5,5000.5);
615 fHmpHistList->Add(fHmpMipCharge1cm);
617 fHmpNumPhots = new TH1F("fHmpNumPhots","HMPID Number of photon clusters on ring;#photon clus.;Entries",51,-0.5,50.5);
618 fHmpHistList->Add(fHmpNumPhots);
620 fHmpTrkFlags = new TH1F("fHmpTrkFlags","HMPID track flags",6,0,6);
621 TString summary[6] = {"NotPerformed","MipDistCut", "MipQdcCut", "NoPhotAccept", "kNoRad", "other"};
622 for(Int_t ibin = 0; ibin < 6; ibin++) fHmpTrkFlags->GetXaxis()->SetBinLabel(ibin+1,Form("%i %s",ibin+1,summary[ibin].Data()));
623 fHmpHistList->Add(fHmpTrkFlags);
625 fPionEff = new TH1F("PionEff","Identified pions",fN1,0,fN1);
626 fKaonEff = new TH1F("KaonEff","Identified kaons",fN1,0,fN1);
627 fProtEff = new TH1F("ProtEff","Identified protons",fN2,0,fN2);
628 fPionTot = new TH1I("PionTot","Total MC pions",fN1,0,fN1);
629 fKaonTot = new TH1I("KaonTot","Total MC kaons",fN1,0,fN1);
630 fProtTot = new TH1I("ProtTot","Total MC protons",fN2,0,fN2);
631 fPionNot = new TH1F("PionNot","Misidentified pions",fN1,0,fN1);
632 fKaonNot = new TH1F("KaonNot","Misidentified kaons",fN1,0,fN1);
633 fProtNot = new TH1F("ProtNot","Misidentified protons",fN2,0,fN2);
634 fPionCon = new TH1I("PionCon","Total not MC pions",fN1,0,fN1);
635 fKaonCon = new TH1I("KaonCon","Total not MC kaons",fN1,0,fN1);
636 fProtCon = new TH1I("ProtCon","Total not MC protons",fN2,0,fN2);
638 fHmpHistList->Add(fPionEff); fHmpHistList->Add(fKaonEff); fHmpHistList->Add(fProtEff);
639 fHmpHistList->Add(fPionTot); fHmpHistList->Add(fKaonTot); fHmpHistList->Add(fProtTot);
640 fHmpHistList->Add(fPionNot); fHmpHistList->Add(fKaonNot); fHmpHistList->Add(fProtNot);
641 fHmpHistList->Add(fPionCon); fHmpHistList->Add(fKaonCon); fHmpHistList->Add(fProtCon);
643 fProtGen = new TH1I("ProtGen","Generated protons",fN2,0,fN2);
644 fPbarGen = new TH1I("PbarGen","Generated antiprotons",fN2,0,fN2);
645 fProtHmp = new TH1I("ProtHmp","HMPID protons",fN2,0,fN2);
646 fPbarHmp = new TH1I("PbarHmp","HMPID antiprotons",fN2,0,fN2);
648 fHmpHistList->Add(fProtGen); fHmpHistList->Add(fPbarGen);
649 fHmpHistList->Add(fProtHmp); fHmpHistList->Add(fPbarHmp);
651 fThetavsPiFromK= new TH2F("ThetavsPiFromK","Theta vs p of pions from K;p_esd (GeV/c);#Theta_C",100,0,10,110,0,1.1);
652 fHmpHistList->Add(fThetavsPiFromK);
654 fThetapivsPesd = new TH2F("ThetapivsPesd","Theta of pions vs p of esd;p_esd (GeV/c);#Theta_C",100,0,10,110,0,1.1);
655 fHmpHistList->Add(fThetapivsPesd);
657 fThetaKvsPesd = new TH2F("ThetaKvsPesd","Theta of kaons vs p of esd;p_esd (GeV/c);#Theta_C",100,0,10,110,0,1.1);
658 fHmpHistList->Add(fThetaKvsPesd);
660 fThetaPvsPesd = new TH2F("ThetaPvsPesd","Theta of protons vs p of esd;p_esd (GeV/c);#Theta_C",100,0,10,110,0,1.1);
661 fHmpHistList->Add(fThetaPvsPesd);
664 fTree = new TTree("Tree","Tree with data");
665 fTree->Branch("Chamber",&fVar[0]);
666 fTree->Branch("pHmp3",&fVar[1]);
667 fTree->Branch("P",&fVar[2]);
668 fTree->Branch("Xpc",&fVar[3]);
669 fTree->Branch("Ypc",&fVar[4]);
670 fTree->Branch("X",&fVar[5]);
671 fTree->Branch("Y",&fVar[6]);
672 fTree->Branch("HMPIDsignal",&fVar[7]);
673 fTree->Branch("Charge",&fVar[8]);
674 fTree->Branch("Theta",&fVar[9]);
675 fTree->Branch("Phi",&fVar[10]);
676 fTree->Branch("Sign",&fVar[11]);
677 fTree->Branch("NumPhotons",&fVar[12]);
678 fTree->Branch("NumTPCclust",&fVar[13]);
679 fTree->Branch("Prob0",&fVar[14]);
680 fTree->Branch("Prob1",&fVar[15]);
681 fTree->Branch("Prob2",&fVar[16]);
682 fTree->Branch("Prob3",&fVar[17]);
683 fTree->Branch("Prob4",&fVar[18]);
684 fTree->Branch("TOFsignal",&fVar[19]);
685 fTree->Branch("KinkIndex",&fVar[20]);
686 fTree->Branch("Xv",&fVar[21]);
687 fTree->Branch("Yv",&fVar[22]);
688 fTree->Branch("Zv",&fVar[23]);
689 fTree->Branch("TPCchi2",&fVar[24]);
690 fTree->Branch("b0",&fVar[25]);
691 fTree->Branch("b1",&fVar[26]);
692 fTree->Branch("ClustSize",&fVar[27]);
694 PostData(1,fHmpHistList);
698 //____________________________________________________________________________________________________________________________________
699 Bool_t AliHMPIDAnalysisTask::Equal(Double_t x, Double_t y, Double_t tolerance)
701 return abs(x - y) <= tolerance ;