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