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