]> 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
change systematic and resolution plots
[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
98 DefineInput(2, TObjArray::Class()); // v0 list
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
3d19c1b0 147 fV0s = dynamic_cast<TObjArray *>(GetInputData(2));
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"))){
206 h = new TProfile2D("PHT", "",
207 xBins, -0.5, xBins - 0.5,
208 AliTRDtrackerV1::GetNTimeBins(), -0.5, AliTRDtrackerV1::GetNTimeBins() - 0.5);
209 } else h->Reset();
210 fPH->AddAt(h, 0);
211 if(!(h = (TProfile2D*)gROOT->FindObject("PHX"))){
212 h = new TProfile2D("PHX", "",
213 xBins, -0.5, xBins - 0.5,
214 AliTRDtrackerV1::GetNTimeBins(), 0., .5*AliTRDgeometry::CamHght()+AliTRDgeometry::CdrHght());
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();
654 for(Int_t ic = 0; ic < AliTRDtrackerV1::GetNTimeBins(); ic++){
655 if(!(cluster = tracklet->GetClusters(ic))) continue;
656 hPHT -> Fill(iBin, cluster->GetLocalTimeBin(), TMath::Abs(cluster->GetQ()));
657 hPHX -> Fill(iBin, x0 - cluster->GetX(), tracklet->GetdQdl(ic));
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 }
849 AliDebug(2, Form("TRACK[%d] species[%s][%d]\n", fkESD->GetId(), AliPID::ParticleShortName(fkMC->GetPID()), fkMC->GetPDG()));
850
851 TH1 *h=dynamic_cast<TH1F*>(fContainer->At(kV0));
852 Int_t sgn(0), n(0); AliTRDv0Info *v0(NULL);
853 for(Int_t iv0(fV0s->GetEntriesFast()); iv0--;){
854 if(!(v0=(AliTRDv0Info*)fV0s->At(iv0))) continue;
855 if(!(sgn = v0->HasTrack(fkESD->GetId()))) continue;
856 //for(Int_t is=AliPID::kSPECIES; is--;) v0->GetPID(is, track);
857 //v0->Print();
858 n++;
859 //break;
860 }
861 h->Fill(n);
862 return h;
863}
1ee39b3a 864
865//________________________________________________________
866Bool_t AliTRDcheckPID::GetRefFigure(Int_t ifig)
867{
868// Steering function to retrieve performance plots
869
870 Bool_t kFIRST = kTRUE;
b91fdd71 871 TGraphErrors *g = NULL;
872 TAxis *ax = NULL;
873 TObjArray *arr = NULL;
874 TH1 *h1 = NULL, *h=NULL;
875 TH2 *h2 = NULL;
876 TList *content = NULL;
1ee39b3a 877 switch(ifig){
878 case kEfficiency:{
879 gPad->Divide(2, 1, 1.e-5, 1.e-5);
880 TList *l=gPad->GetListOfPrimitives();
881 TVirtualPad *pad = ((TVirtualPad*)l->At(0));pad->cd();
882 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
883
2e6aa247 884 TLegend *legEff = new TLegend(.64, .84, .98, .98);
885 legEff->SetBorderSize(1);legEff->SetTextSize(0.03255879);
1ee39b3a 886 legEff->SetFillColor(0);
887 h=new TH1S("hEff", "", 1, .5, 11.);
888 h->SetLineColor(1);h->SetLineWidth(1);
889 ax = h->GetXaxis();
890 ax->SetTitle("p [GeV/c]");
891 ax->SetRangeUser(.5, 11.);
892 ax->SetMoreLogLabels();
893 ax = h->GetYaxis();
894 ax->SetTitle("#epsilon_{#pi} [%]");
895 ax->CenterTitle();
896 ax->SetRangeUser(1.e-2, 10.);
897 h->Draw();
2e6aa247 898 content = (TList *)fGraph->FindObject(Form("Eff_%s", AliTRDCalPID::GetPartName(AliPID::kPion)));
1ee39b3a 899 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kLQ))) break;
900 if(!g->GetN()) break;
2e6aa247 901 legEff->SetHeader("PID Method [PION]");
1ee39b3a 902 g->Draw("pc"); legEff->AddEntry(g, "LQ 2D", "pl");
903 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kNN))) break;
904 g->Draw("pc"); legEff->AddEntry(g, "NN", "pl");
905 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kESD))) break;
906 g->Draw("p"); legEff->AddEntry(g, "ESD", "pl");
907 legEff->Draw();
908 gPad->SetLogy();
909 gPad->SetLogx();
910 gPad->SetGridy();
911 gPad->SetGridx();
912
2e6aa247 913
1ee39b3a 914 pad = ((TVirtualPad*)l->At(1));pad->cd();
915 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
916 h=new TH1S("hThr", "", 1, .5, 11.);
917 h->SetLineColor(1);h->SetLineWidth(1);
918 ax = h->GetXaxis();
919 ax->SetTitle("p [GeV/c]");
920 ax->SetMoreLogLabels();
921 ax = h->GetYaxis();
922 ax->SetTitle("Threshold [%]");
923 ax->SetRangeUser(5.e-2, 1.);
924 h->Draw();
2e6aa247 925 content = (TList *)fGraph->FindObject("Thres");
1ee39b3a 926 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kLQ))) break;
927 if(!g->GetN()) break;
928 g->Draw("pc");
929 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kNN))) break;
930 g->Draw("pc");
931 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kESD))) break;
932 g->Draw("p");
933 gPad->SetLogx();
934 gPad->SetGridy();
935 gPad->SetGridx();
936 return kTRUE;
937 }
2e6aa247 938 case kEfficiencyKa:{
939 gPad->Divide(2, 1, 1.e-5, 1.e-5);
940 TList *l=gPad->GetListOfPrimitives();
941 TVirtualPad *pad = ((TVirtualPad*)l->At(0));pad->cd();
942 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
943
944 TLegend *legEff = new TLegend(.64, .84, .98, .98);
945 legEff->SetBorderSize(1);legEff->SetTextSize(0.03255879);
946 legEff->SetFillColor(0);
947 h = (TH1S*)gROOT->FindObject("hEff");
948 h=(TH1S*)h->Clone("hEff_K");
949 h->SetYTitle("#epsilon_{K} [%]");
950 h->GetYaxis()->SetRangeUser(1.e-2, 1.e2);
951 h->Draw();
952 content = (TList *)fGraph->FindObject(Form("Eff_%s", AliTRDCalPID::GetPartName(AliPID::kKaon)));
953 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kLQ))) break;
954 if(!g->GetN()) break;
955 legEff->SetHeader("PID Method [KAON]");
956 g->Draw("pc"); legEff->AddEntry(g, "LQ 2D", "pl");
957 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kNN))) break;
958 g->Draw("pc"); legEff->AddEntry(g, "NN", "pl");
959 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kESD))) break;
960 g->Draw("p"); legEff->AddEntry(g, "ESD", "pl");
961 legEff->Draw();
962 gPad->SetLogy();
963 gPad->SetLogx();
964 gPad->SetGridy();
965 gPad->SetGridx();
966
967 TLegend *legEff2 = new TLegend(.64, .84, .98, .98);
968 legEff2->SetBorderSize(1);legEff2->SetTextSize(0.03255879);
969 legEff2->SetFillColor(0);
970 pad = ((TVirtualPad*)l->At(1));pad->cd();
971 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
972 h=(TH1S*)h->Clone("hEff_p");
973 h->SetYTitle("#epsilon_{p} [%]");
974 h->GetYaxis()->SetRangeUser(1.e-2, 1.e2);
975 h->Draw();
976 content = (TList *)fGraph->FindObject(Form("Eff_%s", AliTRDCalPID::GetPartName(AliPID::kProton)));
977 if(!(g = (TGraphErrors*)content->At(AliTRDpidUtil::kLQ))) break;
978 if(!g->GetN()) break;
979 legEff2->SetHeader("PID Method [PROTON]");
980 g->Draw("pc"); legEff2->AddEntry(g, "LQ 2D", "pl");
981 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kNN))) break;
982 g->Draw("pc"); legEff2->AddEntry(g, "NN", "pl");
983 if(! (g = (TGraphErrors*)content->At(AliTRDpidUtil::kESD))) break;
984 g->Draw("p"); legEff2->AddEntry(g, "ESD", "pl");
985 legEff2->Draw();
986 gPad->SetLogy();
987 gPad->SetLogx();
988 gPad->SetGridy();
989 gPad->SetGridx();
990 return kTRUE;
991 }
1ee39b3a 992 case kdEdx:{
993 // save 2.0 GeV projection as reference
994 TLegend *legdEdx = new TLegend(.7, .7, .98, .98);
995 legdEdx->SetBorderSize(1);
996 kFIRST = kTRUE;
997 if(!(h2 = (TH2F*)(fContainer->At(kdEdx)))) break;
998 legdEdx->SetHeader("Particle Species");
999 gPad->SetMargin(0.1, 0.01, 0.1, 0.01);
1000 for(Int_t is = AliPID::kSPECIES-1; is>=0; is--){
1001 Int_t bin = FindBin(is, 2.);
1002 h1 = h2->ProjectionY(Form("px%d", is), bin, bin);
1003 if(!h1->GetEntries()) continue;
1004 h1->Scale(1./h1->Integral());
1005 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
1006 if(kFIRST){
1007 h1->GetXaxis()->SetTitle("dE/dx (a.u.)");
1008 h1->GetYaxis()->SetTitle("<Entries>");
1009 }
1010 h = (TH1F*)h1->DrawClone(kFIRST ? "c" : "samec");
1011 legdEdx->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l");
1012 kFIRST = kFALSE;
1013 }
1014 if(kFIRST) break;
1015 legdEdx->Draw();
1016 gPad->SetLogy();
1017 gPad->SetLogx(0);
1018 gPad->SetGridy();
1019 gPad->SetGridx();
1020 return kTRUE;
1021 }
1022 case kdEdxSlice:
1023 break;
1024 case kPH:{
1025 gPad->Divide(2, 1, 1.e-5, 1.e-5);
1026 TList *l=gPad->GetListOfPrimitives();
1027
1028 // save 2.0 GeV projection as reference
1029 TLegend *legPH = new TLegend(.4, .7, .68, .98);
1030 legPH->SetBorderSize(1);legPH->SetFillColor(0);
1031 legPH->SetHeader("Particle Species");
1032 if(!(arr = (TObjArray*)(fContainer->At(kPH)))) break;
1033 if(!(h2 = (TProfile2D*)(arr->At(0)))) break;
1034
1035 TVirtualPad *pad = ((TVirtualPad*)l->At(0));pad->cd();
1036 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
1037 kFIRST = kTRUE;
1038 for(Int_t is=0; is<AliPID::kSPECIES; is++){
1039 Int_t bin = FindBin(is, 2.);
1040 h1 = h2->ProjectionY(Form("pyt%d", is), bin, bin);
1041 if(!h1->GetEntries()) continue;
1042 h1->SetMarkerStyle(24);
1043 h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
1044 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
1045 if(kFIRST){
1046 h1->GetXaxis()->SetTitle("t_{drift} [100*ns]");
1047 h1->GetYaxis()->SetTitle("<dQ/dt> [a.u.]");
1048 }
1049 h = (TH1F*)h1->DrawClone(kFIRST ? "c" : "samec");
1050 legPH->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "pl");
1051 kFIRST = kFALSE;
1052 }
1053
1054 pad = ((TVirtualPad*)l->At(1));pad->cd();
1055 pad->SetMargin(0.1, 0.01, 0.1, 0.01);
1056 if(!(h2 = (TProfile2D*)(arr->At(1)))) break;
1057 kFIRST = kTRUE;
1058 for(Int_t is=0; is<AliPID::kSPECIES; is++){
1059 Int_t bin = FindBin(is, 2.);
1060 h1 = h2->ProjectionY(Form("pyx%d", is), bin, bin);
1061 if(!h1->GetEntries()) continue;
1062 h1->SetMarkerStyle(24);
1063 h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
1064 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
1065 if(kFIRST){
1066 h1->GetXaxis()->SetTitle("x_{drift} [cm]");
1067 h1->GetYaxis()->SetTitle("<dQ/dl> [a.u./cm]");
1068 }
1069 h = (TH1F*)h1->DrawClone(kFIRST ? "c" : "samec");
1070 kFIRST = kFALSE;
1071 }
1072
1073 if(kFIRST) break;
1074 legPH->Draw();
1075 gPad->SetLogy(0);
1076 gPad->SetLogx(0);
1077 gPad->SetGridy();
1078 gPad->SetGridx();
1079 return kTRUE;
1080 }
1081 case kNClus:{
1082 // save 2.0 GeV projection as reference
1083 TLegend *legNClus = new TLegend(.13, .7, .4, .98);
1084 legNClus->SetBorderSize(1);
1085 legNClus->SetFillColor(0);
1086
1087 kFIRST = kTRUE;
1088 if(!(h2 = (TH2F*)(fContainer->At(kNClus)))) break;
1089 legNClus->SetHeader("Particle Species");
1090 for(Int_t is=0; is<AliPID::kSPECIES; is++){
1091 Int_t bin = FindBin(is, 2.);
1092 h1 = h2->ProjectionY(Form("pyNClus%d", is), bin, bin);
1093 if(!h1->GetEntries()) continue;
1094 h1->Scale(100./h1->Integral());
1095 //h1->SetMarkerStyle(24);
1096 //h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
1097 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
1098 if(kFIRST){
1099 h1->GetXaxis()->SetTitle("N^{cl}/tracklet");
1100 h1->GetYaxis()->SetTitle("Prob. [%]");
1101 h = (TH1F*)h1->DrawClone("c");
94f7dff7 1102 h->SetMaximum(20.);
1ee39b3a 1103 h->GetXaxis()->SetRangeUser(0., 35.);
1104 kFIRST = kFALSE;
1105 } else h = (TH1F*)h1->DrawClone("samec");
1106
1107 legNClus->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l");
1108 }
1109 if(kFIRST) break;
1110 legNClus->Draw();
1111 gPad->SetLogy(0);
1112 gPad->SetLogx(0);
1113 gPad->SetGridy();
1114 gPad->SetGridx();
1115 return kTRUE;
1116 }
1117 case kMomentum:
1118 case kMomentumBin:
1119 break;
1120 case kNTracklets:{
1121 TLegend *legNClus = new TLegend(.4, .7, .68, .98);
1122 legNClus->SetBorderSize(1);
1123 kFIRST = kTRUE;
1124 if(!(h2 = (TH2F*)(fContainer->At(kNTracklets)))) break;
1125 legNClus->SetHeader("Particle Species");
1126 for(Int_t is=0; is<AliPID::kSPECIES; is++){
1127 Int_t bin = FindBin(is, 2.);
1128 h1 = h2->ProjectionY(Form("pyNTracklets%d", is), bin, bin);
1129 if(!h1->GetEntries()) continue;
1130 h1->Scale(100./h1->Integral());
1131 //h1->SetMarkerStyle(24);
1132 //h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
1133 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
1134 if(kFIRST){
1135 h1->GetXaxis()->SetTitle("N^{trklt}/track");
1136 h1->GetXaxis()->SetRangeUser(1.,6.);
1137 h1->GetYaxis()->SetTitle("Prob. [%]");
1138 }
1139 h = (TH1F*)h1->DrawClone(kFIRST ? "c" : "samec");
1140 legNClus->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l");
1141 kFIRST = kFALSE;
1142 }
1143 if(kFIRST) break;
1144 legNClus->Draw();
1145 gPad->SetLogy(0);
1146 gPad->SetLogx(0);
1147 gPad->SetGridy();
1148 gPad->SetGridx();
1149 return kTRUE;
1150 }
1151 }
1152 AliInfo(Form("Reference plot [%d] missing result", ifig));
1153 return kFALSE;
1154}
1155
1156//________________________________________________________________________
1157Bool_t AliTRDcheckPID::PostProcess()
1158{
1159 // Draw result to the screen
1160 // Called once at the end of the query
1161
1162 if (!fContainer) {
1163 Printf("ERROR: list not available");
1164 return kFALSE;
1165 }
2e6aa247 1166
1167 TObjArray *eff(NULL);
1ee39b3a 1168 if(!fGraph){
2e6aa247 1169 fGraph = new TObjArray(2*AliPID::kSPECIES);
1ee39b3a 1170 fGraph->SetOwner();
2e6aa247 1171
1172 if(!(eff = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){
1173 AliError("Efficiency container for Electrons missing.");
1174 return kFALSE;
1175 }
1176 EvaluateEfficiency(eff, fGraph, AliPID::kPion, 0.9);
1177 EvaluateEfficiency(eff, fGraph, AliPID::kKaon, 0.9);
1178 EvaluateEfficiency(eff, fGraph, AliPID::kProton, 0.9);
1ee39b3a 1179 }
2e6aa247 1180 fNRefFigures = 12;
1ee39b3a 1181 return kTRUE;
1182}
1183
1184//________________________________________________________________________
64d57299 1185void AliTRDcheckPID::EvaluateEfficiency(const TObjArray * const histoContainer, TObjArray *results, Int_t species, Float_t electronEfficiency){
1ee39b3a 1186// Process PID information for pion efficiency
1187
1188 fUtil->SetElectronEfficiency(electronEfficiency);
1189
2e6aa247 1190
1191 const Int_t kNmethodsPID=Int_t(sizeof(fgMethod)/sizeof(Char_t*));
1192 Color_t colors[kNmethodsPID] = {kBlue, kGreen+2, kRed};
1193 Int_t markerStyle[kNmethodsPID] = {7, 7, 24};
1ee39b3a 1194 // efficiency graphs
2e6aa247 1195 TGraphErrors *g(NULL);
1196 TObjArray *eff = new TObjArray(kNmethodsPID); eff->SetOwner(); eff->SetName(Form("Eff_%s", AliTRDCalPID::GetPartName(species)));
1197 results->AddAt(eff, species);
1198 for(Int_t iMethod = 0; iMethod < kNmethodsPID; iMethod++){
1199 eff->AddAt(g = new TGraphErrors(), iMethod);
1200 g->SetName(Form("%s", fgMethod[iMethod]));
1ee39b3a 1201 g->SetLineColor(colors[iMethod]);
1202 g->SetMarkerColor(colors[iMethod]);
1203 g->SetMarkerStyle(markerStyle[iMethod]);
1204 }
1205
2e6aa247 1206 // Threshold graphs if not already
1207 TObjArray *thres(NULL);
1208 if(!(results->At(AliPID::kSPECIES))){
1209 thres = new TObjArray(kNmethodsPID); thres->SetOwner();
1210 thres->SetName("Thres");
1211 results->AddAt(thres, AliPID::kSPECIES);
1212 for(Int_t iMethod = 0; iMethod < kNmethodsPID; iMethod++){
1213 thres->AddAt(g = new TGraphErrors(), iMethod);
1214 g->SetName(Form("%s", fgMethod[iMethod]));
1215 g->SetLineColor(colors[iMethod]);
1216 g->SetMarkerColor(colors[iMethod]);
1217 g->SetMarkerStyle(markerStyle[iMethod]);
1218 }
1ee39b3a 1219 }
1ee39b3a 1220
2e6aa247 1221 TH2F *hPtr[kNmethodsPID]={
1222 (TH2F*)histoContainer->At(AliTRDpidUtil::kLQ),
1223 (TH2F*)histoContainer->At(AliTRDpidUtil::kNN),
1224 (TH2F*)histoContainer->At(AliTRDpidUtil::kESD)
1225 };
1ee39b3a 1226 for(Int_t iMom = 0; iMom < fMomentumAxis->GetNbins(); iMom++){
2e6aa247 1227 Float_t mom(fMomentumAxis->GetBinCenter(iMom+1));
1ee39b3a 1228
2e6aa247 1229 Int_t binEl(fMomentumAxis->GetNbins() * AliPID::kElectron + iMom + 1),
1230 binXX(fMomentumAxis->GetNbins() * species + iMom + 1);
1231
1232 for(Int_t iMethod = 0; iMethod < kNmethodsPID; iMethod++){
1233 // Calculate the Species Efficiency at electronEfficiency% electron efficiency for each Method
1234
1235 TH1D *histo1 = hPtr[iMethod] -> ProjectionY(Form("%s_el", fgMethod[iMethod]), binEl, binEl);
1236 TH1D *histo2 = hPtr[iMethod] -> ProjectionY(Form("%s_%s", fgMethod[iMethod], AliTRDCalPID::GetPartName(species)), binXX, binXX);
1ee39b3a 1237
1238 if(!fUtil->CalculatePionEffi(histo1, histo2)) continue;
1239
2e6aa247 1240 g=(TGraphErrors*)eff->At(iMethod);
1241 g->SetPoint(iMom, mom, 1.e2*fUtil->GetPionEfficiency());
1242 g->SetPointError(iMom, 0., 1.e2*fUtil->GetError());
1243 AliDebug(2, Form("%s Efficiency for %s is : %f +/- %f", AliTRDCalPID::GetPartName(species), fgMethod[iMethod], fUtil->GetPionEfficiency(), fUtil->GetError()));
1244
1245 if(!thres) continue;
1246 g=(TGraphErrors*)thres->At(iMethod);
1247 g->SetPoint(iMom, mom, fUtil->GetThreshold());
1248 g->SetPointError(iMom, 0., 0.);
1ee39b3a 1249 }
1250 }
1251}