TRD/qaRec/AliTRDpidChecker.cxx

   1 #include "TROOT.h"
   2 #include "TPDGCode.h"
   3 #include "TCanvas.h"
   4 #include "TF1.h"
   5 #include "TH1F.h"
   6 #include "TH1D.h"
   7 #include "TH2F.h"
   8 #include "TProfile.h"
   9 #include "TProfile2D.h"
  10 #include "TGraph.h"
  11 #include "TGraphErrors.h"
  12 #include "TLegend.h"
  13
  14 #include <TClonesArray.h>
  15 #include <TObjArray.h>
  16 #include <TList.h>
  17
  18 #include "AliESDEvent.h"
  19 #include "AliESDInputHandler.h"
  20 #include "AliTrackReference.h"
  21
  22 #include "AliAnalysisTask.h"
  23
  24 #include "AliTRDtrackerV1.h"
  25 #include "AliTRDtrackV1.h"
  26 #include "AliTRDcluster.h"
  27 #include "AliTRDReconstructor.h"
  28 #include "AliCDBManager.h"
  29 #include "AliTRDpidUtil.h"
  30
  31 #include "AliTRDpidChecker.h"
  32 #include "AliTRDtrackInfo/AliTRDtrackInfo.h"
  33
  34 // calculate pion efficiency at 90% electron efficiency for 11 momentum bins
  35 // this task should be used with simulated data only
  36
  37 ClassImp(AliTRDpidChecker)
  38
  39 //________________________________________________________________________
  40 AliTRDpidChecker::AliTRDpidChecker()
  41   :AliTRDrecoTask("PID", "PID Checker")
  42   ,fReconstructor(0x0)
  43   ,fUtil(0x0)
  44   ,fGraph(0x0)
  45   ,fEfficiency(0x0)
  46 {
  47   //
  48   // Default constructor
  49   //
  50
  51   fReconstructor = new AliTRDReconstructor();
  52   fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
  53
  54   fUtil = new AliTRDpidUtil();
  55
  56   InitFunctorList();
  57 }
  58
  59
  60 //________________________________________________________________________
  61 AliTRDpidChecker::~AliTRDpidChecker()
  62 {
  63  if(fGraph){fGraph->Delete(); delete fGraph;}
  64  if(fReconstructor) delete fReconstructor;
  65   if(fUtil) delete fUtil;
  66 }
  67
  68
  69 //________________________________________________________________________
  70 void AliTRDpidChecker::CreateOutputObjects()
  71 {
  72   // Create histograms
  73   // Called once
  74
  75   OpenFile(0, "RECREATE");
  76   fContainer = Histos();
  77 }
  78
  79
  80 //_______________________________________________________
  81 TObjArray * AliTRDpidChecker::Histos(){
  82
  83   //
  84   // Create QA histograms
  85   //
  86   if(fContainer) return fContainer;
  87
  88   Int_t xBins = AliPID::kSPECIES*AliTRDCalPID::kNMom;
  89   fContainer = new TObjArray(); fContainer->Expand(7);
  90
  91   const Float_t epsilon = 1./(2*(AliTRDpidUtil::kBins-1));     // get nice histos with bin center at 0 and 1
  92
  93   // histos of the electron probability of all 5 particle species and 11 momenta for the 2-dim LQ method
  94   fEfficiency = new TObjArray(); fEfficiency->Expand(3);
  95   fContainer->AddAt(fEfficiency, kEfficiency);
  96
  97   TH1 *h = 0x0;
  98   if(!(h = (TH2F*)gROOT->FindObject("PID_LQ"))){
  99     h = new TH2F("PID_LQ", "", xBins, -0.5, xBins - 0.5,
 100       AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon);
 101   } else h->Reset();
 102   fEfficiency->AddAt(h, kLQ);
 103
 104   // histos of the electron probability of all 5 particle species and 11 momenta for the neural network method
 105   if(!(h = (TH2F*)gROOT->FindObject("PID_NN"))){
 106     h = new TH2F("PID_NN", "",
 107       xBins, -0.5, xBins - 0.5,
 108       AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon);
 109   } else h->Reset();
 110   fEfficiency->AddAt(h, kNN);
 111
 112   // histos of the electron probability of all 5 particle species and 11 momenta for the ESD output
 113   if(!(h = (TH2F*)gROOT->FindObject("PID_ESD"))){
 114     h = new TH2F("PID_ESD", "",
 115       xBins, -0.5, xBins - 0.5,
 116       AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon);
 117   } else h->Reset();
 118   fEfficiency->AddAt(h, kESD);
 119
 120   // histos of the dE/dx distribution for all 5 particle species and 11 momenta
 121   if(!(h = (TH2F*)gROOT->FindObject("dEdx"))){
 122     h = new TH2F("dEdx", "",
 123       xBins, -0.5, xBins - 0.5,
 124       200, 0, 10000);
 125   } else h->Reset();
 126   fContainer->AddAt(h, kdEdx);
 127
 128   // histos of the dE/dx slices for all 5 particle species and 11 momenta
 129   if(!(h = (TH2F*)gROOT->FindObject("dEdxSlice"))){
 130     h = new TH2F("dEdxSlice", "",
 131       xBins*AliTRDReconstructor::kLQslices, -0.5, xBins*AliTRDReconstructor::kLQslices - 0.5,
 132       200, 0, 5000);
 133   } else h->Reset();
 134   fContainer->AddAt(h, kdEdxSlice);
 135
 136   // histos of the pulse height distribution for all 5 particle species and 11 momenta
 137   if(!(h = (TH2F*)gROOT->FindObject("PH"))){
 138     h = new TProfile2D("PH", "",
 139       xBins, -0.5, xBins - 0.5,
 140       AliTRDtrackerV1::GetNTimeBins(), -0.5, AliTRDtrackerV1::GetNTimeBins() - 0.5);
 141   } else h->Reset();
 142   fContainer->AddAt(h, kPH);
 143
 144   // histos of the number of clusters distribution for all 5 particle species and 11 momenta
 145   if(!(h = (TH2F*)gROOT->FindObject("NClus"))){
 146     h = new TH2F("NClus", "",
 147       xBins, -0.5, xBins - 0.5,
 148       AliTRDtrackerV1::GetNTimeBins(), -0.5, AliTRDtrackerV1::GetNTimeBins() - 0.5);
 149   } else h->Reset();
 150   fContainer->AddAt(h, kNClus);
 151
 152
 153   // momentum distributions - absolute and in momentum bins
 154   if(!(h = (TH1F*)gROOT->FindObject("hMom"))){
 155     h = new TH1F("hMom", "momentum distribution", 100, 0., 12.);
 156   } else h->Reset();
 157   fContainer->AddAt(h, kMomentum);
 158
 159   if(!(h = (TH1F*)gROOT->FindObject("hMomBin"))){
 160     h = new TH1F("hMomBin", "momentum distribution in momentum bins", AliTRDCalPID::kNMom, 0.5, 11.5);
 161   } else h->Reset();
 162   fContainer->AddAt(h, kMomentumBin);
 163
 164
 165   return fContainer;
 166 }
 167
 168
 169 //________________________________________________________________________
 170 Bool_t AliTRDpidChecker::CheckTrackQuality(const AliTRDtrackV1* track)
 171 {
 172   //
 173   // Check if the track is ok for PID
 174   //
 175
 176   if(track->GetNumberOfTracklets() == AliTRDgeometry::kNlayer) return 1;
 177 //   if(!fESD)
 178 //     return 0;
 179
 180   return 0;
 181 }
 182
 183
 184 //________________________________________________________________________
 185 Int_t AliTRDpidChecker::CalcPDG(AliTRDtrackV1* track)
 186 {
 187
 188   track -> SetReconstructor(fReconstructor);
 189
 190   fReconstructor -> SetOption("nn");
 191   track -> CookPID();
 192
 193   if(track -> GetPID(AliPID::kElectron) > track -> GetPID(AliPID::kMuon) + track -> GetPID(AliPID::kPion)  + track -> GetPID(AliPID::kKaon) + track -> GetPID(AliPID::kProton)){
 194     return kElectron;
 195   }
 196   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)){
 197     return kProton;
 198   }
 199   else if(track -> GetPID(AliPID::kKaon) > track -> GetPID(AliPID::kMuon)  && track -> GetPID(AliPID::kKaon) > track -> GetPID(AliPID::kPion)){
 200     return kKPlus;
 201   }
 202   else if(track -> GetPID(AliPID::kMuon) > track -> GetPID(AliPID::kPion)){
 203     return kMuonPlus;
 204   }
 205   else{
 206     return kPiPlus;
 207   }
 208 }
 209
 210
 211 //_______________________________________________________
 212 TH1 *AliTRDpidChecker::PlotLQ(const AliTRDtrackV1 *track)
 213 {
 214   //
 215   // Plot the probabilities for electrons using 2-dim LQ
 216   //
 217
 218   if(track) fTrack = track;
 219   if(!fTrack){
 220     AliWarning("No Track defined.");
 221     return 0x0;
 222   }
 223
 224   if(!CheckTrackQuality(fTrack)) return 0x0;
 225
 226   if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){
 227     AliWarning("No Histogram defined.");
 228     return 0x0;
 229   }
 230   TH2F *hPIDLQ = 0x0;
 231   if(!(hPIDLQ = dynamic_cast<TH2F *>(fEfficiency->At(kLQ)))){
 232     AliWarning("No Histogram defined.");
 233     return 0x0;
 234   }
 235
 236   AliTRDtrackV1 cTrack(*fTrack);
 237   cTrack.SetReconstructor(fReconstructor);
 238
 239   Int_t pdg = 0;
 240   Float_t momentum = 0.;
 241
 242   if(fMC){
 243     if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
 244     pdg = fMC->GetPDG();
 245   } else{
 246     //AliWarning("No MC info available!");
 247     momentum = cTrack.GetMomentum(0);
 248     pdg = CalcPDG(&cTrack);
 249   }
 250   if(momentum < 0.4) return 0x0;;
 251   if(momentum > 12.) return 0x0;;
 252
 253   fReconstructor -> SetOption("!nn");
 254   cTrack.CookPID();
 255   Int_t iMomBin = fUtil->GetMomentumBin(momentum);
 256
 257   switch(pdg){
 258   case kElectron:
 259   case kPositron:
 260     hPIDLQ -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron));
 261     break;
 262   case kMuonPlus:
 263   case kMuonMinus:
 264     hPIDLQ -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, cTrack .GetPID(AliPID::kElectron));
 265     break;
 266   case kPiPlus:
 267   case kPiMinus:
 268     hPIDLQ -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, cTrack .GetPID(AliPID::kElectron));
 269     break;
 270   case kKPlus:
 271   case kKMinus:
 272     hPIDLQ -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, cTrack .GetPID(AliPID::kElectron));
 273     break;
 274   case kProton:
 275   case kProtonBar:
 276     hPIDLQ -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron));
 277     break;
 278   }
 279
 280   return hPIDLQ;
 281 }
 282
 283
 284 //_______________________________________________________
 285 TH1 *AliTRDpidChecker::PlotNN(const AliTRDtrackV1 *track)
 286 {
 287   //
 288   // Plot the probabilities for electrons using 2-dim LQ
 289   //
 290
 291   if(track) fTrack = track;
 292   if(!fTrack){
 293     AliWarning("No Track defined.");
 294     return 0x0;
 295   }
 296
 297   if(!CheckTrackQuality(fTrack)) return 0x0;
 298
 299   if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){
 300     AliWarning("No Histogram defined.");
 301     return 0x0;
 302   }
 303   TH2F *hPIDNN;
 304   if(!(hPIDNN = dynamic_cast<TH2F *>(fEfficiency->At(kNN)))){
 305     AliWarning("No Histogram defined.");
 306     return 0x0;
 307   }
 308
 309
 310   AliTRDtrackV1 cTrack(*fTrack);
 311   cTrack.SetReconstructor(fReconstructor);
 312
 313   Int_t pdg = 0;
 314   Float_t momentum = 0.;
 315   if(fMC){
 316     if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
 317     pdg = fMC->GetPDG();
 318   } else {
 319     //AliWarning("No MC info available!");
 320     momentum = cTrack.GetMomentum(0);
 321     pdg = CalcPDG(&cTrack);
 322   }
 323   if(momentum < 0.4) return 0x0;;
 324   if(momentum > 12.) return 0x0;;
 325
 326   fReconstructor -> SetOption("nn");
 327   cTrack.CookPID();
 328   Int_t iMomBin = fUtil -> GetMomentumBin(momentum);
 329
 330   switch(pdg){
 331   case kElectron:
 332   case kPositron:
 333     hPIDNN -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron));
 334     break;
 335   case kMuonPlus:
 336   case kMuonMinus:
 337     hPIDNN -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron));
 338     break;
 339   case kPiPlus:
 340   case kPiMinus:
 341     hPIDNN -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron));
 342     break;
 343   case kKPlus:
 344   case kKMinus:
 345     hPIDNN -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron));
 346     break;
 347   case kProton:
 348   case kProtonBar:
 349     hPIDNN -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron));
 350     break;
 351   }
 352   return hPIDNN;
 353 }
 354
 355
 356 //_______________________________________________________
 357 TH1 *AliTRDpidChecker::PlotESD(const AliTRDtrackV1 *track)
 358 {
 359   //
 360   // Plot the probabilities for electrons using 2-dim LQ
 361   //
 362
 363   if(!fESD){
 364     AliWarning("No ESD info available.");
 365     return 0x0;
 366   }
 367
 368   if(track) fTrack = track;
 369   if(!fTrack){
 370     AliWarning("No Track defined.");
 371     return 0x0;
 372   }
 373
 374   if(!CheckTrackQuality(fTrack)) return 0x0;
 375
 376   if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){
 377     AliWarning("No Histogram defined.");
 378     return 0x0;
 379   }
 380   TH2F *hPIDESD = 0x0;
 381   if(!(hPIDESD = dynamic_cast<TH2F *>(fEfficiency->At(kESD)))){
 382     AliWarning("No Histogram defined.");
 383     return 0x0;
 384   }
 385
 386
 387   Int_t pdg = 0;
 388   Float_t momentum = 0.;
 389   if(fMC){
 390     if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
 391     pdg = fMC->GetPDG();
 392   } else {
 393     //AliWarning("No MC info available!");
 394     AliTRDtrackV1 cTrack(*fTrack);
 395     cTrack.SetReconstructor(fReconstructor);
 396     momentum = cTrack.GetMomentum(0);
 397     pdg = CalcPDG(&cTrack);
 398   }
 399   if(momentum < 0.4) return 0x0;;
 400   if(momentum > 12.) return 0x0;;
 401
 402
 403   Int_t iMomBin = fUtil->GetMomentumBin(momentum);
 404
 405
 406 //   Double32_t pidESD[AliPID::kSPECIES];
 407   const Double32_t *pidESD = fESD->GetResponseIter();
 408
 409   switch(pdg){
 410   case kElectron:
 411   case kPositron:
 412     hPIDESD -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, pidESD[0]);
 413     break;
 414   case kMuonPlus:
 415   case kMuonMinus:
 416     hPIDESD -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, pidESD[0]);
 417     break;
 418   case kPiPlus:
 419   case kPiMinus:
 420     hPIDESD -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, pidESD[0]);
 421     break;
 422   case kKPlus:
 423   case kKMinus:
 424     hPIDESD -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, pidESD[0]);
 425     break;
 426   case kProton:
 427   case kProtonBar:
 428     hPIDESD -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, pidESD[0]);
 429     break;
 430   }
 431   return hPIDESD;
 432 }
 433
 434
 435 //_______________________________________________________
 436 TH1 *AliTRDpidChecker::PlotdEdx(const AliTRDtrackV1 *track)
 437 {
 438   //
 439   // Plot the probabilities for electrons using 2-dim LQ
 440   //
 441
 442   if(track) fTrack = track;
 443   if(!fTrack){
 444     AliWarning("No Track defined.");
 445     return 0x0;
 446   }
 447
 448   if(!CheckTrackQuality(fTrack)) return 0x0;
 449
 450   TH2F *hdEdx;
 451   if(!(hdEdx = dynamic_cast<TH2F *>(fContainer->At(kdEdx)))){
 452     AliWarning("No Histogram defined.");
 453     return 0x0;
 454   }
 455
 456
 457   Int_t pdg = 0;
 458   Float_t momentum = 0.;
 459   if(fMC){
 460     if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
 461     pdg = fMC->GetPDG();
 462   } else {
 463     //AliWarning("No MC info available!");
 464     AliTRDtrackV1 cTrack(*fTrack);
 465     cTrack.SetReconstructor(fReconstructor);
 466     momentum = cTrack.GetMomentum(0);
 467     pdg = CalcPDG(&cTrack);
 468   }
 469   if(momentum < 0.4) return 0x0;;
 470   if(momentum > 12.) return 0x0;;
 471
 472   Int_t iMomBin = fUtil -> GetMomentumBin(momentum);
 473
 474
 475
 476   Float_t SumdEdx[AliTRDgeometry::kNlayer];
 477   AliTRDseedV1 *TRDtracklet[AliTRDgeometry::kNlayer];
 478   for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = 0x0;
 479   for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = fTrack -> GetTracklet(iChamb);
 480
 481   Float_t *fdEdx;
 482   Float_t dEdxSlice[AliTRDgeometry::kNlayer][AliTRDReconstructor::kLQslices];
 483
 484   for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
 485     SumdEdx[iChamb] = 0.;
 486     fdEdx = TRDtracklet[iChamb] -> GetdEdx();
 487     SumdEdx[iChamb] += fdEdx[0] + fdEdx[1] + fdEdx[2];
 488     for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){
 489       dEdxSlice[iChamb][iSlice] = fdEdx[iSlice];
 490     }
 491   }
 492
 493   switch(pdg){
 494   case kElectron:
 495   case kPositron:
 496     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 497       hdEdx -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]);
 498     break;
 499   case kMuonPlus:
 500   case kMuonMinus:
 501     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 502       hdEdx -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]);
 503     break;
 504   case kPiPlus:
 505   case kPiMinus:
 506     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 507       hdEdx -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]);
 508     break;
 509   case kKPlus:
 510   case kKMinus:
 511     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 512       hdEdx -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]);
 513     break;
 514   case kProton:
 515   case kProtonBar:
 516     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 517       hdEdx -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]);
 518     break;
 519   }
 520
 521   return hdEdx;
 522 }
 523
 524
 525 //_______________________________________________________
 526 TH1 *AliTRDpidChecker::PlotdEdxSlice(const AliTRDtrackV1 *track)
 527 {
 528   //
 529   // Plot the probabilities for electrons using 2-dim LQ
 530   //
 531
 532   if(track) fTrack = track;
 533   if(!fTrack){
 534     AliWarning("No Track defined.");
 535     return 0x0;
 536   }
 537
 538   if(!CheckTrackQuality(fTrack)) return 0x0;
 539
 540   TH2F *hdEdxSlice;
 541   if(!(hdEdxSlice = dynamic_cast<TH2F *>(fContainer->At(kdEdxSlice)))){
 542     AliWarning("No Histogram defined.");
 543     return 0x0;
 544   }
 545
 546   Int_t pdg = 0;
 547   Float_t momentum = 0.;
 548   if(fMC){
 549     if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
 550     pdg = fMC->GetPDG();
 551   } else {
 552     //AliWarning("No MC info available!");
 553     AliTRDtrackV1 cTrack(*fTrack);
 554     cTrack.SetReconstructor(fReconstructor);
 555     momentum = cTrack.GetMomentum(0);
 556     pdg = CalcPDG(&cTrack);
 557   }
 558   if(momentum < 0.4) return 0x0;
 559   if(momentum > 12.) return 0x0;;
 560
 561   Int_t iMomBin = fUtil -> GetMomentumBin(momentum);
 562
 563
 564
 565   AliTRDseedV1 *TRDtracklet[AliTRDgeometry::kNlayer];
 566   for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = 0x0;
 567   for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = fTrack -> GetTracklet(iChamb);
 568
 569   Float_t *fdEdx;
 570   Float_t dEdxSlice[AliTRDgeometry::kNlayer][AliTRDReconstructor::kLQslices];
 571
 572   for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
 573     fdEdx = TRDtracklet[iChamb] -> GetdEdx();
 574     for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){
 575       dEdxSlice[iChamb][iSlice] = fdEdx[iSlice];
 576     }
 577   }
 578
 579   switch(pdg){
 580   case kElectron:
 581   case kPositron:
 582     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
 583       for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){
 584         hdEdxSlice -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice, dEdxSlice[iChamb][iSlice]);
 585       }
 586     }
 587     break;
 588   case kMuonPlus:
 589   case kMuonMinus:
 590     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
 591       for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){
 592         hdEdxSlice -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice,
 593         dEdxSlice[iChamb][iSlice]);
 594       }
 595     }
 596     break;
 597   case kPiPlus:
 598   case kPiMinus:
 599     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 600       for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++)
 601         hdEdxSlice -> Fill(AliPID::kPion * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice,
 602                            dEdxSlice[iChamb][iSlice]);
 603     break;
 604   case kKPlus:
 605   case kKMinus:
 606     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 607       for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++)
 608         hdEdxSlice -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice,
 609                            dEdxSlice[iChamb][iSlice]);
 610     break;
 611   case kProton:
 612   case kProtonBar:
 613     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 614       for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++)
 615         hdEdxSlice -> Fill(AliPID::kProton * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice,
 616                            dEdxSlice[iChamb][iSlice]);
 617     break;
 618   }
 619
 620   return hdEdxSlice;
 621
 622 }
 623
 624
 625 //_______________________________________________________
 626 TH1 *AliTRDpidChecker::PlotPH(const AliTRDtrackV1 *track)
 627 {
 628   //
 629   // Plot the probabilities for electrons using 2-dim LQ
 630   //
 631
 632   if(track) fTrack = track;
 633   if(!fTrack){
 634     AliWarning("No Track defined.");
 635     return 0x0;
 636   }
 637
 638   if(!CheckTrackQuality(fTrack)) return 0x0;
 639
 640   TProfile2D *hPH;
 641   if(!(hPH = dynamic_cast<TProfile2D *>(fContainer->At(kPH)))){
 642     AliWarning("No Histogram defined.");
 643     return 0x0;
 644   }
 645
 646   Int_t pdg = 0;
 647   Float_t momentum = 0.;
 648   if(fMC){
 649     if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
 650     pdg = fMC->GetPDG();
 651   } else {
 652     //AliWarning("No MC info available!");
 653     AliTRDtrackV1 cTrack(*fTrack);
 654     cTrack.SetReconstructor(fReconstructor);
 655     momentum = cTrack.GetMomentum(0);
 656     pdg = CalcPDG(&cTrack);
 657   }
 658   if(momentum < 0.4) return 0x0;;
 659   if(momentum > 12.) return 0x0;;
 660
 661   Int_t iMomBin = fUtil -> GetMomentumBin(momentum);
 662
 663   AliTRDseedV1 *TRDtracklet[AliTRDgeometry::kNlayer];
 664   for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = 0x0;
 665   for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = fTrack -> GetTracklet(iChamb);
 666
 667   AliTRDcluster *TRDcluster = 0x0;
 668
 669   switch(pdg){
 670   case kElectron:
 671   case kPositron:
 672     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
 673       for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){
 674         if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus)))
 675           continue;
 676         hPH -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus));
 677       }
 678     }
 679     break;
 680   case kMuonPlus:
 681   case kMuonMinus:
 682     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
 683       for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){
 684         if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus)))
 685           continue;
 686         hPH -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus));
 687       }
 688     }
 689     break;
 690   case kPiPlus:
 691   case kPiMinus:
 692     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
 693       for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){
 694         if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus)))
 695           continue;
 696         hPH -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus));
 697       }
 698     }
 699     break;
 700   case kKPlus:
 701   case kKMinus:
 702     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
 703       for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){
 704         if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus)))
 705           continue;
 706         hPH -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus));
 707       }
 708     }
 709     break;
 710   case kProton:
 711   case kProtonBar:
 712     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
 713       for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){
 714         if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus)))
 715           continue;
 716         hPH -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus));
 717       }
 718     }
 719     break;
 720   }
 721
 722   return hPH;
 723 }
 724
 725
 726 //_______________________________________________________
 727 TH1 *AliTRDpidChecker::PlotNClus(const AliTRDtrackV1 *track)
 728 {
 729   //
 730   // Plot the probabilities for electrons using 2-dim LQ
 731   //
 732
 733   if(track) fTrack = track;
 734   if(!fTrack){
 735     AliWarning("No Track defined.");
 736     return 0x0;
 737   }
 738
 739   if(!CheckTrackQuality(fTrack)) return 0x0;
 740
 741   TH2F *hNClus;
 742   if(!(hNClus = dynamic_cast<TH2F *>(fContainer->At(kNClus)))){
 743     AliWarning("No Histogram defined.");
 744     return 0x0;
 745   }
 746
 747
 748   Int_t pdg = 0;
 749   Float_t momentum = 0.;
 750   if(fMC){
 751     if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
 752     pdg = fMC->GetPDG();
 753   } else {
 754     //AliWarning("No MC info available!");
 755     AliTRDtrackV1 cTrack(*fTrack);
 756     cTrack.SetReconstructor(fReconstructor);
 757     momentum = cTrack.GetMomentum(0);
 758     pdg = CalcPDG(&cTrack);
 759   }
 760   if(momentum < 0.4) return 0x0;;
 761   if(momentum > 12.) return 0x0;;
 762
 763   Int_t iMomBin = fUtil -> GetMomentumBin(momentum);
 764
 765
 766   Int_t iNClus[AliTRDgeometry::kNlayer];
 767   memset(iNClus, 0, sizeof(Int_t) * AliTRDgeometry::kNlayer);
 768
 769   AliTRDseedV1 *TRDtracklet[AliTRDgeometry::kNlayer];
 770   for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = 0x0;
 771   for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
 772     TRDtracklet[iChamb] = fTrack -> GetTracklet(iChamb);
 773     iNClus[iChamb] = TRDtracklet[iChamb] -> GetN();
 774   }
 775
 776   switch(pdg){
 777   case kElectron:
 778   case kPositron:
 779     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 780       hNClus -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]);
 781     break;
 782   case kMuonPlus:
 783   case kMuonMinus:
 784     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 785       hNClus -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]);
 786     break;
 787   case kPiPlus:
 788   case kPiMinus:
 789     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 790       hNClus -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]);
 791     break;
 792   case kKPlus:
 793   case kKMinus:
 794     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 795       hNClus -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]);
 796     break;
 797   case kProton:
 798   case kProtonBar:
 799     for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
 800       hNClus -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]);
 801     break;
 802   }
 803
 804   return hNClus;
 805 }
 806
 807
 808 //_______________________________________________________
 809 TH1 *AliTRDpidChecker::PlotMom(const AliTRDtrackV1 *track)
 810 {
 811   //
 812   // Plot the probabilities for electrons using 2-dim LQ
 813   //
 814
 815   if(track) fTrack = track;
 816   if(!fTrack){
 817     AliWarning("No Track defined.");
 818     return 0x0;
 819   }
 820
 821   if(!CheckTrackQuality(fTrack)) return 0x0;
 822
 823   TH1F *hMom;
 824   if(!(hMom = dynamic_cast<TH1F *>(fContainer->At(kMomentum)))){
 825     AliWarning("No Histogram defined.");
 826     return 0x0;
 827   }
 828
 829
 830   Int_t pdg = 0;
 831   Float_t momentum = 0.;
 832   if(fMC){
 833     if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
 834     pdg = fMC->GetPDG();
 835   } else {
 836     //AliWarning("No MC info available!");
 837     AliTRDtrackV1 cTrack(*fTrack);
 838     cTrack.SetReconstructor(fReconstructor);
 839     momentum = cTrack.GetMomentum(0);
 840     pdg = CalcPDG(&cTrack);
 841   }
 842   if(momentum < 0.4) return 0x0;
 843   if(momentum > 12.) return 0x0;;
 844
 845   hMom -> Fill(momentum);
 846   return hMom;
 847 }
 848
 849
 850 //_______________________________________________________
 851 TH1 *AliTRDpidChecker::PlotMomBin(const AliTRDtrackV1 *track)
 852 {
 853   //
 854   // Plot the probabilities for electrons using 2-dim LQ
 855   //
 856
 857   if(track) fTrack = track;
 858   if(!fTrack){
 859     AliWarning("No Track defined.");
 860     return 0x0;
 861   }
 862
 863   if(!CheckTrackQuality(fTrack)) return 0x0;
 864
 865   TH1F *hMomBin;
 866   if(!(hMomBin = dynamic_cast<TH1F *>(fContainer->At(kMomentumBin)))){
 867     AliWarning("No Histogram defined.");
 868     return 0x0;
 869   }
 870
 871
 872   Int_t pdg = 0;
 873   Float_t momentum = 0.;
 874
 875   if(fMC){
 876     if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P();
 877     pdg = fMC->GetPDG();
 878   } else {
 879     //AliWarning("No MC info available!");
 880     AliTRDtrackV1 cTrack(*fTrack);
 881     cTrack.SetReconstructor(fReconstructor);
 882     momentum = cTrack.GetMomentum(0);
 883     pdg = CalcPDG(&cTrack);
 884   }
 885   if(momentum < 0.4) return 0x0;
 886   if(momentum > 12.) return 0x0;;
 887
 888   Int_t iMomBin = fUtil -> GetMomentumBin(momentum);
 889   hMomBin -> Fill(iMomBin);
 890   return hMomBin;
 891 }
 892
 893
 894 //________________________________________________________
 895 Bool_t AliTRDpidChecker::GetRefFigure(Int_t ifig)
 896 {
 897   Bool_t FIRST = kTRUE;
 898   TLegend *leg = new TLegend(.7, .7, .98, .98);
 899   leg->SetBorderSize(1);
 900   TGraphErrors *g = 0x0;
 901   TAxis *ax = 0x0;
 902   TH1 *h1 = 0x0, *h=0x0;
 903   TH2 *h2 = 0x0;
 904   switch(ifig){
 905   case kEfficiency:
 906     if(!(g = (TGraphErrors*)fGraph->At(kLQ))) break;
 907     if(!g->GetN()) break;
 908     leg->SetHeader("PID Method");
 909     g->Draw("apl");
 910     ax = g->GetHistogram()->GetXaxis();
 911     ax->SetTitle("p [GeV/c]");
 912     ax->SetRangeUser(.6, 10.5);
 913     ax->SetMoreLogLabels();
 914     ax = g->GetHistogram()->GetYaxis();
 915     ax->SetTitle("#epsilon_{#pi} [%]");
 916     ax->SetRangeUser(1.e-3, 1.e-1);
 917     leg->AddEntry(g, "2D LQ", "pl");
 918     if(! (g = (TGraphErrors*)fGraph->At(kNN))) break;
 919     g->Draw("pl");
 920     leg->AddEntry(g, "NN", "pl");
 921     if(! (g = (TGraphErrors*)fGraph->At(kESD))) break;
 922     g->Draw("p");
 923     leg->AddEntry(g, "ESD", "pl");
 924     leg->Draw();
 925     gPad->SetLogy();
 926     gPad->SetLogx();
 927     gPad->SetGridy();
 928     gPad->SetGridx();
 929     return kTRUE;
 930   case kdEdx:
 931     // save 2.0 GeV projection as reference
 932     FIRST = kTRUE;
 933     if(!(h2 = (TH2F*)(fContainer->At(kdEdx)))) break;
 934     leg->SetHeader("Particle Species");
 935     for(Int_t is = AliPID::kSPECIES-1; is>=0; is--){
 936       Int_t bin = is*AliTRDCalPID::kNMom+4;
 937       h1 = h2->ProjectionY("px", bin, bin);
 938       if(!h1->GetEntries()) continue;
 939       h1->Scale(1./h1->Integral());
 940       h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
 941       h = (TH1F*)h1->DrawClone(FIRST ? "c" : "samec");
 942       leg->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l");
 943       FIRST = kFALSE;
 944     }
 945     if(FIRST) break;
 946     leg->Draw();
 947     gPad->SetLogy();
 948     gPad->SetLogx(0);
 949     gPad->SetGridy();
 950     gPad->SetGridx();
 951     return kTRUE;
 952   case kdEdxSlice:
 953     break;
 954   case kPH:
 955     // save 2.0 GeV projection as reference
 956     FIRST = kTRUE;
 957     if(!(h2 = (TH2F*)(fContainer->At(kPH)))) break;;
 958     leg->SetHeader("Particle Species");
 959     for(Int_t is=0; is<AliPID::kSPECIES; is++){
 960       Int_t bin = is*AliTRDCalPID::kNMom+4;
 961       h1 = h2->ProjectionY("py", bin, bin);
 962       if(!h1->GetEntries()) continue;
 963       h1->SetMarkerStyle(24);
 964       h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
 965       h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
 966       if(FIRST){
 967         h1->GetXaxis()->SetTitle("tb[1/100 ns^{-1}]");
 968         h1->GetYaxis()->SetTitle("<PH> [a.u.]");
 969       }
 970       h = (TH1F*)h1->DrawClone(FIRST ? "c" : "samec");
 971       leg->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "pl");
 972       FIRST = kFALSE;
 973     }
 974     if(FIRST) break;
 975     leg->Draw();
 976     gPad->SetLogy(0);
 977     gPad->SetLogx(0);
 978     gPad->SetGridy();
 979     gPad->SetGridx();
 980     return kTRUE;
 981   case kNClus:
 982     // save 2.0 GeV projection as reference
 983     FIRST = kTRUE;
 984     if(!(h2 = (TH2F*)(fContainer->At(kNClus)))) break;
 985     leg->SetHeader("Particle Species");
 986     for(Int_t is=0; is<AliPID::kSPECIES; is++){
 987       Int_t bin = is*AliTRDCalPID::kNMom+4;
 988       h1 = h2->ProjectionY("py", bin, bin);
 989       if(!h1->GetEntries()) continue;
 990       //h1->SetMarkerStyle(24);
 991       //h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
 992       h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
 993       if(FIRST) h1->GetXaxis()->SetTitle("N^{cl}/tracklet");
 994       h = (TH1F*)h1->DrawClone(FIRST ? "c" : "samec");
 995       leg->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l");
 996       FIRST = kFALSE;
 997     }
 998     if(FIRST) break;
 999     leg->Draw();
1000     gPad->SetLogy();
1001     gPad->SetLogx(0);
1002     gPad->SetGridy();
1003     gPad->SetGridx();
1004     return kTRUE;
1005   case kMomentum:
1006   case kMomentumBin:
1007     break;
1008   case kThresh:
1009     if(!(g = (TGraphErrors*)fGraph->At(kLQ+3))) break;
1010     if(!g->GetN()) break;
1011     leg->SetHeader("PID Method");
1012     g->Draw("apl");
1013     ax = g->GetHistogram()->GetXaxis();
1014     ax->SetTitle("p [GeV/c]");
1015     ax->SetRangeUser(.6, 10.5);
1016     ax->SetMoreLogLabels();
1017     ax = g->GetHistogram()->GetYaxis();
1018     ax->SetTitle("threshold");
1019     ax->SetRangeUser(5.e-2, 1.);
1020     leg->AddEntry(g, "2D LQ", "pl");
1021     if(!(g = (TGraphErrors*)fGraph->At(kNN+3))) break;
1022     g->Draw("pl");
1023     leg->AddEntry(g, "NN", "pl");
1024     if(!(g = (TGraphErrors*)fGraph->At(kESD+3))) break;
1025     g->Draw("p");
1026     leg->AddEntry(g, "ESD", "pl");
1027     leg->Draw();
1028     gPad->SetLogx();
1029     gPad->SetGridy();
1030     gPad->SetGridx();
1031     return kTRUE;
1032   }
1033   AliInfo(Form("Reference plot [%d] missing result", ifig));
1034   return kFALSE;
1035 }
1036
1037 //________________________________________________________________________
1038 Bool_t AliTRDpidChecker::PostProcess()
1039 {
1040   // Draw result to the screen
1041   // Called once at the end of the query
1042
1043   if (!fContainer) {
1044     Printf("ERROR: list not available");
1045     return kFALSE;
1046   }
1047   if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){
1048     AliError("Efficiency container missing.");
1049     return 0x0;
1050   }
1051
1052   TGraphErrors *g = 0x0;
1053   if(!fGraph){
1054     fGraph = new TObjArray(6);
1055     fGraph->SetOwner();
1056
1057     fGraph->AddAt(g = new TGraphErrors(), kLQ);
1058     g->SetLineColor(kBlue);
1059     g->SetMarkerColor(kBlue);
1060     g->SetMarkerStyle(7);
1061
1062     fGraph->AddAt(g = new TGraphErrors(), kNN);
1063     g->SetLineColor(kGreen);
1064     g->SetMarkerColor(kGreen);
1065     g->SetMarkerStyle(7);
1066
1067     fGraph -> AddAt(g = new TGraphErrors(), kESD);
1068     g->SetLineColor(kRed);
1069     g->SetMarkerColor(kRed);
1070     g->SetMarkerStyle(24);
1071
1072     fGraph->AddAt(g = new TGraphErrors(), 3+kLQ);
1073     g->SetLineColor(kBlue);
1074     g->SetMarkerColor(kBlue);
1075     g->SetMarkerStyle(7);
1076
1077     fGraph->AddAt(g = new TGraphErrors(), 3+kNN);
1078     g->SetLineColor(kGreen);
1079     g->SetMarkerColor(kGreen);
1080     g->SetMarkerStyle(7);
1081
1082     fGraph -> AddAt(g = new TGraphErrors(), 3+kESD);
1083     g->SetLineColor(kRed);
1084     g->SetMarkerColor(kRed);
1085     g->SetMarkerStyle(24);
1086   }
1087
1088   Float_t mom = 0.;
1089   TH1D *Histo1=0x0, *Histo2=0x0;
1090
1091   // calculate the pion efficiencies and the errors for 90% electron efficiency (2-dim LQ)
1092   TH2F *hPIDLQ = (TH2F*)fEfficiency->At(kLQ);
1093   for(Int_t iMom = 0; iMom < AliTRDCalPID::kNMom; iMom++){
1094     mom = AliTRDCalPID::GetMomentum(iMom);
1095
1096     Histo1 = hPIDLQ -> ProjectionY("LQ_ele",AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1,AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1);
1097     Histo2 = hPIDLQ -> ProjectionY("LQ_pio",AliTRDCalPID::kNMom*AliPID::kPion+iMom+1,AliTRDCalPID::kNMom*AliPID::kPion+iMom+1);
1098
1099     if(!fUtil->CalculatePionEffi(Histo1, Histo2)) continue;
1100
1101     g = (TGraphErrors*)fGraph->At(kLQ);
1102     g->SetPoint(iMom, mom, fUtil->GetPionEfficiency());
1103     g->SetPointError(iMom, 0., fUtil->GetError());
1104     g = (TGraphErrors*)fGraph->At(3 + kLQ);
1105     g->SetPoint(iMom, mom, fUtil->GetThreshold());
1106     g->SetPointError(iMom, 0., 0.);
1107
1108     if(fDebugLevel>=2) Printf("Pion Efficiency for 2-dim LQ is : %f +/- %f\n\n", fUtil->GetPionEfficiency(), fUtil->GetError());
1109   }
1110
1111
1112   // calculate the pion efficiencies and the errors for 90% electron efficiency (NN)
1113   TH2F *hPIDNN = (TH2F*)fEfficiency->At(kNN);
1114   for(Int_t iMom = 0; iMom < AliTRDCalPID::kNMom; iMom++){
1115     mom = AliTRDCalPID::GetMomentum(iMom);
1116
1117     Histo1 = hPIDNN -> ProjectionY("NN_ele",AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1,AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1);
1118     Histo2 = hPIDNN -> ProjectionY("NN_pio",AliTRDCalPID::kNMom*AliPID::kPion+iMom+1,AliTRDCalPID::kNMom*AliPID::kPion+iMom+1);
1119
1120     if(!fUtil -> CalculatePionEffi(Histo1, Histo2)) continue;
1121
1122     g = (TGraphErrors*)fGraph->At(kNN);
1123     g->SetPoint(iMom, mom, fUtil->GetPionEfficiency());
1124     g->SetPointError(iMom, 0., fUtil->GetError());
1125     g = (TGraphErrors*)fGraph->At(3+kNN);
1126     g->SetPoint(iMom, mom, fUtil->GetThreshold());
1127     g->SetPointError(iMom, 0., 0.);
1128
1129     if(fDebugLevel>=2) Printf("Pion Efficiency for NN is : %f +/- %f\n\n", fUtil->GetPionEfficiency(), fUtil->GetError());
1130   }
1131
1132
1133   // calculate the pion efficiencies and the errors for 90% electron efficiency (ESD)
1134   TH2F *hPIDESD = (TH2F*)fEfficiency->At(kESD);
1135   for(Int_t iMom = 0; iMom < AliTRDCalPID::kNMom; iMom++){
1136     mom = AliTRDCalPID::GetMomentum(iMom);
1137
1138     Histo1 = hPIDESD -> ProjectionY("NN_ele",AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1,AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1);
1139     Histo2 = hPIDESD -> ProjectionY("NN_pio",AliTRDCalPID::kNMom*AliPID::kPion+iMom+1,AliTRDCalPID::kNMom*AliPID::kPion+iMom+1);
1140
1141     if(!fUtil->CalculatePionEffi(Histo1, Histo2)) continue;
1142
1143     g = (TGraphErrors*)fGraph->At(kESD);
1144     g->SetPoint(iMom, mom, fUtil->GetPionEfficiency());
1145     g->SetPointError(iMom, 0., fUtil->GetError());
1146     g = (TGraphErrors*)fGraph->At(3+kESD);
1147     g->SetPoint(iMom, mom, fUtil->GetThreshold());
1148     g->SetPointError(iMom, 0., 0.);
1149
1150     if(fDebugLevel>=2) Printf("Pion Efficiency for ESD is : %f +/- %f\n\n", fUtil->GetPionEfficiency(), fUtil->GetError());
1151   }
1152
1153   fNRefFigures = 8;
1154
1155   return kTRUE;
1156 }
1157
1158
1159 //________________________________________________________________________
1160 void AliTRDpidChecker::Terminate(Option_t *)
1161 {
1162   // Draw result to the screen
1163   // Called once at the end of the query
1164
1165   fContainer = dynamic_cast<TObjArray*>(GetOutputData(0));
1166   if (!fContainer) {
1167     Printf("ERROR: list not available");
1168     return;
1169   }
1170 }
1171
1172