10 #include "TProfile2D.h"
12 #include "TGraphErrors.h"
15 #include <TClonesArray.h>
16 #include <TObjArray.h>
19 #include "AliESDEvent.h"
20 #include "AliESDInputHandler.h"
21 #include "AliTrackReference.h"
23 #include "AliAnalysisTask.h"
25 #include "AliTRDtrackerV1.h"
26 #include "AliTRDtrackV1.h"
27 #include "AliTRDcluster.h"
28 #include "AliTRDReconstructor.h"
29 #include "AliCDBManager.h"
30 #include "AliTRDpidUtil.h"
32 #include "AliTRDcheckPID.h"
33 #include "info/AliTRDtrackInfo.h"
35 // calculate pion efficiency at 90% electron efficiency for 11 momentum bins
36 // this task should be used with simulated data only
38 ClassImp(AliTRDcheckPID)
40 //________________________________________________________________________
41 AliTRDcheckPID::AliTRDcheckPID()
42 :AliTRDrecoTask("checkPID", "TRD PID checker")
48 ,fMinNTracklets(AliTRDgeometry::kNlayer)
49 ,fMaxNTracklets(AliTRDgeometry::kNlayer)
52 // Default constructor
55 fReconstructor = new AliTRDReconstructor();
56 fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
58 // Initialize momentum axis with default values
59 Double_t defaultMomenta[AliTRDCalPID::kNMom+1];
60 for(Int_t imom = 0; imom < AliTRDCalPID::kNMom+1; imom++)
61 defaultMomenta[imom] = AliTRDCalPID::GetMomentumBinning(imom);
62 SetMomentumBinning(AliTRDCalPID::kNMom, defaultMomenta);
64 fUtil = new AliTRDpidUtil();
69 //________________________________________________________________________
70 AliTRDcheckPID::~AliTRDcheckPID()
72 if(fGraph){fGraph->Delete(); delete fGraph;}
73 if(fReconstructor) delete fReconstructor;
74 if(fUtil) delete fUtil;
78 //________________________________________________________________________
79 void AliTRDcheckPID::CreateOutputObjects()
84 OpenFile(0, "RECREATE");
85 fContainer = Histos();
89 //_______________________________________________________
90 TObjArray * AliTRDcheckPID::Histos(){
93 // Create QA histograms
95 if(fContainer) return fContainer;
97 Int_t xBins = AliPID::kSPECIES*fMomentumAxis->GetNbins();
98 fContainer = new TObjArray(); fContainer->Expand(kNPlots);
100 const Float_t epsilon = 1./(2*(AliTRDpidUtil::kBins-1)); // get nice histos with bin center at 0 and 1
102 // histos of the electron probability of all 5 particle species and 11 momenta for the 2-dim LQ method
103 fEfficiency = new TObjArray(3);
104 fEfficiency->SetOwner(); fEfficiency->SetName("Efficiency");
105 fContainer->AddAt(fEfficiency, kEfficiency);
108 if(!(h = (TH2F*)gROOT->FindObject("PID_LQ"))){
109 h = new TH2F("PID_LQ", "", xBins, -0.5, xBins - 0.5,
110 AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon);
112 fEfficiency->AddAt(h, AliTRDpidUtil::kLQ);
114 // histos of the electron probability of all 5 particle species and 11 momenta for the neural network method
115 if(!(h = (TH2F*)gROOT->FindObject("PID_NN"))){
116 h = new TH2F("PID_NN", "",
117 xBins, -0.5, xBins - 0.5,
118 AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon);
120 fEfficiency->AddAt(h, AliTRDpidUtil::kNN);
122 // histos of the electron probability of all 5 particle species and 11 momenta for the ESD output
123 if(!(h = (TH2F*)gROOT->FindObject("PID_ESD"))){
124 h = new TH2F("PID_ESD", "",
125 xBins, -0.5, xBins - 0.5,
126 AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon);
128 fEfficiency->AddAt(h, AliTRDpidUtil::kESD);
130 // histos of the dE/dx distribution for all 5 particle species and 11 momenta
131 if(!(h = (TH2F*)gROOT->FindObject("dEdx"))){
132 h = new TH2F("dEdx", "",
133 xBins, -0.5, xBins - 0.5,
136 fContainer->AddAt(h, kdEdx);
138 // histos of the dE/dx slices for all 5 particle species and 11 momenta
139 if(!(h = (TH2F*)gROOT->FindObject("dEdxSlice"))){
140 h = new TH2F("dEdxSlice", "",
141 xBins*AliTRDpidUtil::kLQslices, -0.5, xBins*AliTRDpidUtil::kLQslices - 0.5,
144 fContainer->AddAt(h, kdEdxSlice);
146 // histos of the pulse height distribution for all 5 particle species and 11 momenta
147 if(!(h = (TH2F*)gROOT->FindObject("PH"))){
148 h = new TProfile2D("PH", "",
149 xBins, -0.5, xBins - 0.5,
150 AliTRDtrackerV1::GetNTimeBins(), -0.5, AliTRDtrackerV1::GetNTimeBins() - 0.5);
152 fContainer->AddAt(h, kPH);
154 // histos of the number of clusters distribution for all 5 particle species and 11 momenta
155 if(!(h = (TH2F*)gROOT->FindObject("NClus"))){
156 h = new TH2F("NClus", "",
157 xBins, -0.5, xBins - 0.5,
158 AliTRDtrackerV1::GetNTimeBins(), -0.5, AliTRDtrackerV1::GetNTimeBins() - 0.5);
160 fContainer->AddAt(h, kNClus);
163 // momentum distributions - absolute and in momentum bins
164 if(!(h = (TH1F*)gROOT->FindObject("hMom"))){
165 h = new TH1F("hMom", "momentum distribution", fMomentumAxis->GetNbins(), fMomentumAxis->GetXmin(), fMomentumAxis->GetXmax());
167 fContainer->AddAt(h, kMomentum);
169 if(!(h = (TH1F*)gROOT->FindObject("hMomBin"))){
170 h = new TH1F("hMomBin", "momentum distribution in momentum bins", fMomentumAxis->GetNbins(), fMomentumAxis->GetXmin(), fMomentumAxis->GetXmax());
172 fContainer->AddAt(h, kMomentumBin);
174 // Number of tracklets per track
175 if(!(h = (TH2F*)gROOT->FindObject("nTracklets"))){
176 h = new TH2F("nTracklets", "",
177 xBins, -0.5, xBins - 0.5,
178 AliTRDgeometry::kNlayer, 0.5, AliTRDgeometry::kNlayer+.5);
180 fContainer->AddAt(h, kNTracklets);
186 //________________________________________________________________________
187 Bool_t AliTRDcheckPID::CheckTrackQuality(const AliTRDtrackV1* track)
190 // Check if the track is ok for PID
193 Int_t ntracklets = track->GetNumberOfTracklets();
194 if(ntracklets >= fMinNTracklets && ntracklets <= fMaxNTracklets) return 1;
201 //________________________________________________________________________
202 Int_t AliTRDcheckPID::CalcPDG(AliTRDtrackV1* track)
205 track -> SetReconstructor(fReconstructor);
207 fReconstructor -> SetOption("nn");
210 if(track -> GetPID(AliPID::kElectron) > track -> GetPID(AliPID::kMuon) + track -> GetPID(AliPID::kPion) + track -> GetPID(AliPID::kKaon) + track -> GetPID(AliPID::kProton)){
213 else if(track -> GetPID(kProton) > track -> GetPID(AliPID::kPion) && track -> GetPID(AliPID::kProton) > track -> GetPID(AliPID::kKaon) && track -> GetPID(AliPID::kProton) > track -> GetPID(AliPID::kMuon)){
216 else if(track -> GetPID(AliPID::kKaon) > track -> GetPID(AliPID::kMuon) && track -> GetPID(AliPID::kKaon) > track -> GetPID(AliPID::kPion)){
219 else if(track -> GetPID(AliPID::kMuon) > track -> GetPID(AliPID::kPion)){
228 //_______________________________________________________
229 TH1 *AliTRDcheckPID::PlotLQ(const AliTRDtrackV1 *track)
232 // Plot the probabilities for electrons using 2-dim LQ
236 AliWarning("No ESD info available.");
240 if(track) fTrack = track;
242 AliWarning("No Track defined.");
247 status = fESD -> GetStatus();
248 if(!(status&AliESDtrack::kTRDin)) return 0x0;
250 if(!CheckTrackQuality(fTrack)) return 0x0;
252 if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){
253 AliWarning("No Histogram defined.");
257 if(!(hPIDLQ = dynamic_cast<TH2F *>(fEfficiency->At(AliTRDpidUtil::kLQ)))){
258 AliWarning("No Histogram defined.");
262 AliTRDtrackV1 cTrack(*fTrack);
263 cTrack.SetReconstructor(fReconstructor);
266 Float_t momentum = 0.;
269 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
272 //AliWarning("No MC info available!");
273 momentum = cTrack.GetMomentum(0);
274 pdg = CalcPDG(&cTrack);
276 if(!IsInRange(momentum)) return 0x0;
278 fReconstructor -> SetOption("!nn");
280 if(cTrack.GetNumberOfTrackletsPID() < fMinNTracklets) return 0x0;
282 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
283 hPIDLQ -> Fill(FindBin(species, momentum), cTrack.GetPID(AliPID::kElectron));
288 //_______________________________________________________
289 TH1 *AliTRDcheckPID::PlotNN(const AliTRDtrackV1 *track)
292 // Plot the probabilities for electrons using 2-dim LQ
296 AliWarning("No ESD info available.");
300 if(track) fTrack = track;
302 AliWarning("No Track defined.");
307 status = fESD -> GetStatus();
308 if(!(status&AliESDtrack::kTRDin)) return 0x0;
310 if(!CheckTrackQuality(fTrack)) return 0x0;
312 if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){
313 AliWarning("No Histogram defined.");
317 if(!(hPIDNN = dynamic_cast<TH2F *>(fEfficiency->At(AliTRDpidUtil::kNN)))){
318 AliWarning("No Histogram defined.");
322 AliTRDtrackV1 cTrack(*fTrack);
323 cTrack.SetReconstructor(fReconstructor);
326 Float_t momentum = 0.;
328 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
331 //AliWarning("No MC info available!");
332 momentum = cTrack.GetMomentum(0);
333 pdg = CalcPDG(&cTrack);
335 if(!IsInRange(momentum)) return 0x0;
337 fReconstructor -> SetOption("nn");
339 if(cTrack.GetNumberOfTrackletsPID() < fMinNTracklets) return 0x0;
341 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
342 hPIDNN -> Fill(FindBin(species, momentum), cTrack.GetPID(AliPID::kElectron));
347 //_______________________________________________________
348 TH1 *AliTRDcheckPID::PlotESD(const AliTRDtrackV1 *track)
351 // Plot the probabilities for electrons using 2-dim LQ
355 AliWarning("No ESD info available.");
359 if(track) fTrack = track;
361 AliWarning("No Track defined.");
366 status = fESD -> GetStatus();
367 if(!(status&AliESDtrack::kTRDin)) return 0x0;
369 if(!CheckTrackQuality(fTrack)) return 0x0;
370 if(fTrack->GetNumberOfTrackletsPID() < fMinNTracklets) return 0x0;
372 if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){
373 AliWarning("No Histogram defined.");
377 if(!(hPIDESD = dynamic_cast<TH2F *>(fEfficiency->At(AliTRDpidUtil::kESD)))){
378 AliWarning("No Histogram defined.");
384 Float_t momentum = 0.;
386 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
389 //AliWarning("No MC info available!");
390 AliTRDtrackV1 cTrack(*fTrack);
391 cTrack.SetReconstructor(fReconstructor);
392 momentum = cTrack.GetMomentum(0);
393 pdg = CalcPDG(&cTrack);
395 if(!IsInRange(momentum)) return 0x0;
397 // Double32_t pidESD[AliPID::kSPECIES];
398 const Double32_t *pidESD = fESD->GetResponseIter();
399 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
400 hPIDESD -> Fill(FindBin(species, momentum), pidESD[0]);
405 //_______________________________________________________
406 TH1 *AliTRDcheckPID::PlotdEdx(const AliTRDtrackV1 *track)
409 // Plot the probabilities for electrons using 2-dim LQ
412 if(track) fTrack = track;
414 AliWarning("No Track defined.");
418 if(!CheckTrackQuality(fTrack)) return 0x0;
421 if(!(hdEdx = dynamic_cast<TH2F *>(fContainer->At(kdEdx)))){
422 AliWarning("No Histogram defined.");
426 AliTRDtrackV1 cTrack(*fTrack);
427 cTrack.SetReconstructor(fReconstructor);
429 Float_t momentum = 0.;
431 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
434 //AliWarning("No MC info available!");
435 momentum = cTrack.GetMomentum(0);
436 pdg = CalcPDG(&cTrack);
438 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
439 if(!IsInRange(momentum)) return 0x0;
441 fReconstructor -> SetOption("!nn");
443 Int_t iBin = FindBin(species, momentum);
444 AliTRDseedV1 *tracklet = 0x0;
445 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
447 tracklet = cTrack.GetTracklet(iChamb);
448 if(!tracklet) continue;
449 tracklet -> CookdEdx(AliTRDpidUtil::kLQslices);
450 for(Int_t i = 3; i--;) SumdEdx += tracklet->GetdEdx()[i];
451 hdEdx -> Fill(iBin, SumdEdx);
457 //_______________________________________________________
458 TH1 *AliTRDcheckPID::PlotdEdxSlice(const AliTRDtrackV1 *track)
461 // Plot the probabilities for electrons using 2-dim LQ
464 if(track) fTrack = track;
466 AliWarning("No Track defined.");
470 if(!CheckTrackQuality(fTrack)) return 0x0;
473 if(!(hdEdxSlice = dynamic_cast<TH2F *>(fContainer->At(kdEdxSlice)))){
474 AliWarning("No Histogram defined.");
478 AliTRDtrackV1 cTrack(*fTrack);
479 cTrack.SetReconstructor(fReconstructor);
481 Float_t momentum = 0.;
483 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
486 //AliWarning("No MC info available!");
487 momentum = cTrack.GetMomentum(0);
488 pdg = CalcPDG(&cTrack);
490 if(!IsInRange(momentum)) return 0x0;
492 fReconstructor -> SetOption("!nn");
493 Int_t iMomBin = fMomentumAxis->FindBin(momentum);
494 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
496 AliTRDseedV1 *tracklet = 0x0;
497 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
498 tracklet = cTrack.GetTracklet(iChamb);
499 if(!tracklet) continue;
500 tracklet -> CookdEdx(AliTRDpidUtil::kLQslices);
501 fdEdx = tracklet->GetdEdx();
502 for(Int_t iSlice = 0; iSlice < AliTRDpidUtil::kLQslices; iSlice++){
503 hdEdxSlice -> Fill(species * fMomentumAxis->GetNbins() * AliTRDpidUtil::kLQslices + (iMomBin-1) * AliTRDpidUtil::kLQslices + iSlice, fdEdx[iSlice]);
511 //_______________________________________________________
512 TH1 *AliTRDcheckPID::PlotPH(const AliTRDtrackV1 *track)
515 // Plot the probabilities for electrons using 2-dim LQ
518 if(track) fTrack = track;
520 AliWarning("No Track defined.");
524 if(!CheckTrackQuality(fTrack)) return 0x0;
527 if(!(hPH = dynamic_cast<TProfile2D *>(fContainer->At(kPH)))){
528 AliWarning("No Histogram defined.");
533 Float_t momentum = 0.;
535 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
538 //AliWarning("No MC info available!");
539 AliTRDtrackV1 cTrack(*fTrack);
540 cTrack.SetReconstructor(fReconstructor);
541 momentum = cTrack.GetMomentum(0);
542 pdg = CalcPDG(&cTrack);
544 if(!IsInRange(momentum)) return 0x0;;
546 AliTRDseedV1 *tracklet = 0x0;
547 AliTRDcluster *TRDcluster = 0x0;
548 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
549 Int_t iBin = FindBin(species, momentum);
550 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
551 tracklet = fTrack->GetTracklet(iChamb);
552 if(!tracklet) continue;
553 for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){
554 if(!(TRDcluster = tracklet->GetClusters(iClus))) continue;
555 hPH -> Fill(iBin, TRDcluster->GetLocalTimeBin(), tracklet->GetdQdl(iClus));
562 //_______________________________________________________
563 TH1 *AliTRDcheckPID::PlotNClus(const AliTRDtrackV1 *track)
566 // Plot the probabilities for electrons using 2-dim LQ
569 if(track) fTrack = track;
571 AliWarning("No Track defined.");
575 if(!CheckTrackQuality(fTrack)) return 0x0;
578 if(!(hNClus = dynamic_cast<TH2F *>(fContainer->At(kNClus)))){
579 AliWarning("No Histogram defined.");
585 Float_t momentum = 0.;
587 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
590 //AliWarning("No MC info available!");
591 AliTRDtrackV1 cTrack(*fTrack);
592 cTrack.SetReconstructor(fReconstructor);
593 momentum = cTrack.GetMomentum(0);
594 pdg = CalcPDG(&cTrack);
596 if(!IsInRange(momentum)) return 0x0;
598 Int_t species = AliTRDpidUtil::AliTRDpidUtil::Pdg2Pid(pdg);
599 Int_t iBin = FindBin(species, momentum);
600 AliTRDseedV1 *tracklet = 0x0;
601 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
602 tracklet = fTrack->GetTracklet(iChamb);
603 if(!tracklet) continue;
604 hNClus -> Fill(iBin, tracklet->GetN());
610 //_______________________________________________________
611 TH1 *AliTRDcheckPID::PlotNTracklets(const AliTRDtrackV1 *track)
614 // Plot the probabilities for electrons using 2-dim LQ
617 if(track) fTrack = track;
619 AliWarning("No Track defined.");
624 if(!(hTracklets = dynamic_cast<TH2F *>(fContainer->At(kNTracklets)))){
625 AliWarning("No Histogram defined.");
629 AliTRDtrackV1 cTrack(*fTrack);
630 cTrack.SetReconstructor(fReconstructor);
632 Float_t momentum = 0.;
634 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
637 //AliWarning("No MC info available!");
638 momentum = cTrack.GetMomentum(0);
639 pdg = CalcPDG(&cTrack);
641 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
642 if(!IsInRange(momentum)) return 0x0;
644 Int_t iBin = FindBin(species, momentum);
645 hTracklets -> Fill(iBin, cTrack.GetNumberOfTracklets());
649 //_______________________________________________________
650 TH1 *AliTRDcheckPID::PlotMom(const AliTRDtrackV1 *track)
653 // Plot the probabilities for electrons using 2-dim LQ
656 if(track) fTrack = track;
658 AliWarning("No Track defined.");
662 if(!CheckTrackQuality(fTrack)) return 0x0;
665 if(!(hMom = dynamic_cast<TH1F *>(fContainer->At(kMomentum)))){
666 AliWarning("No Histogram defined.");
672 Float_t momentum = 0.;
674 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
677 //AliWarning("No MC info available!");
678 AliTRDtrackV1 cTrack(*fTrack);
679 cTrack.SetReconstructor(fReconstructor);
680 momentum = cTrack.GetMomentum(0);
681 pdg = CalcPDG(&cTrack);
683 if(IsInRange(momentum)) hMom -> Fill(momentum);
688 //_______________________________________________________
689 TH1 *AliTRDcheckPID::PlotMomBin(const AliTRDtrackV1 *track)
692 // Plot the probabilities for electrons using 2-dim LQ
695 if(track) fTrack = track;
697 AliWarning("No Track defined.");
701 if(!CheckTrackQuality(fTrack)) return 0x0;
704 if(!(hMomBin = dynamic_cast<TH1F *>(fContainer->At(kMomentumBin)))){
705 AliWarning("No Histogram defined.");
711 Float_t momentum = 0.;
714 if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
717 //AliWarning("No MC info available!");
718 AliTRDtrackV1 cTrack(*fTrack);
719 cTrack.SetReconstructor(fReconstructor);
720 momentum = cTrack.GetMomentum(0);
722 if(IsInRange(momentum)) hMomBin -> Fill(fMomentumAxis->FindBin(momentum));
727 //________________________________________________________
728 Bool_t AliTRDcheckPID::GetRefFigure(Int_t ifig)
730 Bool_t FIRST = kTRUE;
731 TGraphErrors *g = 0x0;
733 TH1 *h1 = 0x0, *h=0x0;
735 TList *content = 0x0;
738 gPad->Divide(2, 1, 0.,0.);
739 TList *l=gPad->GetListOfPrimitives();
740 TVirtualPad *pad = ((TVirtualPad*)l->At(0));pad->cd();
741 pad->SetMargin(0.07, 0.07, 0.1, 0.001);
743 TLegend *legEff = new TLegend(.7, .7, .98, .98);
744 legEff->SetBorderSize(1);
745 content = (TList *)fGraph->FindObject("Efficiencies");
746 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kLQ))) break;
747 if(!g->GetN()) break;
748 legEff->SetHeader("PID Method");
750 ax = g->GetHistogram()->GetXaxis();
751 ax->SetTitle("p (GeV/c)");
752 ax->SetRangeUser(.5, 11.);
753 ax->SetMoreLogLabels();
754 ax = g->GetHistogram()->GetYaxis();
755 ax->SetTitle("#epsilon_{#pi}");
756 ax->SetRangeUser(1.e-3, 1.e-1);
757 legEff->AddEntry(g, "2D LQ", "pl");
758 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kNN))) break;
760 legEff->AddEntry(g, "NN", "pl");
761 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kESD))) break;
763 legEff->AddEntry(g, "ESD", "pl");
770 pad = ((TVirtualPad*)l->At(1));pad->cd();
771 pad->SetMargin(0.07, 0.07, 0.1, 0.001);
772 content = (TList *)fGraph->FindObject("Thresholds");
773 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kLQ))) break;
774 if(!g->GetN()) break;
776 ax = g->GetHistogram()->GetXaxis();
777 ax->SetTitle("p (GeV/c)");
778 ax->SetRangeUser(.5, 11.5);
779 ax->SetMoreLogLabels();
780 ax = g->GetHistogram()->GetYaxis();
781 ax->SetTitle("threshold (%)");
782 ax->SetRangeUser(5.e-2, 1.);
783 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kNN))) break;
785 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kESD))) break;
793 // save 2.0 GeV projection as reference
794 TLegend *legdEdx = new TLegend(.7, .7, .98, .98);
795 legdEdx->SetBorderSize(1);
797 if(!(h2 = (TH2F*)(fContainer->At(kdEdx)))) break;
798 legdEdx->SetHeader("Particle Species");
799 for(Int_t is = AliPID::kSPECIES-1; is>=0; is--){
800 Int_t bin = FindBin(is, 2.);
801 h1 = h2->ProjectionY("px", bin, bin);
802 if(!h1->GetEntries()) continue;
803 h1->Scale(1./h1->Integral());
804 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
806 h1->GetXaxis()->SetTitle("dE/dx (a.u.)");
807 h1->GetYaxis()->SetTitle("<Entries>");
809 h = (TH1F*)h1->DrawClone(FIRST ? "c" : "samec");
810 legdEdx->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l");
824 // save 2.0 GeV projection as reference
825 TLegend *legPH = new TLegend(.4, .7, .68, .98);
826 legPH->SetBorderSize(1);
828 if(!(h2 = (TH2F*)(fContainer->At(kPH)))) break;;
829 legPH->SetHeader("Particle Species");
830 for(Int_t is=0; is<AliPID::kSPECIES; is++){
831 Int_t bin = FindBin(is, 2.);
832 h1 = h2->ProjectionY("py", bin, bin);
833 if(!h1->GetEntries()) continue;
834 h1->SetMarkerStyle(24);
835 h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
836 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
838 h1->GetXaxis()->SetTitle("tb(1/100 ns^{-1})");
839 h1->GetYaxis()->SetTitle("<PH> (a.u.)");
841 h = (TH1F*)h1->DrawClone(FIRST ? "c" : "samec");
842 legPH->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "pl");
854 // save 2.0 GeV projection as reference
855 TLegend *legNClus = new TLegend(.4, .7, .68, .98);
856 legNClus->SetBorderSize(1);
858 if(!(h2 = (TH2F*)(fContainer->At(kNClus)))) break;
859 legNClus->SetHeader("Particle Species");
860 for(Int_t is=0; is<AliPID::kSPECIES; is++){
861 Int_t bin = FindBin(is, 2.);
862 h1 = h2->ProjectionY("py", bin, bin);
863 if(!h1->GetEntries()) continue;
864 h1->Scale(1./h1->Integral());
865 //h1->SetMarkerStyle(24);
866 //h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
867 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
868 if(FIRST) h1->GetXaxis()->SetTitle("N^{cl}/tracklet");
869 if(FIRST) h1->GetYaxis()->SetTitle("<Entries>");
870 h = (TH1F*)h1->DrawClone(FIRST ? "c" : "samec");
871 legNClus->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l");
886 TLegend *legNClus = new TLegend(.4, .7, .68, .98);
887 legNClus->SetBorderSize(1);
889 if(!(h2 = (TH2F*)(fContainer->At(kNTracklets)))) break;
890 legNClus->SetHeader("Particle Species");
891 for(Int_t is=0; is<AliPID::kSPECIES; is++){
892 Int_t bin = FindBin(is, 2.);
893 h1 = h2->ProjectionY("py", bin, bin);
894 if(!h1->GetEntries()) continue;
895 h1->Scale(1./h1->Integral());
896 //h1->SetMarkerStyle(24);
897 //h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
898 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
899 if(FIRST) h1->GetXaxis()->SetTitle("N^{tl}/track");
900 if(FIRST) h1->GetYaxis()->SetTitle("<Entries>");
901 h = (TH1F*)h1->DrawClone(FIRST ? "c" : "samec");
902 legNClus->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l");
914 AliInfo(Form("Reference plot [%d] missing result", ifig));
918 //________________________________________________________________________
919 Bool_t AliTRDcheckPID::PostProcess()
921 // Draw result to the screen
922 // Called once at the end of the query
925 Printf("ERROR: list not available");
928 if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){
929 AliError("Efficiency container missing.");
933 fGraph = new TObjArray(6);
935 EvaluatePionEfficiency(fEfficiency, fGraph, 0.9);
941 //________________________________________________________________________
942 void AliTRDcheckPID::EvaluatePionEfficiency(TObjArray *histoContainer, TObjArray *results, Float_t electron_efficiency){
943 fUtil->SetElectronEfficiency(electron_efficiency);
945 Color_t colors[3] = {kBlue, kGreen+2, kRed};
946 Int_t markerStyle[3] = {7, 7, 24};
947 TString MethodName[3] = {"LQ", "NN", "ESD"};
949 TGraphErrors *g, *gPtrEff[3], *gPtrThres[3];
950 TObjArray *eff = new TObjArray(3); eff->SetOwner(); eff->SetName("Efficiencies");
951 results->AddAt(eff, 0);
952 for(Int_t iMethod = 0; iMethod < 3; iMethod++){
953 eff->AddAt(g = gPtrEff[iMethod] = new TGraphErrors(), iMethod);
954 g->SetName(Form("efficiency_%s", MethodName[iMethod].Data()));
955 g->SetLineColor(colors[iMethod]);
956 g->SetMarkerColor(colors[iMethod]);
957 g->SetMarkerStyle(markerStyle[iMethod]);
961 TObjArray *thres = new TObjArray(3); thres->SetOwner(); thres->SetName("Thresholds");
962 results->AddAt(thres, 1);
963 for(Int_t iMethod = 0; iMethod < 3; iMethod++){
964 thres->AddAt(g = gPtrThres[iMethod] = new TGraphErrors(), iMethod);
965 g->SetName(Form("threshold_%s", MethodName[iMethod].Data()));
966 g->SetLineColor(colors[iMethod]);
967 g->SetMarkerColor(colors[iMethod]);
968 g->SetMarkerStyle(markerStyle[iMethod]);
972 TH1D *Histo1=0x0, *Histo2=0x0;
975 hPtr[0] = (TH2F*)histoContainer->At(AliTRDpidUtil::kLQ);
976 hPtr[1] = (TH2F*)histoContainer->At(AliTRDpidUtil::kNN);
977 hPtr[2] = (TH2F*)histoContainer->At(AliTRDpidUtil::kESD);
979 for(Int_t iMom = 0; iMom < fMomentumAxis->GetNbins(); iMom++){
980 mom = fMomentumAxis->GetBinCenter(iMom+1);
982 Int_t binEl = fMomentumAxis->GetNbins() * AliPID::kElectron + iMom + 1,
983 binPi = fMomentumAxis->GetNbins() * AliPID::kPion + iMom + 1;
984 for(Int_t iMethod = 0; iMethod < 3; iMethod++){
985 // Calculate the Pion Efficiency at 90% electron efficiency for each Method
986 Histo1 = hPtr[iMethod] -> ProjectionY(Form("%s_ele", MethodName[iMethod].Data()), binEl, binEl);
987 Histo2 = hPtr[iMethod] -> ProjectionY(Form("%s_pio", MethodName[iMethod].Data()), binPi, binPi);
989 if(!fUtil->CalculatePionEffi(Histo1, Histo2)) continue;
991 gPtrEff[iMethod]->SetPoint(iMom, mom, fUtil->GetPionEfficiency());
992 gPtrEff[iMethod]->SetPointError(iMom, 0., fUtil->GetError());
993 gPtrThres[iMethod]->SetPoint(iMom, mom, fUtil->GetThreshold());
994 gPtrThres[iMethod]->SetPointError(iMom, 0., 0.);
996 if(fDebugLevel>=2) Printf(Form("Pion Efficiency for 2-dim LQ is : %f +/- %f\n\n", MethodName[iMethod].Data()), fUtil->GetPionEfficiency(), fUtil->GetError());
1001 //________________________________________________________________________
1002 void AliTRDcheckPID::Terminate(Option_t *)
1004 // Draw result to the screen
1005 // Called once at the end of the query
1007 fContainer = dynamic_cast<TObjArray*>(GetOutputData(0));
1009 Printf("ERROR: list not available");