PWGPP/TRD/AliTRDcheckPID.cxx

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