]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG1/TRD/AliTRDcheckPID.cxx
fix compilation problem in GSI train
[u/mrichter/AliRoot.git] / PWG1 / TRD / AliTRDcheckPID.cxx
CommitLineData
1ee39b3a 1//////////////////////////////////////////////////////
2//
3// PID performance checker of the TRD
4//
5// Performs checks of ESD information for TRD-PID and recalculate PID based on OCDB information
6// Also provides performance plots on detector based on the PID information - for the moment only
7// MC source is used but V0 is also possible. Here is a list of detector performance checks
8// - Integrated dE/dx per chamber
9// - <PH> as function of time bin and local radial position
10// - number of clusters/tracklet
11// - number of tracklets/track
12//
13// Author : Alex Wilk <wilka@uni-muenster.de>
14// Alex Bercuci <A.Bercuci@gsi.de>
15// Markus Fasel <M.Fasel@gsi.de>
16//
17///////////////////////////////////////////////////////
18
19#include "TAxis.h"
20#include "TROOT.h"
21#include "TPDGCode.h"
22#include "TCanvas.h"
23#include "TF1.h"
24#include "TH1F.h"
25#include "TH1D.h"
26#include "TH2F.h"
27#include "TProfile.h"
28#include "TProfile2D.h"
29#include "TGraph.h"
30#include "TGraphErrors.h"
31#include "TLegend.h"
32
33#include <TClonesArray.h>
34#include <TObjArray.h>
35#include <TList.h>
36
37#include "AliESDEvent.h"
38#include "AliESDInputHandler.h"
39#include "AliTrackReference.h"
40
41#include "AliAnalysisTask.h"
42
43#include "AliTRDtrackerV1.h"
44#include "AliTRDtrackV1.h"
45#include "AliTRDcluster.h"
46#include "AliTRDReconstructor.h"
47#include "AliCDBManager.h"
48#include "AliTRDpidUtil.h"
49
2e6aa247 50#include "Cal/AliTRDCalPID.h"
51#include "Cal/AliTRDCalPIDNN.h"
1ee39b3a 52#include "AliTRDcheckPID.h"
53#include "info/AliTRDtrackInfo.h"
b91fdd71 54#include "info/AliTRDpidInfo.h"
1ee39b3a 55
2e6aa247 56Char_t const * AliTRDcheckPID::fgMethod[3] = {"LQ", "NN", "ESD"};
57
1ee39b3a 58ClassImp(AliTRDcheckPID)
59
60//________________________________________________________________________
61AliTRDcheckPID::AliTRDcheckPID()
f2e89a4c 62 :AliTRDrecoTask()
f8f46e4d 63 ,fReconstructor(NULL)
64 ,fUtil(NULL)
65 ,fGraph(NULL)
66 ,fPID(NULL)
67 ,fMomentumAxis(NULL)
68 ,fMinNTracklets(AliTRDgeometry::kNlayer)
69 ,fMaxNTracklets(AliTRDgeometry::kNlayer)
70 {
71 //
72 // Default constructor
73 //
f2e89a4c 74 SetNameTitle("checkPID", "Check TRD PID");
2c5222ab 75 LocalInit();
f8f46e4d 76}
77
f2e89a4c 78//________________________________________________________________________
f8f46e4d 79AliTRDcheckPID::AliTRDcheckPID(char* name )
f2e89a4c 80 :AliTRDrecoTask(name, "Check TRD PID")
b91fdd71 81 ,fReconstructor(NULL)
82 ,fUtil(NULL)
83 ,fGraph(NULL)
84 ,fPID(NULL)
85 ,fMomentumAxis(NULL)
1ee39b3a 86 ,fMinNTracklets(AliTRDgeometry::kNlayer)
87 ,fMaxNTracklets(AliTRDgeometry::kNlayer)
88 {
89 //
90 // Default constructor
91 //
92
2c5222ab 93 LocalInit();
94 InitFunctorList();
95
96 DefineInput(2, TObjArray::Class()); // v0 list
97 DefineOutput(2, TObjArray::Class()); // pid info list
98}
99
100
101//________________________________________________________________________
102void AliTRDcheckPID::LocalInit()
103{
104// Initialize working data
1ee39b3a 105 fReconstructor = new AliTRDReconstructor();
106 fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
107
108 // Initialize momentum axis with default values
109 Double_t defaultMomenta[AliTRDCalPID::kNMom+1];
110 for(Int_t imom = 0; imom < AliTRDCalPID::kNMom+1; imom++)
111 defaultMomenta[imom] = AliTRDCalPID::GetMomentumBinning(imom);
112 SetMomentumBinning(AliTRDCalPID::kNMom, defaultMomenta);
113
2e6aa247 114 memset(fEfficiency, 0, AliPID::kSPECIES*sizeof(TObjArray*));
115
1ee39b3a 116 fUtil = new AliTRDpidUtil();
1ee39b3a 117}
118
1ee39b3a 119//________________________________________________________________________
120AliTRDcheckPID::~AliTRDcheckPID()
121{
b91fdd71 122 if(fPID){fPID->Delete(); delete fPID;}
123 if(fGraph){fGraph->Delete(); delete fGraph;}
124 if(fReconstructor) delete fReconstructor;
125 if(fUtil) delete fUtil;
1ee39b3a 126}
127
128
129//________________________________________________________________________
f8f46e4d 130void AliTRDcheckPID::UserCreateOutputObjects()
1ee39b3a 131{
132 // Create histograms
133 // Called once
134
f8f46e4d 135 OpenFile(1, "RECREATE");
1ee39b3a 136 fContainer = Histos();
b91fdd71 137
138 fPID = new TObjArray();
139 fPID->SetOwner(kTRUE);
140}
141
142//________________________________________________________
f8f46e4d 143void AliTRDcheckPID::UserExec(Option_t *opt)
b91fdd71 144{
145 //
146 // Execution part
147 //
148
149 fPID->Delete();
150
b4414720 151 AliTRDrecoTask::UserExec(opt);
b91fdd71 152
f8f46e4d 153 PostData(2, fPID);
1ee39b3a 154}
155
156
157//_______________________________________________________
158TObjArray * AliTRDcheckPID::Histos(){
159
160 //
161 // Create QA histograms
162 //
163 if(fContainer) return fContainer;
164
165 Int_t xBins = AliPID::kSPECIES*fMomentumAxis->GetNbins();
166 fContainer = new TObjArray(); fContainer->Expand(kNPlots);
167
168 const Float_t epsilon = 1./(2*(AliTRDpidUtil::kBins-1)); // get nice histos with bin center at 0 and 1
b91fdd71 169 TH1 *h = NULL;
1ee39b3a 170
2e6aa247 171 const Int_t kNmethodsPID=Int_t(sizeof(fgMethod)/sizeof(Char_t*));
172 // histos with posterior probabilities for all particle species
173 for(Int_t is=AliPID::kSPECIES; is--;){
174 fEfficiency[is] = new TObjArray(kNmethodsPID);
175 fEfficiency[is]->SetOwner();
176 fEfficiency[is]->SetName(Form("Eff_%s", AliPID::ParticleShortName(is)));
177 fContainer->AddAt(fEfficiency[is], is?kEfficiencyMu+is-1:kEfficiency);
178 for(Int_t im=kNmethodsPID; im--;){
179 if(!(h = (TH2F*)gROOT->FindObject(Form("PID_%s_%s", fgMethod[im], AliPID::ParticleShortName(is))))) {
180 h = new TH2F(Form("PID_%s_%s", fgMethod[im], AliPID::ParticleShortName(is)), "", xBins, -0.5, xBins - 0.5,
181 AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon);
182 } else h->Reset();
183 fEfficiency[is]->AddAt(h, im);
184 }
185 }
1ee39b3a 186
187 // histos of the dE/dx distribution for all 5 particle species and 11 momenta
188 if(!(h = (TH2F*)gROOT->FindObject("dEdx"))){
189 h = new TH2F("dEdx", "",
190 xBins, -0.5, xBins - 0.5,
2e6aa247 191 200, 0, 15);
192// 200, 0, 10000);
1ee39b3a 193 } else h->Reset();
194 fContainer->AddAt(h, kdEdx);
195
196 // histos of the dE/dx slices for all 5 particle species and 11 momenta
197 if(!(h = (TH2F*)gROOT->FindObject("dEdxSlice"))){
198 h = new TH2F("dEdxSlice", "",
199 xBins*AliTRDpidUtil::kLQslices, -0.5, xBins*AliTRDpidUtil::kLQslices - 0.5,
200 200, 0, 5000);
201 } else h->Reset();
202 fContainer->AddAt(h, kdEdxSlice);
203
204 // histos of the pulse height distribution for all 5 particle species and 11 momenta
205 TObjArray *fPH = new TObjArray(2);
206 fPH->SetOwner(); fPH->SetName("PH");
207 fContainer->AddAt(fPH, kPH);
208 if(!(h = (TProfile2D*)gROOT->FindObject("PHT"))){
209 h = new TProfile2D("PHT", "",
210 xBins, -0.5, xBins - 0.5,
211 AliTRDtrackerV1::GetNTimeBins(), -0.5, AliTRDtrackerV1::GetNTimeBins() - 0.5);
212 } else h->Reset();
213 fPH->AddAt(h, 0);
214 if(!(h = (TProfile2D*)gROOT->FindObject("PHX"))){
215 h = new TProfile2D("PHX", "",
216 xBins, -0.5, xBins - 0.5,
217 AliTRDtrackerV1::GetNTimeBins(), 0., .5*AliTRDgeometry::CamHght()+AliTRDgeometry::CdrHght());
218 } else h->Reset();
219 fPH->AddAt(h, 1);
220
221 // histos of the number of clusters distribution for all 5 particle species and 11 momenta
222 if(!(h = (TH2F*)gROOT->FindObject("NClus"))){
223 h = new TH2F("NClus", "",
224 xBins, -0.5, xBins - 0.5,
225 50, -0.5, 49.5);
226 } else h->Reset();
227 fContainer->AddAt(h, kNClus);
228
229
230 // momentum distributions - absolute and in momentum bins
231 if(!(h = (TH1F*)gROOT->FindObject("hMom"))){
232 h = new TH1F("hMom", "momentum distribution", fMomentumAxis->GetNbins(), fMomentumAxis->GetXmin(), fMomentumAxis->GetXmax());
233 } else h->Reset();
234 fContainer->AddAt(h, kMomentum);
235
236 if(!(h = (TH1F*)gROOT->FindObject("hMomBin"))){
237 h = new TH1F("hMomBin", "momentum distribution in momentum bins", fMomentumAxis->GetNbins(), fMomentumAxis->GetXmin(), fMomentumAxis->GetXmax());
238 } else h->Reset();
239 fContainer->AddAt(h, kMomentumBin);
240
241 // Number of tracklets per track
242 if(!(h = (TH2F*)gROOT->FindObject("nTracklets"))){
243 h = new TH2F("nTracklets", "",
244 xBins, -0.5, xBins - 0.5,
245 AliTRDgeometry::kNlayer, 0.5, AliTRDgeometry::kNlayer+.5);
246 } else h->Reset();
247 fContainer->AddAt(h, kNTracklets);
248
249 return fContainer;
250}
251
252
253//________________________________________________________________________
254Bool_t AliTRDcheckPID::CheckTrackQuality(const AliTRDtrackV1* track) const
255{
256 //
257 // Check if the track is ok for PID
258 //
259
260 Int_t ntracklets = track->GetNumberOfTracklets();
261 if(ntracklets >= fMinNTracklets && ntracklets <= fMaxNTracklets) return 1;
262// if(!fESD)
263// return 0;
264
265 return 0;
266}
267
268//________________________________________________________________________
269Int_t AliTRDcheckPID::CalcPDG(AliTRDtrackV1* track)
270{
271// Documentation to come
272
273 track -> SetReconstructor(fReconstructor);
274 (const_cast<AliTRDrecoParam*>(fReconstructor->GetRecoParam()))->SetPIDNeuralNetwork();
275 track -> CookPID();
276
277 if(track -> GetPID(AliPID::kElectron) > track -> GetPID(AliPID::kMuon) + track -> GetPID(AliPID::kPion) + track -> GetPID(AliPID::kKaon) + track -> GetPID(AliPID::kProton)){
278 return kElectron;
279 }
280 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)){
281 return kProton;
282 }
283 else if(track -> GetPID(AliPID::kKaon) > track -> GetPID(AliPID::kMuon) && track -> GetPID(AliPID::kKaon) > track -> GetPID(AliPID::kPion)){
284 return kKPlus;
285 }
286 else if(track -> GetPID(AliPID::kMuon) > track -> GetPID(AliPID::kPion)){
287 return kMuonPlus;
288 }
289 else{
290 return kPiPlus;
291 }
292}
293
294
295//_______________________________________________________
296TH1 *AliTRDcheckPID::PlotLQ(const AliTRDtrackV1 *track)
297{
298 //
299 // Plot the probabilities for electrons using 2-dim LQ
300 //
301
302 if(!fkESD){
f232df0d 303 AliDebug(2, "No ESD info available.");
2e6aa247 304 return NULL;
1ee39b3a 305 }
306
307 if(track) fkTrack = track;
308 if(!fkTrack){
f232df0d 309 AliDebug(2, "No Track defined.");
2e6aa247 310 return NULL;
1ee39b3a 311 }
312
313 ULong_t status;
314 status = fkESD -> GetStatus();
2e6aa247 315 if(!(status&AliESDtrack::kTRDin)) return NULL;
1ee39b3a 316
2e6aa247 317 if(!CheckTrackQuality(fkTrack)) return NULL;
1ee39b3a 318
319 AliTRDtrackV1 cTrack(*fkTrack);
320 cTrack.SetReconstructor(fReconstructor);
321
322 Int_t pdg = 0;
323 Float_t momentum = 0.;
324
325 if(fkMC){
326 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
327 pdg = fkMC->GetPDG();
328 } else{
329 //AliWarning("No MC info available!");
330 momentum = cTrack.GetMomentum(0);
331 pdg = CalcPDG(&cTrack);
332 }
2e6aa247 333 if(!IsInRange(momentum)) return NULL;
1ee39b3a 334
335 (const_cast<AliTRDrecoParam*>(fReconstructor->GetRecoParam()))->SetPIDNeuralNetwork(kFALSE);
336 cTrack.CookPID();
2e6aa247 337 if(cTrack.GetNumberOfTrackletsPID() < fMinNTracklets) return NULL;
1ee39b3a 338 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
2e6aa247 339
340 TH2F *hPIDLQ(NULL);
341 TObjArray *eff(NULL);
342 for(Int_t is=AliPID::kSPECIES; is--;){
343 if(!(eff = dynamic_cast<TObjArray *>(fContainer->At(is?kEfficiencyMu+is-1:kEfficiency)))){
344 AliWarning("No Histogram List defined.");
345 return NULL;
346 }
347 if(!(hPIDLQ = dynamic_cast<TH2F *>(eff->At(AliTRDpidUtil::kLQ)))){
348 AliWarning("No Histogram defined.");
349 return NULL;
350 }
351
352 hPIDLQ -> Fill(FindBin(species, momentum), cTrack.GetPID(is));
353 }
1ee39b3a 354 return hPIDLQ;
355}
356
357
2e6aa247 358
359
1ee39b3a 360//_______________________________________________________
361TH1 *AliTRDcheckPID::PlotNN(const AliTRDtrackV1 *track)
362{
363 //
364 // Plot the probabilities for electrons using 2-dim LQ
365 //
366
367 if(!fkESD){
f232df0d 368 AliDebug(2, "No ESD info available.");
2e6aa247 369 return NULL;
1ee39b3a 370 }
371
372 if(track) fkTrack = track;
373 if(!fkTrack){
f232df0d 374 AliDebug(2, "No Track defined.");
2e6aa247 375 return NULL;
1ee39b3a 376 }
377
378 ULong_t status;
379 status = fkESD -> GetStatus();
2e6aa247 380 if(!(status&AliESDtrack::kTRDin)) return NULL;
1ee39b3a 381
2e6aa247 382 if(!CheckTrackQuality(fkTrack)) return NULL;
1ee39b3a 383
1ee39b3a 384 AliTRDtrackV1 cTrack(*fkTrack);
385 cTrack.SetReconstructor(fReconstructor);
386
387 Int_t pdg = 0;
388 Float_t momentum = 0.;
389 if(fkMC){
390 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
391 pdg = fkMC->GetPDG();
392 } else {
393 //AliWarning("No MC info available!");
394 momentum = cTrack.GetMomentum(0);
395 pdg = CalcPDG(&cTrack);
396 }
2e6aa247 397 if(!IsInRange(momentum)) return NULL;
1ee39b3a 398
399 (const_cast<AliTRDrecoParam*>(fReconstructor->GetRecoParam()))->SetPIDNeuralNetwork();
400 cTrack.CookPID();
2e6aa247 401 if(cTrack.GetNumberOfTrackletsPID() < fMinNTracklets) return NULL;
1ee39b3a 402
403 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
2e6aa247 404
405
406 TH2F *hPIDNN(NULL);
407 TObjArray *eff(NULL);
408 for(Int_t is=AliPID::kSPECIES; is--;){
409 if(!(eff = dynamic_cast<TObjArray *>(fContainer->At(is?kEfficiencyMu+is-1:kEfficiency)))){
410 AliWarning("No Histogram List defined.");
411 return NULL;
412 }
413 if(!(hPIDNN = dynamic_cast<TH2F *>(eff->At(AliTRDpidUtil::kNN)))){
414 AliWarning("No Histogram defined.");
415 return NULL;
416 }
417
418 hPIDNN->Fill(FindBin(species, momentum), cTrack.GetPID(is));
419 }
1ee39b3a 420 return hPIDNN;
421}
422
423
2e6aa247 424
1ee39b3a 425//_______________________________________________________
426TH1 *AliTRDcheckPID::PlotESD(const AliTRDtrackV1 *track)
427{
428 //
429 // Plot the probabilities for electrons using 2-dim LQ
430 //
431
432 if(!fkESD){
f232df0d 433 AliDebug(2, "No ESD info available.");
2e6aa247 434 return NULL;
1ee39b3a 435 }
436
437 if(track) fkTrack = track;
438 if(!fkTrack){
f232df0d 439 AliDebug(2, "No Track defined.");
2e6aa247 440 return NULL;
1ee39b3a 441 }
442
443 ULong_t status;
444 status = fkESD -> GetStatus();
2e6aa247 445 if(!(status&AliESDtrack::kTRDin)) return NULL;
1ee39b3a 446
2e6aa247 447 if(!CheckTrackQuality(fkTrack)) return NULL;
448 if(fkTrack->GetNumberOfTrackletsPID() < fMinNTracklets) return NULL;
1ee39b3a 449
1ee39b3a 450 Int_t pdg = 0;
451 Float_t momentum = 0.;
452 if(fkMC){
453 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
454 pdg = fkMC->GetPDG();
455 } else {
456 //AliWarning("No MC info available!");
457 AliTRDtrackV1 cTrack(*fkTrack);
458 cTrack.SetReconstructor(fReconstructor);
459 momentum = cTrack.GetMomentum(0);
460 pdg = CalcPDG(&cTrack);
461 }
2e6aa247 462 if(!IsInRange(momentum)) return NULL;
1ee39b3a 463
1ee39b3a 464 const Double32_t *pidESD = fkESD->GetResponseIter();
465 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
2e6aa247 466
467 TH2F *hPID(NULL);
468 TObjArray *eff(NULL);
469 for(Int_t is=AliPID::kSPECIES; is--;){
470 if(!(eff = dynamic_cast<TObjArray *>(fContainer->At(is?kEfficiencyMu+is-1:kEfficiency)))){
471 AliWarning("No Histogram List defined.");
472 return NULL;
473 }
474 if(!(hPID = dynamic_cast<TH2F *>(eff->At(AliTRDpidUtil::kESD)))){
475 AliWarning("No Histogram defined.");
476 return NULL;
477 }
b91fdd71 478
479 AliTRDpidInfo *pid = new AliTRDpidInfo();
480 fPID->Add(pid);
2e6aa247 481 hPID->Fill(FindBin(species, momentum), pidESD[is]);
482 }
483 return hPID;
1ee39b3a 484}
485
486
2e6aa247 487
1ee39b3a 488//_______________________________________________________
489TH1 *AliTRDcheckPID::PlotdEdx(const AliTRDtrackV1 *track)
490{
491 //
492 // Plot the probabilities for electrons using 2-dim LQ
493 //
494
495 if(track) fkTrack = track;
496 if(!fkTrack){
f232df0d 497 AliDebug(2, "No Track defined.");
2e6aa247 498 return NULL;
1ee39b3a 499 }
500
2e6aa247 501 if(!CheckTrackQuality(fkTrack)) return NULL;
1ee39b3a 502
503 TH2F *hdEdx;
504 if(!(hdEdx = dynamic_cast<TH2F *>(fContainer->At(kdEdx)))){
505 AliWarning("No Histogram defined.");
2e6aa247 506 return NULL;
1ee39b3a 507 }
508
509 AliTRDtrackV1 cTrack(*fkTrack);
510 cTrack.SetReconstructor(fReconstructor);
511 Int_t pdg = 0;
512 Float_t momentum = 0.;
513 if(fkMC){
514 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
515 pdg = fkMC->GetPDG();
516 } else {
517 //AliWarning("No MC info available!");
518 momentum = cTrack.GetMomentum(0);
519 pdg = CalcPDG(&cTrack);
520 }
521 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
2e6aa247 522 if(!IsInRange(momentum)) return NULL;
1ee39b3a 523
2e6aa247 524 (const_cast<AliTRDrecoParam*>(fReconstructor->GetRecoParam()))->SetPIDNeuralNetwork(kTRUE);
525// (const_cast<AliTRDrecoParam*>(fReconstructor->GetRecoParam()))->SetPIDNeuralNetwork(kFALSE);
1ee39b3a 526 Float_t sumdEdx = 0;
527 Int_t iBin = FindBin(species, momentum);
b91fdd71 528 AliTRDseedV1 *tracklet = NULL;
1ee39b3a 529 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
530 sumdEdx = 0;
531 tracklet = cTrack.GetTracklet(iChamb);
532 if(!tracklet) continue;
2e6aa247 533 tracklet -> CookdEdx(AliTRDpidUtil::kNNslices);
534 for(Int_t i = AliTRDpidUtil::kNNslices; i--;) sumdEdx += tracklet->GetdEdx()[i];
535// tracklet -> CookdEdx(AliTRDpidUtil::kLQslices);
536// for(Int_t i = 3; i--;) sumdEdx += tracklet->GetdEdx()[i];
537 sumdEdx /= AliTRDCalPIDNN::kMLPscale;
1ee39b3a 538 hdEdx -> Fill(iBin, sumdEdx);
539 }
540 return hdEdx;
541}
542
543
544//_______________________________________________________
545TH1 *AliTRDcheckPID::PlotdEdxSlice(const AliTRDtrackV1 *track)
546{
547 //
548 // Plot the probabilities for electrons using 2-dim LQ
549 //
550
551 if(track) fkTrack = track;
552 if(!fkTrack){
f232df0d 553 AliDebug(2, "No Track defined.");
2e6aa247 554 return NULL;
1ee39b3a 555 }
556
2e6aa247 557 if(!CheckTrackQuality(fkTrack)) return NULL;
1ee39b3a 558
559 TH2F *hdEdxSlice;
560 if(!(hdEdxSlice = dynamic_cast<TH2F *>(fContainer->At(kdEdxSlice)))){
561 AliWarning("No Histogram defined.");
2e6aa247 562 return NULL;
1ee39b3a 563 }
564
565 AliTRDtrackV1 cTrack(*fkTrack);
566 cTrack.SetReconstructor(fReconstructor);
567 Int_t pdg = 0;
568 Float_t momentum = 0.;
569 if(fkMC){
570 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
571 pdg = fkMC->GetPDG();
572 } else {
573 //AliWarning("No MC info available!");
574 momentum = cTrack.GetMomentum(0);
575 pdg = CalcPDG(&cTrack);
576 }
2e6aa247 577 if(!IsInRange(momentum)) return NULL;
1ee39b3a 578
579 (const_cast<AliTRDrecoParam*>(fReconstructor->GetRecoParam()))->SetPIDNeuralNetwork(kFALSE);
580 Int_t iMomBin = fMomentumAxis->FindBin(momentum);
581 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
582 Float_t *fdEdx;
b91fdd71 583 AliTRDseedV1 *tracklet = NULL;
1ee39b3a 584 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
585 tracklet = cTrack.GetTracklet(iChamb);
586 if(!tracklet) continue;
587 tracklet -> CookdEdx(AliTRDpidUtil::kLQslices);
588 fdEdx = const_cast<Float_t *>(tracklet->GetdEdx());
589 for(Int_t iSlice = 0; iSlice < AliTRDpidUtil::kLQslices; iSlice++){
590 hdEdxSlice -> Fill(species * fMomentumAxis->GetNbins() * AliTRDpidUtil::kLQslices + (iMomBin-1) * AliTRDpidUtil::kLQslices + iSlice, fdEdx[iSlice]);
591 }
592 }
593
594 return hdEdxSlice;
595}
596
597
598//_______________________________________________________
599TH1 *AliTRDcheckPID::PlotPH(const AliTRDtrackV1 *track)
600{
601 //
602 // Plot the probabilities for electrons using 2-dim LQ
603 //
604
605 if(track) fkTrack = track;
606 if(!fkTrack){
f232df0d 607 AliDebug(2, "No Track defined.");
2e6aa247 608 return NULL;
1ee39b3a 609 }
610
2e6aa247 611 if(!CheckTrackQuality(fkTrack)) return NULL;
1ee39b3a 612
b91fdd71 613 TObjArray *arr = NULL;
1ee39b3a 614 TProfile2D *hPHX, *hPHT;
615 if(!(arr = dynamic_cast<TObjArray *>(fContainer->At(kPH)))){
616 AliWarning("No Histogram defined.");
2e6aa247 617 return NULL;
1ee39b3a 618 }
619 hPHT = (TProfile2D*)arr->At(0);
620 hPHX = (TProfile2D*)arr->At(1);
621
622 Int_t pdg = 0;
623 Float_t momentum = 0.;
624 if(fkMC){
625 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
626 pdg = fkMC->GetPDG();
627 } else {
628 //AliWarning("No MC info available!");
629 AliTRDtrackV1 cTrack(*fkTrack);
630 cTrack.SetReconstructor(fReconstructor);
631 momentum = cTrack.GetMomentum(0);
632 pdg = CalcPDG(&cTrack);
633 }
2e6aa247 634 if(!IsInRange(momentum)) return NULL;;
1ee39b3a 635
b91fdd71 636 AliTRDseedV1 *tracklet = NULL;
637 AliTRDcluster *cluster = NULL;
1ee39b3a 638 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
639 Int_t iBin = FindBin(species, momentum);
640 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
641 tracklet = fkTrack->GetTracklet(iChamb);
642 if(!tracklet) continue;
643 Float_t x0 = tracklet->GetX0();
644 for(Int_t ic = 0; ic < AliTRDtrackerV1::GetNTimeBins(); ic++){
645 if(!(cluster = tracklet->GetClusters(ic))) continue;
646 hPHT -> Fill(iBin, cluster->GetLocalTimeBin(), TMath::Abs(cluster->GetQ()));
647 hPHX -> Fill(iBin, x0 - cluster->GetX(), tracklet->GetdQdl(ic));
648 }
649 }
650 return hPHT;
651}
652
653
654//_______________________________________________________
655TH1 *AliTRDcheckPID::PlotNClus(const AliTRDtrackV1 *track)
656{
657 //
658 // Plot the probabilities for electrons using 2-dim LQ
659 //
660
661 if(track) fkTrack = track;
662 if(!fkTrack){
f232df0d 663 AliDebug(2, "No Track defined.");
2e6aa247 664 return NULL;
1ee39b3a 665 }
666
2e6aa247 667 if(!CheckTrackQuality(fkTrack)) return NULL;
1ee39b3a 668
669 TH2F *hNClus;
670 if(!(hNClus = dynamic_cast<TH2F *>(fContainer->At(kNClus)))){
671 AliWarning("No Histogram defined.");
2e6aa247 672 return NULL;
1ee39b3a 673 }
674
675
676 Int_t pdg = 0;
677 Float_t momentum = 0.;
678 if(fkMC){
679 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
680 pdg = fkMC->GetPDG();
681 } else {
682 //AliWarning("No MC info available!");
683 AliTRDtrackV1 cTrack(*fkTrack);
684 cTrack.SetReconstructor(fReconstructor);
685 momentum = cTrack.GetMomentum(0);
686 pdg = CalcPDG(&cTrack);
687 }
2e6aa247 688 if(!IsInRange(momentum)) return NULL;
1ee39b3a 689
690 Int_t species = AliTRDpidUtil::AliTRDpidUtil::Pdg2Pid(pdg);
691 Int_t iBin = FindBin(species, momentum);
b91fdd71 692 AliTRDseedV1 *tracklet = NULL;
1ee39b3a 693 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
694 tracklet = fkTrack->GetTracklet(iChamb);
695 if(!tracklet) continue;
696 hNClus -> Fill(iBin, tracklet->GetN());
697 }
698
699 return hNClus;
700}
701
702//_______________________________________________________
703TH1 *AliTRDcheckPID::PlotNTracklets(const AliTRDtrackV1 *track)
704{
705 //
706 // Plot the probabilities for electrons using 2-dim LQ
707 //
708
709 if(track) fkTrack = track;
710 if(!fkTrack){
f232df0d 711 AliDebug(2, "No Track defined.");
2e6aa247 712 return NULL;
1ee39b3a 713 }
714
715 TH2F *hTracklets;
716 if(!(hTracklets = dynamic_cast<TH2F *>(fContainer->At(kNTracklets)))){
717 AliWarning("No Histogram defined.");
2e6aa247 718 return NULL;
1ee39b3a 719 }
720
721 AliTRDtrackV1 cTrack(*fkTrack);
722 cTrack.SetReconstructor(fReconstructor);
723 Int_t pdg = 0;
724 Float_t momentum = 0.;
725 if(fkMC){
726 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
727 pdg = fkMC->GetPDG();
728 } else {
729 //AliWarning("No MC info available!");
730 momentum = cTrack.GetMomentum(0);
731 pdg = CalcPDG(&cTrack);
732 }
733 Int_t species = AliTRDpidUtil::Pdg2Pid(pdg);
2e6aa247 734 if(!IsInRange(momentum)) return NULL;
1ee39b3a 735
736 Int_t iBin = FindBin(species, momentum);
737 hTracklets -> Fill(iBin, cTrack.GetNumberOfTracklets());
738 return hTracklets;
739}
740
741//_______________________________________________________
742TH1 *AliTRDcheckPID::PlotMom(const AliTRDtrackV1 *track)
743{
744 //
745 // Plot the probabilities for electrons using 2-dim LQ
746 //
747
748 if(track) fkTrack = track;
749 if(!fkTrack){
f232df0d 750 AliDebug(2, "No Track defined.");
2e6aa247 751 return NULL;
1ee39b3a 752 }
753
2e6aa247 754 if(!CheckTrackQuality(fkTrack)) return NULL;
1ee39b3a 755
756 TH1F *hMom;
757 if(!(hMom = dynamic_cast<TH1F *>(fContainer->At(kMomentum)))){
758 AliWarning("No Histogram defined.");
2e6aa247 759 return NULL;
1ee39b3a 760 }
761
762
763 Int_t pdg = 0;
764 Float_t momentum = 0.;
765 if(fkMC){
766 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
767 pdg = fkMC->GetPDG();
768 } else {
769 //AliWarning("No MC info available!");
770 AliTRDtrackV1 cTrack(*fkTrack);
771 cTrack.SetReconstructor(fReconstructor);
772 momentum = cTrack.GetMomentum(0);
773 pdg = CalcPDG(&cTrack);
774 }
775 if(IsInRange(momentum)) hMom -> Fill(momentum);
776 return hMom;
777}
778
779
780//_______________________________________________________
781TH1 *AliTRDcheckPID::PlotMomBin(const AliTRDtrackV1 *track)
782{
783 //
784 // Plot the probabilities for electrons using 2-dim LQ
785 //
786
787 if(track) fkTrack = track;
788 if(!fkTrack){
f232df0d 789 AliDebug(2, "No Track defined.");
2e6aa247 790 return NULL;
1ee39b3a 791 }
792
2e6aa247 793 if(!CheckTrackQuality(fkTrack)) return NULL;
1ee39b3a 794
795 TH1F *hMomBin;
796 if(!(hMomBin = dynamic_cast<TH1F *>(fContainer->At(kMomentumBin)))){
797 AliWarning("No Histogram defined.");
2e6aa247 798 return NULL;
1ee39b3a 799 }
800
801
802 Int_t pdg = 0;
803 Float_t momentum = 0.;
804
805 if(fkMC){
806 if(fkMC->GetTrackRef()) momentum = fkMC->GetTrackRef()->P();
807 pdg = fkMC->GetPDG();
808 } else {
809 //AliWarning("No MC info available!");
810 AliTRDtrackV1 cTrack(*fkTrack);
811 cTrack.SetReconstructor(fReconstructor);
812 momentum = cTrack.GetMomentum(0);
813 }
814 if(IsInRange(momentum)) hMomBin -> Fill(fMomentumAxis->FindBin(momentum));
815 return hMomBin;
816}
817
818
819//________________________________________________________
820Bool_t AliTRDcheckPID::GetRefFigure(Int_t ifig)
821{
822// Steering function to retrieve performance plots
823
824 Bool_t kFIRST = kTRUE;
b91fdd71 825 TGraphErrors *g = NULL;
826 TAxis *ax = NULL;
827 TObjArray *arr = NULL;
828 TH1 *h1 = NULL, *h=NULL;
829 TH2 *h2 = NULL;
830 TList *content = NULL;
1ee39b3a 831 switch(ifig){
832 case kEfficiency:{
833 gPad->Divide(2, 1, 1.e-5, 1.e-5);
834 TList *l=gPad->GetListOfPrimitives();
835 TVirtualPad *pad = ((TVirtualPad*)l->At(0));pad->cd();
836 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
837
2e6aa247 838 TLegend *legEff = new TLegend(.64, .84, .98, .98);
839 legEff->SetBorderSize(1);legEff->SetTextSize(0.03255879);
1ee39b3a 840 legEff->SetFillColor(0);
841 h=new TH1S("hEff", "", 1, .5, 11.);
842 h->SetLineColor(1);h->SetLineWidth(1);
843 ax = h->GetXaxis();
844 ax->SetTitle("p [GeV/c]");
845 ax->SetRangeUser(.5, 11.);
846 ax->SetMoreLogLabels();
847 ax = h->GetYaxis();
848 ax->SetTitle("#epsilon_{#pi} [%]");
849 ax->CenterTitle();
850 ax->SetRangeUser(1.e-2, 10.);
851 h->Draw();
2e6aa247 852 content = (TList *)fGraph->FindObject(Form("Eff_%s", AliTRDCalPID::GetPartName(AliPID::kPion)));
1ee39b3a 853 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kLQ))) break;
854 if(!g->GetN()) break;
2e6aa247 855 legEff->SetHeader("PID Method [PION]");
1ee39b3a 856 g->Draw("pc"); legEff->AddEntry(g, "LQ 2D", "pl");
857 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kNN))) break;
858 g->Draw("pc"); legEff->AddEntry(g, "NN", "pl");
859 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kESD))) break;
860 g->Draw("p"); legEff->AddEntry(g, "ESD", "pl");
861 legEff->Draw();
862 gPad->SetLogy();
863 gPad->SetLogx();
864 gPad->SetGridy();
865 gPad->SetGridx();
866
2e6aa247 867
1ee39b3a 868 pad = ((TVirtualPad*)l->At(1));pad->cd();
869 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
870 h=new TH1S("hThr", "", 1, .5, 11.);
871 h->SetLineColor(1);h->SetLineWidth(1);
872 ax = h->GetXaxis();
873 ax->SetTitle("p [GeV/c]");
874 ax->SetMoreLogLabels();
875 ax = h->GetYaxis();
876 ax->SetTitle("Threshold [%]");
877 ax->SetRangeUser(5.e-2, 1.);
878 h->Draw();
2e6aa247 879 content = (TList *)fGraph->FindObject("Thres");
1ee39b3a 880 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kLQ))) break;
881 if(!g->GetN()) break;
882 g->Draw("pc");
883 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kNN))) break;
884 g->Draw("pc");
885 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kESD))) break;
886 g->Draw("p");
887 gPad->SetLogx();
888 gPad->SetGridy();
889 gPad->SetGridx();
890 return kTRUE;
891 }
2e6aa247 892 case kEfficiencyKa:{
893 gPad->Divide(2, 1, 1.e-5, 1.e-5);
894 TList *l=gPad->GetListOfPrimitives();
895 TVirtualPad *pad = ((TVirtualPad*)l->At(0));pad->cd();
896 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
897
898 TLegend *legEff = new TLegend(.64, .84, .98, .98);
899 legEff->SetBorderSize(1);legEff->SetTextSize(0.03255879);
900 legEff->SetFillColor(0);
901 h = (TH1S*)gROOT->FindObject("hEff");
902 h=(TH1S*)h->Clone("hEff_K");
903 h->SetYTitle("#epsilon_{K} [%]");
904 h->GetYaxis()->SetRangeUser(1.e-2, 1.e2);
905 h->Draw();
906 content = (TList *)fGraph->FindObject(Form("Eff_%s", AliTRDCalPID::GetPartName(AliPID::kKaon)));
907 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kLQ))) break;
908 if(!g->GetN()) break;
909 legEff->SetHeader("PID Method [KAON]");
910 g->Draw("pc"); legEff->AddEntry(g, "LQ 2D", "pl");
911 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kNN))) break;
912 g->Draw("pc"); legEff->AddEntry(g, "NN", "pl");
913 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kESD))) break;
914 g->Draw("p"); legEff->AddEntry(g, "ESD", "pl");
915 legEff->Draw();
916 gPad->SetLogy();
917 gPad->SetLogx();
918 gPad->SetGridy();
919 gPad->SetGridx();
920
921 TLegend *legEff2 = new TLegend(.64, .84, .98, .98);
922 legEff2->SetBorderSize(1);legEff2->SetTextSize(0.03255879);
923 legEff2->SetFillColor(0);
924 pad = ((TVirtualPad*)l->At(1));pad->cd();
925 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
926 h=(TH1S*)h->Clone("hEff_p");
927 h->SetYTitle("#epsilon_{p} [%]");
928 h->GetYaxis()->SetRangeUser(1.e-2, 1.e2);
929 h->Draw();
930 content = (TList *)fGraph->FindObject(Form("Eff_%s", AliTRDCalPID::GetPartName(AliPID::kProton)));
931 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kLQ))) break;
932 if(!g->GetN()) break;
933 legEff2->SetHeader("PID Method [PROTON]");
934 g->Draw("pc"); legEff2->AddEntry(g, "LQ 2D", "pl");
935 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kNN))) break;
936 g->Draw("pc"); legEff2->AddEntry(g, "NN", "pl");
937 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kESD))) break;
938 g->Draw("p"); legEff2->AddEntry(g, "ESD", "pl");
939 legEff2->Draw();
940 gPad->SetLogy();
941 gPad->SetLogx();
942 gPad->SetGridy();
943 gPad->SetGridx();
944 return kTRUE;
945 }
1ee39b3a 946 case kdEdx:{
947 // save 2.0 GeV projection as reference
948 TLegend *legdEdx = new TLegend(.7, .7, .98, .98);
949 legdEdx->SetBorderSize(1);
950 kFIRST = kTRUE;
951 if(!(h2 = (TH2F*)(fContainer->At(kdEdx)))) break;
952 legdEdx->SetHeader("Particle Species");
953 gPad->SetMargin(0.1, 0.01, 0.1, 0.01);
954 for(Int_t is = AliPID::kSPECIES-1; is>=0; is--){
955 Int_t bin = FindBin(is, 2.);
956 h1 = h2->ProjectionY(Form("px%d", is), bin, bin);
957 if(!h1->GetEntries()) continue;
958 h1->Scale(1./h1->Integral());
959 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
960 if(kFIRST){
961 h1->GetXaxis()->SetTitle("dE/dx (a.u.)");
962 h1->GetYaxis()->SetTitle("<Entries>");
963 }
964 h = (TH1F*)h1->DrawClone(kFIRST ? "c" : "samec");
965 legdEdx->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l");
966 kFIRST = kFALSE;
967 }
968 if(kFIRST) break;
969 legdEdx->Draw();
970 gPad->SetLogy();
971 gPad->SetLogx(0);
972 gPad->SetGridy();
973 gPad->SetGridx();
974 return kTRUE;
975 }
976 case kdEdxSlice:
977 break;
978 case kPH:{
979 gPad->Divide(2, 1, 1.e-5, 1.e-5);
980 TList *l=gPad->GetListOfPrimitives();
981
982 // save 2.0 GeV projection as reference
983 TLegend *legPH = new TLegend(.4, .7, .68, .98);
984 legPH->SetBorderSize(1);legPH->SetFillColor(0);
985 legPH->SetHeader("Particle Species");
986 if(!(arr = (TObjArray*)(fContainer->At(kPH)))) break;
987 if(!(h2 = (TProfile2D*)(arr->At(0)))) break;
988
989 TVirtualPad *pad = ((TVirtualPad*)l->At(0));pad->cd();
990 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
991 kFIRST = kTRUE;
992 for(Int_t is=0; is<AliPID::kSPECIES; is++){
993 Int_t bin = FindBin(is, 2.);
994 h1 = h2->ProjectionY(Form("pyt%d", is), bin, bin);
995 if(!h1->GetEntries()) continue;
996 h1->SetMarkerStyle(24);
997 h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
998 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
999 if(kFIRST){
1000 h1->GetXaxis()->SetTitle("t_{drift} [100*ns]");
1001 h1->GetYaxis()->SetTitle("<dQ/dt> [a.u.]");
1002 }
1003 h = (TH1F*)h1->DrawClone(kFIRST ? "c" : "samec");
1004 legPH->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "pl");
1005 kFIRST = kFALSE;
1006 }
1007
1008 pad = ((TVirtualPad*)l->At(1));pad->cd();
1009 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
1010 if(!(h2 = (TProfile2D*)(arr->At(1)))) break;
1011 kFIRST = kTRUE;
1012 for(Int_t is=0; is<AliPID::kSPECIES; is++){
1013 Int_t bin = FindBin(is, 2.);
1014 h1 = h2->ProjectionY(Form("pyx%d", is), bin, bin);
1015 if(!h1->GetEntries()) continue;
1016 h1->SetMarkerStyle(24);
1017 h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
1018 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
1019 if(kFIRST){
1020 h1->GetXaxis()->SetTitle("x_{drift} [cm]");
1021 h1->GetYaxis()->SetTitle("<dQ/dl> [a.u./cm]");
1022 }
1023 h = (TH1F*)h1->DrawClone(kFIRST ? "c" : "samec");
1024 kFIRST = kFALSE;
1025 }
1026
1027 if(kFIRST) break;
1028 legPH->Draw();
1029 gPad->SetLogy(0);
1030 gPad->SetLogx(0);
1031 gPad->SetGridy();
1032 gPad->SetGridx();
1033 return kTRUE;
1034 }
1035 case kNClus:{
1036 // save 2.0 GeV projection as reference
1037 TLegend *legNClus = new TLegend(.13, .7, .4, .98);
1038 legNClus->SetBorderSize(1);
1039 legNClus->SetFillColor(0);
1040
1041 kFIRST = kTRUE;
1042 if(!(h2 = (TH2F*)(fContainer->At(kNClus)))) break;
1043 legNClus->SetHeader("Particle Species");
1044 for(Int_t is=0; is<AliPID::kSPECIES; is++){
1045 Int_t bin = FindBin(is, 2.);
1046 h1 = h2->ProjectionY(Form("pyNClus%d", is), bin, bin);
1047 if(!h1->GetEntries()) continue;
1048 h1->Scale(100./h1->Integral());
1049 //h1->SetMarkerStyle(24);
1050 //h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
1051 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
1052 if(kFIRST){
1053 h1->GetXaxis()->SetTitle("N^{cl}/tracklet");
1054 h1->GetYaxis()->SetTitle("Prob. [%]");
1055 h = (TH1F*)h1->DrawClone("c");
1056 h->SetMaximum(55.);
1057 h->GetXaxis()->SetRangeUser(0., 35.);
1058 kFIRST = kFALSE;
1059 } else h = (TH1F*)h1->DrawClone("samec");
1060
1061 legNClus->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l");
1062 }
1063 if(kFIRST) break;
1064 legNClus->Draw();
1065 gPad->SetLogy(0);
1066 gPad->SetLogx(0);
1067 gPad->SetGridy();
1068 gPad->SetGridx();
1069 return kTRUE;
1070 }
1071 case kMomentum:
1072 case kMomentumBin:
1073 break;
1074 case kNTracklets:{
1075 TLegend *legNClus = new TLegend(.4, .7, .68, .98);
1076 legNClus->SetBorderSize(1);
1077 kFIRST = kTRUE;
1078 if(!(h2 = (TH2F*)(fContainer->At(kNTracklets)))) break;
1079 legNClus->SetHeader("Particle Species");
1080 for(Int_t is=0; is<AliPID::kSPECIES; is++){
1081 Int_t bin = FindBin(is, 2.);
1082 h1 = h2->ProjectionY(Form("pyNTracklets%d", is), bin, bin);
1083 if(!h1->GetEntries()) continue;
1084 h1->Scale(100./h1->Integral());
1085 //h1->SetMarkerStyle(24);
1086 //h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
1087 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
1088 if(kFIRST){
1089 h1->GetXaxis()->SetTitle("N^{trklt}/track");
1090 h1->GetXaxis()->SetRangeUser(1.,6.);
1091 h1->GetYaxis()->SetTitle("Prob. [%]");
1092 }
1093 h = (TH1F*)h1->DrawClone(kFIRST ? "c" : "samec");
1094 legNClus->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l");
1095 kFIRST = kFALSE;
1096 }
1097 if(kFIRST) break;
1098 legNClus->Draw();
1099 gPad->SetLogy(0);
1100 gPad->SetLogx(0);
1101 gPad->SetGridy();
1102 gPad->SetGridx();
1103 return kTRUE;
1104 }
1105 }
1106 AliInfo(Form("Reference plot [%d] missing result", ifig));
1107 return kFALSE;
1108}
1109
1110//________________________________________________________________________
1111Bool_t AliTRDcheckPID::PostProcess()
1112{
1113 // Draw result to the screen
1114 // Called once at the end of the query
1115
1116 if (!fContainer) {
1117 Printf("ERROR: list not available");
1118 return kFALSE;
1119 }
2e6aa247 1120
1121 TObjArray *eff(NULL);
1ee39b3a 1122 if(!fGraph){
2e6aa247 1123 fGraph = new TObjArray(2*AliPID::kSPECIES);
1ee39b3a 1124 fGraph->SetOwner();
2e6aa247 1125
1126 if(!(eff = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){
1127 AliError("Efficiency container for Electrons missing.");
1128 return kFALSE;
1129 }
1130 EvaluateEfficiency(eff, fGraph, AliPID::kPion, 0.9);
1131 EvaluateEfficiency(eff, fGraph, AliPID::kKaon, 0.9);
1132 EvaluateEfficiency(eff, fGraph, AliPID::kProton, 0.9);
1ee39b3a 1133 }
2e6aa247 1134 fNRefFigures = 12;
1ee39b3a 1135 return kTRUE;
1136}
1137
1138//________________________________________________________________________
2e6aa247 1139void AliTRDcheckPID::EvaluateEfficiency(TObjArray * const histoContainer, TObjArray *results, Int_t species, Float_t electronEfficiency){
1ee39b3a 1140// Process PID information for pion efficiency
1141
1142 fUtil->SetElectronEfficiency(electronEfficiency);
1143
2e6aa247 1144
1145 const Int_t kNmethodsPID=Int_t(sizeof(fgMethod)/sizeof(Char_t*));
1146 Color_t colors[kNmethodsPID] = {kBlue, kGreen+2, kRed};
1147 Int_t markerStyle[kNmethodsPID] = {7, 7, 24};
1ee39b3a 1148 // efficiency graphs
2e6aa247 1149 TGraphErrors *g(NULL);
1150 TObjArray *eff = new TObjArray(kNmethodsPID); eff->SetOwner(); eff->SetName(Form("Eff_%s", AliTRDCalPID::GetPartName(species)));
1151 results->AddAt(eff, species);
1152 for(Int_t iMethod = 0; iMethod < kNmethodsPID; iMethod++){
1153 eff->AddAt(g = new TGraphErrors(), iMethod);
1154 g->SetName(Form("%s", fgMethod[iMethod]));
1ee39b3a 1155 g->SetLineColor(colors[iMethod]);
1156 g->SetMarkerColor(colors[iMethod]);
1157 g->SetMarkerStyle(markerStyle[iMethod]);
1158 }
1159
2e6aa247 1160 // Threshold graphs if not already
1161 TObjArray *thres(NULL);
1162 if(!(results->At(AliPID::kSPECIES))){
1163 thres = new TObjArray(kNmethodsPID); thres->SetOwner();
1164 thres->SetName("Thres");
1165 results->AddAt(thres, AliPID::kSPECIES);
1166 for(Int_t iMethod = 0; iMethod < kNmethodsPID; iMethod++){
1167 thres->AddAt(g = new TGraphErrors(), iMethod);
1168 g->SetName(Form("%s", fgMethod[iMethod]));
1169 g->SetLineColor(colors[iMethod]);
1170 g->SetMarkerColor(colors[iMethod]);
1171 g->SetMarkerStyle(markerStyle[iMethod]);
1172 }
1ee39b3a 1173 }
1ee39b3a 1174
2e6aa247 1175 TH2F *hPtr[kNmethodsPID]={
1176 (TH2F*)histoContainer->At(AliTRDpidUtil::kLQ),
1177 (TH2F*)histoContainer->At(AliTRDpidUtil::kNN),
1178 (TH2F*)histoContainer->At(AliTRDpidUtil::kESD)
1179 };
1ee39b3a 1180 for(Int_t iMom = 0; iMom < fMomentumAxis->GetNbins(); iMom++){
2e6aa247 1181 Float_t mom(fMomentumAxis->GetBinCenter(iMom+1));
1ee39b3a 1182
2e6aa247 1183 Int_t binEl(fMomentumAxis->GetNbins() * AliPID::kElectron + iMom + 1),
1184 binXX(fMomentumAxis->GetNbins() * species + iMom + 1);
1185
1186 for(Int_t iMethod = 0; iMethod < kNmethodsPID; iMethod++){
1187 // Calculate the Species Efficiency at electronEfficiency% electron efficiency for each Method
1188
1189 TH1D *histo1 = hPtr[iMethod] -> ProjectionY(Form("%s_el", fgMethod[iMethod]), binEl, binEl);
1190 TH1D *histo2 = hPtr[iMethod] -> ProjectionY(Form("%s_%s", fgMethod[iMethod], AliTRDCalPID::GetPartName(species)), binXX, binXX);
1ee39b3a 1191
1192 if(!fUtil->CalculatePionEffi(histo1, histo2)) continue;
1193
2e6aa247 1194 g=(TGraphErrors*)eff->At(iMethod);
1195 g->SetPoint(iMom, mom, 1.e2*fUtil->GetPionEfficiency());
1196 g->SetPointError(iMom, 0., 1.e2*fUtil->GetError());
1197 AliDebug(2, Form("%s Efficiency for %s is : %f +/- %f", AliTRDCalPID::GetPartName(species), fgMethod[iMethod], fUtil->GetPionEfficiency(), fUtil->GetError()));
1198
1199 if(!thres) continue;
1200 g=(TGraphErrors*)thres->At(iMethod);
1201 g->SetPoint(iMom, mom, fUtil->GetThreshold());
1202 g->SetPointError(iMom, 0., 0.);
1ee39b3a 1203 }
1204 }
1205}