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