PWG1/TRD/AliTRDcheckESD.cxx

   1 /**************************************************************************
   2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3 *                                                                        *
   4 * Author: The ALICE Off-line Project.                                    *
   5 * Contributors are mentioned in the code where appropriate.              *
   6 *                                                                        *
   7 * Permission to use, copy, modify and distribute this software and its   *
   8 * documentation strictly for non-commercial purposes is hereby granted   *
   9 * without fee, provided that the above copyright notice appears in all   *
  10 * copies and that both the copyright notice and this permission notice   *
  11 * appear in the supporting documentation. The authors make no claims     *
  12 * about the suitability of this software for any purpose. It is          *
  13 * provided "as is" without express or implied warranty.                  *
  14 **************************************************************************/
  15
  16 /////////////////////////////////////////////////////
  17 //
  18 // Check basic detector results at ESD level
  19 //   - Geometrical efficiency
  20 //   - Tracking efficiency
  21 //   - PID efficiency
  22 //   - Refit efficiency
  23 //
  24 // Author
  25 //   Alex Bercuci <A.Bercuci@gsi.de>
  26 //
  27 //////////////////////////////////////////////////////
  28
  29 #include <TClonesArray.h>
  30 #include <TCanvas.h>
  31 #include <TObjArray.h>
  32 #include <TPad.h>
  33 #include <TLegend.h>
  34 #include <TLatex.h>
  35 #include <TLine.h>
  36 #include <TF1.h>
  37 #include <TH1D.h>
  38 #include <TH2D.h>
  39 #include <TH3D.h>
  40 #include <TH2I.h>
  41 #include <TH2F.h>
  42 #include <TH3S.h>
  43 #include <TH3F.h>
  44 #include <TProfile2D.h>
  45 #include <TProfile.h>
  46 #include <TGraphErrors.h>
  47 #include <TGraphAsymmErrors.h>
  48 #include <TFile.h>
  49 #include <TTree.h>
  50 #include <TROOT.h>
  51 #include <TChain.h>
  52 #include <TParticle.h>
  53 #include <TTimeStamp.h>
  54 #include <TRandom.h>
  55
  56 #include "AliLog.h"
  57 #include "AliAnalysisManager.h"
  58 #include "AliAnalysisCuts.h"
  59 #include "AliESDEvent.h"
  60 #include "AliESDkink.h"
  61 #include "AliMCEvent.h"
  62 #include "AliESDInputHandler.h"
  63 #include "AliMCEventHandler.h"
  64 #include "AliESDpid.h"
  65
  66 #include "AliESDtrack.h"
  67 #include "AliMCParticle.h"
  68 #include "AliPID.h"
  69 #include "AliStack.h"
  70 #include "AliTrackReference.h"
  71 //#include "AliESDCentrality.h"
  72 #include "AliMultiplicity.h"
  73 #include "AliCFContainer.h"
  74
  75 #include "AliTRDcheckESD.h"
  76 #include <iostream>
  77 using namespace std;
  78
  79 ClassImp(AliTRDcheckESD)
  80
  81 const Float_t AliTRDcheckESD::fgkxTPC = 290.;
  82 const Float_t AliTRDcheckESD::fgkxTOF = 365.;
  83 const UChar_t AliTRDcheckESD::fgkNgraph[AliTRDcheckESD::kNrefs] ={
  84 8, 4, 2, 20};
  85 FILE* AliTRDcheckESD::fgFile = NULL;
  86
  87 const Float_t AliTRDcheckESD::fgkEvVertexZ = 15.;
  88 const Int_t   AliTRDcheckESD::fgkEvVertexN = 1;
  89 const Float_t AliTRDcheckESD::fgkTrkDCAxy  = 40.;
  90 const Float_t AliTRDcheckESD::fgkTrkDCAz   = 15.;
  91 const Int_t   AliTRDcheckESD::fgkNclTPC    = 100;
  92 const Float_t AliTRDcheckESD::fgkPt        = 0.2;
  93 const Float_t AliTRDcheckESD::fgkEta       = 0.9;
  94 const Float_t AliTRDcheckESD::fgkQs        = 0.002;
  95
  96 //____________________________________________________________________
  97 AliTRDcheckESD::AliTRDcheckESD():
  98   AliAnalysisTaskSE()
  99   ,fStatus(0)
 100   ,fNRefFigures(0)
 101   ,fESD(NULL)
 102   ,fMC(NULL)
 103   ,fESDpid(new AliESDpid)
 104   ,fHistos(NULL)
 105   ,fResults(NULL)
 106   ,fCfContainer(NULL)
 107   ,fReferenceTrackFilter(NULL)
 108 {
 109   //
 110   // Default constructor
 111   //
 112   SetNameTitle("TRDcheckESD", "Check TRD @ ESD level");
 113   SetMC(kTRUE);
 114 }
 115
 116 //____________________________________________________________________
 117 AliTRDcheckESD::AliTRDcheckESD(char* name):
 118   AliAnalysisTaskSE(name)
 119   ,fStatus(0)
 120   ,fNRefFigures(0)
 121   ,fESD(NULL)
 122   ,fMC(NULL)
 123   ,fESDpid(new AliESDpid)
 124   ,fHistos(NULL)
 125   ,fResults(NULL)
 126   ,fCfContainer(NULL)
 127   ,fReferenceTrackFilter(NULL)
 128 {
 129   //
 130   // Default constructor
 131   //
 132   SetMC(kTRUE);
 133   SetTitle("Check TRD @ ESD level");
 134   DefineOutput(1, TObjArray::Class());
 135 }
 136
 137 //____________________________________________________________________
 138 AliTRDcheckESD::~AliTRDcheckESD()
 139 {
 140 // Destructor
 141   if(fHistos){
 142     //fHistos->Delete();
 143     delete fHistos;
 144   }
 145   if(fResults){
 146     fResults->Delete();
 147     delete fResults;
 148   }
 149 }
 150
 151 //____________________________________________________________________
 152 void AliTRDcheckESD::UserCreateOutputObjects()
 153 {
 154   //
 155   // Create Output Containers (TObjectArray containing 1D histograms)
 156   //
 157   Histos();
 158   PostData(1, fHistos);
 159 }
 160
 161 //____________________________________________________________________
 162 void AliTRDcheckESD::MakeSummaryFromCF(Double_t* trendValues, Bool_t useIsolatedBC, Bool_t cutTOFbc){
 163   //
 164   // Draw summary plots for the ESDcheck task using the CF container
 165   //
 166
 167   cout << "Make summary from CF" << endl;
 168   TCanvas *cOut=0x0;
 169   if(gROOT->FindObject("trackingSummary")) delete gROOT->FindObject("trackingSummary");
 170   cOut = new TCanvas("trackingSummary", "Tracking summary for the ESD task", 1600, 1200);
 171   cOut->cd();
 172   PlotTrackingSummaryFromCF(0, trendValues, useIsolatedBC, cutTOFbc);
 173   cOut->SaveAs("trackingSummary.gif");
 174
 175   if(gROOT->FindObject("pidSummary")) delete gROOT->FindObject("pidSummary");
 176   cOut = new TCanvas("pidSummary", "PID summary for the ESD task", 1600, 1200);
 177   cOut->cd();
 178   //GetRefFigure(6);
 179   PlotPidSummaryFromCF(0, trendValues, useIsolatedBC, cutTOFbc);
 180   cOut->SaveAs("pidSummary.gif");
 181
 182   if(gROOT->FindObject("centSummary")) delete gROOT->FindObject("centSummary");
 183   cOut = new TCanvas("centSummary", "Centrality summary for the ESD task", 1600, 1200);
 184   cOut->cd();
 185   //GetRefFigure(7);
 186   PlotCentSummaryFromCF(trendValues, useIsolatedBC, cutTOFbc);
 187   cOut->SaveAs("centSummary.gif");
 188
 189 }
 190
 191
 192 //____________________________________________________________________
 193 void AliTRDcheckESD::MakeSummary(Double_t* trendValues){
 194   //
 195   // Draw summary plots for the ESDcheck task
 196   //
 197   // Old method to draw summary pictures from histograms. Use the MakeSummaryFromCF() when CF container is present
 198
 199   cout << "Make summary" << endl;
 200   TCanvas *cTracking=0x0;
 201   if(gROOT->FindObject("trackingSummary")) delete gROOT->FindObject("trackingSummary");
 202   cTracking = new TCanvas("trackingSummary", "Tracking summary for the ESD task", 1600, 1200);
 203   cTracking->cd();
 204   //GetRefFigure(5);
 205   if(PlotTrackingSummary(0, trendValues))
 206     cTracking->SaveAs("trackingSummary.gif");
 207
 208   TCanvas* cPid = 0x0;
 209   if(gROOT->FindObject("pidSummary")) delete gROOT->FindObject("pidSummary");
 210   cPid = new TCanvas("pidSummary", "PID summary for the ESD task", 1600, 1200);
 211   cPid->cd();
 212   //GetRefFigure(6);
 213   if(PlotPidSummary(0, trendValues))
 214     cPid->SaveAs("pidSummary.gif");
 215
 216   TCanvas* cCent=0x0;
 217   if(gROOT->FindObject("centSummary")) delete gROOT->FindObject("centSummary");
 218   cCent = new TCanvas("centSummary", "Centrality summary for the ESD task", 1600, 1200);
 219   cCent->cd();
 220   //GetRefFigure(7);
 221   if(PlotCentSummary(trendValues))
 222     cCent->SaveAs("centSummary.gif");
 223 }
 224
 225 //____________________________________________________________________
 226 Bool_t AliTRDcheckESD::GetRefFigure(Int_t ifig)
 227 {
 228   //
 229   // Produce reference Plots during PostProcessing
 230   //
 231   if(ifig>=fNRefFigures){
 232     AliWarning(Form("Ref plot %d not available. Valid only up to %d", ifig, fNRefFigures));
 233     return kFALSE;
 234   }
 235   if(!gPad){
 236     AliWarning("Please provide a canvas to draw results.");
 237     return kFALSE;
 238   } else {
 239     gPad->SetLogx(0);gPad->SetLogy(0);
 240     gPad->SetMargin(0.125, 0.015, 0.1, 0.015);
 241   }
 242
 243   const Char_t *title[20];
 244   TH1 *hF(NULL);
 245   if((hF=(TH1S*)gROOT->FindObject("hFcheckESD"))) delete hF;
 246   TLegend *leg(NULL);
 247   TList *l(NULL); TVirtualPad *pad(NULL);
 248   TGraphErrors *g(NULL);TGraphAsymmErrors *ga(NULL);
 249   TObjArray *arr(NULL);
 250   TLatex *lat=new TLatex();
 251   lat->SetTextSize(0.07);
 252   lat->SetTextColor(2);
 253   TLine line;
 254   TTimeStamp now;
 255   switch(ifig){
 256   case kNCl: // number of clusters/track
 257     if(!(arr = (TObjArray*)fResults->At(kNCl))) return kFALSE;
 258
 259     leg = new TLegend(.83, .7, .99, .96);
 260     leg->SetHeader("Species");
 261     leg->SetBorderSize(0); leg->SetFillStyle(0);
 262     for(Int_t ig(0); ig<fgkNgraph[kNCl-1]; ig++){
 263       if(!(g = (TGraphErrors*)arr->At(ig))) return kFALSE;
 264       if(!g->GetN()) continue;
 265       g->Draw(ig?"pc":"apc"); leg->AddEntry(g, g->GetTitle(), "pl");
 266       if(ig) continue;
 267       hF=g->GetHistogram();
 268       hF->SetXTitle("no of clusters");
 269       hF->SetYTitle("entries");
 270       hF->GetYaxis()->CenterTitle(1);
 271       hF->GetYaxis()->SetTitleOffset(1.2);
 272       hF->SetMinimum(5);
 273     }
 274     leg->Draw(); gPad->SetLogy();
 275     break;
 276   case kTRDstat: // Efficiency
 277     if(!(arr = (TObjArray*)fResults->At(kTRDstat))) return kFALSE;
 278     leg = new TLegend(.62, .77, .98, .98);
 279     leg->SetHeader("TRD Efficiency");
 280     leg->SetBorderSize(0); leg->SetFillStyle(0);
 281     title[0] = "Geometrical (TRDin/TPCout)";
 282     title[1] = "Tracking (TRDout/TRDin)";
 283     title[2] = "PID (TRDpid/TRDin)";
 284     title[3] = "Refit (TRDrefit/TRDin)";
 285     hF = new TH1S("hFcheckESD", ";p [GeV/c];Efficiency", 10, 0.1, 10.);
 286     hF->SetMaximum(1.4);
 287     hF->GetXaxis()->SetMoreLogLabels();
 288     hF->GetYaxis()->CenterTitle(1);
 289     hF->Draw("p");
 290     for(Int_t ig(0); ig<fgkNgraph[kTRDstat-1]; ig++){
 291       if(!(g = (TGraphErrors*)arr->At(ig))) return kFALSE;
 292       g->Draw("pl"); leg->AddEntry(g, title[ig], "pl");
 293       //PutTrendValue(name[id], g->GetMean(2));
 294       //PutTrendValue(Form("%sRMS", name[id]), g->GetRMS(2));
 295     }
 296     leg->Draw(); gPad->SetLogx();
 297     break;
 298   case kTRDmom: // Energy loss
 299     if(!(arr = (TObjArray*)fResults->At(kTRDmom))) return kFALSE;
 300     leg = new TLegend(.65, .7, .95, .99);
 301     leg->SetHeader("Energy Loss");
 302     leg->SetBorderSize(1); leg->SetFillColor(0);
 303     title[0] = "Max & 90% quantile";
 304     title[1] = "Mean & 60% quantile";
 305     hF = new TH1S("hFcheckESD", ";layer;#Delta E", 6, -0.5, 5.5);
 306     hF->SetMaximum(1.3);hF->SetMinimum(-.3);
 307     hF->Draw("p");
 308     for(Int_t ig(0); ig<fgkNgraph[kTRDmom-1]; ig++){
 309       if(!(ga = (TGraphAsymmErrors*)arr->At(ig))) return kFALSE;
 310       ga->Draw("pl"); leg->AddEntry(ga, title[ig], "pl");
 311       //PutTrendValue(name[id], g->GetMean(2));
 312       //PutTrendValue(Form("%sRMS", name[id]), g->GetRMS(2));
 313     }
 314     leg->Draw();gPad->SetLogx(kFALSE);
 315     break;
 316   case kPtRes: // Pt resolution @ vertex
 317     if(!(arr = (TObjArray*)fResults->At(kPtRes))) return kFALSE;
 318     gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
 319     pad = ((TVirtualPad*)l->At(0)); pad->cd(); pad->SetLogx();
 320     pad->SetMargin(0.1, 0.022, 0.1, 0.023);
 321     hF = new TH1S("hFcheckESD", "ITS+TPC+TRD;p_{t} [GeV/c];#Delta p_{t} / p_{t} [%]", 10, 0.2, 10.);
 322     hF->SetMaximum(10.);hF->SetMinimum(-3.);
 323     hF->GetXaxis()->SetMoreLogLabels();
 324     hF->GetXaxis()->SetTitleOffset(1.2);
 325     hF->GetYaxis()->CenterTitle();
 326     hF->Draw("p");
 327     //for(Int_t ig(0); ig<fgkNgraph[kPtRes-1]/2; ig++){
 328     for(Int_t ig(2); ig<6; ig++){
 329       if(!(g = (TGraphErrors*)arr->At(ig))) continue;
 330       if(!g->GetN()) continue;
 331       g->Draw("pl");
 332       //PutTrendValue(name[id], g->GetMean(2));
 333       //PutTrendValue(Form("%sRMS", name[id]), g->GetRMS(2));
 334     }
 335     pad = ((TVirtualPad*)l->At(1)); pad->cd(); pad->SetLogx();
 336     pad->SetMargin(0.1, 0.22, 0.1, 0.023);
 337     hF = (TH1*)hF->Clone("hFcheckESD1");
 338     hF->SetTitle("ITS+TPC");
 339     hF->SetMaximum(10.);hF->SetMinimum(-3.);
 340     hF->Draw("p");
 341     leg = new TLegend(.78, .1, .99, .98);
 342     leg->SetHeader("P_{t} @ DCA");
 343     leg->SetBorderSize(1); leg->SetFillColor(0);
 344     leg->SetTextAlign(22);
 345     leg->SetTextFont(12);
 346     leg->SetTextSize(0.03813559);
 347     {
 348       Int_t nPlots(0);
 349       //for(Int_t ig(fgkNgraph[kPtRes-1]/2); ig<fgkNgraph[kPtRes-1]; ig++){
 350       for(Int_t ig(12); ig<16; ig++){
 351         if(!(g = (TGraphErrors*)arr->At(ig))) continue;
 352         if(!g->GetN()) continue;
 353         nPlots++;
 354         g->Draw("pl"); leg->AddEntry(g, g->GetTitle(), "pl");
 355         //PutTrendValue(name[id], g->GetMean(2));
 356         //PutTrendValue(Form("%sRMS", name[id]), g->GetRMS(2));
 357       }
 358       if(nPlots) leg->Draw();
 359     }
 360     break;
 361   case 5: // plot a 3x3 canvas with tracking related histograms
 362     PlotTrackingSummary(0);
 363     break;
 364
 365   case 6: // plot a 3x3 canvas with PID related histograms
 366     PlotPidSummary(0);
 367     break;
 368
 369   case 7: // plot a 3x3 canvas with centrality dependence histograms
 370     PlotCentSummary();
 371     break;
 372
 373   }
 374   return kTRUE;
 375 }
 376
 377 //____________________________________________________________________
 378 void AliTRDcheckESD::UserExec(Option_t *){
 379   //
 380   // Run the Analysis
 381   //
 382   fESD = dynamic_cast<AliESDEvent*>(InputEvent());
 383   fMC = MCEvent();
 384
 385   if(!fESD){
 386     AliError("ESD event missing.");
 387     return;
 388   }
 389
 390   // Get MC information if available
 391   AliStack * fStack = NULL;
 392   if(HasMC()){
 393     if(!fMC){
 394       AliWarning("MC event missing");
 395       SetMC(kFALSE);
 396     } else {
 397       if(!(fStack = fMC->Stack())){
 398         AliWarning("MC stack missing");
 399         SetMC(kFALSE);
 400       }
 401     }
 402   }
 403   TH1 *h(NULL);
 404
 405   Double_t values[kNTrdCfVariables];      // array where the CF container variables are stored
 406   values[kEventVtxZ] = fESD->GetPrimaryVertex()->GetZv();
 407   values[kEventBC] = fESD->GetBunchCrossNumber();
 408
 409   const AliMultiplicity* mult=fESD->GetMultiplicity();
 410   Double_t itsNTracklets = mult->GetNumberOfTracklets();
 411   if(itsNTracklets<1) return;
 412   Int_t multLimits[6] = {0, 700, 1400, 2100, 2800, 3500};
 413   Int_t centralityClass = 0;
 414   for(Int_t iCent=0; iCent<5; ++iCent) {
 415     if(itsNTracklets>=multLimits[iCent] && itsNTracklets<multLimits[iCent+1])
 416       centralityClass=iCent+1;
 417   }
 418   values[kEventMult] = itsNTracklets;
 419   if(centralityClass == 0) return;
 420
 421   AliESDtrack *esdTrack(NULL);
 422   for(Int_t itrk = 0; itrk < fESD->GetNumberOfTracks(); itrk++){
 423     esdTrack = fESD->GetTrack(itrk);
 424     if(!fReferenceTrackFilter->IsSelected(esdTrack)) continue;
 425
 426     ULong_t status = esdTrack->GetStatus(); //PrintStatus(status);
 427
 428     // pid quality
 429     Bool_t kBarrel = Bool_t(status & AliESDtrack::kTRDin);
 430
 431     // find position and momentum of the track at entrance in TRD
 432     Double_t localCoord[3] = {0., 0., 0.};
 433     Bool_t localCoordGood = esdTrack->GetXYZAt(298., fESD->GetMagneticField(), localCoord);
 434     Double_t localMom[3] = {0., 0., 0.};
 435     Bool_t localMomGood = esdTrack->GetPxPyPzAt(298., fESD->GetMagneticField(), localMom);
 436     //Double_t localPhi = (localMomGood ? TMath::ATan2(localMom[1], localMom[0]) : 0.0);
 437     Double_t localSagitaPhi = (localCoordGood ? TMath::ATan2(localCoord[1], localCoord[0]) : 0.0);
 438
 439     values[kTrackTOFdeltaBC] = esdTrack->GetTOFDeltaBC();
 440     values[kTrackTOFBC]      = esdTrack->GetTOFBunchCrossing(fESD->GetMagneticField());
 441     Float_t dcaXY=0.0; Float_t dcaZ=0.0;
 442     esdTrack->GetImpactParameters(dcaXY, dcaZ);
 443     values[kTrackDCAxy]  = dcaXY;
 444     values[kTrackDCAz]   = dcaZ;
 445     values[kTrackCharge] = esdTrack->Charge();
 446     values[kTrackPhi]    = localSagitaPhi;
 447     values[kTrackEta]    = esdTrack->Eta();
 448     values[kTrackPt]     = esdTrack->Pt();
 449     values[kTrackP]      = esdTrack->P();
 450     values[kTrackTrdTracklets] = esdTrack->GetTRDntracklets();
 451     values[kTrackTrdClusters] = esdTrack->GetTRDncls();
 452     for(Int_t i=0; i<6; ++i) values[kTrackQtot+i] = 0.0;
 453
 454     if(localCoordGood && localMomGood) fCfContainer->Fill(values, 0);      // fill the TPC reference step
 455
 456     // TRD reference tracks
 457     if(esdTrack->GetTRDntracklets()>=1) {
 458       // (slicePH,sliceNo) distribution and Qtot from slices
 459       for(Int_t iPlane=0; iPlane<6; iPlane++) {
 460         Float_t qtot=esdTrack->GetTRDslice(iPlane, 0);
 461         for(Int_t iSlice=0; iSlice<8; iSlice++) {
 462           if(esdTrack->GetTRDslice(iPlane, iSlice)>20.) {
 463             h = (TH2F*)fHistos->At(kPHSlice); h->Fill(iSlice, esdTrack->GetTRDslice(iPlane, iSlice));
 464             h = (TH2F*)fHistos->At(kPHSlice+centralityClass); h->Fill(iSlice, esdTrack->GetTRDslice(iPlane, iSlice));
 465           }
 466         }
 467         values[kTrackQtot+iPlane] = fgkQs*qtot;
 468       }
 469
 470       if(localCoordGood && localMomGood) {
 471         fCfContainer->Fill(values, 1);
 472         if(Bool_t(status & AliESDtrack::kTOFpid)) fCfContainer->Fill(values, 2);
 473       }
 474     }  // end if nTRDtrkl>=1
 475
 476     // look at external track param
 477     const AliExternalTrackParam *op = esdTrack->GetOuterParam();
 478     const AliExternalTrackParam *ip = esdTrack->GetInnerParam();
 479
 480     Double_t pt(0.), pt0(0.), ptTRD(0.);
 481     // read MC info if available
 482     Bool_t kFOUND(kFALSE), kPhysPrim(kFALSE);
 483     AliMCParticle *mcParticle(NULL);
 484     if(HasMC()){
 485       AliTrackReference *ref(NULL);
 486       Int_t fLabel(esdTrack->GetLabel());
 487       Int_t fIdx(TMath::Abs(fLabel));
 488       if(!fStack || fIdx > fStack->GetNtrack()) continue;
 489
 490       // read MC particle
 491       if(!(mcParticle = (AliMCParticle*) fMC->GetTrack(fIdx))) {
 492         AliWarning(Form("MC particle missing. Label[ %d].", fLabel));
 493         continue;
 494       }
 495
 496       pt   = esdTrack->Pt();
 497       pt0  = mcParticle->Pt();
 498       //Double_t eta0 = mcParticle->Eta();
 499       //Double_t phi0 = mcParticle->Phi();
 500       kPhysPrim = fMC->IsPhysicalPrimary(fIdx);
 501
 502       // read track references
 503       Int_t nRefs = mcParticle->GetNumberOfTrackReferences();
 504       if(!nRefs){
 505         AliWarning(Form("No TR found for track @ Label[%d].", fLabel));
 506         continue;
 507       }
 508       Int_t iref = 0;
 509       while(iref<nRefs){
 510         ref = mcParticle->GetTrackReference(iref);
 511         if(ref->LocalX() > fgkxTPC) break;
 512         ref=NULL; iref++;
 513       }
 514       if(ref){
 515         if(ref->LocalX() > fgkxTOF){ // track skipping TRD fiducial volume
 516           ref = mcParticle->GetTrackReference(TMath::Max(iref-1, 0));
 517         }
 518       } else { // track stopped in TPC
 519         ref = mcParticle->GetTrackReference(TMath::Max(iref-1, 0));
 520       }
 521       ptTRD = ref->Pt();kFOUND=kTRUE;
 522     } else { // use reconstructed values
 523       if(op){
 524         Double_t x(op->GetX());
 525         if(x<fgkxTOF && x>fgkxTPC){
 526           ptTRD=op->Pt();
 527           kFOUND=kTRUE;
 528         }
 529       }
 530
 531       if(!kFOUND && ip){
 532         ptTRD=ip->Pt();
 533         kFOUND=kTRUE;
 534       }
 535     }     // end if(HasMC())
 536
 537     if(kFOUND){
 538       h = (TH2I*)fHistos->At(kTRDstat);
 539       if(status & AliESDtrack::kTPCout) h->Fill(ptTRD, kTPCout);
 540       if(status & AliESDtrack::kTRDin) h->Fill(ptTRD, kTRDin);
 541       if(kBarrel && (status & AliESDtrack::kTRDout)) h->Fill(ptTRD, kTRDout);
 542       if(kBarrel && (status & AliESDtrack::kTRDpid)) h->Fill(ptTRD, kTRDpid);
 543       if(kBarrel && (status & AliESDtrack::kTRDrefit)) h->Fill(ptTRD, kTRDref);
 544     }
 545     Int_t idx(HasMC() ? Pdg2Idx(TMath::Abs(mcParticle->PdgCode())): 0)
 546          ,sgn(esdTrack->Charge()<0?0:1);
 547     if(kBarrel && kPhysPrim) {
 548       TH3 *h3 = (TH3S*)fHistos->At(kPtRes);
 549       Int_t offset = (status & AliESDtrack::kTRDrefit) ? 0 : 10;
 550       h3->Fill(pt0, 1.e2*(pt/pt0-1.),
 551         offset + 2*idx + sgn);
 552     }
 553     ((TH1*)fHistos->At(kNCl))->Fill(esdTrack->GetTRDncls(), 2*idx + sgn);
 554     if(ip){
 555       h = (TH2I*)fHistos->At(kTRDmom);
 556       Float_t pTRD(0.);
 557       for(Int_t ily=6; ily--;){
 558         if((pTRD=esdTrack->GetTRDmomentum(ily))<0.) continue;
 559         h->Fill(ip->GetP()-pTRD, ily);
 560       }
 561     }
 562   }  // end loop over tracks
 563
 564   // fill the number of tracks histograms
 565   //h = (TH1I*)fHistos->At(kNTracksAcc);
 566   //h->Fill(nTracksAcc);
 567   //h = (TH1I*)fHistos->At(kNTracksTPC);
 568   //h->Fill(nTracksTPC);
 569   PostData(1, fHistos);
 570 }
 571
 572 //____________________________________________________________________
 573 TObjArray* AliTRDcheckESD::Histos()
 574 {
 575 // Retrieve histograms array if already build or build it
 576
 577   if(fHistos) return fHistos;
 578
 579   fHistos = new TObjArray(kNhistos+1);
 580   fHistos->SetOwner(kTRUE);
 581
 582   TH1 *h = NULL;
 583
 584   // clusters per track
 585   const Int_t kNpt(30);
 586   Float_t pt(0.2);
 587   Float_t binsPt[kNpt+1];
 588   for(Int_t i=0;i<kNpt+1; i++,pt+=(TMath::Exp(i*i*.001)-1.)) binsPt[i]=pt;
 589   if(!(h = (TH2I*)gROOT->FindObject("hNCl"))){
 590     h = new TH2I("hNCl", "Clusters per TRD track;N_{cl}^{TRD};SPECIES;entries", 60, 0., 180., 10, -0.5, 9.5);
 591     TAxis *ay(h->GetYaxis());
 592     ay->SetLabelOffset(0.015);
 593     for(Int_t i(0); i<AliPID::kSPECIES; i++){
 594       ay->SetBinLabel(2*i+1, Form("%s^{-}", AliPID::ParticleLatexName(i)));
 595       ay->SetBinLabel(2*i+2, Form("%s^{+}", AliPID::ParticleLatexName(i)));
 596     }
 597   } else h->Reset();
 598   fHistos->AddAt(h, kNCl); fNRefFigures++;
 599
 600   // status bits histogram
 601   const Int_t kNbits(5);
 602   Float_t bits(.5);
 603   Float_t binsBits[kNbits+1];
 604   for(Int_t i=0; i<kNbits+1; i++,bits+=1.) binsBits[i]=bits;
 605   if(!(h = (TH2I*)gROOT->FindObject("hTRDstat"))){
 606     h = new TH2I("hTRDstat", "TRD status bits;p_{t} @ TRD [GeV/c];status;entries", kNpt, binsPt, kNbits, binsBits);
 607     TAxis *ay(h->GetYaxis());
 608     ay->SetBinLabel(1, "kTPCout");
 609     ay->SetBinLabel(2, "kTRDin");
 610     ay->SetBinLabel(3, "kTRDout");
 611     ay->SetBinLabel(4, "kTRDpid");
 612     ay->SetBinLabel(5, "kTRDrefit");
 613   } else h->Reset();
 614   fHistos->AddAt(h, kTRDstat);
 615
 616   // energy loss
 617   if(!(h = (TH2I*)gROOT->FindObject("hTRDmom"))){
 618     h = new TH2I("hTRDmom", "TRD energy loss;p_{inner} - p_{ly} [GeV/c];ly;entries", 100, -1., 2., 6, -0.5, 5.5);
 619   } else h->Reset();
 620   fHistos->AddAt(h, kTRDmom);
 621   //if(!HasMC()) return fHistos;
 622
 623   // pt resolution
 624   const Int_t kNdpt(100), kNspec(4*AliPID::kSPECIES);
 625   Float_t dpt(-3.), spec(-0.5);
 626   Float_t binsDPt[kNdpt+1], binsSpec[kNspec+1];
 627   for(Int_t i=0; i<kNdpt+1; i++,dpt+=6.e-2) binsDPt[i]=dpt;
 628   for(Int_t i=0; i<kNspec+1; i++,spec+=1.) binsSpec[i]=spec;
 629   if(!(h = (TH3S*)gROOT->FindObject("hPtRes"))){
 630     h = new TH3S("hPtRes", "P_{t} resolution @ DCA;p_{t}^{MC} [GeV/c];#Delta p_{t}/p_{t}^{MC} [%];SPECIES", kNpt, binsPt, kNdpt, binsDPt, kNspec, binsSpec);
 631     TAxis *az(h->GetZaxis());
 632     az->SetLabelOffset(0.015);
 633     for(Int_t i(0); i<AliPID::kSPECIES; i++){
 634       az->SetBinLabel(2*i+1, Form("%s^{-}", AliPID::ParticleLatexName(i)));
 635       az->SetBinLabel(2*i+2, Form("%s^{+}", AliPID::ParticleLatexName(i)));
 636       az->SetBinLabel(10+2*i+1, Form("%s^{-}", AliPID::ParticleLatexName(i)));
 637       az->SetBinLabel(10+2*i+2, Form("%s^{+}", AliPID::ParticleLatexName(i)));
 638     }
 639   } else h->Reset();
 640   fHistos->AddAt(h, kPtRes);
 641
 642   // TPC event vertex distribution
 643   if(!(h = (TH1F*)gROOT->FindObject("hTPCVertex"))){
 644     h = new TH1F("hTPCVertex", "Event vertex Z coord. from TPC tracks", 100, -25., 25.);
 645   } else h->Reset();
 646   fHistos->AddAt(h, kTPCVertex);
 647
 648   // Event vertex
 649   if(!(h = (TH1F*)gROOT->FindObject("hEventVertex"))){
 650     h = new TH1F("hEventVertex", "Event vertex Z coord.", 100, -25., 25.);
 651   } else h->Reset();
 652   fHistos->AddAt(h, kEventVertex);
 653
 654   // Number of all tracks
 655   if(!(h = (TH1I*)gROOT->FindObject("hNTracksAll"))){
 656     h = new TH1I("hNTracksAll", "Number of tracks per event, event vertex cuts", 5000, 0, 5000);
 657   } else h->Reset();
 658   fHistos->AddAt(h, kNTracksAll);
 659
 660   // Number of tracks in acceptance and DCA cut
 661   if(!(h = (TH1I*)gROOT->FindObject("hNTracksAcc"))){
 662     h = new TH1I("hNTracksAcc", Form("Number of tracks per event, |#eta|<%.1f, |DCAxy|<%.1f, |DCAz|<%.1f",
 663                                      fgkEta, fgkTrkDCAxy, fgkTrkDCAz), 5000, 0, 5000);
 664   } else h->Reset();
 665   fHistos->AddAt(h, kNTracksAcc);
 666
 667   // Number of tracks in TPC (Ncls>10)
 668   if(!(h = (TH1I*)gROOT->FindObject("hNTracksTPC"))){
 669     h = new TH1I("hNTracksTPC", Form("Number of tracks per event, |#eta|<%.1f, pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 670                                      fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 5000, 0, 5000);
 671   } else h->Reset();
 672   fHistos->AddAt(h, kNTracksTPC);
 673
 674   // Distribution of DCA-xy
 675   if(!(h = (TH1F*)gROOT->FindObject("hDCAxy"))){
 676     h = new TH1F("hDCAxy", "Distribution of transverse DCA", 100, -100., 100.);
 677   } else h->Reset();
 678   fHistos->AddAt(h, kDCAxy);
 679
 680   // Distribution of DCA-z
 681   if(!(h = (TH1F*)gROOT->FindObject("hDCAz"))){
 682     h = new TH1F("hDCAz", "Distribution of longitudinal DCA", 100, -100., 100.);
 683   } else h->Reset();
 684   fHistos->AddAt(h, kDCAz);
 685
 686   Double_t binPtLimits[33] = {0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9,
 687                              1.0, 1.1, 1.2, 1.3, 1.4,
 688                              1.6, 1.8, 2.0, 2.2, 2.4, 2.6, 2.8, 3.0,
 689                              3.4, 3.8, 4.2, 4.6, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0};
 690   // Pt distributions
 691   if(!(h = (TH1F*)gROOT->FindObject("hPt1"))){
 692     h = new TH1F("hPt1", Form("dN/dpt, |#eta|<%.1f and pt>%.1f", fgkEta, fgkPt), 32, binPtLimits);
 693   } else h->Reset();
 694   fHistos->AddAt(h, kPt1);
 695
 696   if(!(h = (TH1F*)gROOT->FindObject("hPt2"))){
 697     h = new TH1F("hPt2", Form("dN/dpt, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f",
 698                               fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz), 32, binPtLimits);
 699   } else h->Reset();
 700   fHistos->AddAt(h, kPt2);
 701
 702   if(!(h = (TH1F*)gROOT->FindObject("hPt3pos"))){
 703     h = new TH1F("hPt3pos", Form("dN/dpt (positives), |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 704                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 32, binPtLimits);
 705   } else h->Reset();
 706   fHistos->AddAt(h, kPt3pos);
 707
 708   if(!(h = (TH1F*)gROOT->FindObject("hPt3neg"))){
 709     h = new TH1F("hPt3neg", Form("dN/dpt (negatives), |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 710                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 32, binPtLimits);
 711   } else h->Reset();
 712   fHistos->AddAt(h, kPt3neg);
 713
 714   if(!(h = (TH1F*)gROOT->FindObject("hPt4pos"))){
 715     h = new TH1F("hPt4pos", Form("dN/dpt (positives), |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq 1",
 716                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 32, binPtLimits);
 717   } else h->Reset();
 718   fHistos->AddAt(h, kPt4pos);
 719
 720   if(!(h = (TH1F*)gROOT->FindObject("hPt4neg"))){
 721     h = new TH1F("hPt4pos", Form("dN/dpt (negatives), |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq 1",
 722                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 32, binPtLimits);
 723   } else h->Reset();
 724   fHistos->AddAt(h, kPt4neg);
 725
 726   // theta distribution of TRD tracks
 727   if(!(h = (TH1F*)gROOT->FindObject("hTheta"))){
 728     h = new TH1F("hTheta", Form("dN/d#theta, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq 1",
 729                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 220,.5,2.7);
 730   } else h->Reset();
 731   fHistos->AddAt(h, kTheta);
 732
 733   // phi distribution of TRD tracks
 734   if(!(h = (TH1F*)gROOT->FindObject("hPhi"))){
 735     h = new TH1F("hPhi", Form("dN/d#varphi, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq 1",
 736                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 157,0,6.28);
 737   } else h->Reset();
 738   fHistos->AddAt(h, kPhi);
 739
 740   // TPC cluster distribution
 741   if(!(h = (TH1F*)gROOT->FindObject("hNTPCCl"))){
 742     h = new TH1I("hNTPCCl", Form("Number of TPC clusters/track, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f",
 743                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz), 160, 0, 160);
 744   } else h->Reset();
 745   fHistos->AddAt(h, kNTPCCl);
 746
 747   if(!(h = (TH1I*)gROOT->FindObject("hNTPCCl2"))){
 748     h = new TH1F("hNTPCCl2", Form("Number of TPC clusters/track, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, pt>1.0 GeV/c",
 749                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz), 160, 0, 160);
 750   } else h->Reset();
 751   fHistos->AddAt(h, kNTPCCl2);
 752
 753   // dE/dx vs P for TPC reference tracks
 754   if(!(h = (TH2F*)gROOT->FindObject("hTPCDedx"))){
 755     h = new TH2F("hTPCDedx", Form("TPC dE/dx vs P, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 756                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, 0.1,10.1, 150, 0, 150.);
 757   } else h->Reset();
 758   fHistos->AddAt(h, kTPCDedx);
 759
 760   // eta,phi distribution of TPC reference tracks
 761   if(!(h = (TH2F*)gROOT->FindObject("hEtaPhi"))){
 762     h = new TH2F("hEtaPhi", Form("TPC (#eta,#varphi), |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 763                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 50, -1, 1, 157, 0, 6.28);
 764   } else h->Reset();
 765   fHistos->AddAt(h, kEtaPhi);
 766
 767   // Nclusters vs eta distribution for TPC tracks
 768   if(!(h = (TH2F*)gROOT->FindObject("hEtaNclsTPC"))){
 769     h = new TH2F("hEtaNclsTPC", Form("TPC Nclusters vs. #eta, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f",
 770                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz), 50, -1, 1, 160, 0, 160.);
 771   } else h->Reset();
 772   fHistos->AddAt(h, kEtaNclsTPC);
 773
 774   // Nclusters vs phi distribution for TPC reference tracks
 775   if(!(h = (TH2F*)gROOT->FindObject("hPhiNclsTPC"))){
 776     h = new TH2F("hPhiNclsTPC", Form("TPC Nclusters vs. #varphi, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f",
 777                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz), 157, 0, 6.28, 160, 0, 160.);
 778   } else h->Reset();
 779   fHistos->AddAt(h, kPhiNclsTPC);
 780
 781   // SPD multiplicity distribution
 782   if(!(h = (TH1F*)gROOT->FindObject("hSPDMult"))){
 783     h = new TH1F("hSPDMult", "SPD multiplicity", 10000, -0.5, 9999.5);
 784   } else h->Reset();
 785   fHistos->AddAt(h, kSPDMult);
 786
 787   // Ntracklets/track vs P for TRD reference tracks
 788   Double_t binsP[19] = {0.0, 0.2, 0.4, 0.6, 0.8, 1.0, 1.2, 1.4, 1.7, 2.0,
 789                         2.5, 3.0, 3.5, 4.0, 5.0, 6.0, 7.0, 9.0, 12.0};
 790   for(Int_t iCent=0; iCent<=5; ++iCent) {
 791     if(!(h = (TH2F*)gROOT->FindObject(Form("hNTrackletsTRD_cent%d",iCent+1)))){
 792       h = new TH2F(Form("hNTrackletsTRD_cent%d",iCent+1), Form("TRD Ntracklets/track vs. P, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 793                                                                iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 18, binsP, 7, -0.5, 6.5);
 794     } else h->Reset();
 795     fHistos->AddAt(h, kNTrackletsTRD+iCent);
 796   }
 797
 798   // Nclusters/track vs P for TRD reference tracks
 799   for(Int_t iCent=0; iCent<=5; ++iCent) {
 800     if(!(h = (TH2F*)gROOT->FindObject(Form("hNClsTrackTRD_cent%d",iCent+1)))){
 801       h = new TH2F(Form("hNClsTrackTRD_cent%d",iCent+1), Form("TRD Nclusters/track vs. P, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 802                                                               iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 18, binsP, 300, 0., 300.);
 803     } else h->Reset();
 804     fHistos->AddAt(h, kNClsTrackTRD+iCent);
 805   }
 806
 807   // <PH> vs slice number for TRD reference tracklets
 808   for(Int_t iCent=0; iCent<=5; ++iCent) {
 809     if(!(h = (TH2F*)gROOT->FindObject(Form("hPHSlice_cent%d",iCent+1)))){
 810       h = new TH2F(Form("hPHSlice_cent%d",iCent+1), Form("<PH> vs sliceNo, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 811                                     iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 8, -0.5, 7.5, 200, 0., 2000.);
 812     } else h->Reset();
 813     fHistos->AddAt(h, kPHSlice+iCent);
 814   }
 815
 816   // <PH> vs slice number for TRD reference tracklets, from TPC pions
 817   for(Int_t iCent=0; iCent<=5; ++iCent) {
 818     if(!(h = (TH2F*)gROOT->FindObject(Form("hPHSliceTPCpions_cent%d",iCent+1)))){
 819       h = new TH2F(Form("hPHSliceTPCpions_cent%d",iCent+1), Form("<PH> vs sliceNo, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, TPC pions",
 820                                                                  iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 8, -0.5, 7.5, 100, 0., 2000.);
 821     } else h->Reset();
 822     fHistos->AddAt(h, kPHSliceTPCpions+iCent);
 823   }
 824
 825   // TPC dE/dx vs P for TRD reference tracks, pions
 826   for(Int_t iCent=0; iCent<=5; ++iCent) {
 827     if(!(h = (TH2F*)gROOT->FindObject(Form("hTPCdedxPions_cent%d",iCent+1)))){
 828       h = new TH2F(Form("hTPCdedxPions_cent%d",iCent+1), Form("TPC dE/dx vs P, TPC pions, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 829                                          iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, 0.1,10.1, 100, 0,100.);
 830     } else h->Reset();
 831     fHistos->AddAt(h, kTPCdedxPions+iCent);
 832   }
 833
 834   // <PH> vs slice number for TRD reference tracklets, from TPC electrons
 835   for(Int_t iCent=0; iCent<=5; ++iCent) {
 836     if(!(h = (TH2F*)gROOT->FindObject(Form("hPHSliceTPCelectrons_cent%d",iCent+1)))){
 837       h = new TH2F(Form("hPHSliceTPCelectrons_cent%d",iCent+1), Form("<PH> vs sliceNo, centrality %d,|#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, TPC electrons",
 838                                                 iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 8, -0.5, 7.5, 100, 0., 2000.);
 839     } else h->Reset();
 840     fHistos->AddAt(h, kPHSliceTPCelectrons+iCent);
 841   }
 842
 843   // TPC dE/dx vs P for TRD reference tracks, electrons
 844   for(Int_t iCent=0; iCent<=5; ++iCent) {
 845     if(!(h = (TH2F*)gROOT->FindObject(Form("hTPCdedxElectrons_cent%d",iCent+1)))){
 846       h = new TH2F(Form("hTPCdedxElectrons_cent%d",iCent+1), Form("TPC dE/dx vs P, TPC electrons, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 847                                              iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, 0.1,10.1, 100, 0,100.);
 848     } else h->Reset();
 849     fHistos->AddAt(h, kTPCdedxElectrons+iCent);
 850   }
 851
 852   // Qtot vs P for TRD reference tracklets
 853   for(Int_t iCent=0; iCent<=5; ++iCent) {
 854     if(!(h = (TH2F*)gROOT->FindObject(Form("hQtotP_cent%d",iCent+1)))){
 855       h = new TH2F(Form("hQtotP_cent%d",iCent+1), Form("Qtot(from slices) vs P, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 856                                   iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 18, binsP, 400, 0., 20);
 857     } else h->Reset();
 858     fHistos->AddAt(h, kQtotP+iCent);
 859   }
 860
 861   // (X,Y,momentum) distribution after AliESDtrack::PropagateTo(r=300.)
 862   if(!(h = (TH3F*)gROOT->FindObject("hPropagXYvsP"))){
 863     h = new TH3F("hPropagXYvsP", Form("(x,y) vs P after AliESDtrack::PropagateTo(r=300.), |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 864                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100,-500,500, 100,-500,500, 10, 0.,10.);
 865   } else h->Reset();
 866   fHistos->AddAt(h, kPropagXYvsP);
 867
 868   // (R,Z,momentum) distribution after AliESDtrack::PropagateTo(r=300.)
 869   if(!(h = (TH3F*)gROOT->FindObject("hPropagRZvsP"))){
 870     h = new TH3F("hPropagRZvsP", Form("(r,z) vs P after AliESDtrack::PropagateTo(r=300.), |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 871                  fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100,-350., 350., 100,0.,500., 10, 0.,10.);
 872   } else h->Reset();
 873   fHistos->AddAt(h, kPropagRZvsP);
 874
 875   Double_t etaBinLimits[101];
 876   for(Int_t i=0; i<101; i++) etaBinLimits[i] = -1.0 + i*2.0/100.;
 877   Double_t phiBinLimits[151];
 878   for(Int_t i=0; i<151; i++) phiBinLimits[i] = -1.1*TMath::Pi() + i*2.2*TMath::Pi()/150.;
 879   // (eta,detector phi,P) distribution of reference TPC positive tracks
 880   for(Int_t iCent=0; iCent<=5; ++iCent) {
 881     if(!(h = (TH3F*)gROOT->FindObject(Form("hTPCRefTracksPos_cent%d",iCent+1)))){
 882       h = new TH3F(Form("hTPCRefTracksPos_cent%d",iCent+1), Form("(#eta,detector #varphi,p) for TPC positive reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 883                                             iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, etaBinLimits, 150, phiBinLimits, 32, binPtLimits);
 884     } else h->Reset();
 885     fHistos->AddAt(h, kTPCRefTracksPos+iCent);
 886   }
 887
 888   // (eta,detector phi,P) distribution of reference TPC negative tracks
 889   for(Int_t iCent=0; iCent<=5; ++iCent) {
 890     if(!(h = (TH3F*)gROOT->FindObject(Form("hTPCRefTracksNeg_cent%d",iCent+1)))){
 891       h = new TH3F(Form("hTPCRefTracksNeg_cent%d",iCent+1), Form("(#eta,detector #varphi,p) for TPC negative reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d",
 892                                             iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, etaBinLimits, 150, phiBinLimits, 32, binPtLimits);
 893     } else h->Reset();
 894     fHistos->AddAt(h, kTPCRefTracksNeg+iCent);
 895   }
 896
 897   // (eta,detector phi,P) distribution of reference TRD positive tracks
 898   for(Int_t iCent=0; iCent<=5; ++iCent) {
 899     if(!(h = (TH3F*)gROOT->FindObject(Form("hTRDRefTracksPos_cent%d",iCent+1)))){
 900       h = new TH3F(Form("hTRDRefTracksPos_cent%d",iCent+1), Form("(#eta,detector #varphi,p) for TRD positive reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq1",
 901                                             iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, etaBinLimits, 150, phiBinLimits, 32, binPtLimits);
 902     } else h->Reset();
 903     fHistos->AddAt(h, kTRDRefTracksPos+iCent);
 904   }
 905
 906   // (eta,detector phi,P) distribution of reference TRD negative tracks
 907   for(Int_t iCent=0; iCent<=5; ++iCent) {
 908     if(!(h = (TH3F*)gROOT->FindObject(Form("hTRDRefTracksNeg_cent%d",iCent+1)))){
 909       h = new TH3F(Form("hTRDRefTracksNeg_cent%d",iCent+1), Form("(#eta,detector #varphi,p) for TRD negative reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq1",
 910                                             iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, etaBinLimits, 150, phiBinLimits, 32, binPtLimits);
 911     } else h->Reset();
 912     fHistos->AddAt(h, kTRDRefTracksNeg+iCent);
 913   }
 914
 915   // (eta,detector phi,P) distribution of reference TRD positive tracks with 4 tracklets
 916   for(Int_t iCent=0; iCent<=5; ++iCent) {
 917     if(!(h = (TH3F*)gROOT->FindObject(Form("hTRDRefTracksPos4_cent%d",iCent+1)))){
 918       h = new TH3F(Form("hTRDRefTracksPos4_cent%d",iCent+1), Form("(#eta,detector #varphi,p) for TRD positive reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq1",
 919                                                                  iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, etaBinLimits, 150, phiBinLimits, 32, binPtLimits);
 920     } else h->Reset();
 921     fHistos->AddAt(h, kTRDRefTracksPos4+iCent);
 922   }
 923
 924   // (eta,detector phi,P) distribution of reference TRD positive tracks with 5 tracklets
 925   for(Int_t iCent=0; iCent<=5; ++iCent) {
 926     if(!(h = (TH3F*)gROOT->FindObject(Form("hTRDRefTracksPos5_cent%d",iCent+1)))){
 927       h = new TH3F(Form("hTRDRefTracksPos5_cent%d",iCent+1), Form("(#eta,detector #varphi,p) for TRD positive reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq1",
 928                                                                   iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, etaBinLimits, 150, phiBinLimits, 32, binPtLimits);
 929     } else h->Reset();
 930     fHistos->AddAt(h, kTRDRefTracksPos5+iCent);
 931   }
 932
 933   // (eta,detector phi,P) distribution of reference TRD positive tracks with 6 tracklets
 934   for(Int_t iCent=0; iCent<=5; ++iCent) {
 935     if(!(h = (TH3F*)gROOT->FindObject(Form("hTRDRefTracksPos6_cent%d",iCent+1)))){
 936       h = new TH3F(Form("hTRDRefTracksPos6_cent%d",iCent+1), Form("(#eta,detector #varphi,p) for TRD positive reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq1",
 937                                                                   iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, etaBinLimits, 150, phiBinLimits, 32, binPtLimits);
 938     } else h->Reset();
 939     fHistos->AddAt(h, kTRDRefTracksPos6+iCent);
 940   }
 941
 942   // (eta,detector phi,P) distribution of reference TRD negative tracks with 4 tracklets
 943   for(Int_t iCent=0; iCent<=5; ++iCent) {
 944     if(!(h = (TH3F*)gROOT->FindObject(Form("hTRDRefTracksNeg4_cent%d",iCent+1)))){
 945       h = new TH3F(Form("hTRDRefTracksNeg4_cent%d",iCent+1), Form("(#eta,detector #varphi,p) for TRD negative reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq1",
 946                                                                  iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, etaBinLimits, 150, phiBinLimits, 32, binPtLimits);
 947     } else h->Reset();
 948     fHistos->AddAt(h, kTRDRefTracksNeg4+iCent);
 949   }
 950
 951   // (eta,detector phi,P) distribution of reference TRD negative tracks with 5 tracklets
 952   for(Int_t iCent=0; iCent<=5; ++iCent) {
 953     if(!(h = (TH3F*)gROOT->FindObject(Form("hTRDRefTracksNeg5_cent%d",iCent+1)))){
 954       h = new TH3F(Form("hTRDRefTracksNeg5_cent%d",iCent+1), Form("(#eta,detector #varphi,p) for TRD negative reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq1",
 955                                                                   iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, etaBinLimits, 150, phiBinLimits, 32, binPtLimits);
 956     } else h->Reset();
 957     fHistos->AddAt(h, kTRDRefTracksNeg5+iCent);
 958   }
 959
 960   // (eta,detector phi,P) distribution of reference TRD negative tracks with 6 tracklets
 961   for(Int_t iCent=0; iCent<=5; ++iCent) {
 962     if(!(h = (TH3F*)gROOT->FindObject(Form("hTRDRefTracksNeg6_cent%d",iCent+1)))){
 963       h = new TH3F(Form("hTRDRefTracksNeg6_cent%d",iCent+1), Form("(#eta,detector #varphi,p) for TRD negative reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq1",
 964                                                                   iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, etaBinLimits, 150, phiBinLimits, 32, binPtLimits);
 965     } else h->Reset();
 966     fHistos->AddAt(h, kTRDRefTracksNeg6+iCent);
 967   }
 968
 969
 970   // (eta,detector phi) profile of average number of TRD tracklets/track
 971   for(Int_t iCent=0; iCent<=5; ++iCent) {
 972     if(!(h = (TProfile2D*)gROOT->FindObject(Form("hTRDEtaPhiAvNtrkl_cent%d",iCent+1)))){
 973       h = new TProfile2D(Form("hTRDEtaPhiAvNtrkl_cent%d",iCent+1), Form("<Ntracklets/track> vs (#eta,detector #varphi) for TRD reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq1",
 974                                                    iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, -1.0, 1.0, 150, -1.1*TMath::Pi(), 1.1*TMath::Pi());
 975     } else h->Reset();
 976     fHistos->AddAt(h, kTRDEtaPhiAvNtrkl+iCent);
 977   }
 978
 979   // (eta,delta phi) profile of average number of TRD tracklets/track
 980   for(Int_t iCent=0; iCent<=5; ++iCent) {
 981     if(!(h = (TProfile2D*)gROOT->FindObject(Form("hTRDEtaDeltaPhiAvNtrkl_cent%d",iCent+1)))){
 982       h = new TProfile2D(Form("hTRDEtaDeltaPhiAvNtrkl_cent%d",iCent+1), Form("<Ntracklets/track> vs (#eta, #Delta#varphi) for TRD reference tracks, centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq1",
 983                                                         iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, -1.0, 1.0, 50, -0.4*TMath::Pi(), 0.4*TMath::Pi());
 984     } else h->Reset();
 985     fHistos->AddAt(h, kTRDEtaDeltaPhiAvNtrkl+iCent);
 986   }
 987
 988   // (eta, detector phi) profile of average tracklet Qtot from slices
 989   for(Int_t iCent=0; iCent<=5; ++iCent) {
 990     for(Int_t iLayer=0;iLayer<6;iLayer++) {
 991       if(!(h = (TProfile2D*)gROOT->FindObject(Form("hTRDEtaPhiAvQtot_Layer%d_cent%d",iLayer,iCent+1)))) {
 992         h = new TProfile2D(Form("hTRDEtaPhiAvQtot_Layer%d_cent%d",iLayer,iCent+1),
 993                            Form("<Q_{tot}> vs (#eta, detector #varphi) for TRD reference tracks (layer %d), centrality %d, |#eta|<%.1f and pt>%.1f, |DCAxy|<%.1f, |DCAz|<%.1f, TPC nclusters>%d, nTRDtracklets#geq1",
 994                                 iLayer, iCent+1, fgkEta, fgkPt, fgkTrkDCAxy, fgkTrkDCAz, fgkNclTPC), 100, -1.0, 1.0, 150, -1.1*TMath::Pi(), 1.1*TMath::Pi());
 995       } else h->Reset();
 996       fHistos->AddAt(h, kTRDEtaPhiAvQtot+iCent*6+iLayer);
 997     }
 998   }
 999
1000   // create a CF container and add it to the list of histograms
1001   Int_t nbinsCf[kNTrdCfVariables];
1002   for(Int_t i=0;i<kNTrdCfVariables;++i) nbinsCf[i]=0;
1003   nbinsCf[kEventVtxZ]         =   12;
1004   nbinsCf[kEventMult]         =    5;
1005   nbinsCf[kEventBC]           = 3500;
1006   nbinsCf[kTrackTOFdeltaBC]   = 2001;
1007   nbinsCf[kTrackTOFBC]        =  201;
1008   nbinsCf[kTrackDCAxy]        =   26;
1009   nbinsCf[kTrackDCAz]         =  101;
1010   nbinsCf[kTrackCharge]       =    2;
1011   nbinsCf[kTrackPhi]          =  150;
1012   nbinsCf[kTrackEta]          =  100;
1013   nbinsCf[kTrackPt]           =   32;
1014   nbinsCf[kTrackP]            =   18;
1015   nbinsCf[kTrackTrdTracklets] =    7;
1016   nbinsCf[kTrackTrdClusters]  =  200;
1017   for(Int_t i=0;i<6;++i) nbinsCf[kTrackQtot+i] = 100;
1018   Double_t evVtxLims[2]      = {-12.,+12.};
1019   Double_t evMultLims[6]     = {0.0, 700., 1400., 2100., 2800., 3500.};
1020   Double_t evBCLims[2]       = {-0.5, +3499.5};
1021   Double_t trkTOFdeltaBClims[2] = {-1000.5, 1000.5};
1022   Double_t trkTOFBClims[2]   = {-100.5, 100.5};
1023   Double_t trkDCAxyLims[27]  = {-50.0, -40.0, -30.0, -25.0, -20.0,
1024                                 -15.0, -10.0,  -8.0,  -6.0,  -4.0,
1025                                  -3.0,  -2.0,  -1.0,   0.0,  +1.0,
1026                                  +2.0,  +3.0,  +4.0,  +6.0,  +8.0,
1027                                 +10.0, +15.0, +20.0, +25.0, +30.0,
1028                                 +40.0, +50.0};
1029   Double_t trkDCAzLims[2]    = {-50.5, +50.5};
1030   Double_t trkChargeLims[2]  = {-1.5, +1.5};
1031   Double_t trkPhiLims[2]     = {-1.1*TMath::Pi(), +1.1*TMath::Pi()};
1032   Double_t trkEtaLims[2]     = {-1.0, +1.0};
1033   Double_t trkPtLims[33]     = {0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9,
1034                                 1.0, 1.1, 1.2, 1.3, 1.4, 1.6, 1.8, 2.0, 2.2, 2.4,
1035                                 2.6, 2.8, 3.0, 3.4, 3.8, 4.2, 4.6, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0};
1036   Double_t trkPLims[19]      = {0.0, 0.2, 0.4, 0.6, 0.8, 1.0, 1.2, 1.4, 1.7, 2.0,
1037                                 2.5, 3.0, 3.5, 4.0, 5.0, 6.0, 7.0, 9.0, 12.0};
1038   Double_t trkTrdNLims[2]    = {-0.5, 6.5};
1039   Double_t trkTrdNclsLims[2] = {-0.5, 199.5};
1040   Double_t trkQtotLims[2]    = {0.0, 20.};
1041   fCfContainer = new AliCFContainer("TrdCfContainer", "TRD ESD CF container", 3, kNTrdCfVariables, nbinsCf);
1042   fCfContainer->SetBinLimits(kEventVtxZ, evVtxLims[0], evVtxLims[1]);
1043   fCfContainer->SetBinLimits(kEventMult, evMultLims);
1044   fCfContainer->SetBinLimits(kEventBC, evBCLims[0], evBCLims[1]);
1045   fCfContainer->SetBinLimits(kTrackTOFdeltaBC, trkTOFdeltaBClims[0], trkTOFdeltaBClims[1]);
1046   fCfContainer->SetBinLimits(kTrackTOFBC, trkTOFBClims[0], trkTOFBClims[1]);
1047   fCfContainer->SetBinLimits(kTrackDCAxy, trkDCAxyLims);
1048   fCfContainer->SetBinLimits(kTrackDCAz, trkDCAzLims[0], trkDCAzLims[1]);
1049   fCfContainer->SetBinLimits(kTrackCharge, trkChargeLims[0], trkChargeLims[1]);
1050   fCfContainer->SetBinLimits(kTrackPhi, trkPhiLims[0], trkPhiLims[1]);
1051   fCfContainer->SetBinLimits(kTrackEta, trkEtaLims[0], trkEtaLims[1]);
1052   fCfContainer->SetBinLimits(kTrackPt, trkPtLims);
1053   fCfContainer->SetBinLimits(kTrackP, trkPLims);
1054   fCfContainer->SetBinLimits(kTrackTrdTracklets, trkTrdNLims[0], trkTrdNLims[1]);
1055   fCfContainer->SetBinLimits(kTrackTrdClusters, trkTrdNclsLims[0], trkTrdNclsLims[1]);
1056   for(Int_t i=0; i<6; ++i) fCfContainer->SetBinLimits(kTrackQtot+i, trkQtotLims[0], trkQtotLims[1]);
1057   fCfContainer->SetVarTitle(kEventVtxZ, "vtxZ");
1058   fCfContainer->SetVarTitle(kEventMult, "multiplicity");
1059   fCfContainer->SetVarTitle(kEventBC, "BC");
1060   fCfContainer->SetVarTitle(kTrackTOFdeltaBC, "TOFdeltaBC");
1061   fCfContainer->SetVarTitle(kTrackTOFBC, "TOFBC");
1062   fCfContainer->SetVarTitle(kTrackDCAxy, "DCAxy");
1063   fCfContainer->SetVarTitle(kTrackDCAz, "DCAz");
1064   fCfContainer->SetVarTitle(kTrackCharge, "charge");
1065   fCfContainer->SetVarTitle(kTrackPhi, "phi");
1066   fCfContainer->SetVarTitle(kTrackEta, "eta");
1067   fCfContainer->SetVarTitle(kTrackPt, "pt");
1068   fCfContainer->SetVarTitle(kTrackP, "P");
1069   fCfContainer->SetVarTitle(kTrackTrdTracklets, "tracklets");
1070   fCfContainer->SetVarTitle(kTrackTrdClusters, "clusters");
1071   for(Int_t i=0; i<6; ++i) fCfContainer->SetVarTitle(kTrackQtot+i, Form("Qtot%d",i));
1072
1073   fCfContainer->SetStepTitle(0, "TPC reference");
1074   fCfContainer->SetStepTitle(1, "TRD");
1075   fCfContainer->SetStepTitle(2, "TOF");
1076
1077   //  fHistos->AddAt(fCfContainer, kNhistos);
1078   fHistos->AddAt(fCfContainer, 0);
1079
1080   return fHistos;
1081 }
1082
1083 //____________________________________________________________________
1084 Bool_t AliTRDcheckESD::Load(const Char_t *file, const Char_t *dir, const Char_t *name)
1085 {
1086 // Load data from performance file
1087
1088   if(!TFile::Open(file)){
1089     AliWarning(Form("Couldn't open file %s.", file));
1090     return kFALSE;
1091   }
1092   if(dir){
1093     if(!gFile->cd(dir)){
1094       AliWarning(Form("Couldn't cd to %s in %s.", dir, file));
1095       return kFALSE;
1096     }
1097   }
1098   TObjArray *o(NULL);
1099   const Char_t *tn=(name ? name : GetName());
1100   if(!(o = (TObjArray*)gDirectory->Get(tn))){
1101     AliWarning(Form("Missing histogram container %s.", tn));
1102     return kFALSE;
1103   }
1104   fHistos = (TObjArray*)o->Clone(GetName());
1105   fCfContainer = (AliCFContainer*)fHistos->At(0);
1106   gFile->Close();
1107   return kTRUE;
1108 }
1109
1110 //_______________________________________________________
1111 Bool_t AliTRDcheckESD::PutTrendValue(const Char_t *name, Double_t val)
1112 {
1113 // Dump trending value to default file
1114
1115   if(!fgFile){
1116     fgFile = fopen("TRD.Performance.txt", "at");
1117   }
1118   fprintf(fgFile, "%s_%s %f\n", GetName(), name, val);
1119   return kTRUE;
1120 }
1121
1122 //____________________________________________________________________
1123 void AliTRDcheckESD::Terminate(Option_t *)
1124 {
1125   // Steer post-processing
1126   if(!fHistos){
1127     fHistos = dynamic_cast<TObjArray *>(GetOutputData(1));
1128     if(!fHistos){
1129       AliError("Histogram container not found in output");
1130       return;
1131     }
1132   }
1133
1134   //  fNRefFigures = 15;
1135   //  return;
1136
1137   const Char_t *name[kNrefs] = {
1138     "Ncl", "Eff", "Eloss", "PtResDCA"
1139   };
1140
1141   TObjArray *arr(NULL); TGraph *g(NULL);
1142   if(!fResults){
1143     fResults = new TObjArray(kNrefs);
1144     fResults->SetOwner();
1145     fResults->SetName("results");
1146     for(Int_t iref(0); iref<kNrefs; iref++){
1147       fResults->AddAt(arr = new TObjArray(fgkNgraph[iref]), iref);
1148       arr->SetName(name[iref]);  arr->SetOwner();
1149       switch(iref+1){
1150       case kNCl:
1151         for(Int_t ig(0); ig<fgkNgraph[iref]; ig++){
1152           arr->AddAt(g = new TGraphErrors(), ig);
1153           g->SetLineColor(ig+1);
1154           g->SetMarkerColor(ig+1);
1155           g->SetMarkerStyle(ig+20);
1156           g->SetName(Form("s%d", ig));
1157           switch(ig){
1158           case 0: g->SetTitle("ALL"); break;
1159           case 1: g->SetTitle("NEG"); break;
1160           case 2: g->SetTitle("POS"); break;
1161           default: g->SetTitle(AliPID::ParticleLatexName(ig-3)); break;
1162           };
1163         }
1164         break;
1165       case kTRDmom:
1166         for(Int_t ig(0); ig<fgkNgraph[iref]; ig++){
1167           arr->AddAt(g = new TGraphAsymmErrors(), ig);
1168           g->SetLineColor(ig+1);
1169           g->SetMarkerColor(ig+1);
1170           g->SetMarkerStyle(ig+20);
1171         }
1172         break;
1173       case kPtRes:
1174         for(Int_t idx(0); idx<AliPID::kSPECIES; idx++){
1175           Int_t ig(2*idx);
1176           arr->AddAt(g = new TGraphErrors(), ig);
1177           g->SetLineColor(kRed-idx);
1178           g->SetMarkerColor(kRed-idx);
1179           g->SetMarkerStyle(20+idx);
1180           g->SetNameTitle(Form("s%d", ig), Form("res %s", AliPID::ParticleLatexName(idx)));
1181           arr->AddAt(g = new TGraphErrors(), ig+1);
1182           g->SetLineColor(kBlue-idx);
1183           g->SetMarkerColor(kBlue-idx);
1184           g->SetMarkerStyle(20+idx);
1185           g->SetNameTitle(Form("m%d", ig+1), Form("sys %s", AliPID::ParticleLatexName(idx)));
1186
1187           ig+=10;
1188           arr->AddAt(g = new TGraphErrors(), ig);
1189           g->SetLineColor(kRed-idx);
1190           g->SetMarkerColor(kRed-idx);
1191           g->SetMarkerStyle(20+idx);
1192           g->SetNameTitle(Form("s%d", ig), Form("sigma %s", AliPID::ParticleLatexName(idx)));
1193           arr->AddAt(g = new TGraphErrors(), ig+1);
1194           g->SetLineColor(kBlue-idx);
1195           g->SetMarkerColor(kBlue-idx);
1196           g->SetMarkerStyle(20+idx);
1197           g->SetNameTitle(Form("m%d", ig+1), Form("mean %s", AliPID::ParticleLatexName(idx)));
1198         }
1199         break;
1200       default:
1201         for(Int_t ig(0); ig<fgkNgraph[iref]; ig++){
1202           arr->AddAt(g = new TGraphErrors(), ig);
1203           g->SetLineColor(ig+1);
1204           g->SetMarkerColor(ig+1);
1205           g->SetMarkerStyle(ig+20);
1206         }
1207         break;
1208       }
1209     }
1210   }
1211   TH1 *h1[2] = {NULL, NULL};
1212   TH2I *h2(NULL);
1213   TAxis *ax(NULL);
1214
1215   // No of clusters
1216   if(!(h2 = (TH2I*)fHistos->At(kNCl))) return;
1217   ax = h2->GetXaxis();
1218   arr = (TObjArray*)fResults->At(kNCl);
1219
1220   // All tracks
1221   h1[0] = h2->ProjectionX("Ncl_px");
1222   TGraphErrors *ge=(TGraphErrors*)arr->At(0);
1223   for(Int_t ib=2; ib<=ax->GetNbins(); ib++){
1224     ge->SetPoint(ib-2, ax->GetBinCenter(ib), h1[0]->GetBinContent(ib));
1225   }
1226
1227   // All charged tracks
1228   TH1 *hNclCh[2] = {(TH1D*)h1[0]->Clone("NEG"), (TH1D*)h1[0]->Clone("POS")};
1229   hNclCh[0]->Reset();hNclCh[1]->Reset();
1230   for(Int_t is(1); is<=AliPID::kSPECIES; is++){
1231     hNclCh[0]->Add(h2->ProjectionX("Ncl_px", 2*is-1, 2*is-1)); // neg
1232     hNclCh[1]->Add(h2->ProjectionX("Ncl_px", 2*is, 2*is));     // pos
1233   }
1234   if(Int_t(hNclCh[0]->GetEntries())){
1235     ge=(TGraphErrors*)arr->At(1);
1236     for(Int_t ib=2; ib<=ax->GetNbins(); ib++){
1237       ge->SetPoint(ib-2, ax->GetBinCenter(ib), hNclCh[0]->GetBinContent(ib));
1238     }
1239   }
1240
1241   if(Int_t(hNclCh[1]->GetEntries())){
1242     ge=(TGraphErrors*)arr->At(2);
1243     for(Int_t ib=2; ib<=ax->GetNbins(); ib++){
1244       ge->SetPoint(ib-2, ax->GetBinCenter(ib), hNclCh[1]->GetBinContent(ib));
1245     }
1246   }
1247   // Species wise
1248   for(Int_t is(1); is<=AliPID::kSPECIES; is++){
1249     h1[0] = h2->ProjectionX("Ncl_px", 2*is-1, 2*is);
1250     if(!Int_t(h1[0]->GetEntries())) continue;
1251     ge=(TGraphErrors*)arr->At(2+is);
1252     for(Int_t ib=2; ib<=ax->GetNbins(); ib++){
1253       ge->SetPoint(ib-2, ax->GetBinCenter(ib), h1[0]->GetBinContent(ib));
1254     }
1255   }
1256
1257   fNRefFigures = 1;
1258
1259   // EFFICIENCY
1260   // geometrical efficiency
1261   if(!(h2 = (TH2I*)fHistos->At(kTRDstat))) return;
1262   arr = (TObjArray*)fResults->At(kTRDstat-1);
1263   h1[0] = h2->ProjectionX("checkESDx0", kTPCout, kTPCout);
1264   h1[1] = h2->ProjectionX("checkESDx1", kTRDin, kTRDin);
1265   Process(h1, (TGraphErrors*)arr->At(0));
1266   delete h1[0];delete h1[1];
1267   // tracking efficiency
1268   h1[0] = h2->ProjectionX("checkESDx0", kTRDin, kTRDin);
1269   h1[1] = h2->ProjectionX("checkESDx1", kTRDout, kTRDout);
1270   Process(h1, (TGraphErrors*)arr->At(1));
1271   delete h1[1];
1272   // PID efficiency
1273   h1[1] = h2->ProjectionX("checkESDx1", kTRDpid, kTRDpid);
1274   Process(h1, (TGraphErrors*)arr->At(2));
1275   delete h1[1];
1276   // Refit efficiency
1277   h1[1] = h2->ProjectionX("checkESDx1", kTRDref, kTRDref);
1278   Process(h1, (TGraphErrors*)arr->At(3));
1279   delete h1[1];
1280   fNRefFigures++;
1281
1282   // ENERGY LOSS
1283   if(!(h2 = dynamic_cast<TH2I*>(fHistos->At(kTRDmom)))) return;
1284   arr = (TObjArray*)fResults->At(kTRDmom-1);
1285   TGraphAsymmErrors *g06 = (TGraphAsymmErrors*)arr->At(0), *g09 = (TGraphAsymmErrors*)arr->At(1);
1286   ax=h2->GetXaxis();
1287   const Int_t nq(4);
1288   const Double_t xq[nq] = {0.05, 0.2, 0.8, 0.95};
1289   Double_t yq[nq];
1290   for(Int_t ily=6; ily--;){
1291     h1[0] = h2->ProjectionX("checkESDp0", ily+1, ily+1);
1292     h1[0]->GetQuantiles(nq,yq,xq);
1293     g06->SetPoint(ily, Float_t(ily), ax->GetBinCenter(h1[0]->GetMaximumBin()));
1294     g06->SetPointError(ily, 0., 0., TMath::Abs(yq[0]), yq[3]);
1295     g09->SetPoint(ily, Float_t(ily), h1[0]->GetMean());
1296     g09->SetPointError(ily, 0., 0., TMath::Abs(yq[1]), yq[2]);
1297
1298     //printf(" max[%f] mean[%f] q[%f %f %f %f]\n", ax->GetBinCenter(h1[0]->GetMaximumBin()), h1[0]->GetMean(), yq[0], yq[1], yq[2], yq[3]);
1299     delete h1[0];
1300   }
1301   fNRefFigures++;
1302 //  if(!HasMC()) return;
1303
1304   // Pt RESOLUTION @ DCA
1305   TH3S* h3(NULL); TGraphErrors *gg[2] = {NULL,NULL};
1306   if(!(h3 = dynamic_cast<TH3S*>(fHistos->At(kPtRes)))) return;
1307   arr = (TObjArray*)fResults->At(kPtRes-1);
1308   TAxis *az(h3->GetZaxis());
1309   for(Int_t i(0); i<AliPID::kSPECIES; i++){
1310     Int_t idx(2*i);
1311     az->SetRange(idx+1, idx+2);
1312     gg[1] = (TGraphErrors*)arr->At(idx);
1313     gg[0] = (TGraphErrors*)arr->At(idx+1);
1314     Process2D((TH2*)h3->Project3D("yx"), gg);
1315
1316     idx+=10;
1317     az->SetRange(idx+1, idx+2);
1318     gg[1] = (TGraphErrors*)arr->At(idx);
1319     gg[0] = (TGraphErrors*)arr->At(idx+1);
1320     Process2D((TH2*)h3->Project3D("yx"), gg);
1321   }
1322   fNRefFigures++;
1323
1324   fNRefFigures++;
1325   // 3x3 tracking summary canvases for every centrality class
1326   fNRefFigures++;
1327   // 3x3 PID summary canvases for every centrality class
1328   fNRefFigures++;
1329   // 3x3 for centrality dependent pictures
1330   fNRefFigures++;
1331
1332   //DoTrending();
1333 }
1334
1335 //____________________________________________________________________
1336 Int_t AliTRDcheckESD::Pdg2Idx(Int_t pdg) const
1337 {
1338   //
1339   // Helper function converting PDG code into AliPID index
1340   //
1341   switch(pdg){
1342   case kElectron:
1343   case kPositron: return AliPID::kElectron;
1344   case kMuonPlus:
1345   case kMuonMinus: return AliPID::kMuon;
1346   case kPiPlus:
1347   case kPiMinus: return AliPID::kPion;
1348   case kKPlus:
1349   case kKMinus: return AliPID::kKaon;
1350   case kProton:
1351   case kProtonBar: return AliPID::kProton;
1352   }
1353   return -1;
1354 }
1355
1356 //____________________________________________________________________
1357 void AliTRDcheckESD::Process(TH1 **h1, TGraphErrors *g)
1358 {
1359 // Generic function to process one reference plot
1360
1361   Int_t n1 = 0, n2 = 0, ip=0;
1362   Double_t eff = 0.;
1363
1364   TAxis *ax = h1[0]->GetXaxis();
1365   for(Int_t ib=1; ib<=ax->GetNbins(); ib++){
1366     if(!(n1 = (Int_t)h1[0]->GetBinContent(ib))) continue;
1367     n2 = (Int_t)h1[1]->GetBinContent(ib);
1368     eff = n2/Float_t(n1);
1369
1370     ip=g->GetN();
1371     g->SetPoint(ip, ax->GetBinCenter(ib), eff);
1372     g->SetPointError(ip, 0., n2 ? eff*TMath::Sqrt(1./n1+1./n2) : 0.);
1373   }
1374 }
1375 //________________________________________________________
1376 void AliTRDcheckESD::Process2D(TH2 * const h2, TGraphErrors **g)
1377 {
1378   //
1379   // Do the processing
1380   //
1381
1382   Int_t n = 0;
1383   if((n=g[0]->GetN())) for(;n--;) g[0]->RemovePoint(n);
1384   if((n=g[1]->GetN())) for(;n--;) g[1]->RemovePoint(n);
1385   TF1 f("fg", "gaus", -3.,3.);
1386   for(Int_t ibin = 1; ibin <= h2->GetNbinsX(); ibin++){
1387     Double_t x = h2->GetXaxis()->GetBinCenter(ibin);
1388     TH1D *h = h2->ProjectionY("py", ibin, ibin);
1389     if(h->GetEntries()<100) continue;
1390     //AdjustF1(h, f);
1391
1392     h->Fit(&f, "QN");
1393     Int_t ip = g[0]->GetN();
1394     g[0]->SetPoint(ip, x, f.GetParameter(1));
1395     g[0]->SetPointError(ip, 0., f.GetParError(1));
1396     g[1]->SetPoint(ip, x, f.GetParameter(2));
1397     g[1]->SetPointError(ip, 0., f.GetParError(2));
1398   }
1399   return;
1400 }
1401 //____________________________________________________________________
1402 void AliTRDcheckESD::PrintStatus(ULong_t status)
1403 {
1404 // Dump track status to stdout
1405
1406   printf("ITS[i(%d) o(%d) r(%d)] TPC[i(%d) o(%d) r(%d) p(%d)] TRD[i(%d) o(%d) r(%d) p(%d) s(%d)] HMPID[o(%d) p(%d)]\n"
1407     ,Bool_t(status & AliESDtrack::kITSin)
1408     ,Bool_t(status & AliESDtrack::kITSout)
1409     ,Bool_t(status & AliESDtrack::kITSrefit)
1410     ,Bool_t(status & AliESDtrack::kTPCin)
1411     ,Bool_t(status & AliESDtrack::kTPCout)
1412     ,Bool_t(status & AliESDtrack::kTPCrefit)
1413     ,Bool_t(status & AliESDtrack::kTPCpid)
1414     ,Bool_t(status & AliESDtrack::kTRDin)
1415     ,Bool_t(status & AliESDtrack::kTRDout)
1416     ,Bool_t(status & AliESDtrack::kTRDrefit)
1417     ,Bool_t(status & AliESDtrack::kTRDpid)
1418     ,Bool_t(status & AliESDtrack::kTRDStop)
1419     ,Bool_t(status & AliESDtrack::kHMPIDout)
1420     ,Bool_t(status & AliESDtrack::kHMPIDpid)
1421   );
1422 }
1423
1424 //____________________________________________________________________
1425 TH1D* AliTRDcheckESD::Proj2D(TH2* hist) {
1426   //
1427   // project the PH vs Slice 2D-histo into a 1D histo
1428   //
1429   /*TH1D* hProjection = new TH1F("hProjection","", hist->GetXaxis()->GetXbins()->GetSize()-1,
1430                                hist->GetXaxis()->GetXbins()->GetArray());*/
1431   TH1D* hProjection = (TH1D*)hist->ProjectionX(Form("hProjection_%f", gRandom->Rndm()));
1432   hProjection->Reset();
1433   //cout << "Proj2D: nbins = " << hist->GetXaxis()->GetXbins()->GetSize()-1 << endl;
1434   TF1* fitLandau = new TF1("landauFunc","landau",0.,2000.);
1435   TH1D *hD;
1436   for(Int_t iBin=1;iBin<=hist->GetXaxis()->GetNbins();iBin++) {
1437     if(gROOT->FindObject("projection"))
1438       delete gROOT->FindObject("projection");
1439     hD = (TH1D*)hist->ProjectionY("projection",iBin,iBin);
1440     hD->Rebin(4);
1441     if(hD->Integral()>10) {
1442       fitLandau->SetParameter(1, hD->GetBinCenter(hD->GetMaximumBin()));
1443       fitLandau->SetParLimits(1, 0.2*hD->GetBinCenter(hD->GetMaximumBin()), 3.0*hD->GetBinCenter(hD->GetMaximumBin()));
1444       fitLandau->SetParameter(0, 1000.);
1445       fitLandau->SetParLimits(0, 1., 10000000.);
1446       fitLandau->SetParameter(2, 0.5*hD->GetBinCenter(hD->GetMaximumBin()));
1447       fitLandau->SetParLimits(2, 0.01*hD->GetBinCenter(hD->GetMaximumBin()), 1.0*hD->GetBinCenter(hD->GetMaximumBin()));
1448       hD->Fit(fitLandau, "Q0", "", hD->GetXaxis()->GetXmin(), hD->GetXaxis()->GetXmax());
1449       hD->Fit(fitLandau, "Q0", "", hD->GetXaxis()->GetXmin(), hD->GetXaxis()->GetXmax());
1450       hProjection->SetBinContent(iBin, fitLandau->GetParameter(1));
1451       hProjection->SetBinError(iBin, fitLandau->GetParameter(2));
1452     }
1453     else{
1454       hProjection->SetBinContent(iBin, 0);
1455       hProjection->SetBinError(iBin, 0);
1456     }
1457   }
1458   return hProjection;
1459 }
1460
1461 //____________________________________________________________________
1462 TH2F* AliTRDcheckESD::Proj3D(TH3* hist, TH2* accMap, Int_t zbinLow, Int_t zbinHigh, Float_t &entries) {
1463   //
1464   //  Project a 3D histogram to a 2D histogram in the Z axis interval [zbinLow,zbinHigh]
1465   //  Return the 2D histogram and also the number of entries into this projection (entries)
1466
1467   Int_t nBinsX = hist->GetXaxis()->GetNbins();   // X and Y axis bins are assumed to be all equal
1468   Float_t minX = hist->GetXaxis()->GetXmin();
1469   Float_t maxX = hist->GetXaxis()->GetXmax();
1470   Int_t nBinsY = hist->GetYaxis()->GetNbins();
1471   Float_t minY = hist->GetYaxis()->GetXmin();
1472   Float_t maxY = hist->GetYaxis()->GetXmax();
1473   Int_t nBinsZ = hist->GetZaxis()->GetNbins();  // Z axis bins (pt) might have different widths
1474
1475   TH2F* projHisto = (TH2F*)gROOT->FindObject("projHisto");
1476   if(projHisto)
1477     projHisto->Reset();
1478   else
1479     projHisto = new TH2F("projHisto", "projection", nBinsX, minX, maxX, nBinsY, minY, maxY);
1480
1481   for(Int_t iZ=1; iZ<=nBinsZ; iZ++) {
1482     if(iZ<zbinLow) continue;
1483     if(iZ>zbinHigh) continue;
1484     for(Int_t iX=1; iX<=nBinsX; iX++) {
1485       for(Int_t iY=1; iY<=nBinsY; iY++) {
1486         if(accMap) {
1487           if(accMap->GetBinContent(iX,iY)>0.1)
1488             projHisto->SetBinContent(iX, iY, projHisto->GetBinContent(iX, iY)+hist->GetBinContent(iX,iY,iZ));
1489         }
1490         else    // no acc. cut
1491           projHisto->SetBinContent(iX, iY, projHisto->GetBinContent(iX, iY)+hist->GetBinContent(iX,iY,iZ));
1492         // count only the entries which are inside the acceptance map
1493         if(accMap) {
1494           if(accMap->GetBinContent(iX,iY)>0.1)
1495             entries+=hist->GetBinContent(iX,iY,iZ);
1496         }
1497         else    // no acc. cut
1498           entries+=hist->GetBinContent(iX,iY,iZ);
1499       }
1500     }
1501   }
1502   return projHisto;
1503 }
1504
1505 //____________________________________________________________________
1506 void AliTRDcheckESD::CheckActiveSM(TH1D* phiProj, Bool_t activeSM[18]) {
1507   //
1508   // Check the active super-modules
1509   //
1510   Double_t entries[18] = {0.0};
1511   Double_t smPhiLimits[19];
1512   for(Int_t ism=0; ism<=18; ++ism) smPhiLimits[ism] = -TMath::Pi() + (2.0*TMath::Pi()/18.0)*ism;
1513   for(Int_t phiBin=1; phiBin<=phiProj->GetXaxis()->GetNbins(); ++phiBin) {
1514     Double_t phi = phiProj->GetBinCenter(phiBin);
1515     Int_t sm = -1;
1516     for(Int_t ism=0; ism<18; ++ism)
1517       if(phi>=smPhiLimits[ism] && phi<smPhiLimits[ism+1]) sm = ism;
1518     if(sm==-1) continue;
1519     entries[sm] += phiProj->GetBinContent(phiBin);
1520   }
1521   Double_t avEntries = Double_t(phiProj->Integral())/18.0;
1522   for(Int_t ism=0; ism<18; ++ism)
1523     if(entries[ism]>0.5*avEntries) activeSM[ism] = kTRUE;
1524 }
1525
1526 //____________________________________________________________________
1527 TH1F* AliTRDcheckESD::EfficiencyTRD(TH3* tpc3D, TH3* trd3D, Bool_t useAcceptance) {
1528   //
1529   // Calculate the TRD-TPC matching efficiency as function of pt
1530   //
1531
1532   if(!tpc3D || !trd3D) return NULL;
1533   Int_t nBinsZ = trd3D->GetZaxis()->GetNbins();
1534   // project everything on the eta-phi map to obtain an acceptance map
1535   Float_t nada = 0.;
1536   TH2F *trdAcc = (useAcceptance ? (TH2F*)Proj3D(trd3D, 0x0, 1, nBinsZ, nada)->Clone(Form("trdAcc%f", gRandom->Rndm())) : 0x0);
1537   TH1D *phiProj = (trdAcc ? trdAcc->ProjectionY(Form("phiProj%f", gRandom->Rndm())) : 0x0);
1538
1539   // prepare the acceptance map
1540   Bool_t activeSM[18] = {kFALSE};
1541   Double_t smPhiLimits[19];
1542   for(Int_t ism=0; ism<=18; ++ism) smPhiLimits[ism] = -TMath::Pi() + (2.0*TMath::Pi()/18.0)*ism;
1543   if(phiProj) {
1544     CheckActiveSM(phiProj, activeSM);   // get the active SMs
1545     trdAcc->Reset();
1546     // Put 1 entry in every bin which belongs to an active SM
1547     for(Int_t iY=1; iY<=trdAcc->GetYaxis()->GetNbins(); ++iY) {
1548       Double_t phi = trdAcc->GetYaxis()->GetBinCenter(iY);
1549       Bool_t isActive = kFALSE;
1550       for(Int_t ism=0; ism<18; ++ism) {
1551         if(phi>=smPhiLimits[ism] && phi<smPhiLimits[ism+1] && activeSM[ism]) {
1552           isActive = kTRUE;
1553         }
1554       }
1555       if(!isActive) continue;
1556       for(Int_t iX=1; iX<=trdAcc->GetXaxis()->GetNbins(); ++iX)
1557         if(trdAcc->GetXaxis()->GetBinCenter(iX)>=-0.85 && trdAcc->GetXaxis()->GetBinCenter(iX)<=0.85) trdAcc->SetBinContent(iX, iY, 1.0);
1558     }  // end for over Y(phi) bins
1559   }  // end if phiProj
1560
1561   // get the bin limits from the Z axis of 3D histos
1562   Float_t *ptBinLimits = new Float_t[nBinsZ+1];
1563   for(Int_t i=1; i<=nBinsZ; i++) {
1564     ptBinLimits[i-1] = trd3D->GetZaxis()->GetBinLowEdge(i);
1565   }
1566   ptBinLimits[nBinsZ] = trd3D->GetZaxis()->GetBinUpEdge(nBinsZ);
1567
1568   TH1F *efficiency = new TH1F(Form("eff%d", Int_t(1000000.0*gRandom->Rndm())), "TRD-TPC matching efficiency", nBinsZ, ptBinLimits);
1569
1570   // loop over Z bins
1571   Bool_t effGood = kFALSE;
1572   for(Int_t i=1; i<=nBinsZ; i++) {
1573     Float_t tpcEntries = 0.0; Float_t trdEntries = 0.0;
1574     Proj3D(tpc3D, trdAcc, i, i, tpcEntries);
1575     Proj3D(trd3D, trdAcc, i, i, trdEntries);
1576     Float_t ratio = 0;
1577     if(tpcEntries>0) ratio = trdEntries/tpcEntries;
1578     Float_t error = 0;
1579     if(tpcEntries>0 && trdEntries>0 && (tpcEntries-trdEntries)>=0.0)
1580       error = TMath::Sqrt(trdEntries*(tpcEntries-trdEntries)/tpcEntries/tpcEntries/tpcEntries);
1581     if(ratio>0.001) {
1582       efficiency->SetBinContent(i,ratio);
1583       efficiency->SetBinError(i,error);
1584       effGood = kTRUE;
1585     }
1586   }     // end loop over Z bins
1587   if(!effGood) return 0x0;
1588
1589   return efficiency;
1590 }
1591
1592 //__________________________________________________________________________________________________
1593 void AliTRDcheckESD::PlotCentSummaryFromCF(Double_t* trendValues, Bool_t useIsolatedBC, Bool_t cutTOFbc) {
1594   //
1595   // Make the centrality summary figure from the CF container
1596   //
1597   if(!fCfContainer) return;
1598
1599   trendValues = trendValues;
1600
1601   TLatex* lat=new TLatex();
1602   lat->SetTextSize(0.06);
1603   lat->SetTextColor(2);
1604
1605   gPad->SetTopMargin(0.05); gPad->SetBottomMargin(0.001); gPad->SetLeftMargin(0.001); gPad->SetRightMargin(0.001);
1606   gPad->Divide(3,3,0.,0.);
1607   TList* l=gPad->GetListOfPrimitives();
1608   TVirtualPad* pad=0x0;
1609
1610   fCfContainer->SetRangeUser(kTrackDCAxy, -0.999, +0.999);
1611   fCfContainer->SetRangeUser(kTrackDCAz, -3.0, +3.0);
1612   if(cutTOFbc) fCfContainer->SetRangeUser(kTrackTOFBC, 0.0, 0.0);
1613
1614   // find all the isolated bunch crossings with entries
1615   Bool_t isIsolated[3500];
1616   TH1D* tempTH1D = (TH1D*)fCfContainer->Project(0, kEventBC);
1617   FindIsolatedBCs(tempTH1D, isIsolated); delete tempTH1D;
1618   Int_t nIsolatedBC = 0;
1619   for(Int_t ibc=0; ibc<3500; ++ibc)
1620     if(isIsolated[ibc]) nIsolatedBC++;
1621
1622   if(gROOT->FindObject("rangeEffPt")) delete gROOT->FindObject("rangeEffPt");
1623   TH2F* rangeEffPt=new TH2F("rangeEffPt", "",10,0.,10.,10,0.,1.3);
1624   rangeEffPt->SetStats(kFALSE);
1625   SetStyle(rangeEffPt->GetXaxis(), "p_{T} [GeV/c]", 0.07, 0.8, kTRUE, 0.05);
1626   SetStyle(rangeEffPt->GetYaxis(), "efficiency", 0.07, 0.8, kTRUE, 0.05);
1627
1628   Int_t padsForEffs[5] = {0,3,6,1,4};
1629   for(Int_t iCent=1; iCent<6; ++iCent) {
1630     pad = ((TVirtualPad*)l->At(padsForEffs[iCent-1])); pad->cd();
1631     pad->SetLeftMargin(0.15); pad->SetRightMargin(0.02); pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
1632     pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
1633     rangeEffPt->Draw();
1634     TLine line;
1635     line.SetLineStyle(2);
1636     line.SetLineWidth(2);
1637     line.DrawLine(rangeEffPt->GetXaxis()->GetXmin(), 0.7, rangeEffPt->GetXaxis()->GetXmax(), 0.7);
1638     line.DrawLine(rangeEffPt->GetXaxis()->GetXmin(), 0.9, rangeEffPt->GetXaxis()->GetXmax(), 0.9);
1639
1640     fCfContainer->SetRangeUser(kEventMult, Double_t(iCent), Double_t(iCent), kTRUE);
1641     fCfContainer->SetRangeUser(kTrackCharge, +1.0, +1.0);        // positive charges
1642     fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
1643
1644     TH3D* h3PosTPC = (TH3D*)fCfContainer->Project(0, kTrackEta, kTrackPhi, kTrackPt);
1645     if(h3PosTPC->GetEntries()<10) {
1646       delete h3PosTPC;
1647       continue;
1648     }
1649     TH3D* h3PosTRDall = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1650     fCfContainer->SetRangeUser(kTrackTrdTracklets, 4.0, 4.0);        // >= 4 TRD tracklets
1651     TH3D* h3PosTRDtrk4 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1652     fCfContainer->SetRangeUser(kTrackTrdTracklets, 5.0, 5.0);        // >= 5 TRD tracklets
1653     TH3D* h3PosTRDtrk5 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1654     fCfContainer->SetRangeUser(kTrackTrdTracklets, 6.0, 6.0);        // >= 6 TRD tracklets
1655     TH3D* h3PosTRDtrk6 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1656
1657     fCfContainer->SetRangeUser(kTrackCharge, -1.0, -1.0);        // negative charges
1658     fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);    // >= 0 TRD tracklets
1659     TH3D* h3NegTPC = (TH3D*)fCfContainer->Project(0, kTrackEta, kTrackPhi, kTrackPt);
1660     TH3D* h3NegTRDall = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1661     fCfContainer->SetRangeUser(kTrackTrdTracklets, 4.0, 4.0);        // 4 TRD tracklets
1662     TH3D* h3NegTRDtrk4 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1663     fCfContainer->SetRangeUser(kTrackTrdTracklets, 5.0, 5.0);        // 5 TRD tracklets
1664     TH3D* h3NegTRDtrk5 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1665     fCfContainer->SetRangeUser(kTrackTrdTracklets, 6.0, 6.0);        // 6 TRD tracklets
1666     TH3D* h3NegTRDtrk6 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1667
1668     TH1F* hEffPosAll = EfficiencyTRD(h3PosTPC, h3PosTRDall, kTRUE);
1669     TH1F* hEffPosTrk4 = EfficiencyTRD(h3PosTPC, h3PosTRDtrk4, kTRUE);
1670     TH1F* hEffPosTrk5 = EfficiencyTRD(h3PosTPC, h3PosTRDtrk5, kTRUE);
1671     TH1F* hEffPosTrk6 = EfficiencyTRD(h3PosTPC, h3PosTRDtrk6, kTRUE);
1672     TH1F* hEffNegAll = EfficiencyTRD(h3NegTPC, h3NegTRDall, kTRUE);
1673     TH1F* hEffNegTrk4 = EfficiencyTRD(h3NegTPC, h3NegTRDtrk4, kTRUE);
1674     TH1F* hEffNegTrk5 = EfficiencyTRD(h3NegTPC, h3NegTRDtrk5, kTRUE);
1675     TH1F* hEffNegTrk6 = EfficiencyTRD(h3NegTPC, h3NegTRDtrk6, kTRUE);
1676
1677     delete h3PosTPC; delete h3NegTPC;
1678     delete h3PosTRDall; delete h3PosTRDtrk4; delete h3PosTRDtrk5; delete h3PosTRDtrk6;
1679     delete h3NegTRDall; delete h3NegTRDtrk4; delete h3NegTRDtrk5; delete h3NegTRDtrk6;
1680
1681     // get matching efficiencies for isolated bunches
1682     TH3D* h3TPCrefPos_IsolatedBC=0x0; TH3D* h3TPCrefNeg_IsolatedBC=0x0;
1683     TH3D* h3TRDrefPosAll_IsolatedBC=0x0; TH3D* h3TRDrefPosTrk4_IsolatedBC=0x0; TH3D* h3TRDrefPosTrk5_IsolatedBC=0x0; TH3D* h3TRDrefPosTrk6_IsolatedBC=0x0;
1684     TH3D* h3TRDrefNegAll_IsolatedBC=0x0; TH3D* h3TRDrefNegTrk4_IsolatedBC=0x0; TH3D* h3TRDrefNegTrk5_IsolatedBC=0x0; TH3D* h3TRDrefNegTrk6_IsolatedBC=0x0;
1685     for(Int_t ibc=0; ibc<3500; ++ibc) {
1686       if(!isIsolated[ibc]) continue;
1687       fCfContainer->SetRangeUser(kEventBC, Double_t(ibc), Double_t(ibc));
1688
1689       TH3D* tempTH3D;
1690       fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
1691       fCfContainer->SetRangeUser(kTrackCharge, +1.0, +1.0);    // positive charges
1692       tempTH3D = (TH3D*)fCfContainer->Project(0, kTrackEta, kTrackPhi, kTrackPt);
1693       if(!h3TPCrefPos_IsolatedBC) h3TPCrefPos_IsolatedBC = tempTH3D;
1694       else {h3TPCrefPos_IsolatedBC->Add(tempTH3D); delete tempTH3D;};
1695       tempTH3D = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1696       if(!h3TRDrefPosAll_IsolatedBC) h3TRDrefPosAll_IsolatedBC = tempTH3D;
1697       else {h3TRDrefPosAll_IsolatedBC->Add(tempTH3D); delete tempTH3D;}
1698       fCfContainer->SetRangeUser(kTrackTrdTracklets, 4.0, 4.0);
1699       tempTH3D = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1700       if(!h3TRDrefPosTrk4_IsolatedBC) h3TRDrefPosTrk4_IsolatedBC = tempTH3D;
1701       else {h3TRDrefPosTrk4_IsolatedBC->Add(tempTH3D); delete tempTH3D;}
1702       fCfContainer->SetRangeUser(kTrackTrdTracklets, 5.0, 5.0);
1703       tempTH3D = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1704       if(!h3TRDrefPosTrk5_IsolatedBC) h3TRDrefPosTrk5_IsolatedBC = tempTH3D;
1705       else {h3TRDrefPosTrk5_IsolatedBC->Add(tempTH3D); delete tempTH3D;}
1706       fCfContainer->SetRangeUser(kTrackTrdTracklets, 6.0, 6.0);
1707       tempTH3D = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1708       if(!h3TRDrefPosTrk6_IsolatedBC) h3TRDrefPosTrk6_IsolatedBC = tempTH3D;
1709       else {h3TRDrefPosTrk6_IsolatedBC->Add(tempTH3D); delete tempTH3D;}
1710
1711       fCfContainer->SetRangeUser(kTrackCharge, -1.0, -1.0);   // negative charges
1712       fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
1713       tempTH3D = (TH3D*)fCfContainer->Project(0, kTrackEta, kTrackPhi, kTrackPt);
1714       if(!h3TPCrefNeg_IsolatedBC) h3TPCrefNeg_IsolatedBC = tempTH3D;
1715       else {h3TPCrefNeg_IsolatedBC->Add(tempTH3D); delete tempTH3D;}
1716       tempTH3D = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1717       if(!h3TRDrefNegAll_IsolatedBC) h3TRDrefNegAll_IsolatedBC = tempTH3D;
1718       else {h3TRDrefNegAll_IsolatedBC->Add(tempTH3D); delete tempTH3D;}
1719       fCfContainer->SetRangeUser(kTrackTrdTracklets, 4.0, 4.0);
1720       tempTH3D = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1721       if(!h3TRDrefNegTrk4_IsolatedBC) h3TRDrefNegTrk4_IsolatedBC = tempTH3D;
1722       else {h3TRDrefNegTrk4_IsolatedBC->Add(tempTH3D); delete tempTH3D;}
1723       fCfContainer->SetRangeUser(kTrackTrdTracklets, 5.0, 5.0);
1724       tempTH3D = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1725       if(!h3TRDrefNegTrk5_IsolatedBC) h3TRDrefNegTrk5_IsolatedBC = tempTH3D;
1726       else {h3TRDrefNegTrk5_IsolatedBC->Add(tempTH3D); delete tempTH3D;}
1727       fCfContainer->SetRangeUser(kTrackTrdTracklets, 6.0, 6.0);
1728       tempTH3D = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
1729       if(!h3TRDrefNegTrk6_IsolatedBC) h3TRDrefNegTrk6_IsolatedBC = tempTH3D;
1730       else {h3TRDrefNegTrk6_IsolatedBC->Add(tempTH3D); delete tempTH3D;}
1731     }
1732     fCfContainer->SetRangeUser(kEventBC, 0.0, 3500.0);
1733     fCfContainer->SetRangeUser(kTrackCharge, -1.0, +1.0);
1734     fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
1735
1736     TH1F* hEffPosAll_IsolatedBC  = EfficiencyTRD(h3TPCrefPos_IsolatedBC, h3TRDrefPosAll_IsolatedBC, kTRUE);
1737     TH1F* hEffPosTrk4_IsolatedBC = EfficiencyTRD(h3TPCrefPos_IsolatedBC, h3TRDrefPosTrk4_IsolatedBC, kTRUE);
1738     TH1F* hEffPosTrk5_IsolatedBC = EfficiencyTRD(h3TPCrefPos_IsolatedBC, h3TRDrefPosTrk5_IsolatedBC, kTRUE);
1739     TH1F* hEffPosTrk6_IsolatedBC = EfficiencyTRD(h3TPCrefPos_IsolatedBC, h3TRDrefPosTrk6_IsolatedBC, kTRUE);
1740     TH1F* hEffNegAll_IsolatedBC  = EfficiencyTRD(h3TPCrefNeg_IsolatedBC, h3TRDrefNegAll_IsolatedBC, kTRUE);
1741     TH1F* hEffNegTrk4_IsolatedBC = EfficiencyTRD(h3TPCrefNeg_IsolatedBC, h3TRDrefNegTrk4_IsolatedBC, kTRUE);
1742     TH1F* hEffNegTrk5_IsolatedBC = EfficiencyTRD(h3TPCrefNeg_IsolatedBC, h3TRDrefNegTrk5_IsolatedBC, kTRUE);
1743     TH1F* hEffNegTrk6_IsolatedBC = EfficiencyTRD(h3TPCrefNeg_IsolatedBC, h3TRDrefNegTrk6_IsolatedBC, kTRUE);
1744
1745     if(h3TPCrefPos_IsolatedBC) delete h3TPCrefPos_IsolatedBC;
1746     if(h3TPCrefNeg_IsolatedBC) delete h3TPCrefNeg_IsolatedBC;
1747     if(h3TRDrefPosAll_IsolatedBC) delete h3TRDrefPosAll_IsolatedBC;
1748     if(h3TRDrefPosTrk4_IsolatedBC) delete h3TRDrefPosTrk4_IsolatedBC;
1749     if(h3TRDrefPosTrk5_IsolatedBC) delete h3TRDrefPosTrk5_IsolatedBC;
1750     if(h3TRDrefPosTrk6_IsolatedBC) delete h3TRDrefPosTrk6_IsolatedBC;
1751     if(h3TRDrefNegAll_IsolatedBC) delete h3TRDrefNegAll_IsolatedBC;
1752     if(h3TRDrefNegTrk4_IsolatedBC) delete h3TRDrefNegTrk4_IsolatedBC;
1753     if(h3TRDrefNegTrk5_IsolatedBC) delete h3TRDrefNegTrk5_IsolatedBC;
1754     if(h3TRDrefNegTrk6_IsolatedBC) delete h3TRDrefNegTrk6_IsolatedBC;
1755
1756     if(!useIsolatedBC) {
1757       SetStyle(hEffPosAll,  1, kRed, 1, 24, kRed, 1);
1758       SetStyle(hEffPosTrk4, 1, kRed, 1, 25, kRed, 1);
1759       SetStyle(hEffPosTrk5, 1, kRed, 1, 26, kRed, 1);
1760       SetStyle(hEffPosTrk6, 1, kRed, 1, 27, kRed, 1);
1761       SetStyle(hEffNegAll,  1, kBlue, 1, 24, kBlue, 1);
1762       SetStyle(hEffNegTrk4, 1, kBlue, 1, 25, kBlue, 1);
1763       SetStyle(hEffNegTrk5, 1, kBlue, 1, 26, kBlue, 1);
1764       SetStyle(hEffNegTrk6, 1, kBlue, 1, 27, kBlue, 1);
1765       hEffPosAll->Draw("same");
1766       hEffNegAll->Draw("same");
1767       hEffPosTrk4->Draw("same");
1768       hEffNegTrk4->Draw("same");
1769       hEffPosTrk5->Draw("same");
1770       hEffNegTrk5->Draw("same");
1771       hEffPosTrk6->Draw("same");
1772       hEffNegTrk6->Draw("same");
1773     }
1774     else {
1775       if(nIsolatedBC>0) {
1776         SetStyle(hEffPosAll_IsolatedBC,  1, kRed, 1, 24, kRed, 1);
1777         SetStyle(hEffPosTrk4_IsolatedBC, 1, kRed, 1, 25, kRed, 1);
1778         SetStyle(hEffPosTrk5_IsolatedBC, 1, kRed, 1, 26, kRed, 1);
1779         SetStyle(hEffPosTrk6_IsolatedBC, 1, kRed, 1, 27, kRed, 1);
1780         SetStyle(hEffNegAll_IsolatedBC,  1, kBlue, 1, 24, kBlue, 1);
1781         SetStyle(hEffNegTrk4_IsolatedBC, 1, kBlue, 1, 25, kBlue, 1);
1782         SetStyle(hEffNegTrk5_IsolatedBC, 1, kBlue, 1, 26, kBlue, 1);
1783         SetStyle(hEffNegTrk6_IsolatedBC, 1, kBlue, 1, 27, kBlue, 1);
1784         hEffPosAll_IsolatedBC->Draw("same");
1785         hEffNegAll_IsolatedBC->Draw("same");
1786         hEffPosTrk4_IsolatedBC->Draw("same");
1787         hEffNegTrk4_IsolatedBC->Draw("same");
1788         hEffPosTrk5_IsolatedBC->Draw("same");
1789         hEffNegTrk5_IsolatedBC->Draw("same");
1790         hEffPosTrk6_IsolatedBC->Draw("same");
1791         hEffNegTrk6_IsolatedBC->Draw("same");
1792       }
1793     }
1794
1795     TLegend* leg=new TLegend(0.18, 0.7, 0.77, 0.89);
1796     if(iCent==1) {
1797       leg->SetFillColor(0);
1798       leg->SetNColumns(2);
1799       leg->SetMargin(0.1);
1800       leg->SetBorderSize(0);
1801       if(useIsolatedBC) {
1802         leg->SetHeader("Isolated bunch crossings");
1803         if(nIsolatedBC) {
1804           leg->AddEntry(hEffPosAll_IsolatedBC,  "pos. (#geq 1 tracklet)", "p");
1805           leg->AddEntry(hEffNegAll_IsolatedBC,  "neg. (#geq 1 tracklet)", "p");
1806           leg->AddEntry(hEffPosTrk4_IsolatedBC, "pos. (4 tracklets)", "p");
1807           leg->AddEntry(hEffNegTrk4_IsolatedBC, "neg. (4 tracklets)", "p");
1808           leg->AddEntry(hEffPosTrk5_IsolatedBC, "pos. (5 tracklets)", "p");
1809           leg->AddEntry(hEffNegTrk5_IsolatedBC, "neg. (5 tracklets)", "p");
1810           leg->AddEntry(hEffPosTrk6_IsolatedBC, "pos. (6 tracklets)", "p");
1811           leg->AddEntry(hEffNegTrk6_IsolatedBC, "neg. (6 tracklets)", "p");
1812         }
1813       }
1814       else {
1815         leg->SetHeader("All bunch crossings");
1816         leg->AddEntry(hEffPosAll,  "pos. (#geq 1 tracklet)", "p");
1817         leg->AddEntry(hEffNegAll,  "neg. (#geq 1 tracklet)", "p");
1818         leg->AddEntry(hEffPosTrk4, "pos. (4 tracklets)", "p");
1819         leg->AddEntry(hEffNegTrk4, "neg. (4 tracklets)", "p");
1820         leg->AddEntry(hEffPosTrk5, "pos. (5 tracklets)", "p");
1821         leg->AddEntry(hEffNegTrk5, "neg. (5 tracklets)", "p");
1822         leg->AddEntry(hEffPosTrk6, "pos. (6 tracklets)", "p");
1823         leg->AddEntry(hEffNegTrk6, "neg. (6 tracklets)", "p");
1824       }
1825       leg->Draw();
1826     }
1827     lat->DrawLatex(0.2, 1.32, Form("Centrality class %d", iCent));
1828   }   // end for loop over multiplicity classes
1829
1830   // Reset the modified user ranges of the CF container
1831   fCfContainer->SetRangeUser(kEventMult, 0, 6, kTRUE);
1832   fCfContainer->SetRangeUser(kTrackCharge, -1.0, +1.0);
1833   fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
1834
1835   // Cluster distributions in all multiplicity classes
1836   pad = ((TVirtualPad*)l->At(2)); pad->cd();
1837   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.02);
1838   pad->SetTopMargin(0.02); pad->SetBottomMargin(0.15);
1839   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
1840   if(gROOT->FindObject("rangeNcls")) delete gROOT->FindObject("rangeNcls");
1841   TH2F* rangeNcls = new TH2F("rangeNcls", "", 10, 0.0, 199.9, 10, 0.0, 1.199);
1842   SetStyle(rangeNcls->GetXaxis(), "# TRD clusters", 0.07, 0.8, kTRUE, 0.05);
1843   SetStyle(rangeNcls->GetYaxis(), "entries (a.u.)", 0.07, 0.8, kTRUE, 0.05);
1844   rangeNcls->SetStats(kFALSE);
1845   rangeNcls->Draw();
1846
1847   TH1D* hNcls[6]={0x0};
1848   TLegend* legCls=new TLegend(0.7, 0.75, 0.97, 0.97);
1849   legCls->SetBorderSize(0);
1850   legCls->SetFillColor(0);
1851   legCls->SetMargin(0.15);
1852   for(Int_t iCent=0; iCent<6; ++iCent) {
1853     if(iCent>0)
1854       fCfContainer->SetRangeUser(kEventMult, Double_t(iCent), Double_t(iCent), kTRUE);
1855     if(!useIsolatedBC) hNcls[iCent] = (TH1D*)fCfContainer->Project(1, kTrackTrdClusters);
1856     else{
1857       for(Int_t ibc=0; ibc<3500; ++ibc) {
1858         if(!isIsolated[ibc]) continue;
1859         fCfContainer->SetRangeUser(kEventBC, Double_t(ibc), Double_t(ibc));
1860         tempTH1D = (TH1D*)fCfContainer->Project(1, kTrackTrdClusters);
1861         if(!hNcls[iCent]) hNcls[iCent] = tempTH1D;
1862         else {hNcls[iCent]->Add(tempTH1D); delete tempTH1D;}
1863       }
1864       fCfContainer->SetRangeUser(kEventBC, 0.0, 3500.0);
1865     }
1866     if(!hNcls[iCent]) continue;
1867     hNcls[iCent]->SetLineColor(iCent<4 ? iCent+1 : iCent+2);
1868     Double_t maximum = hNcls[iCent]->GetMaximum();
1869     if(maximum>1.0)
1870       hNcls[iCent]->Scale(1.0/maximum);
1871     hNcls[iCent]->SetStats(kFALSE);
1872     hNcls[iCent]->SetTitle("");
1873     hNcls[iCent]->SetLineWidth(2);
1874
1875     if(hNcls[iCent]->Integral()>0.01) {
1876       hNcls[iCent]->Draw("same");
1877       legCls->AddEntry(hNcls[iCent], (iCent==0 ? "all centralities" : Form("centrality class %d", iCent)), "l");
1878     }
1879   }
1880   if(useIsolatedBC) legCls->SetHeader("Isolated bunch crossings");
1881   else legCls->SetHeader("All bunch crossings");
1882   legCls->Draw();
1883   fCfContainer->SetRangeUser(kEventMult, 0.0, 6.0, kTRUE);
1884
1885   // Qtot vs P
1886   pad = ((TVirtualPad*)l->At(5)); pad->cd();
1887   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.02);
1888   pad->SetTopMargin(0.02); pad->SetBottomMargin(0.15);
1889   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
1890   if(gROOT->FindObject("rangeQtot")) delete gROOT->FindObject("rangeQtot");
1891   TH2F* rangeQtot = new TH2F("rangeQtot", "", 10, 0.0, 9.999, 10, 0.0, 1.199);
1892   SetStyle(rangeQtot->GetXaxis(), "Q_{tot} (a.u.)", 0.07, 0.8, kTRUE, 0.05);
1893   SetStyle(rangeQtot->GetYaxis(), "entries (a.u.)", 0.07, 0.8, kTRUE, 0.05);
1894   rangeQtot->SetStats(kFALSE);
1895   rangeQtot->Draw();
1896
1897   TH1D* hQtot[6]={0x0};
1898   TLegend* leg2=new TLegend(0.6, 0.7, 0.9, 0.97);
1899   leg2->SetFillColor(0);
1900   leg2->SetBorderSize(0);
1901   for(Int_t iCent=0; iCent<6; ++iCent) {
1902     if(iCent>0)
1903       fCfContainer->SetRangeUser(kEventMult, Double_t(iCent), Double_t(iCent), kTRUE);
1904     if(useIsolatedBC) {
1905       for(Int_t ibc=0; ibc<3500; ++ibc) {
1906         if(!isIsolated[ibc]) continue;
1907         fCfContainer->SetRangeUser(kEventBC, Double_t(ibc), Double_t(ibc));
1908         for(Int_t il=0; il<6; ++il) {
1909           tempTH1D = (TH1D*)fCfContainer->Project(1, kTrackQtot+il);
1910           if(!hQtot[iCent] && il==0) hQtot[iCent] = tempTH1D;
1911           else {hQtot[iCent]->Add(tempTH1D); delete tempTH1D;}
1912         }
1913       }
1914       fCfContainer->SetRangeUser(kEventBC, 0.0, 3500.0);
1915     }
1916     else {
1917       for(Int_t il=0; il<6; ++il) {
1918         tempTH1D = (TH1D*)fCfContainer->Project(1, kTrackQtot+il);
1919         if(il==0) hQtot[iCent] = tempTH1D;
1920         else hQtot[iCent]->Add(tempTH1D);
1921       }
1922     } // end if(useIsolatedBC)
1923
1924     if(!hQtot[iCent]) continue;
1925     hQtot[iCent]->SetBinContent(1, 0);
1926     Double_t maximum = hQtot[iCent]->GetMaximum();
1927     if(maximum>1.0)
1928       hQtot[iCent]->Scale(1.0/maximum);
1929     hQtot[iCent]->SetLineColor(iCent<4 ? iCent+1 : iCent+2);
1930     hQtot[iCent]->SetStats(kFALSE);
1931     hQtot[iCent]->SetTitle("");
1932     hQtot[iCent]->SetLineWidth(2);
1933     if(hQtot[iCent]->Integral()>0.01) {
1934       hQtot[iCent]->Draw(iCent==0 ? "" : "same");
1935       leg2->AddEntry(hQtot[iCent], (iCent==0 ? "all centralities" : Form("centrality class %d", iCent)), "l");
1936     }
1937   }
1938   if(useIsolatedBC) leg2->SetHeader("Isolated bunch crossings");
1939   else leg2->SetHeader("All bunch crossings");
1940   leg2->Draw();
1941   fCfContainer->SetRangeUser(kEventMult, 0.0, 5.0, kTRUE);
1942   if(cutTOFbc) fCfContainer->SetRangeUser(kTrackTOFBC, -1000.0, +1000.0);  // reset the cut on TOFbc
1943 }
1944
1945
1946 //_________________________________________________________________
1947 void AliTRDcheckESD::PlotTrackingSummaryFromCF(Int_t centralityClass, Double_t* trendValues, Bool_t useIsolatedBC, Bool_t cutTOFbc) {
1948
1949   if(!fCfContainer) return;
1950
1951   TLatex *lat=new TLatex();
1952   lat->SetTextSize(0.06);
1953   lat->SetTextColor(2);
1954
1955   gPad->SetTopMargin(0.05); gPad->SetBottomMargin(0.001);
1956   gPad->SetLeftMargin(0.001); gPad->SetRightMargin(0.001);
1957   gPad->Divide(3,3,0.,0.);
1958   TList* l=gPad->GetListOfPrimitives();
1959
1960   // eta-phi distr. for positive TPC tracks
1961   TVirtualPad* pad = ((TVirtualPad*)l->At(0)); pad->cd();
1962   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
1963   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
1964   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
1965
1966   fCfContainer->SetRangeUser(kTrackDCAxy, -0.999, +0.999);
1967   fCfContainer->SetRangeUser(kTrackDCAz, -3.0, +3.0);
1968   if(cutTOFbc) fCfContainer->SetRangeUser(kTrackTOFBC, 0.0, 0.0);
1969
1970   // find all the isolated bunch crossings with entries
1971   Bool_t isIsolated[3500];
1972   TH1D* tempTH1D = (TH1D*)fCfContainer->Project(0, kEventBC);
1973   FindIsolatedBCs(tempTH1D, isIsolated); delete tempTH1D;
1974   Int_t nIsolatedBC = 0;
1975   for(Int_t ibc=0; ibc<3500; ++ibc)
1976     if(isIsolated[ibc]) nIsolatedBC++;
1977
1978   if(centralityClass>0) // select the multiplicity class
1979     fCfContainer->SetRangeUser(kEventMult, Double_t(centralityClass), Double_t(centralityClass), kTRUE);
1980
1981   TH2D* hTPCrefPos = 0x0; TH2D* hTRDrefPos = 0x0; TH2D* hTOFrefPos = 0x0;
1982   TH2D* hTPCrefNeg = 0x0; TH2D* hTRDrefNeg = 0x0; TH2D* hTOFrefNeg = 0x0;
1983   if(!useIsolatedBC) {
1984     fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
1985     fCfContainer->SetRangeUser(kTrackCharge, +1.0, +1.0);      // positive charges
1986     hTPCrefPos = (TH2D*)fCfContainer->Project(0, kTrackEta, kTrackPhi);
1987     hTRDrefPos = (TH2D*)fCfContainer->Project(1, kTrackEta, kTrackPhi);
1988     hTOFrefPos = (TH2D*)fCfContainer->Project(2, kTrackEta, kTrackPhi);
1989     fCfContainer->SetRangeUser(kTrackCharge, -1.0, -1.0);      // negative charges
1990     hTPCrefNeg = (TH2D*)fCfContainer->Project(0, kTrackEta, kTrackPhi);
1991     hTRDrefNeg = (TH2D*)fCfContainer->Project(1, kTrackEta, kTrackPhi);
1992     hTOFrefNeg = (TH2D*)fCfContainer->Project(2, kTrackEta, kTrackPhi);
1993     fCfContainer->SetRangeUser(kTrackCharge, -1.0, +1.0);      // reset charge cut
1994   }
1995
1996   TH2D* hTPCrefPos_IsolatedBC=0x0; TH2D* hTPCrefNeg_IsolatedBC=0x0;
1997   TH2D* hTRDrefPos_IsolatedBC=0x0; TH2D* hTRDrefNeg_IsolatedBC=0x0;
1998   TH2D* hTOFrefPos_IsolatedBC=0x0; TH2D* hTOFrefNeg_IsolatedBC=0x0;
1999   if(useIsolatedBC) {
2000     fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
2001     for(Int_t ibc=0; ibc<3500; ++ibc) {
2002       if(!isIsolated[ibc]) continue;
2003       fCfContainer->SetRangeUser(kEventBC, Double_t(ibc), Double_t(ibc));
2004       fCfContainer->SetRangeUser(kTrackCharge, +1.0, +1.0);      // positive charges
2005       TH2D* tempTH2D;
2006       tempTH2D = (TH2D*)fCfContainer->Project(0, kTrackEta, kTrackPhi);
2007       if(!hTPCrefPos_IsolatedBC) hTPCrefPos_IsolatedBC = tempTH2D;
2008       else {hTPCrefPos_IsolatedBC->Add(tempTH2D); delete tempTH2D;}
2009       tempTH2D = (TH2D*)fCfContainer->Project(1, kTrackEta, kTrackPhi);
2010       if(!hTRDrefPos_IsolatedBC) hTRDrefPos_IsolatedBC = tempTH2D;
2011       else {hTRDrefPos_IsolatedBC->Add(tempTH2D); delete tempTH2D;}
2012       tempTH2D = (TH2D*)fCfContainer->Project(2, kTrackEta, kTrackPhi);
2013       if(!hTOFrefPos_IsolatedBC) hTOFrefPos_IsolatedBC = tempTH2D;
2014       else {hTOFrefPos_IsolatedBC->Add(tempTH2D); delete tempTH2D;}
2015       fCfContainer->SetRangeUser(kTrackCharge, -1.0, -1.0);      // negative charges
2016       tempTH2D = (TH2D*)fCfContainer->Project(0, kTrackEta, kTrackPhi);
2017       if(!hTPCrefNeg_IsolatedBC) hTPCrefNeg_IsolatedBC = tempTH2D;
2018       else {hTPCrefNeg_IsolatedBC->Add(tempTH2D); delete tempTH2D;}
2019       tempTH2D = (TH2D*)fCfContainer->Project(1, kTrackEta, kTrackPhi);
2020       if(!hTRDrefNeg_IsolatedBC) hTRDrefNeg_IsolatedBC = tempTH2D;
2021       else {hTRDrefNeg_IsolatedBC->Add(tempTH2D); delete tempTH2D;}
2022       tempTH2D = (TH2D*)fCfContainer->Project(2, kTrackEta, kTrackPhi);
2023       if(!hTOFrefNeg_IsolatedBC) hTOFrefNeg_IsolatedBC = tempTH2D;
2024       else {hTOFrefNeg_IsolatedBC->Add(tempTH2D); delete tempTH2D;}
2025     }
2026     fCfContainer->SetRangeUser(kEventBC, 0.0, 3500.);   // reset the BC range
2027     fCfContainer->SetRangeUser(kTrackCharge, -1.0, +1.0);      // reset charge cut
2028   }
2029
2030   if(gROOT->FindObject("rangeEtaPhi")) delete gROOT->FindObject("rangeEtaPhi");
2031   TH2F* rangeEtaPhi = new TH2F("rangeEtaPhi", "", 10, -0.99, +0.99, 10, -3.4, +3.4);
2032   SetStyle(rangeEtaPhi->GetXaxis(), "#eta", 0.07, 0.8, kTRUE, 0.05);
2033   SetStyle(rangeEtaPhi->GetYaxis(), "detector #varphi", 0.07, 0.8, kTRUE, 0.05);
2034   rangeEtaPhi->SetStats(kFALSE);
2035
2036   //----------------------------------------------
2037   // eta-phi efficiency for positive TRD tracks
2038   pad = ((TVirtualPad*)l->At(0)); pad->cd();
2039   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
2040   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
2041   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2042   rangeEtaPhi->Draw();
2043
2044   if(useIsolatedBC) {
2045     if(nIsolatedBC>0) {
2046       TH2D* hTRDeffPos_IsolatedBC = (TH2D*)hTRDrefPos_IsolatedBC->Clone("hTRDeffPos_IsolatedBC");
2047       hTRDeffPos_IsolatedBC->Reset();
2048       hTRDeffPos_IsolatedBC->SetStats(kFALSE);
2049       hTRDeffPos_IsolatedBC->Divide(hTRDrefPos_IsolatedBC, hTPCrefPos_IsolatedBC);
2050       hTRDeffPos_IsolatedBC->SetMaximum(1.0);
2051       hTRDeffPos_IsolatedBC->Draw("samecolz");
2052     }
2053   }
2054   else {
2055     TH2D* hTRDeffPos = (TH2D*)hTRDrefPos->Clone("hTRDeffPos");
2056     hTRDeffPos->Reset();
2057     hTRDeffPos->SetStats(kFALSE);
2058     hTRDeffPos->Divide(hTRDrefPos, hTPCrefPos);
2059     hTRDeffPos->SetMaximum(1.0);
2060     hTRDeffPos->Draw("samecolz");
2061   }
2062   lat->DrawLatex(-0.9, 3.6, "TPC-TRD matching for positive tracks");
2063   DrawTRDGrid();
2064
2065   //----------------------------------------------
2066   // eta-phi efficiency for negative TRD tracks
2067   pad = ((TVirtualPad*)l->At(3)); pad->cd();
2068   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
2069   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
2070   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2071   rangeEtaPhi->Draw();
2072
2073   if(useIsolatedBC) {
2074     if(nIsolatedBC>0) {
2075       TH2D* hTRDeffNeg_IsolatedBC = (TH2D*)hTRDrefNeg_IsolatedBC->Clone("hTRDeffNeg_IsolatedBC");
2076       hTRDeffNeg_IsolatedBC->Reset();
2077       hTRDeffNeg_IsolatedBC->SetStats(kFALSE);
2078       hTRDeffNeg_IsolatedBC->Divide(hTRDrefNeg_IsolatedBC, hTPCrefNeg_IsolatedBC);
2079       hTRDeffNeg_IsolatedBC->SetMaximum(1.0);
2080       hTRDeffNeg_IsolatedBC->Draw("samecolz");
2081     }
2082   }
2083   else {
2084     TH2D* hTRDeffNeg = (TH2D*)hTRDrefNeg->Clone("hTRDeffNeg");
2085     hTRDeffNeg->Reset();
2086     hTRDeffNeg->SetStats(kFALSE);
2087     hTRDeffNeg->Divide(hTRDrefNeg, hTPCrefNeg);
2088     hTRDeffNeg->SetMaximum(1.0);
2089     hTRDeffNeg->Draw("samecolz");
2090   }
2091   lat->DrawLatex(-0.9, 3.6, "TPC-TRD matching for negative tracks");
2092   DrawTRDGrid();
2093
2094   //----------------------------------------------
2095   // eta-phi TRD-TOF matching efficiency for positive tracks
2096   pad = ((TVirtualPad*)l->At(1)); pad->cd();
2097   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
2098   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
2099   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2100   rangeEtaPhi->Draw();
2101
2102   if(useIsolatedBC) {
2103     if(nIsolatedBC>0) {
2104       TH2D* hTOFeffPos_IsolatedBC = (TH2D*)hTOFrefPos_IsolatedBC->Clone("hTOFeffPos_IsolatedBC");
2105       hTOFeffPos_IsolatedBC->Reset();
2106       hTOFeffPos_IsolatedBC->SetStats(kFALSE);
2107       hTOFeffPos_IsolatedBC->Divide(hTOFrefPos_IsolatedBC, hTRDrefPos_IsolatedBC);
2108       hTOFeffPos_IsolatedBC->SetMaximum(1.0);
2109       hTOFeffPos_IsolatedBC->Draw("samecolz");
2110     }
2111   }
2112   else {
2113     TH2D* hTOFeffPos = (TH2D*)hTOFrefPos->Clone("hTOFeffPos");
2114     hTOFeffPos->Reset();
2115     hTOFeffPos->SetStats(kFALSE);
2116     hTOFeffPos->Divide(hTOFrefPos, hTRDrefPos);
2117     hTOFeffPos->SetMaximum(1.0);
2118     hTOFeffPos->Draw("samecolz");
2119   }
2120   lat->DrawLatex(-0.9, 3.6, "TRD-TOF matching for positive tracks");
2121   DrawTRDGrid();
2122
2123   //----------------------------------------------
2124   // eta-phi TRD-TOF matching efficiency for negative tracks
2125   pad = ((TVirtualPad*)l->At(4)); pad->cd();
2126   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
2127   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
2128   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2129   rangeEtaPhi->Draw();
2130
2131   if(useIsolatedBC) {
2132     if(nIsolatedBC) {
2133       TH2D* hTOFeffNeg_IsolatedBC = (TH2D*)hTOFrefNeg_IsolatedBC->Clone("hTOFeffNeg_IsolatedBC");
2134       hTOFeffNeg_IsolatedBC->Reset();
2135       hTOFeffNeg_IsolatedBC->SetStats(kFALSE);
2136       hTOFeffNeg_IsolatedBC->Divide(hTOFrefNeg_IsolatedBC, hTRDrefNeg_IsolatedBC);
2137       hTOFeffNeg_IsolatedBC->SetMaximum(1.0);
2138       hTOFeffNeg_IsolatedBC->Draw("samecolz");
2139     }
2140   }
2141   else {
2142     TH2D* hTOFeffNeg = (TH2D*)hTOFrefNeg->Clone("hTOFeffNeg");
2143     hTOFeffNeg->Reset();
2144     hTOFeffNeg->SetStats(kFALSE);
2145     hTOFeffNeg->Divide(hTOFrefNeg, hTRDrefNeg);
2146     hTOFeffNeg->SetMaximum(1.0);
2147     hTOFeffNeg->Draw("samecolz");
2148   }
2149   lat->DrawLatex(-0.9, 3.6, "TRD-TOF matching for negative tracks");
2150   DrawTRDGrid();
2151
2152   if(hTRDrefPos) delete hTRDrefPos; if(hTPCrefPos) delete hTPCrefPos; if(hTOFrefPos) delete hTOFrefPos;
2153   if(hTRDrefNeg) delete hTRDrefNeg; if(hTPCrefNeg) delete hTPCrefNeg; if(hTOFrefNeg) delete hTOFrefNeg;
2154   if(hTRDrefPos_IsolatedBC) delete hTRDrefPos_IsolatedBC; if(hTPCrefPos_IsolatedBC) delete hTPCrefPos_IsolatedBC; if(hTOFrefPos_IsolatedBC) delete hTOFrefPos_IsolatedBC;
2155   if(hTRDrefNeg_IsolatedBC) delete hTRDrefNeg_IsolatedBC; if(hTPCrefNeg_IsolatedBC) delete hTPCrefNeg_IsolatedBC; if(hTOFrefNeg_IsolatedBC) delete hTOFrefNeg_IsolatedBC;
2156
2157   fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
2158   fCfContainer->SetRangeUser(kTrackCharge, +1.0, +1.0);    // positive charges
2159   TH3D* h3TPCrefPos = (TH3D*)fCfContainer->Project(0, kTrackEta, kTrackPhi, kTrackPt);
2160   TH3D* h3TRDrefPosAll = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2161   TH3D* h3TOFrefPosAll = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2162   fCfContainer->SetRangeUser(kTrackTrdTracklets, 4.0, 4.0);
2163   TH3D* h3TRDrefPosTrk4 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2164   TH3D* h3TOFrefPosTrk4 = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2165   fCfContainer->SetRangeUser(kTrackTrdTracklets, 5.0, 5.0);
2166   TH3D* h3TRDrefPosTrk5 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2167   TH3D* h3TOFrefPosTrk5 = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2168   fCfContainer->SetRangeUser(kTrackTrdTracklets, 6.0, 6.0);
2169   TH3D* h3TRDrefPosTrk6 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2170   TH3D* h3TOFrefPosTrk6 = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2171
2172   fCfContainer->SetRangeUser(kTrackCharge, -1.0, -1.0);   // negative charges
2173   fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
2174   TH3D* h3TPCrefNeg = (TH3D*)fCfContainer->Project(0, kTrackEta, kTrackPhi, kTrackPt);
2175   TH3D* h3TRDrefNegAll = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2176   TH3D* h3TOFrefNegAll = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2177   fCfContainer->SetRangeUser(kTrackTrdTracklets, 4.0, 4.0);
2178   TH3D* h3TRDrefNegTrk4 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2179   TH3D* h3TOFrefNegTrk4 = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2180   fCfContainer->SetRangeUser(kTrackTrdTracklets, 5.0, 5.0);
2181   TH3D* h3TRDrefNegTrk5 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2182   TH3D* h3TOFrefNegTrk5 = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2183   fCfContainer->SetRangeUser(kTrackTrdTracklets, 6.0, 6.0);
2184   TH3D* h3TRDrefNegTrk6 = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2185   TH3D* h3TOFrefNegTrk6 = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2186   fCfContainer->SetRangeUser(kTrackCharge, -1.0, +1.0);
2187   fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
2188
2189   TH1F* hTRDEffPtPosAll = EfficiencyTRD(h3TPCrefPos, h3TRDrefPosAll, kTRUE);
2190   TH1F* hTRDEffPtNegAll = EfficiencyTRD(h3TPCrefNeg, h3TRDrefNegAll, kTRUE);
2191   TH1F* hTRDEffPtPosTrk4 = EfficiencyTRD(h3TPCrefPos, h3TRDrefPosTrk4, kTRUE);
2192   TH1F* hTRDEffPtNegTrk4 = EfficiencyTRD(h3TPCrefNeg, h3TRDrefNegTrk4, kTRUE);
2193   TH1F* hTRDEffPtPosTrk5 = EfficiencyTRD(h3TPCrefPos, h3TRDrefPosTrk5, kTRUE);
2194   TH1F* hTRDEffPtNegTrk5 = EfficiencyTRD(h3TPCrefNeg, h3TRDrefNegTrk5, kTRUE);
2195   TH1F* hTRDEffPtPosTrk6 = EfficiencyTRD(h3TPCrefPos, h3TRDrefPosTrk6, kTRUE);
2196   TH1F* hTRDEffPtNegTrk6 = EfficiencyTRD(h3TPCrefNeg, h3TRDrefNegTrk6, kTRUE);
2197
2198   TH1F* hTOFEffPtPosAll = EfficiencyTRD(h3TRDrefPosAll, h3TOFrefPosAll, kFALSE);
2199   TH1F* hTOFEffPtNegAll = EfficiencyTRD(h3TRDrefNegAll, h3TOFrefNegAll, kFALSE);
2200   TH1F* hTOFEffPtPosTrk4 = EfficiencyTRD(h3TRDrefPosTrk4, h3TOFrefPosTrk4, kFALSE);
2201   TH1F* hTOFEffPtNegTrk4 = EfficiencyTRD(h3TRDrefNegTrk4, h3TOFrefNegTrk4, kFALSE);
2202   TH1F* hTOFEffPtPosTrk5 = EfficiencyTRD(h3TRDrefPosTrk5, h3TOFrefPosTrk5, kFALSE);
2203   TH1F* hTOFEffPtNegTrk5 = EfficiencyTRD(h3TRDrefNegTrk5, h3TOFrefNegTrk5, kFALSE);
2204   TH1F* hTOFEffPtPosTrk6 = EfficiencyTRD(h3TRDrefPosTrk6, h3TOFrefPosTrk6, kFALSE);
2205   TH1F* hTOFEffPtNegTrk6 = EfficiencyTRD(h3TRDrefNegTrk6, h3TOFrefNegTrk6, kFALSE);
2206
2207   delete h3TPCrefPos; delete h3TPCrefNeg;
2208   delete h3TRDrefPosAll; delete h3TRDrefPosTrk4; delete h3TRDrefPosTrk5; delete h3TRDrefPosTrk6;
2209   delete h3TRDrefNegAll; delete h3TRDrefNegTrk4; delete h3TRDrefNegTrk5; delete h3TRDrefNegTrk6;
2210   delete h3TOFrefPosAll; delete h3TOFrefPosTrk4; delete h3TOFrefPosTrk5; delete h3TOFrefPosTrk6;
2211   delete h3TOFrefNegAll; delete h3TOFrefNegTrk4; delete h3TOFrefNegTrk5; delete h3TOFrefNegTrk6;
2212
2213   TF1* funcConst = new TF1("constFunc", "[0]", 1.0, 3.0);
2214   if(trendValues && hTRDEffPtPosAll){
2215     if(hTRDEffPtPosAll->Integral()>0.1) {
2216       hTRDEffPtPosAll->Fit(funcConst, "Q0ME", "goff", 1.0, 3.0);
2217       trendValues[0] = funcConst->GetParameter(0);
2218       trendValues[1] = funcConst->GetParError(0);
2219     }
2220   }
2221   if(trendValues && hTRDEffPtNegAll){
2222     if(hTRDEffPtNegAll->Integral()>0.1) {
2223       hTRDEffPtNegAll->Fit(funcConst, "Q0ME", "goff", 1.0, 3.0);
2224       trendValues[2] = funcConst->GetParameter(0);
2225       trendValues[3] = funcConst->GetParError(0);
2226     }
2227   }
2228   if(trendValues && hTOFEffPtPosAll){
2229     if(hTOFEffPtPosAll->Integral()>0.1) {
2230       hTOFEffPtPosAll->Fit(funcConst, "Q0ME", "goff", 1.0, 3.0);
2231       trendValues[4] = funcConst->GetParameter(0);
2232       trendValues[5] = funcConst->GetParError(0);
2233     }
2234   }
2235   if(trendValues && hTOFEffPtNegAll){
2236     if(hTOFEffPtNegAll->Integral()>0.1) {
2237       hTOFEffPtNegAll->Fit(funcConst, "Q0ME", "goff", 1.0, 3.0);
2238       trendValues[6] = funcConst->GetParameter(0);
2239       trendValues[7] = funcConst->GetParError(0);
2240     }
2241   }
2242   // get matching efficiencies for isolated bunches
2243   TH3D* h3TPCrefPos_IsolatedBC=0x0; TH3D* h3TPCrefNeg_IsolatedBC=0x0;
2244   TH3D* h3TRDrefPosAll_IsolatedBC=0x0; TH3D* h3TRDrefPosTrk4_IsolatedBC=0x0; TH3D* h3TRDrefPosTrk5_IsolatedBC=0x0; TH3D* h3TRDrefPosTrk6_IsolatedBC=0x0;
2245   TH3D* h3TRDrefNegAll_IsolatedBC=0x0; TH3D* h3TRDrefNegTrk4_IsolatedBC=0x0; TH3D* h3TRDrefNegTrk5_IsolatedBC=0x0; TH3D* h3TRDrefNegTrk6_IsolatedBC=0x0;
2246   TH3D* h3TOFrefPosAll_IsolatedBC=0x0; TH3D* h3TOFrefPosTrk4_IsolatedBC=0x0; TH3D* h3TOFrefPosTrk5_IsolatedBC=0x0; TH3D* h3TOFrefPosTrk6_IsolatedBC=0x0;
2247   TH3D* h3TOFrefNegAll_IsolatedBC=0x0; TH3D* h3TOFrefNegTrk4_IsolatedBC=0x0; TH3D* h3TOFrefNegTrk5_IsolatedBC=0x0; TH3D* h3TOFrefNegTrk6_IsolatedBC=0x0;
2248   for(Int_t ibc=0; ibc<3500; ++ibc) {
2249     if(!isIsolated[ibc]) continue;
2250     fCfContainer->SetRangeUser(kEventBC, Double_t(ibc), Double_t(ibc));
2251
2252     TH3D* temp;
2253     fCfContainer->SetRangeUser(kTrackCharge, +1.0, +1.0);    // positive charges
2254     fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
2255     temp = (TH3D*)fCfContainer->Project(0, kTrackEta, kTrackPhi, kTrackPt);
2256     if(!h3TPCrefPos_IsolatedBC) h3TPCrefPos_IsolatedBC = temp;
2257     else {h3TPCrefPos_IsolatedBC->Add(temp); delete temp;}
2258     temp = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2259     if(!h3TRDrefPosAll_IsolatedBC) h3TRDrefPosAll_IsolatedBC = temp;
2260     else {h3TRDrefPosAll_IsolatedBC->Add(temp); delete temp;}
2261     temp = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2262     if(!h3TOFrefPosAll_IsolatedBC) h3TOFrefPosAll_IsolatedBC = temp;
2263     else {h3TOFrefPosAll_IsolatedBC->Add(temp); delete temp;}
2264     fCfContainer->SetRangeUser(kTrackTrdTracklets, 4.0, 4.0);
2265     temp = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2266     if(!h3TRDrefPosTrk4_IsolatedBC) h3TRDrefPosTrk4_IsolatedBC = temp;
2267     else {h3TRDrefPosTrk4_IsolatedBC->Add(temp); delete temp;}
2268     temp = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2269     if(!h3TOFrefPosTrk4_IsolatedBC) h3TOFrefPosTrk4_IsolatedBC = temp;
2270     else {h3TOFrefPosTrk4_IsolatedBC->Add(temp); delete temp;}
2271     fCfContainer->SetRangeUser(kTrackTrdTracklets, 5.0, 5.0);
2272     temp = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2273     if(!h3TRDrefPosTrk5_IsolatedBC) h3TRDrefPosTrk5_IsolatedBC = temp;
2274     else {h3TRDrefPosTrk5_IsolatedBC->Add(temp); delete temp;}
2275     temp = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2276     if(!h3TOFrefPosTrk5_IsolatedBC) h3TOFrefPosTrk5_IsolatedBC = temp;
2277     else {h3TOFrefPosTrk5_IsolatedBC->Add(temp); delete temp;}
2278     fCfContainer->SetRangeUser(kTrackTrdTracklets, 6.0, 6.0);
2279     temp = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2280     if(!h3TRDrefPosTrk6_IsolatedBC) h3TRDrefPosTrk6_IsolatedBC = temp;
2281     else {h3TRDrefPosTrk6_IsolatedBC->Add(temp); delete temp;}
2282     temp = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2283     if(!h3TOFrefPosTrk6_IsolatedBC) h3TOFrefPosTrk6_IsolatedBC = temp;
2284     else {h3TOFrefPosTrk6_IsolatedBC->Add(temp); delete temp;}
2285
2286     fCfContainer->SetRangeUser(kTrackCharge, -1.0, -1.0);   // negative charges
2287     fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
2288     temp = (TH3D*)fCfContainer->Project(0, kTrackEta, kTrackPhi, kTrackPt);
2289     if(!h3TPCrefNeg_IsolatedBC) h3TPCrefNeg_IsolatedBC = temp;
2290     else {h3TPCrefNeg_IsolatedBC->Add(temp); delete temp;}
2291     temp = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2292     if(!h3TRDrefNegAll_IsolatedBC) h3TRDrefNegAll_IsolatedBC = temp;
2293     else {h3TRDrefNegAll_IsolatedBC->Add(temp); delete temp;}
2294     temp = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2295     if(!h3TOFrefNegAll_IsolatedBC) h3TOFrefNegAll_IsolatedBC = temp;
2296     else {h3TOFrefNegAll_IsolatedBC->Add(temp); delete temp;}
2297     fCfContainer->SetRangeUser(kTrackTrdTracklets, 4.0, 4.0);
2298     temp = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2299     if(!h3TRDrefNegTrk4_IsolatedBC) h3TRDrefNegTrk4_IsolatedBC = temp;
2300     else {h3TRDrefNegTrk4_IsolatedBC->Add(temp); delete temp;}
2301     temp = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2302     if(!h3TOFrefNegTrk4_IsolatedBC) h3TOFrefNegTrk4_IsolatedBC = temp;
2303     else {h3TOFrefNegTrk4_IsolatedBC->Add(temp); delete temp;}
2304     fCfContainer->SetRangeUser(kTrackTrdTracklets, 5.0, 5.0);
2305     temp = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2306     if(!h3TRDrefNegTrk5_IsolatedBC) h3TRDrefNegTrk5_IsolatedBC = temp;
2307     else {h3TRDrefNegTrk5_IsolatedBC->Add(temp); delete temp;}
2308     temp = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2309     if(!h3TOFrefNegTrk5_IsolatedBC) h3TOFrefNegTrk5_IsolatedBC = temp;
2310     else {h3TOFrefNegTrk5_IsolatedBC->Add(temp); delete temp;}
2311     fCfContainer->SetRangeUser(kTrackTrdTracklets, 6.0, 6.0);
2312     temp = (TH3D*)fCfContainer->Project(1, kTrackEta, kTrackPhi, kTrackPt);
2313     if(!h3TRDrefNegTrk6_IsolatedBC) h3TRDrefNegTrk6_IsolatedBC = temp;
2314     else {h3TRDrefNegTrk6_IsolatedBC->Add(temp); delete temp;}
2315     temp = (TH3D*)fCfContainer->Project(2, kTrackEta, kTrackPhi, kTrackPt);
2316     if(!h3TOFrefNegTrk6_IsolatedBC) h3TOFrefNegTrk6_IsolatedBC = temp;
2317     else {h3TOFrefNegTrk6_IsolatedBC->Add(temp); delete temp;}
2318   }
2319   fCfContainer->SetRangeUser(kEventBC, 0.0, 3500.0);
2320   fCfContainer->SetRangeUser(kTrackCharge, -1.0, +1.0);
2321   fCfContainer->SetRangeUser(kTrackTrdTracklets, 0.0, 6.0);
2322
2323   TH1F* hTRDEffPtPosAll_IsolatedBC  = EfficiencyTRD(h3TPCrefPos_IsolatedBC, h3TRDrefPosAll_IsolatedBC, kTRUE);
2324   TH1F* hTRDEffPtNegAll_IsolatedBC  = EfficiencyTRD(h3TPCrefNeg_IsolatedBC, h3TRDrefNegAll_IsolatedBC, kTRUE);
2325   TH1F* hTRDEffPtPosTrk4_IsolatedBC = EfficiencyTRD(h3TPCrefPos_IsolatedBC, h3TRDrefPosTrk4_IsolatedBC, kTRUE);
2326   TH1F* hTRDEffPtNegTrk4_IsolatedBC = EfficiencyTRD(h3TPCrefNeg_IsolatedBC, h3TRDrefNegTrk4_IsolatedBC, kTRUE);
2327   TH1F* hTRDEffPtPosTrk5_IsolatedBC = EfficiencyTRD(h3TPCrefPos_IsolatedBC, h3TRDrefPosTrk5_IsolatedBC, kTRUE);
2328   TH1F* hTRDEffPtNegTrk5_IsolatedBC = EfficiencyTRD(h3TPCrefNeg_IsolatedBC, h3TRDrefNegTrk5_IsolatedBC, kTRUE);
2329   TH1F* hTRDEffPtPosTrk6_IsolatedBC = EfficiencyTRD(h3TPCrefPos_IsolatedBC, h3TRDrefPosTrk6_IsolatedBC, kTRUE);
2330   TH1F* hTRDEffPtNegTrk6_IsolatedBC = EfficiencyTRD(h3TPCrefNeg_IsolatedBC, h3TRDrefNegTrk6_IsolatedBC, kTRUE);
2331
2332   TH1F* hTOFEffPtPosAll_IsolatedBC  = EfficiencyTRD(h3TRDrefPosAll_IsolatedBC, h3TOFrefPosAll_IsolatedBC, kFALSE);
2333   TH1F* hTOFEffPtNegAll_IsolatedBC  = EfficiencyTRD(h3TRDrefNegAll_IsolatedBC, h3TOFrefNegAll_IsolatedBC, kFALSE);
2334   TH1F* hTOFEffPtPosTrk4_IsolatedBC = EfficiencyTRD(h3TRDrefPosTrk4_IsolatedBC, h3TOFrefPosTrk4_IsolatedBC, kFALSE);
2335   TH1F* hTOFEffPtNegTrk4_IsolatedBC = EfficiencyTRD(h3TRDrefNegTrk4_IsolatedBC, h3TOFrefNegTrk4_IsolatedBC, kFALSE);
2336   TH1F* hTOFEffPtPosTrk5_IsolatedBC = EfficiencyTRD(h3TRDrefPosTrk5_IsolatedBC, h3TOFrefPosTrk5_IsolatedBC, kFALSE);
2337   TH1F* hTOFEffPtNegTrk5_IsolatedBC = EfficiencyTRD(h3TRDrefNegTrk5_IsolatedBC, h3TOFrefNegTrk5_IsolatedBC, kFALSE);
2338   TH1F* hTOFEffPtPosTrk6_IsolatedBC = EfficiencyTRD(h3TRDrefPosTrk6_IsolatedBC, h3TOFrefPosTrk6_IsolatedBC, kFALSE);
2339   TH1F* hTOFEffPtNegTrk6_IsolatedBC = EfficiencyTRD(h3TRDrefNegTrk6_IsolatedBC, h3TOFrefNegTrk6_IsolatedBC, kFALSE);
2340
2341   if(h3TPCrefPos_IsolatedBC) delete h3TPCrefPos_IsolatedBC; if(h3TPCrefNeg_IsolatedBC) delete h3TPCrefNeg_IsolatedBC;
2342   if(h3TRDrefPosAll_IsolatedBC) delete h3TRDrefPosAll_IsolatedBC; if(h3TRDrefPosTrk4_IsolatedBC) delete h3TRDrefPosTrk4_IsolatedBC;
2343   if(h3TRDrefPosTrk5_IsolatedBC) delete h3TRDrefPosTrk5_IsolatedBC; if(h3TRDrefPosTrk6_IsolatedBC) delete h3TRDrefPosTrk6_IsolatedBC;
2344   if(h3TRDrefNegAll_IsolatedBC) delete h3TRDrefNegAll_IsolatedBC; if(h3TRDrefNegTrk4_IsolatedBC) delete h3TRDrefNegTrk4_IsolatedBC;
2345   if(h3TRDrefNegTrk5_IsolatedBC) delete h3TRDrefNegTrk5_IsolatedBC; if(h3TRDrefNegTrk6_IsolatedBC) delete h3TRDrefNegTrk6_IsolatedBC;
2346   if(h3TOFrefPosAll_IsolatedBC) delete h3TOFrefPosAll_IsolatedBC; if(h3TOFrefPosTrk4_IsolatedBC) delete h3TOFrefPosTrk4_IsolatedBC;
2347   if(h3TOFrefPosTrk5_IsolatedBC) delete h3TOFrefPosTrk5_IsolatedBC; if(h3TOFrefPosTrk6_IsolatedBC) delete h3TOFrefPosTrk6_IsolatedBC;
2348   if(h3TOFrefNegAll_IsolatedBC) delete h3TOFrefNegAll_IsolatedBC; if(h3TOFrefNegTrk4_IsolatedBC) delete h3TOFrefNegTrk4_IsolatedBC;
2349   if(h3TOFrefNegTrk5_IsolatedBC) delete h3TOFrefNegTrk5_IsolatedBC; if(h3TOFrefNegTrk6_IsolatedBC) delete h3TOFrefNegTrk6_IsolatedBC;
2350
2351   if(trendValues &&  hTRDEffPtPosAll_IsolatedBC && hTRDEffPtPosAll_IsolatedBC->Integral()>0.1) {
2352     hTRDEffPtPosAll_IsolatedBC->Fit(funcConst, "Q0ME", "goff", 1.0, 3.0);
2353     trendValues[18] = funcConst->GetParameter(0);
2354     trendValues[19] = funcConst->GetParError(0);
2355   }
2356   if(trendValues &&  hTRDEffPtNegAll_IsolatedBC && hTRDEffPtNegAll_IsolatedBC->Integral()>0.1) {
2357     hTRDEffPtNegAll_IsolatedBC->Fit(funcConst, "Q0ME", "goff", 1.0, 3.0);
2358     trendValues[20] = funcConst->GetParameter(0);
2359     trendValues[21] = funcConst->GetParError(0);
2360   }
2361   if(trendValues &&  hTOFEffPtPosAll_IsolatedBC && hTOFEffPtPosAll_IsolatedBC->Integral()>0.1) {
2362     hTOFEffPtPosAll_IsolatedBC->Fit(funcConst, "Q0ME", "goff", 1.0, 3.0);
2363     trendValues[22] = funcConst->GetParameter(0);
2364     trendValues[23] = funcConst->GetParError(0);
2365   }
2366   if(trendValues &&  hTOFEffPtNegAll_IsolatedBC && hTOFEffPtNegAll_IsolatedBC->Integral()>0.1) {
2367     hTOFEffPtNegAll_IsolatedBC->Fit(funcConst, "Q0ME", "goff", 1.0, 3.0);
2368     trendValues[24] = funcConst->GetParameter(0);
2369     trendValues[25] = funcConst->GetParError(0);
2370   }
2371
2372   //---------------------------------------------------------
2373   // TPC-TRD matching efficiency vs pt
2374   pad = ((TVirtualPad*)l->At(6)); pad->cd();
2375   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.02);
2376   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
2377   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2378
2379   if(gROOT->FindObject("rangeEffPt2")) delete gROOT->FindObject("rangeEffPt2");
2380   TH2F* rangeEffPt=new TH2F("rangeEffPt2", "",10,0.,10.,10,0.,1.4);
2381   rangeEffPt->SetStats(kFALSE);
2382   SetStyle(rangeEffPt->GetXaxis(), "p_{T} [GeV/c]", 0.07, 0.8, kTRUE, 0.05);
2383   SetStyle(rangeEffPt->GetYaxis(), "efficiency", 0.07, 0.8, kTRUE, 0.05);
2384   rangeEffPt->Draw();
2385   lat->DrawLatex(0.2, 1.42, "TPC-TRD matching efficiency");
2386   //++++++++++++++++++
2387   TLine line;
2388   line.SetLineStyle(2);
2389   line.SetLineWidth(2);
2390   line.DrawLine(rangeEffPt->GetXaxis()->GetXmin(), 0.7, rangeEffPt->GetXaxis()->GetXmax(), 0.7);
2391   line.DrawLine(rangeEffPt->GetXaxis()->GetXmin(), 0.9, rangeEffPt->GetXaxis()->GetXmax(), 0.9);
2392   TLegend* leg=new TLegend(0.2, 0.7, 0.7, 0.89);
2393   leg->SetNColumns(2);
2394   leg->SetMargin(0.15);
2395   leg->SetBorderSize(0);
2396   leg->SetFillColor(0);
2397
2398   if(!useIsolatedBC) {
2399     SetStyle(hTRDEffPtPosAll, 1, kRed, 1, 24, kRed, 1);
2400     SetStyle(hTRDEffPtNegAll, 1, kBlue, 1, 24, kBlue, 1);
2401     SetStyle(hTRDEffPtPosTrk4, 1, kRed, 1, 25, kRed, 1);
2402     SetStyle(hTRDEffPtNegTrk4, 1, kBlue, 1, 25, kBlue, 1);
2403     SetStyle(hTRDEffPtPosTrk5, 1, kRed, 1, 26, kRed, 1);
2404     SetStyle(hTRDEffPtNegTrk5, 1, kBlue, 1, 26, kBlue, 1);
2405     SetStyle(hTRDEffPtPosTrk6, 1, kRed, 1, 27, kRed, 1);
2406     SetStyle(hTRDEffPtNegTrk6, 1, kBlue, 1, 27, kBlue, 1);
2407     leg->SetHeader("All bunch crossings");
2408     hTRDEffPtPosAll->Draw("same"); leg->AddEntry(hTRDEffPtPosAll, "pos. (#geq 1 tracklet)", "p");
2409     hTRDEffPtNegAll->Draw("same"); leg->AddEntry(hTRDEffPtNegAll, "neg. (#geq 1 tracklet)", "p");
2410     hTRDEffPtPosTrk4->Draw("same"); leg->AddEntry(hTRDEffPtPosTrk4, "pos. (4 tracklets)", "p");
2411     hTRDEffPtNegTrk4->Draw("same"); leg->AddEntry(hTRDEffPtNegTrk4, "neg. (4 tracklets)", "p");
2412     hTRDEffPtPosTrk5->Draw("same"); leg->AddEntry(hTRDEffPtPosTrk5, "pos. (5 tracklets)", "p");
2413     hTRDEffPtNegTrk5->Draw("same"); leg->AddEntry(hTRDEffPtNegTrk5, "neg. (5 tracklets)", "p");
2414     hTRDEffPtPosTrk6->Draw("same"); leg->AddEntry(hTRDEffPtPosTrk6, "pos. (6 tracklets)", "p");
2415     hTRDEffPtNegTrk6->Draw("same"); leg->AddEntry(hTRDEffPtNegTrk6, "neg. (6 tracklets)", "p");
2416   }
2417   else {
2418     if(nIsolatedBC>0) {
2419       SetStyle(hTRDEffPtPosAll_IsolatedBC, 1, kRed, 1, 24, kRed, 1);
2420       SetStyle(hTRDEffPtNegAll_IsolatedBC, 1, kBlue, 1, 24, kBlue, 1);
2421       SetStyle(hTRDEffPtPosTrk4_IsolatedBC, 1, kRed, 1, 25, kRed, 1);
2422       SetStyle(hTRDEffPtNegTrk4_IsolatedBC, 1, kBlue, 1, 25, kBlue, 1);
2423       SetStyle(hTRDEffPtPosTrk5_IsolatedBC, 1, kRed, 1, 26, kRed, 1);
2424       SetStyle(hTRDEffPtNegTrk5_IsolatedBC, 1, kBlue, 1, 26, kBlue, 1);
2425       SetStyle(hTRDEffPtPosTrk6_IsolatedBC, 1, kRed, 1, 27, kRed, 1);
2426       SetStyle(hTRDEffPtNegTrk6_IsolatedBC, 1, kBlue, 1, 27, kBlue, 1);
2427       leg->SetHeader("Isolated bunch crossings");
2428       hTRDEffPtPosAll_IsolatedBC->Draw("same"); leg->AddEntry(hTRDEffPtPosAll_IsolatedBC, "pos. (#geq 1 tracklet)", "p");
2429       hTRDEffPtNegAll_IsolatedBC->Draw("same"); leg->AddEntry(hTRDEffPtNegAll_IsolatedBC, "neg. (#geq 1 tracklet)", "p");
2430       hTRDEffPtPosTrk4_IsolatedBC->Draw("same"); leg->AddEntry(hTRDEffPtPosTrk4_IsolatedBC, "pos. (4 tracklets)", "p");
2431       hTRDEffPtNegTrk4_IsolatedBC->Draw("same"); leg->AddEntry(hTRDEffPtNegTrk4_IsolatedBC, "neg. (4 tracklets)", "p");
2432       hTRDEffPtPosTrk5_IsolatedBC->Draw("same"); leg->AddEntry(hTRDEffPtPosTrk5_IsolatedBC, "pos. (5 tracklets)", "p");
2433       hTRDEffPtNegTrk5_IsolatedBC->Draw("same"); leg->AddEntry(hTRDEffPtNegTrk5_IsolatedBC, "neg. (5 tracklets)", "p");
2434       hTRDEffPtPosTrk6_IsolatedBC->Draw("same"); leg->AddEntry(hTRDEffPtPosTrk6_IsolatedBC, "pos. (6 tracklets)", "p");
2435       hTRDEffPtNegTrk6_IsolatedBC->Draw("same"); leg->AddEntry(hTRDEffPtNegTrk6_IsolatedBC, "neg. (6 tracklets)", "p");
2436     }
2437   }
2438
2439   leg->Draw();
2440
2441   //---------------------------------------------------------
2442   // TRD-TOF matching efficiency vs pt
2443   pad = ((TVirtualPad*)l->At(7)); pad->cd();
2444   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.02);
2445   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
2446   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2447
2448   rangeEffPt->Draw();
2449   lat->DrawLatex(0.2, 1.42, "TRD-TOF matching efficiency");
2450
2451   if(!useIsolatedBC) {
2452     SetStyle(hTOFEffPtPosAll, 1, kRed, 1, 24, kRed, 1);
2453     SetStyle(hTOFEffPtPosTrk4, 1, kRed, 1, 25, kRed, 1);
2454     SetStyle(hTOFEffPtPosTrk5, 1, kRed, 1, 26, kRed, 1);
2455     SetStyle(hTOFEffPtPosTrk6, 1, kRed, 1, 27, kRed, 1);
2456     SetStyle(hTOFEffPtNegAll, 1, kBlue, 1, 24, kBlue, 1);
2457     SetStyle(hTOFEffPtNegTrk4, 1, kBlue, 1, 25, kBlue, 1);
2458     SetStyle(hTOFEffPtNegTrk5, 1, kBlue, 1, 26, kBlue, 1);
2459     SetStyle(hTOFEffPtNegTrk6, 1, kBlue, 1, 27, kBlue, 1);
2460     hTOFEffPtPosAll->Draw("same");
2461     hTOFEffPtPosTrk4->Draw("same");
2462     hTOFEffPtPosTrk5->Draw("same");
2463     hTOFEffPtPosTrk6->Draw("same");
2464     hTOFEffPtNegAll->Draw("same");
2465     hTOFEffPtNegTrk4->Draw("same");
2466     hTOFEffPtNegTrk5->Draw("same");
2467     hTOFEffPtNegTrk6->Draw("same");
2468   }
2469   else {
2470     if(nIsolatedBC>0) {
2471       SetStyle(hTOFEffPtPosAll_IsolatedBC, 1, kRed, 1, 24, kRed, 1);
2472       SetStyle(hTOFEffPtPosTrk4_IsolatedBC, 1, kRed, 1, 25, kRed, 1);
2473       SetStyle(hTOFEffPtPosTrk5_IsolatedBC, 1, kRed, 1, 26, kRed, 1);
2474       SetStyle(hTOFEffPtPosTrk6_IsolatedBC, 1, kRed, 1, 27, kRed, 1);
2475       SetStyle(hTOFEffPtNegAll_IsolatedBC, 1, kBlue, 1, 24, kBlue, 1);
2476       SetStyle(hTOFEffPtNegTrk4_IsolatedBC, 1, kBlue, 1, 25, kBlue, 1);
2477       SetStyle(hTOFEffPtNegTrk5_IsolatedBC, 1, kBlue, 1, 26, kBlue, 1);
2478       SetStyle(hTOFEffPtNegTrk6_IsolatedBC, 1, kBlue, 1, 27, kBlue, 1);
2479       hTOFEffPtPosAll_IsolatedBC->Draw("same");
2480       hTOFEffPtPosTrk4_IsolatedBC->Draw("same");
2481       hTOFEffPtPosTrk5_IsolatedBC->Draw("same");
2482       hTOFEffPtPosTrk6_IsolatedBC->Draw("same");
2483       hTOFEffPtNegAll_IsolatedBC->Draw("same");
2484       hTOFEffPtNegTrk4_IsolatedBC->Draw("same");
2485       hTOFEffPtNegTrk5_IsolatedBC->Draw("same");
2486       hTOFEffPtNegTrk6_IsolatedBC->Draw("same");
2487     }
2488   }
2489
2490   //-----------------------------------------------------
2491   // <ntracklets> vs (phi,eta)
2492   pad = ((TVirtualPad*)l->At(2)); pad->cd();
2493   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
2494   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
2495   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2496
2497   rangeEtaPhi->Draw();
2498   lat->DrawLatex(-0.9, 3.6, "TRD <N_{tracklets}>");
2499
2500   TH3D* hNtracklets = (TH3D*)fCfContainer->Project(1, kTrackPhi, kTrackEta, kTrackTrdTracklets);
2501   TH3D* hNtracklets_IsolatedBC=0x0;
2502   for(Int_t ibc=0; ibc<3500; ++ibc) {
2503     if(!isIsolated[ibc]) continue;
2504     fCfContainer->SetRangeUser(kEventBC, Double_t(ibc), Double_t(ibc));
2505     TH3D* temp = (TH3D*)fCfContainer->Project(1, kTrackPhi, kTrackEta, kTrackTrdTracklets);
2506     if(!hNtracklets_IsolatedBC) hNtracklets_IsolatedBC = temp;
2507     else {hNtracklets_IsolatedBC->Add(temp); delete temp;}
2508   }
2509   fCfContainer->SetRangeUser(kEventBC, 0.0, 3500.0);   // reset the BC range
2510
2511   TProfile2D* hNtrackletsProf=0x0;
2512   if(useIsolatedBC) {
2513     if(nIsolatedBC>0)
2514       hNtrackletsProf = hNtracklets_IsolatedBC->Project3DProfile();
2515   }
2516   else hNtrackletsProf = hNtracklets->Project3DProfile();
2517   delete hNtracklets; if(hNtracklets_IsolatedBC) delete hNtracklets_IsolatedBC;
2518   if(hNtrackletsProf) {
2519     hNtrackletsProf->SetStats(kFALSE);
2520     hNtrackletsProf->SetMinimum(0.);
2521     hNtrackletsProf->SetMaximum(6.);
2522     hNtrackletsProf->Draw("samecolz");
2523     DrawTRDGrid();
2524   }
2525
2526   // calculate the trend value for tracklets/track
2527   TH2D* hNtrackletsVsP = (TH2D*)fCfContainer->Project(1, kTrackP, kTrackTrdTracklets);
2528   if(trendValues &&  hNtrackletsVsP && hNtrackletsVsP->GetEntries()>0.1) {
2529     TProfile* hNtrackletsVsPprof = hNtrackletsVsP->ProfileX("hNtrackletsVsPprof");
2530     hNtrackletsVsPprof->Fit(funcConst, "QME0", "goff", 1.0, 3.0);
2531     trendValues[8] = funcConst->GetParameter(0);
2532     trendValues[9] = funcConst->GetParError(0);
2533     delete hNtrackletsVsP;
2534   }
2535   TH2D* hNtrackletsVsP_IsolatedBC=0x0;
2536   for(Int_t ibc=0; ibc<3500; ++ibc) {
2537     if(!isIsolated[ibc]) continue;
2538     fCfContainer->SetRangeUser(kEventBC, Double_t(ibc), Double_t(ibc));
2539     TH2D* temp = (TH2D*)fCfContainer->Project(1, kTrackP, kTrackTrdTracklets);
2540     if(!hNtrackletsVsP_IsolatedBC) hNtrackletsVsP_IsolatedBC = temp;
2541     else {hNtrackletsVsP_IsolatedBC->Add(temp); delete temp;}
2542   }
2543   fCfContainer->SetRangeUser(kEventBC, 0.0, 3500.0);   // reset the BC range
2544   if(trendValues &&  hNtrackletsVsP_IsolatedBC && hNtrackletsVsP_IsolatedBC->GetEntries()>0.1) {
2545     TProfile* hNtrackletsVsPprof_IsolatedBC = hNtrackletsVsP_IsolatedBC->ProfileX("hNtrackletsVsPprof_IsolatedBC");
2546     hNtrackletsVsPprof_IsolatedBC->Fit(funcConst, "QME0", "goff", 1.0, 3.0);
2547     trendValues[26] = funcConst->GetParameter(0);
2548     trendValues[27] = funcConst->GetParError(0);
2549     delete hNtrackletsVsPprof_IsolatedBC;
2550   }
2551
2552   //--------------------------------------------------------------
2553   // Nclusters per TRD track vs momentum
2554   pad = ((TVirtualPad*)l->At(5)); pad->cd();
2555   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.12);
2556   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
2557   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2558   pad->SetLogz();
2559
2560   if(gROOT->FindObject("rangeNclsP")) delete gROOT->FindObject("rangeNclsP");
2561   TH2F* rangeNclsP = new TH2F("rangeNclsP", "", 10, 0.0, 11.99, 10, 0.0, 199.0);
2562   SetStyle(rangeNclsP->GetXaxis(), "p [GeV/c]", 0.07, 0.8, kTRUE, 0.05);
2563   SetStyle(rangeNclsP->GetYaxis(), "#clusters", 0.07, 0.8, kTRUE, 0.05);
2564   rangeNclsP->SetStats(kFALSE);
2565   rangeNclsP->Draw();
2566   lat->DrawLatex(1.0, 205., "TRD Clusters / track");
2567
2568   TH2D* hNclsVsP = (TH2D*)fCfContainer->Project(1, kTrackP, kTrackTrdClusters);
2569   TH2D* hNclsVsP_IsolatedBC=0x0;
2570   for(Int_t ibc=0; ibc<3500; ++ibc) {
2571     if(!isIsolated[ibc]) continue;
2572     fCfContainer->SetRangeUser(kEventBC, Double_t(ibc), Double_t(ibc));
2573     TH2D* temp = (TH2D*)fCfContainer->Project(1, kTrackP, kTrackTrdClusters);
2574     if(!hNclsVsP_IsolatedBC) hNclsVsP_IsolatedBC = temp;
2575     else {hNclsVsP_IsolatedBC->Add(temp); delete temp;}
2576   }
2577   fCfContainer->SetRangeUser(kEventBC, 0.0, 3500.0);   // reset the BC range
2578
2579   if(useIsolatedBC && hNclsVsP_IsolatedBC) {
2580     hNclsVsP_IsolatedBC->SetStats(kFALSE);
2581     hNclsVsP_IsolatedBC->Draw("samecolz");
2582   }
2583   if(!useIsolatedBC && hNclsVsP) {
2584     hNclsVsP->SetStats(kFALSE);
2585     hNclsVsP->Draw("samecolz");
2586   }
2587
2588   if(trendValues && hNclsVsP && hNclsVsP->GetEntries()>10) {
2589     TProfile* hNclsVsPprof = hNclsVsP->ProfileX("hNclsVsPprof");
2590     hNclsVsPprof->Fit(funcConst, "QME0", "goff", 1.0, 3.0);
2591     trendValues[10] = funcConst->GetParameter(0);
2592     trendValues[11] = funcConst->GetParError(0);
2593   }
2594   if(trendValues && hNclsVsP_IsolatedBC && hNclsVsP_IsolatedBC->GetEntries()>10) {
2595     TProfile* hNclsVsPprof_IsolatedBC = hNclsVsP_IsolatedBC->ProfileX("hNclsVsPprof_IsolateBC");
2596     hNclsVsPprof_IsolatedBC->Fit(funcConst, "QME0", "goff", 1.0, 3.0);
2597     trendValues[28] = funcConst->GetParameter(0);
2598     trendValues[29] = funcConst->GetParError(0);
2599   }
2600
2601   //--------------------------------------------------------------
2602   // TRD-TPC and TOF-TRD matching efficiency vs bunch crossing
2603   pad = ((TVirtualPad*)l->At(8)); pad->cd();
2604   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.02);
2605   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
2606   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2607
2608   TH1D* phiProj = (TH1D*)fCfContainer->Project(1, kTrackPhi);
2609   Double_t smPhiLimits[19];
2610   Bool_t activeSM[18] = {kFALSE};
2611   for(Int_t ism=0; ism<=18; ++ism) smPhiLimits[ism] = -TMath::Pi() + (2.0*TMath::Pi()/18.0)*ism;
2612   CheckActiveSM(phiProj, activeSM);
2613   for(Int_t ism=0; ism<18; ++ism) cout << "sm " << ism << " is active : " << (activeSM[ism] ? "yes" : "no") << endl;
2614
2615   fCfContainer->SetRangeUser(kTrackPt, 1.01, 2.99);   // 1.0 < pt < 3.0 GeV/c
2616   TH2D* hTPCPhiBC = (TH2D*)fCfContainer->Project(0, kEventBC, kTrackPhi);
2617   TH2D* hTRDPhiBC = (TH2D*)fCfContainer->Project(1, kEventBC, kTrackPhi);
2618   TH2D* hTOFPhiBC = (TH2D*)fCfContainer->Project(2, kEventBC, kTrackPhi);
2619   TH1D* projectionBC = (TH1D*)fCfContainer->Project(0, kEventBC);
2620   fCfContainer->SetRangeUser(kTrackPt, 0.0, 100.0);   // reset the pt range
2621   TH1D* hTRDEffBC = new TH1D("hTRDEffBC", "", hTPCPhiBC->GetXaxis()->GetNbins(), hTPCPhiBC->GetXaxis()->GetXmin(), hTPCPhiBC->GetXaxis()->GetXmax());
2622   TH1D* hTOFEffBC = new TH1D("hTOFEffBC", "", hTPCPhiBC->GetXaxis()->GetNbins(), hTPCPhiBC->GetXaxis()->GetXmin(), hTPCPhiBC->GetXaxis()->GetXmax());
2623
2624   for(Int_t bcBin=1; bcBin<=hTPCPhiBC->GetXaxis()->GetNbins(); ++bcBin) {
2625     if(projectionBC->GetBinContent(bcBin)<0.1) continue;
2626     Double_t tpcEntries = 0.0; Double_t trdEntries = 0.0; Double_t tofEntries = 0.0;
2627     for(Int_t phiBin=1; phiBin<=hTPCPhiBC->GetYaxis()->GetNbins(); ++phiBin) {
2628       Double_t phi = hTPCPhiBC->GetYaxis()->GetBinCenter(phiBin);
2629       for(Int_t ism=0; ism<18; ++ism) {
2630         if(phi>=smPhiLimits[ism] && phi<smPhiLimits[ism+1] && activeSM[ism]) {
2631           tpcEntries += hTPCPhiBC->GetBinContent(bcBin, phiBin);
2632           trdEntries += hTRDPhiBC->GetBinContent(bcBin, phiBin);
2633           tofEntries += hTOFPhiBC->GetBinContent(bcBin, phiBin);
2634         }
2635       }  // end loop over super-modules
2636     }  // end loop over phi bins
2637     hTRDEffBC->SetBinContent(bcBin, (tpcEntries>0.01 ? trdEntries/tpcEntries : 0.0));
2638     if(tpcEntries>0.01 && trdEntries>0.01 && (tpcEntries-trdEntries)>=0.01)
2639     hTRDEffBC->SetBinError(bcBin, TMath::Sqrt(trdEntries*(tpcEntries-trdEntries)/tpcEntries/tpcEntries/tpcEntries));
2640     hTOFEffBC->SetBinContent(bcBin, (trdEntries>0.01 ? tofEntries/trdEntries : 0.0));
2641     if(trdEntries>0.01 && tofEntries>0.01 && (trdEntries-tofEntries)>=0.01)
2642     hTOFEffBC->SetBinError(bcBin, TMath::Sqrt(tofEntries*(trdEntries-tofEntries)/trdEntries/trdEntries/trdEntries));
2643   }  // end loop over BC bins
2644   delete hTPCPhiBC; delete hTRDPhiBC; delete hTOFPhiBC;
2645   delete projectionBC;
2646
2647   if(gROOT->FindObject("rangeBC")) delete gROOT->FindObject("rangeBC");
2648   TH2F* rangeBC = new TH2F("rangeBC", "", 10, -0.5, 3499.5, 10, 0.0, 1.4);
2649   rangeBC->SetStats(kFALSE);
2650   SetStyle(rangeBC->GetXaxis(), "Bunch crossing", 0.07, 0.8, kTRUE, 0.05);
2651   SetStyle(rangeBC->GetYaxis(), "efficiency", 0.07, 0.8, kTRUE, 0.05);
2652   rangeBC->Draw();
2653
2654   TLegend* legBC=new TLegend(0.8, 0.7, 0.95, 0.89);
2655   legBC->SetBorderSize(0);
2656   legBC->SetMargin(0.15);
2657   legBC->SetFillColor(0);
2658   if(hTRDEffBC) {
2659     hTRDEffBC->SetStats(kFALSE);
2660     SetStyle(hTRDEffBC, 1, kRed, 2, 24, kRed, 1); legBC->AddEntry(hTRDEffBC, "TPC-TRD", "p");
2661     SetStyle(hTOFEffBC, 1, kBlue, 2, 24, kBlue, 1); legBC->AddEntry(hTOFEffBC, "TRD-TOF", "p");
2662     hTRDEffBC->Draw("same");
2663     hTOFEffBC->Draw("same");
2664     legBC->Draw();
2665     lat->DrawLatex(200., 1.42, "Matching efficiency at 1<p_{T}<3 GeV/c");
2666   }
2667
2668   // reset the user range on the event multiplicity
2669   fCfContainer->SetRangeUser(kEventMult, 0.0, 6.0, kTRUE);
2670   if(cutTOFbc) fCfContainer->SetRangeUser(kTrackTOFBC, -1000.0, +1000.0);  // reset the cut on TOFbc
2671
2672   delete funcConst;
2673 }
2674
2675
2676 //_________________________________________________________________
2677 void AliTRDcheckESD::PlotPidSummaryFromCF(Int_t centralityClass, Double_t* trendValues, Bool_t useIsolatedBC, Bool_t cutTOFbc) {
2678
2679   if(!fCfContainer) return;
2680
2681   TLatex *lat=new TLatex();
2682   lat->SetTextSize(0.07);
2683   lat->SetTextColor(2);
2684   gPad->SetTopMargin(0.05); gPad->SetBottomMargin(0.001);
2685   gPad->SetLeftMargin(0.001); gPad->SetRightMargin(0.001);
2686   gPad->Divide(3,3,0.,0.);
2687   TList* l=gPad->GetListOfPrimitives();
2688
2689   fCfContainer->SetRangeUser(kTrackDCAxy, -0.999, +0.999);
2690   fCfContainer->SetRangeUser(kTrackDCAz, -3.0, +3.0);
2691   if(cutTOFbc) fCfContainer->SetRangeUser(kTrackTOFBC, 0.0, 0.0);
2692
2693   // find all the isolated bunch crossings with entries
2694   Bool_t isIsolated[3500];
2695   TH1D* tempTH1D = (TH1D*)fCfContainer->Project(0, kEventBC);
2696   FindIsolatedBCs(tempTH1D, isIsolated); delete tempTH1D;
2697   Int_t nIsolatedBC = 0;
2698   for(Int_t ibc=0; ibc<3500; ++ibc)
2699     if(isIsolated[ibc]) nIsolatedBC++;
2700
2701   if(centralityClass>0) // select the multiplicity class
2702     fCfContainer->SetRangeUser(kEventMult, Double_t(centralityClass), Double_t(centralityClass), kTRUE);
2703
2704   if(gROOT->FindObject("rangeEtaPhi2")) delete gROOT->FindObject("rangeEtaPhi2");
2705   TH2F* rangeEtaPhi = new TH2F("rangeEtaPhi2", "", 10, -0.99, +0.99, 10, -3.4, +3.4);
2706   SetStyle(rangeEtaPhi->GetXaxis(), "#eta", 0.07, 0.8, kTRUE, 0.05);
2707   SetStyle(rangeEtaPhi->GetYaxis(), "detector #varphi", 0.07, 0.8, kTRUE, 0.05);
2708   rangeEtaPhi->SetStats(kFALSE);
2709
2710   // eta-phi distr. for <Qtot> in layer 0
2711   TVirtualPad* pad;
2712   TProfile2D* hProf2D;
2713   for(Int_t iLayer=0; iLayer<6; ++iLayer) {
2714     pad = ((TVirtualPad*)l->At((iLayer<3 ? iLayer*3 : (iLayer-3)*3+1))); pad->cd();
2715     pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
2716     pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
2717     pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2718     rangeEtaPhi->Draw();
2719
2720     TH3D* hQtotEtaPhi=0x0;
2721     if(useIsolatedBC) {
2722       for(Int_t ibc=0; ibc<3500; ++ibc) {
2723         if(!isIsolated[ibc]) continue;
2724         fCfContainer->SetRangeUser(kEventBC, Double_t(ibc), Double_t(ibc));
2725         TH3D* temp = (TH3D*)fCfContainer->Project(1, kTrackPhi, kTrackEta, kTrackQtot+iLayer);
2726         if(!hQtotEtaPhi) hQtotEtaPhi = temp;
2727         else {hQtotEtaPhi->Add(temp); delete temp;}
2728       }
2729       fCfContainer->SetRangeUser(kEventBC, 0.0, 3500.0);
2730     }
2731     else hQtotEtaPhi = (TH3D*)fCfContainer->Project(1, kTrackPhi, kTrackEta, kTrackQtot+iLayer);
2732     hProf2D = (hQtotEtaPhi ? hQtotEtaPhi->Project3DProfile() : 0x0);
2733     if(hQtotEtaPhi) delete hQtotEtaPhi;
2734
2735     if(hProf2D) {
2736       hProf2D->SetStats(kFALSE);
2737       hProf2D->SetMinimum(0.);
2738       hProf2D->SetMaximum(4.);
2739       hProf2D->Draw("samecolz");
2740     }
2741     lat->DrawLatex(-0.9, 3.6, Form("TRD <Q_{tot}> Layer %d", iLayer));
2742     DrawTRDGrid();
2743   }
2744
2745   // PH versus slice number
2746   pad = ((TVirtualPad*)l->At(2)); pad->cd();
2747   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
2748   pad->SetTopMargin(0.03); pad->SetBottomMargin(0.15);
2749   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2750
2751   if(gROOT->FindObject("rangePHslice")) delete gROOT->FindObject("rangePHslice");
2752   TH2F* rangePHslice=new TH2F("rangePHslice", "", 8, -0.5, 7.5, 10, 0.0, 2000.);
2753   rangePHslice->SetStats(kFALSE);
2754   SetStyle(rangePHslice->GetXaxis(), "slice", 0.07, 0.8, kTRUE, 0.05);
2755   SetStyle(rangePHslice->GetYaxis(), "PH", 0.07, 0.8, kTRUE, 0.05);
2756   rangePHslice->Draw();
2757
2758   TF1* funcPol1 = new TF1("funcPol1", "[0]+[1]*x", 2.9, 6.4);
2759
2760   TH2F* h2F;
2761   TH1D* hF;
2762   if((h2F = dynamic_cast<TH2F*>(fHistos->At(kPHSlice+centralityClass)))) {
2763     hF = Proj2D(h2F);
2764     h2F->SetStats(kFALSE);
2765     h2F->Draw("samecolz");
2766     if(trendValues) {
2767       hF->Fit(funcPol1, "QME0", "goff", 2.9, 6.4);
2768       trendValues[12] = funcPol1->GetParameter(0);  // PH plateau
2769       trendValues[13] = funcPol1->GetParError(0);   // PH plateau
2770       trendValues[14] = funcPol1->GetParameter(1);  // PH slope
2771       trendValues[15] = funcPol1->GetParError(1);   // PH slope
2772     }
2773     hF->SetLineWidth(2);
2774     hF->SetLineStyle(2);
2775     hF->Draw("same");
2776   }
2777   delete funcPol1;
2778
2779   // Qtot vs P
2780   pad = ((TVirtualPad*)l->At(5)); pad->cd();
2781   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
2782   pad->SetTopMargin(0.03); pad->SetBottomMargin(0.15);
2783   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2784   pad->SetLogz();
2785
2786   if(gROOT->FindObject("rangeQtotP")) delete gROOT->FindObject("rangeQtotP");
2787   TH2F* rangeQtotP = new TH2F("rangeQtotP", "", 10, 0.0, 11.99, 10, 0.0, 11.99);
2788   SetStyle(rangeQtotP->GetXaxis(), "P [GeV/c]", 0.07, 0.8, kTRUE, 0.05);
2789   SetStyle(rangeQtotP->GetYaxis(), "Q_{tot}", 0.07, 0.8, kTRUE, 0.05);
2790   rangeQtotP->SetStats(kFALSE);
2791   rangeQtotP->Draw();
2792
2793   TH2D* hQtotP_IsolatedBC = 0x0;
2794   for(Int_t ibc=0; ibc<3500; ++ibc) {
2795     if(!isIsolated[ibc]) continue;
2796     fCfContainer->SetRangeUser(kEventBC, Double_t(ibc), Double_t(ibc));
2797     TH2D* temp = (TH2D*)fCfContainer->Project(1, kTrackP, kTrackQtot);
2798     if(!hQtotP_IsolatedBC) hQtotP_IsolatedBC = temp;
2799     else {hQtotP_IsolatedBC->Add(temp); delete temp;}
2800     for(Int_t il=1; il<6; ++il) {
2801       temp = (TH2D*)fCfContainer->Project(1, kTrackP, kTrackQtot+il);
2802       hQtotP_IsolatedBC->Add(temp); delete temp;
2803     }
2804   }
2805   fCfContainer->SetRangeUser(kEventBC, 0.0, 3500.0);
2806
2807   TH2D* temp = (TH2D*)fCfContainer->Project(1, kTrackP, kTrackQtot);
2808   TH2D* hQtotP = temp;
2809   for(Int_t il=1; il<6; ++il) {
2810     temp = (TH2D*)fCfContainer->Project(1, kTrackP, kTrackQtot+il);
2811     hQtotP->Add(temp); delete temp;
2812   }
2813
2814   if(hQtotP_IsolatedBC)
2815     for(Int_t i=1; i<=hQtotP_IsolatedBC->GetXaxis()->GetNbins(); ++i)
2816       hQtotP_IsolatedBC->SetBinContent(i, 1, 0.0);
2817   TH1D* hQtotProj_IsolatedBC = (hQtotP_IsolatedBC ? Proj2D(hQtotP_IsolatedBC) : 0x0);
2818   if(hQtotProj_IsolatedBC) SetStyle(hQtotProj_IsolatedBC, 2, kBlue, 2, 1, kBlue, 1);
2819   if(trendValues && hQtotProj_IsolatedBC && hQtotProj_IsolatedBC->GetEntries()>2) {
2820     trendValues[30] = hQtotProj_IsolatedBC->GetBinContent(hQtotProj_IsolatedBC->FindBin(1.0));   // Landau MPV at 1GeV/c
2821     trendValues[31] = hQtotProj_IsolatedBC->GetBinError(hQtotProj_IsolatedBC->FindBin(1.0));     // Landau width at 1 GeV/c
2822   }
2823   if(hQtotP)
2824     for(Int_t i=1; i<=hQtotP->GetXaxis()->GetNbins(); ++i)
2825       hQtotP->SetBinContent(i, 1, 0.0);
2826   TH1D* hQtotProj = (hQtotP ? Proj2D(hQtotP) : 0x0);
2827   if(hQtotProj) SetStyle(hQtotProj, 2, kBlue, 2, 1, kBlue, 1);
2828   if(trendValues && hQtotProj && hQtotProj->GetEntries()>2) {
2829     trendValues[16] = hQtotProj->GetBinContent(hQtotProj->FindBin(1.0));   // Landau MPV at 1GeV/c
2830     trendValues[17] = hQtotProj->GetBinError(hQtotProj->FindBin(1.0));     // Landau width at 1 GeV/c
2831   }
2832   if(useIsolatedBC && hQtotP_IsolatedBC) {
2833     hQtotP_IsolatedBC->SetStats(kFALSE);
2834     for(Int_t i=1; i<=hQtotP_IsolatedBC->GetXaxis()->GetNbins(); ++i) hQtotP_IsolatedBC->SetBinContent(i, 1, 0.0);
2835     hQtotP_IsolatedBC->Draw("samecolz");
2836     hQtotProj_IsolatedBC->Draw("same");
2837   }
2838   if(!useIsolatedBC && hQtotP) {
2839     hQtotP->SetStats(kFALSE);
2840     hQtotP->Draw("samecolz");
2841     hQtotProj->Draw("same");
2842   }
2843
2844   // reset the user range on the event multiplicity
2845   fCfContainer->SetRangeUser(kEventMult, 0.0, 6.0, kTRUE);
2846   if(cutTOFbc) fCfContainer->SetRangeUser(kTrackTOFBC, -1000.0, +1000.0);  // reset the cut on TOFbc
2847 }
2848
2849
2850 //_________________________________________________________________
2851 Bool_t AliTRDcheckESD::PlotCentSummary(Double_t* trendValues) {
2852
2853   Bool_t isGoodForSaving=kFALSE;
2854
2855   trendValues = trendValues;
2856
2857   TLatex* lat=new TLatex();
2858   lat->SetTextSize(0.06);
2859   lat->SetTextColor(2);
2860
2861   gPad->SetTopMargin(0.05); gPad->SetBottomMargin(0.001);
2862   gPad->SetLeftMargin(0.001); gPad->SetRightMargin(0.001);
2863   gPad->Divide(3,3,0.,0.);
2864   TList* l=gPad->GetListOfPrimitives();
2865
2866   TPad* pad=0x0;
2867
2868   if(gROOT->FindObject("rangeEffPt")) delete gROOT->FindObject("rangeEffPt");
2869   TH2F* rangeEffPt=new TH2F("rangeEffPt", "",10,0.,10.,10,0.,1.4);
2870   rangeEffPt->SetStats(kFALSE);
2871   SetStyle(rangeEffPt->GetXaxis(), "p_{T} [GeV/c]", 0.07, 0.8, kTRUE, 0.05);
2872   SetStyle(rangeEffPt->GetYaxis(), "efficiency", 0.07, 0.8, kTRUE, 0.05);
2873
2874   TH3F *h3(NULL), *h3p(NULL), *h3n(NULL);
2875   Int_t padsForEffs[5] = {0,3,6,1,4};
2876   for(Int_t iCent=1; iCent<6; ++iCent) {
2877     // TPC-TRD matching efficiencies
2878     pad = ((TPad*)l->At(padsForEffs[iCent-1])); pad->cd();
2879     pad->SetLeftMargin(0.15); pad->SetRightMargin(0.02); pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
2880     pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2881
2882     if(!(h3p = dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksPos+iCent)))) continue;
2883     if(!(h3n = dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksNeg+iCent)))) continue;
2884     // =============================================
2885     if(!(h3 = dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksPos+iCent)))) continue;
2886     TH1F* hFeffP = EfficiencyTRD(h3p, h3, kTRUE);
2887     //
2888     if(!(h3 = dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksNeg+iCent)))) continue;
2889     TH1F* hFeffN = EfficiencyTRD(h3n, h3, kTRUE);
2890     // =============================================
2891     if(!(h3 = dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksPos4+iCent)))) continue;
2892     TH1F* hFeffP4 = EfficiencyTRD(h3p, h3, kTRUE);
2893     //
2894     if(!(h3 = dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksNeg4+iCent)))) continue;
2895     TH1F* hFeffN4 = EfficiencyTRD(h3n, h3, kTRUE);
2896     // =============================================
2897     if(!(h3 = dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksPos5+iCent)))) continue;
2898     TH1F* hFeffP5 = EfficiencyTRD(h3p, h3, kTRUE);
2899     //
2900     if(!(h3 = dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksNeg5+iCent)))) continue;
2901     TH1F* hFeffN5 = EfficiencyTRD(h3n, h3, kTRUE);
2902     // =============================================
2903     if(!(h3 = dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksPos6+iCent)))) continue;
2904     TH1F* hFeffP6 = EfficiencyTRD(h3p, h3, kTRUE);
2905     //
2906     if(!(h3 = dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksNeg6+iCent)))) continue;
2907     TH1F* hFeffN6 = EfficiencyTRD(h3n, h3, kTRUE);
2908
2909     rangeEffPt->Draw();
2910
2911     TLine line;
2912     line.SetLineStyle(2);
2913     line.SetLineWidth(2);
2914     line.DrawLine(rangeEffPt->GetXaxis()->GetXmin(), 0.7, rangeEffPt->GetXaxis()->GetXmax(), 0.7);
2915     line.DrawLine(rangeEffPt->GetXaxis()->GetXmin(), 0.9, rangeEffPt->GetXaxis()->GetXmax(), 0.9);
2916     line.SetLineStyle(1);
2917     line.SetLineWidth(1);
2918     line.DrawLine(rangeEffPt->GetXaxis()->GetXmin(), 1.0, rangeEffPt->GetXaxis()->GetXmax(), 1.0);
2919     if(hFeffP) SetStyle(hFeffP, 1, kRed, 1, 24, kRed, 1);
2920     if(hFeffP4) SetStyle(hFeffP4, 1, kRed, 1, 25, kRed, 1);
2921     if(hFeffP5) SetStyle(hFeffP5, 1, kRed, 1, 26, kRed, 1);
2922     if(hFeffP6) SetStyle(hFeffP6, 1, kRed, 1, 27, kRed, 1);
2923     if(hFeffN) SetStyle(hFeffN, 1, kBlue, 1, 24, kBlue, 1);
2924     if(hFeffN4) SetStyle(hFeffN4, 1, kBlue, 1, 25, kBlue, 1);
2925     if(hFeffN5) SetStyle(hFeffN5, 1, kBlue, 1, 26, kBlue, 1);
2926     if(hFeffN6) SetStyle(hFeffN6, 1, kBlue, 1, 27, kBlue, 1);
2927
2928     TLegend* leg=new TLegend(0.16, 0.7, 0.61, 0.89);
2929     leg->SetFillColor(0);
2930     leg->SetNColumns(2);
2931     leg->SetTextSize(0.039);
2932     leg->SetMargin(0.1);
2933     if(hFeffP && hFeffP->Integral()>0.001) {
2934       isGoodForSaving = kTRUE;
2935       hFeffP->Draw("same");
2936       leg->AddEntry(hFeffP, "pos. (#geq 1 trcklt)", "p");
2937     }
2938     if(hFeffN && hFeffN->Integral()>0.001) {
2939       isGoodForSaving = kTRUE;
2940       hFeffN->Draw("same");
2941       leg->AddEntry(hFeffN, "neg. (#geq 1 trcklt)", "p");
2942     }
2943     if(hFeffP4 && hFeffP4->Integral()>0.001) {
2944       isGoodForSaving = kTRUE;
2945       hFeffP4->Draw("same");
2946       leg->AddEntry(hFeffP4, "pos. (4 trcklts)", "p");
2947     }
2948     if(hFeffN4 && hFeffN4->Integral()>0.001) {
2949       isGoodForSaving = kTRUE;
2950       hFeffN4->Draw("same");
2951       leg->AddEntry(hFeffN4, "neg. (4 trcklts)", "p");
2952     }
2953     if(hFeffP5 && hFeffP5->Integral()>0.001) {
2954       isGoodForSaving = kTRUE;
2955       hFeffP5->Draw("same");
2956       leg->AddEntry(hFeffP5, "pos. (5 trcklts)", "p");
2957     }
2958     if(hFeffN5 && hFeffN5->Integral()>0.001) {
2959       isGoodForSaving = kTRUE;
2960       hFeffN5->Draw("same");
2961       leg->AddEntry(hFeffN5, "neg. (5 trcklts)", "p");
2962     }
2963     if(hFeffP6 && hFeffP6->Integral()>0.001) {
2964       isGoodForSaving = kTRUE;
2965       hFeffP6->Draw("same");
2966       leg->AddEntry(hFeffP6, "pos. (6 trcklts)", "p");
2967     }
2968     if(hFeffN6 && hFeffN6->Integral()>0.001) {
2969       isGoodForSaving = kTRUE;
2970       hFeffN6->Draw("same");
2971       leg->AddEntry(hFeffN6, "neg. (6 trklts)", "p");
2972     }
2973
2974     if(isGoodForSaving) {
2975       if(iCent==1) leg->Draw();
2976       lat->DrawLatex(5.6, 1.3, Form("Centrality class %d", iCent));
2977       lat->DrawLatex(0.5, 1.42, "TPC-TRD matching efficiency");
2978     }
2979   }   // end loop over multiplicity intervals
2980
2981   // Number of clusters per TRD track
2982   pad = ((TPad*)l->At(2)); pad->cd();
2983   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.02);
2984   pad->SetTopMargin(0.02); pad->SetBottomMargin(0.15);
2985   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
2986   if(gROOT->FindObject("rangeNcls")) delete gROOT->FindObject("rangeNcls");
2987   TH2F* rangeNcls = new TH2F("rangeNcls", "", 10, 0.0, 199.9, 10, 0.0, 1.199);
2988   SetStyle(rangeNcls->GetXaxis(), "# TRD clusters", 0.07, 0.8, kTRUE, 0.05);
2989   SetStyle(rangeNcls->GetYaxis(), "entries (a.u.)", 0.07, 0.8, kTRUE, 0.05);
2990   rangeNcls->SetStats(kFALSE);
2991   rangeNcls->Draw();
2992
2993   TH2F* h2F[6]; TH1D* proj[6];
2994   TLegend* leg=new TLegend(0.2, 0.7, 0.5, 0.95);
2995   leg->SetFillColor(0);
2996   Bool_t isGood=kFALSE;
2997   for(Int_t iCent=0; iCent<6; ++iCent) {
2998     h2F[iCent] = dynamic_cast<TH2F*>(fHistos->At(kNClsTrackTRD+iCent));
2999     proj[iCent] = (h2F[iCent] && h2F[iCent]->GetEntries()>10 ? h2F[iCent]->ProjectionY(Form("projCent%d",iCent)) : 0x0);
3000     if(proj[iCent]) {
3001       proj[iCent]->SetLineColor(iCent<4 ? iCent+1 : iCent+2);
3002       Double_t maximum = proj[iCent]->GetMaximum();
3003       if(maximum>1.0)
3004         proj[iCent]->Scale(1.0/maximum);
3005       proj[iCent]->SetStats(kFALSE);
3006       proj[iCent]->Draw("same");
3007       leg->AddEntry(proj[iCent], (iCent==0 ? "all centralities" : Form("centrality class %d", iCent)), "l");
3008       isGood = kTRUE;
3009     }
3010   }
3011   if(isGood) leg->Draw();
3012   isGoodForSaving = isGoodForSaving || isGood;
3013
3014   // Qtot vs P
3015   pad = ((TPad*)l->At(5)); pad->cd();
3016   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.02);
3017   pad->SetTopMargin(0.02); pad->SetBottomMargin(0.15);
3018   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3019   if(gROOT->FindObject("rangeQtot")) delete gROOT->FindObject("rangeQtot");
3020   TH2F* rangeQtot = new TH2F("rangeQtot", "", 10, 0.0, 9.999, 10, 0.0, 1.199);
3021   SetStyle(rangeQtot->GetXaxis(), "Q_{tot} (a.u.)", 0.07, 0.8, kTRUE, 0.05);
3022   SetStyle(rangeQtot->GetYaxis(), "entries (a.u.)", 0.07, 0.8, kTRUE, 0.05);
3023   rangeQtot->SetStats(kFALSE);
3024   rangeQtot->Draw();
3025
3026   TH1D* projQ[6];
3027   TLegend* leg2=new TLegend(0.6, 0.7, 0.9, 0.95);
3028   leg2->SetFillColor(0);
3029   isGood = kFALSE;
3030   for(Int_t iCent=0; iCent<6; ++iCent) {
3031     h2F[iCent] = dynamic_cast<TH2F*>(fHistos->At(kQtotP+iCent));
3032     projQ[iCent] = (h2F[iCent] && h2F[iCent]->GetEntries()>10 ? h2F[iCent]->ProjectionY(Form("projQCent%d",iCent)) : 0x0);
3033     if(projQ[iCent]) {
3034       projQ[iCent]->SetLineColor(iCent<4 ? iCent+1 : iCent+2);
3035       Double_t maximum = projQ[iCent]->GetMaximum();
3036       if(maximum>1.0)
3037         projQ[iCent]->Scale(1.0/maximum);
3038       projQ[iCent]->SetStats(kFALSE);
3039       projQ[iCent]->Draw("same");
3040       leg2->AddEntry(projQ[iCent], (iCent==0 ? "all centralities" : Form("centrality class %d", iCent)), "l");
3041       isGood = kTRUE;
3042     }
3043   }
3044   if(isGood) leg2->Draw();
3045   isGoodForSaving = isGoodForSaving || isGood;
3046   return isGoodForSaving;
3047 }
3048
3049
3050 //_________________________________________________________________
3051 Bool_t AliTRDcheckESD::PlotTrackingSummary(Int_t centralityClass, Double_t* trendValues) {
3052
3053   Bool_t isGoodForSaving=kFALSE;
3054
3055   TLatex *lat=new TLatex();
3056   lat->SetTextSize(0.07);
3057   lat->SetTextColor(2);
3058   gPad->SetTopMargin(0.05); gPad->SetBottomMargin(0.001);
3059   gPad->SetLeftMargin(0.001); gPad->SetRightMargin(0.001);
3060   gPad->Divide(3,3,0.,0.);
3061   TList* l=gPad->GetListOfPrimitives();
3062   // eta-phi distr. for positive TPC tracks
3063   TVirtualPad* pad = ((TVirtualPad*)l->At(0)); pad->cd();
3064   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
3065   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3066   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3067   if(gROOT->FindObject("rangeEtaPhi")) delete gROOT->FindObject("rangeEtaPhi");
3068   TH2F* rangeEtaPhi = new TH2F("rangeEtaPhi", "", 10, -0.99, +0.99, 10, -3.4, +3.4);
3069   SetStyle(rangeEtaPhi->GetXaxis(), "#eta", 0.07, 0.8, kTRUE, 0.05);
3070   SetStyle(rangeEtaPhi->GetYaxis(), "detector #varphi", 0.07, 0.8, kTRUE, 0.05);
3071   rangeEtaPhi->SetStats(kFALSE);
3072   rangeEtaPhi->Draw();
3073   lat->DrawLatex(-0.9, 3.6, "TPC positive ref. tracks");
3074
3075   TH3F* h3F = dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksPos+centralityClass));
3076   TH2F* h2FtpcP = 0x0;
3077   Float_t nada=0.0;
3078   if(h3F && h3F->GetEntries()>10) {
3079     h2FtpcP = (TH2F*)Proj3D((TH3F*)h3F, 0x0, 1, h3F->GetZaxis()->GetNbins(), nada)->Clone();
3080     h2FtpcP->SetStats(kFALSE);
3081     h2FtpcP->Draw("samecolz"); isGoodForSaving = kTRUE;
3082     isGoodForSaving = kTRUE;
3083   }
3084   //-----------------
3085   // eta-phi distr. for negative TPC tracks
3086   pad = ((TVirtualPad*)l->At(1)); pad->cd();
3087   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
3088   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3089   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3090   rangeEtaPhi->Draw();
3091   lat->DrawLatex(-0.9, 3.6, "TPC negative ref. tracks");
3092
3093   h3F = dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksNeg+centralityClass));
3094   TH2F* h2FtpcN = 0x0;
3095   if(h3F && h3F->GetEntries()>10) {
3096     h2FtpcN = (TH2F*)Proj3D((TH3F*)h3F, 0x0, 1, h3F->GetZaxis()->GetNbins(), nada)->Clone();
3097     h2FtpcN->SetStats(kFALSE);
3098     h2FtpcN->Draw("samecolz");
3099     isGoodForSaving = kTRUE;
3100   }
3101   //----------------------------------------------
3102   // eta-phi distr. for positive TRD tracks
3103   pad = ((TVirtualPad*)l->At(3)); pad->cd();
3104   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
3105   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3106   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3107   rangeEtaPhi->Draw();
3108   lat->DrawLatex(-0.9, 3.6, "TRD positive ref. tracks");
3109
3110   h3F = dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksPos+centralityClass));
3111   TH2F* h2FtrdP = 0x0;
3112   if(h3F && h3F->GetEntries()>10) {
3113     h2FtrdP = (TH2F*)Proj3D((TH3F*)h3F, 0x0, 1, h3F->GetZaxis()->GetNbins(), nada)->Clone();
3114     h2FtrdP->SetStats(kFALSE);
3115     h2FtrdP->SetMaximum((h2FtpcP ? h2FtpcP->GetMaximum() : h2FtrdP->GetMaximum()));
3116     h2FtrdP->Draw("samecolz");
3117     isGoodForSaving=kTRUE;
3118   }
3119   //--------------------------------------------
3120   // eta-phi distr. for negative TRD tracks
3121   pad = ((TVirtualPad*)l->At(4)); pad->cd();
3122   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
3123   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3124   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3125   rangeEtaPhi->Draw();
3126   lat->DrawLatex(-0.9, 3.6, "TRD negative ref. tracks");
3127
3128   h3F = dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksNeg+centralityClass));
3129   TH2F* h2FtrdN = 0x0;
3130   if(h3F && h3F->GetEntries()>10) {
3131     h2FtrdN = (TH2F*)Proj3D((TH3F*)h3F, 0x0, 1, h3F->GetZaxis()->GetNbins(), nada)->Clone();
3132     h2FtrdN->SetStats(kFALSE);
3133     h2FtrdN->SetMaximum(h2FtpcN ? h2FtpcN->GetMaximum() : h2FtrdN->GetMaximum());
3134     h2FtrdN->Draw("samecolz");
3135     isGoodForSaving=kTRUE;
3136   }
3137   //----------------------------------------------
3138   // eta-phi efficiency for positive TRD tracks
3139   pad = ((TVirtualPad*)l->At(6)); pad->cd();
3140   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
3141   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3142   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3143   rangeEtaPhi->Draw();
3144   lat->DrawLatex(-0.9, 3.6, "Efficiency positive tracks");
3145
3146   TH2F* h2Feff = (h2FtrdP ? (TH2F*)h2FtrdP->Clone("h2FeffPos") : 0x0);
3147   if(h2Feff) {
3148     h2Feff->Reset();
3149     h2Feff->SetStats(kFALSE);
3150     h2Feff->Divide(h2FtrdP, h2FtpcP);
3151     h2Feff->SetMaximum(1.0);
3152     if(h2Feff->GetEntries()>1) { h2Feff->Draw("samecolz"); isGoodForSaving=kTRUE; }
3153   }
3154   //-------------------------------------------------
3155   // eta-phi efficiency for negative TRD tracks
3156   pad = ((TVirtualPad*)l->At(7)); pad->cd();
3157   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
3158   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3159   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3160   rangeEtaPhi->Draw();
3161   lat->DrawLatex(-0.9, 3.6, "Efficiency negative tracks");
3162
3163   h2Feff = (h2FtrdN ? (TH2F*)h2FtrdN->Clone("h2FeffNeg") : 0x0);
3164   if(h2Feff) {
3165     h2Feff->Reset();
3166     h2Feff->SetStats(kFALSE);
3167     h2Feff->Divide(h2FtrdN, h2FtpcN);
3168     h2Feff->SetMaximum(1.0);
3169     if(h2Feff->GetEntries()>0.1) { h2Feff->Draw("samecolz"); isGoodForSaving=kTRUE; }
3170   }
3171   //-----------------------------------------------------
3172   // <ntracklets> vs (phi,eta)
3173   pad = ((TVirtualPad*)l->At(2)); pad->cd();
3174   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
3175   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3176   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3177   rangeEtaPhi->Draw();
3178   lat->DrawLatex(-0.9, 3.6, "TRD <N_{tracklets}>");
3179
3180   TProfile2D* hProf2D;
3181   if((hProf2D = dynamic_cast<TProfile2D*>(fHistos->At(kTRDEtaPhiAvNtrkl+centralityClass)))) {
3182     if(hProf2D->GetEntries()>10) {
3183       hProf2D->SetStats(kFALSE);
3184       hProf2D->SetMinimum(0.);
3185       hProf2D->SetMaximum(6.);
3186       if(hProf2D->GetEntries()>1) { hProf2D->Draw("samecolz"); isGoodForSaving = kTRUE; }
3187     }
3188   }
3189   //---------------------------------------------------------
3190   // TPC-TRD matching efficiency vs pt
3191   pad = ((TVirtualPad*)l->At(5)); pad->cd();
3192   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.02);
3193   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3194   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3195   TH1F* hFeffP = EfficiencyTRD(dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksPos+centralityClass)),
3196               dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksPos+centralityClass)), kTRUE);
3197   TH1F* hFeffN = EfficiencyTRD(dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksNeg+centralityClass)),
3198               dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksNeg+centralityClass)), kTRUE);
3199   TH1F* hFeffP4 = EfficiencyTRD(dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksPos+centralityClass)),
3200         dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksPos4+centralityClass)), kTRUE);
3201   TH1F* hFeffN4 = EfficiencyTRD(dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksNeg+centralityClass)),
3202         dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksNeg4+centralityClass)), kTRUE);
3203   TH1F* hFeffP5 = EfficiencyTRD(dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksPos+centralityClass)),
3204         dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksPos5+centralityClass)), kTRUE);
3205   TH1F* hFeffN5 = EfficiencyTRD(dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksNeg+centralityClass)),
3206         dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksNeg5+centralityClass)), kTRUE);
3207   TH1F* hFeffP6 = EfficiencyTRD(dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksPos+centralityClass)),
3208         dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksPos6+centralityClass)), kTRUE);
3209   TH1F* hFeffN6 = EfficiencyTRD(dynamic_cast<TH3F*>(fHistos->At(kTPCRefTracksNeg+centralityClass)),
3210         dynamic_cast<TH3F*>(fHistos->At(kTRDRefTracksNeg6+centralityClass)), kTRUE);
3211
3212   TF1* funcConst = new TF1("funcConst", "[0]", 1.0, 3.0);
3213
3214   if(gROOT->FindObject("rangeEffPt2")) delete gROOT->FindObject("rangeEffPt2");
3215   TH2F* rangeEffPt2=new TH2F("rangeEffPt2", "",10,0.,10.,10,0.,1.4);
3216   rangeEffPt2->SetStats(kFALSE);
3217   SetStyle(rangeEffPt2->GetXaxis(), "p_{T} [GeV/c]", 0.07, 0.8, kTRUE, 0.05);
3218   SetStyle(rangeEffPt2->GetYaxis(), "efficiency", 0.07, 0.8, kTRUE, 0.05);
3219   rangeEffPt2->Draw();
3220   lat->DrawLatex(0.5, 1.42, "TRD-TPC matching efficiency");
3221   //++++++++++++++++++
3222   TLine line;
3223   line.SetLineStyle(2);
3224   line.SetLineWidth(2);
3225   line.DrawLine(rangeEffPt2->GetXaxis()->GetXmin(), 0.7, rangeEffPt2->GetXaxis()->GetXmax(), 0.7);
3226   line.DrawLine(rangeEffPt2->GetXaxis()->GetXmin(), 0.9, rangeEffPt2->GetXaxis()->GetXmax(), 0.9);
3227   line.SetLineStyle(1);
3228   line.SetLineWidth(1);
3229   line.DrawLine(rangeEffPt2->GetXaxis()->GetXmin(), 1.0, rangeEffPt2->GetXaxis()->GetXmax(), 1.0);
3230   TLegend* leg=new TLegend(0.2, 0.7, 0.6, 0.89);
3231   leg->SetNColumns(2);
3232   leg->SetFillColor(0);
3233   if(hFeffP){
3234     hFeffP->SetMarkerStyle(24);
3235     hFeffP->SetMarkerColor(2);
3236     hFeffP->SetLineColor(2);
3237     if(trendValues && hFeffP->GetEntries()>1) {
3238       hFeffP->Fit(funcConst, "QME0", "goff", 1.0, 3.0);
3239       trendValues[0] = funcConst->GetParameter(0);
3240       trendValues[1] = funcConst->GetParError(0);
3241     }
3242     if(hFeffP->Integral()>0.001) {
3243       hFeffP->Draw("same");
3244       leg->AddEntry(hFeffP, "positives (#geq 1 tracklet)", "p");
3245     }
3246   }
3247   if(hFeffN){
3248     hFeffN->SetMarkerStyle(24);
3249     hFeffN->SetMarkerColor(4);
3250     hFeffN->SetLineColor(4);
3251     if(trendValues && hFeffN->GetEntries()>1) {
3252       hFeffN->Fit(funcConst, "QME0", "goff", 1.0, 3.0);
3253       trendValues[2] = funcConst->GetParameter(0);
3254       trendValues[3] = funcConst->GetParError(0);
3255     }
3256     if(hFeffN->Integral()>0.001) {
3257       hFeffN->Draw("same");
3258       leg->AddEntry(hFeffN, "negatives (#geq 1 tracklet)", "p");
3259     }
3260   }
3261   if(hFeffP4){
3262     hFeffP4->SetMarkerStyle(25);
3263     hFeffP4->SetMarkerColor(2);
3264     hFeffP4->SetLineColor(2);
3265     if(hFeffP4->Integral()>0.001) {
3266       hFeffP4->Draw("same");
3267       leg->AddEntry(hFeffP4, "positives (4 tracklets)", "p");
3268     }
3269   }
3270   if(hFeffN4){
3271     hFeffN4->SetMarkerStyle(25);
3272     hFeffN4->SetMarkerColor(4);
3273     hFeffN4->SetLineColor(4);
3274     if(hFeffN4->Integral()>0.001) {
3275       hFeffN4->Draw("same");
3276       leg->AddEntry(hFeffN4, "negatives (4 tracklets)", "p");
3277     }
3278   }
3279   if(hFeffP5){
3280     hFeffP5->SetMarkerStyle(26);
3281     hFeffP5->SetMarkerColor(2);
3282     hFeffP5->SetLineColor(2);
3283     if(hFeffP5->Integral()>0.001) {
3284       hFeffP5->Draw("same");
3285       leg->AddEntry(hFeffP5, "positives (5 tracklets)", "p");
3286     }
3287   }
3288   if(hFeffN5){
3289     hFeffN5->SetMarkerStyle(26);
3290     hFeffN5->SetMarkerColor(4);
3291     hFeffN5->SetLineColor(4);
3292     if(hFeffN5->Integral()>0.001) {
3293       hFeffN5->Draw("same");
3294       leg->AddEntry(hFeffN5, "negatives (5 tracklets)", "p");
3295     }
3296   }
3297   if(hFeffP6){
3298     hFeffP6->SetMarkerStyle(27);
3299     hFeffP6->SetMarkerColor(2);
3300     hFeffP6->SetLineColor(2);
3301     if(hFeffP6->Integral()>0.001) {
3302       hFeffP6->Draw("same");
3303       leg->AddEntry(hFeffP6, "positives (6 tracklets)", "p");
3304     }
3305   }
3306   if(hFeffN6){
3307     hFeffN6->SetMarkerStyle(27);
3308     hFeffN6->SetMarkerColor(4);
3309     hFeffN6->SetLineColor(4);
3310     if(hFeffN6->Integral()>0.001) {
3311       hFeffN6->Draw("same");
3312       leg->AddEntry(hFeffN6, "negatives (6 tracklets)", "p");
3313     }
3314   }
3315   leg->Draw();
3316
3317   //--------------------------------------------------------------
3318   // Nclusters per TRD track
3319   pad = ((TVirtualPad*)l->At(8)); pad->cd();
3320   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.12);
3321   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3322   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3323   pad->SetLogz();
3324   if(gROOT->FindObject("rangeNclsP")) delete gROOT->FindObject("rangeNclsP");
3325   TH2F* rangeNclsP = new TH2F("rangeNclsP", "", 10, 0.0, 11.99, 10, 0.0, 199.0);
3326   SetStyle(rangeNclsP->GetXaxis(), "p [GeV/c]", 0.07, 0.8, kTRUE, 0.05);
3327   SetStyle(rangeNclsP->GetYaxis(), "#clusters", 0.07, 0.8, kTRUE, 0.05);
3328   rangeNclsP->SetStats(kFALSE);
3329   rangeNclsP->Draw();
3330   lat->DrawLatex(1.0, 205., "TRD Clusters / track");
3331
3332   TH2F* hNclsVsP=0x0;
3333   if((hNclsVsP = dynamic_cast<TH2F*>(fHistos->At(kNClsTrackTRD+centralityClass)))) {
3334     hNclsVsP->SetStats(kFALSE);
3335     if(hNclsVsP->GetEntries()>10) {
3336       hNclsVsP->Draw("samecolz"); isGoodForSaving=kTRUE;
3337       if(trendValues) {
3338         TProfile* h2FProf = hNclsVsP->ProfileX("nclsVsPprof");
3339         h2FProf->Fit(funcConst, "QME0", "goff", 1.0, 3.0);
3340         trendValues[4] = funcConst->GetParameter(0);
3341         trendValues[5] = funcConst->GetParError(0);
3342       }
3343     }
3344   }
3345
3346   delete funcConst;
3347   return isGoodForSaving;
3348 }
3349
3350
3351 //_________________________________________________________________
3352 Bool_t AliTRDcheckESD::PlotPidSummary(Int_t centralityClass, Double_t* trendValues) {
3353
3354   Bool_t isGoodForSaving=kFALSE;
3355
3356   TLatex *lat=new TLatex();
3357   lat->SetTextSize(0.07);
3358   lat->SetTextColor(2);
3359   gPad->SetTopMargin(0.05); gPad->SetBottomMargin(0.001);
3360   gPad->SetLeftMargin(0.001); gPad->SetRightMargin(0.001);
3361   gPad->Divide(3,3,0.,0.);
3362   TList* l=gPad->GetListOfPrimitives();
3363   // eta-phi distr. for <Qtot> in layer 0
3364   TVirtualPad* pad;
3365   TProfile2D* hProf2D;
3366   if(gROOT->FindObject("rangeEtaPhi2")) delete gROOT->FindObject("rangeEtaPhi2");
3367   TH2F* rangeEtaPhi = new TH2F("rangeEtaPhi2", "", 10, -0.99, +0.99, 10, -3.4, +3.4);
3368   SetStyle(rangeEtaPhi->GetXaxis(), "#eta", 0.07, 0.8, kTRUE, 0.05);
3369   SetStyle(rangeEtaPhi->GetYaxis(), "detector #varphi", 0.07, 0.8, kTRUE, 0.05);
3370   rangeEtaPhi->SetStats(kFALSE);
3371
3372   for(Int_t iLayer=0; iLayer<6; ++iLayer) {
3373     pad = ((TVirtualPad*)l->At((iLayer<3 ? iLayer*3 : (iLayer-3)*3+1))); pad->cd();
3374     pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
3375     pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3376     pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3377     rangeEtaPhi->Draw();
3378     lat->DrawLatex(-0.9, 3.6, Form("TRD <Q_{tot}> Layer %d", iLayer));
3379
3380     if(!(hProf2D = dynamic_cast<TProfile2D*>(fHistos->At(kTRDEtaPhiAvQtot+6*centralityClass+iLayer)))) continue;
3381     if(hProf2D && hProf2D->GetEntries()>10) {
3382       hProf2D->SetStats(kFALSE);
3383       hProf2D->SetMinimum(0.);
3384       hProf2D->SetMaximum(4.);
3385       if(hProf2D->GetEntries()>10) { hProf2D->Draw("samecolz"); isGoodForSaving=kTRUE; }
3386     }
3387   }
3388
3389   // PH versus slice number
3390   pad = ((TVirtualPad*)l->At(2)); pad->cd();
3391   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
3392   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3393   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3394   if(gROOT->FindObject("rangePHslice")) delete gROOT->FindObject("rangePHslice");
3395   TH2F* rangePHslice = new TH2F("rangePHslice", "", 10, -0.5, 7.5, 10, 0.0, 2000.0);
3396   SetStyle(rangePHslice->GetXaxis(), "slice", 0.07, 0.8, kTRUE, 0.05);
3397   SetStyle(rangePHslice->GetYaxis(), "PH", 0.07, 0.8, kTRUE, 0.05);
3398   rangePHslice->SetStats(kFALSE);
3399   rangePHslice->Draw();
3400
3401   TF1* funcPol1 = new TF1("funcPol1", "[0]+[1]*x", 2.9, 6.4);
3402
3403   TH2F* h2F;
3404   TH1D* hF;
3405   if((h2F = dynamic_cast<TH2F*>(fHistos->At(kPHSlice+centralityClass)))) {
3406     if(h2F && h2F->GetEntries()>10) {
3407       hF = Proj2D(h2F);
3408       h2F->SetStats(kFALSE);
3409       h2F->Draw("samecolz");
3410       isGoodForSaving=kTRUE;
3411       if(trendValues) {
3412         hF->Fit(funcPol1, "QME0", "goff", 2.9, 6.4);
3413         trendValues[6] = funcPol1->GetParameter(0);
3414         trendValues[7] = funcPol1->GetParError(0);
3415         trendValues[8] = funcPol1->GetParameter(1);
3416         trendValues[9] = funcPol1->GetParError(1);
3417       }
3418       hF->SetLineWidth(2);
3419       hF->Draw("same");
3420     }
3421   }
3422   delete funcPol1;
3423
3424   // Qtot vs P
3425   pad = ((TVirtualPad*)l->At(5)); pad->cd();
3426   pad->SetLeftMargin(0.15); pad->SetRightMargin(0.1);
3427   pad->SetTopMargin(0.1); pad->SetBottomMargin(0.15);
3428   pad->SetGridx(kFALSE); pad->SetGridy(kFALSE);
3429   pad->SetLogz();
3430   if(gROOT->FindObject("rangeQtotP")) delete gROOT->FindObject("rangeQtotP");
3431   TH2F* rangeQtotP = new TH2F("rangeQtotP", "", 10, 0.0, 11.99, 10, 0.0, 11.99);
3432   SetStyle(rangeQtotP->GetXaxis(), "P [GeV/c]", 0.07, 0.8, kTRUE, 0.05);
3433   SetStyle(rangeQtotP->GetYaxis(), "Q_{tot}", 0.07, 0.8, kTRUE, 0.05);
3434   rangeQtotP->SetStats(kFALSE);
3435   rangeQtotP->Draw();
3436
3437   if((h2F = dynamic_cast<TH2F*>(fHistos->At(kQtotP+centralityClass)))) {
3438     if(h2F && h2F->GetEntries()>10) {
3439       h2F->SetStats(kFALSE);
3440       h2F->Draw("samecolz");
3441       isGoodForSaving=kTRUE;
3442       hF = Proj2D(h2F);
3443       hF->SetLineWidth(2);
3444       hF->Draw("same");
3445       if(trendValues) {
3446         trendValues[10] = hF->GetBinContent(hF->FindBin(1.0));
3447         trendValues[11] = hF->GetBinError(hF->FindBin(1.0));
3448       }
3449     }
3450   }
3451   return isGoodForSaving;
3452 }
3453
3454 //__________________________________________________________________________________________________
3455 void AliTRDcheckESD::DrawTRDGrid() {
3456   //
3457   //   Draw a grid of lines showing the TRD supermodule and stack structure in (eta,phi) coordinates.
3458   //   The canvas on which to draw must already exist.
3459   //
3460   TLine line;
3461   line.SetLineColor(2);
3462   line.SetLineWidth(1);
3463   line.SetLineStyle(2);
3464   for(Int_t i=0; i<=9; ++i) {
3465     line.DrawLine(-1.0, 2.0*TMath::Pi()/18.0*i, +1.0, 2.0*TMath::Pi()/18.0*i);
3466     line.DrawLine(-1.0, -2.0*TMath::Pi()/18.0*i, +1.0, -2.0*TMath::Pi()/18.0*i);
3467   }
3468   line.DrawLine(-0.85, -3.2, -0.85, +3.2);
3469   line.DrawLine(-0.54, -3.2, -0.54, +3.2);
3470   line.DrawLine(-0.16, -3.2, -0.16, +3.2);
3471   line.DrawLine(+0.16, -3.2, +0.16, +3.2);
3472   line.DrawLine(+0.54, -3.2, +0.54, +3.2);
3473   line.DrawLine(+0.85, -3.2, +0.85, +3.2);
3474 }
3475
3476 //_________________________________________________________________
3477 void AliTRDcheckESD::SetStyle(TH1* hist,
3478                               Int_t lineStyle, Int_t lineColor, Int_t lineWidth,
3479                               Int_t markerStyle, Int_t markerColor, Int_t markerSize) {
3480   //
3481   // Set style settings for histograms
3482   //
3483   hist->SetLineStyle(lineStyle);
3484   hist->SetLineColor(lineColor);
3485   hist->SetLineWidth(lineWidth);
3486   hist->SetMarkerStyle(markerStyle);
3487   hist->SetMarkerColor(markerColor);
3488   hist->SetMarkerSize(markerSize);
3489 }
3490
3491 //____________________________________________________________________
3492 void AliTRDcheckESD::SetStyle(TAxis* axis, const Char_t* title, Float_t titleSize, Float_t titleOffset, Bool_t centerTitle,
3493                               Float_t labelSize) {
3494   //
3495   // Set style settings for axes
3496   //
3497   axis->SetTitle(title);
3498   axis->SetTitleSize(titleSize);
3499   axis->SetTitleOffset(titleOffset);
3500   axis->CenterTitle(centerTitle);
3501   axis->SetLabelSize(labelSize);
3502 }
3503
3504 //____________________________________________________________________
3505 void AliTRDcheckESD::FindIsolatedBCs(TH1D* bcHist, Bool_t isIsolated[3500]) {
3506   //
3507   // Find the isolated bunch crossings
3508   //
3509   Int_t isolationSize = 10;      // number of free bunches in both directions
3510   for(Int_t bcBin=1; bcBin<=bcHist->GetXaxis()->GetNbins(); ++bcBin) {
3511     Int_t bc = TMath::Nint(bcHist->GetBinCenter(bcBin));
3512     if(bc<-0.001 || bc>3499.01) {
3513       isIsolated[bc] = kFALSE;
3514       continue;
3515     }
3516     Double_t entries = bcHist->GetBinContent(bcBin);
3517     if(entries<0.001) {
3518       isIsolated[bc] = kFALSE;
3519       continue;     // no entries
3520     }
3521
3522     // check isolation
3523     isIsolated[bc] = kTRUE;
3524     for(Int_t ibc = TMath::Max(1,bcBin-isolationSize); ibc <= TMath::Min(3499, bcBin+isolationSize); ++ibc) {
3525       if(ibc==bcBin) continue;
3526       if(bcHist->GetBinContent(ibc)>0.01) {
3527         isIsolated[bc] = kFALSE;
3528         break;
3529       }
3530     }
3531   }   // end loop over BC bins
3532
3533   cout << "Isolated bunches: " << endl;
3534   for(Int_t ibc=0; ibc<3500; ++ibc)
3535     if(isIsolated[ibc]) cout << "BC #" << ibc << endl;
3536 }