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