PWGDQ/dielectron/macrosJPSI/ConfigBJpsi_ff_PbPb.C

   1 void InitHistograms(AliDielectron *die, Int_t cutDefinition);
   2 void InitCF(AliDielectron* die, Int_t cutDefinition);
   3
   4 void SetupTrackCuts(Bool_t isESD, AliDielectron *die, Int_t cutDefinition);
   5 void SetupPairCuts(AliDielectron *die, Int_t cutDefinition);
   6
   7 void AddMCSignals(AliDielectron *die);
   8 void SetEtaCorrection();
   9 TVectorD *GetRunNumbers();
  10
  11 TString names=("TOFTRDany;TOFTRDfirst");
  12 enum { kTOFTRD, kTOFTRD2};
  13
  14 TObjArray *arrNames=names.Tokenize(";");
  15 const Int_t nDie=arrNames->GetEntries();
  16
  17 Bool_t hasMC=kFALSE;
  18
  19 AliDielectron* ConfigBJpsi_ff_PbPb(Int_t cutDefinition, Bool_t isMC=kFALSE)
  20 {
  21   //
  22   // Setup the instance of AliDielectron
  23   //
  24
  25
  26   // MC event handler?
  27   hasMC=isMC;
  28     //(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()!=0x0);
  29
  30   //ESD handler?
  31   Bool_t isESD=(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->IsA()==AliESDInputHandler::Class());
  32
  33   // switch off some configurations
  34   //switch(cutDefinition) {
  35    //case kTOFTRD:
  36    //case kTOFTRD2:
  37   //default:            return 0x0;      break;
  38   //}
  39
  40   // create the actual framework object
  41   TString name=Form("%02d",cutDefinition);
  42   if (cutDefinition<arrNames->GetEntriesFast()){
  43     name=arrNames->At(cutDefinition)->GetName();
  44   }
  45   AliDielectron *die = new AliDielectron(Form("%s",name.Data()),
  46                                          Form("Track cuts: %s",name.Data()));
  47   // debugTree
  48   AliDielectronDebugTree *tree = new AliDielectronDebugTree(Form("%s",name.Data()),Form("Track cuts: %s",name.Data()));
  49   tree->SetOutputFileName(Form("jpsi_debug_tree_%s.root",name.Data()));
  50
  51   if(tree){
  52            tree->AddLegVariable(AliDielectronVarManager::kPt);
  53            tree->AddLegVariable(AliDielectronVarManager::kPx);
  54            tree->AddLegVariable(AliDielectronVarManager::kPy);
  55            tree->AddLegVariable(AliDielectronVarManager::kEta);
  56            tree->AddLegVariable(AliDielectronVarManager::kXv);
  57            tree->AddLegVariable(AliDielectronVarManager::kYv);
  58            tree->AddLegVariable(AliDielectronVarManager::kE);
  59            tree->AddLegVariable(AliDielectronVarManager::kTPCnSigmaEle);
  60            tree->AddLegVariable(AliDielectronVarManager::kTPCsignal);
  61            tree->AddLegVariable(AliDielectronVarManager::kTRDntracklets);
  62            tree->AddLegVariable(AliDielectronVarManager::kNclsTRD);
  63            tree->AddLegVariable(AliDielectronVarManager::kTOFnSigmaEle);
  64            tree->AddLegVariable(AliDielectronVarManager::kITSLayerFirstCls);
  65            tree->AddPairVariable(AliDielectronVarManager::kM);
  66            tree->AddPairVariable(AliDielectronVarManager::kE);
  67            tree->AddPairVariable(AliDielectronVarManager::kP);
  68            tree->AddPairVariable(AliDielectronVarManager::kPt);
  69            tree->AddPairVariable(AliDielectronVarManager::kY);
  70            tree->AddPairVariable(AliDielectronVarManager::kEta);
  71            tree->AddPairVariable(AliDielectronVarManager::kPairType);
  72            tree->AddPairVariable(AliDielectronVarManager::kPseudoProperTime);
  73            tree->AddPairVariable(AliDielectronVarManager::kPseudoProperTimeErr);
  74            if(hasMC) tree->AddPairVariable(AliDielectronVarManager::kPseudoProperTimeResolution);
  75            if(hasMC) tree->AddPairVariable(AliDielectronVarManager::kPseudoProperTimePull);
  76            if(hasMC) tree->AddPairVariable(AliDielectronVarManager::kPdgCode);
  77            if(hasMC) tree->AddPairVariable(AliDielectronVarManager::kPdgCodeMother);
  78            if(hasMC) tree->AddPairVariable(AliDielectronVarManager::kPdgCodeGrandMother);
  79            if(hasMC) tree->AddPairVariable(AliDielectronVarManager::kIsJpsiPrimary);
  80            tree->AddPairVariable(AliDielectronVarManager::kCentralitySPD);
  81            tree->AddPairVariable(AliDielectronVarManager::kCentrality);
  82            tree->AddPairVariable(AliDielectronVarManager::kNevents);
  83            tree->AddPairVariable(AliDielectronVarManager::kXvPrim);
  84            tree->AddPairVariable(AliDielectronVarManager::kYvPrim);
  85            tree->AddPairVariable(AliDielectronVarManager::kZvPrim);
  86            tree->AddPairVariable(AliDielectronVarManager::kXRes);
  87            tree->AddPairVariable(AliDielectronVarManager::kYRes);
  88            tree->AddPairVariable(AliDielectronVarManager::kZRes);
  89
  90             }
  91   //
  92   // QA histogram setup
  93   //
  94   die->SetDebugTree(tree);
  95
  96
  97   // Monte Carlo Signals and TRD efficiency tables
  98   if(hasMC) {
  99     AddMCSignals(die);
 100
 101     // trd tables
 102     TString pidTab="$TRAIN_ROOT/util/dielectron/dielectron/TRDpidEff_eleProb07_TRDntr4_6.root";
 103     TString trainRoot=gSystem->Getenv("TRAIN_ROOT");
 104     if (trainRoot.IsNull()) pidTab="$ALICE_ROOT/PWGDQ/dielectron/files/TRDpidEff_eleProb07_TRDntr4_6.root";
 105
 106     if (gSystem->AccessPathName(gSystem->ExpandPathName(pidTab.Data())))
 107       Error("ConfigPbPb","PID table not found: %s",pidTab.Data());
 108     else
 109       die->SetTRDcorrectionFilename(pidTab.Data());
 110   }
 111
 112   // cut setup
 113   SetupTrackCuts(isESD,die,cutDefinition);
 114   SetupPairCuts(die,cutDefinition);
 115
 116   // histogram setup
 117   if(cutDefinition == kTOFTRD  ||
 118         cutDefinition == kTOFTRD2   )
 119     InitHistograms(die,cutDefinition);
 120
 121   // CF container setup
 122   //  InitCF(die,cutDefinition);
 123
 124   // bgrd estimators
 125   if(!hasMC) {
 126
 127     if(cutDefinition == kTOFTRD) {
 128       // rotations
 129       AliDielectronTrackRotator *rot=new AliDielectronTrackRotator;
 130       rot->SetIterations(10);
 131       rot->SetConeAnglePhi(TMath::Pi());
 132       rot->SetStartAnglePhi(TMath::Pi());
 133       die->SetTrackRotator(rot);
 134       // mixing
 135       AliDielectronMixingHandler *mix=new AliDielectronMixingHandler;
 136       mix->AddVariable(AliDielectronVarManager::kZvPrim,20,-10.,10.);
 137       mix->AddVariable(AliDielectronVarManager::kCentrality,"0,5,10,20,50,80");
 138       mix->SetMixType(AliDielectronMixingHandler::kAll);
 139       mix->SetDepth(50);
 140       die->SetMixingHandler(mix);
 141     }
 142
 143
 144     }
 145
 146
 147   // prefilter settings
 148   //if(cutDefinition == kTOFTRD2)
 149   //  die->SetPreFilterAllSigns();
 150   //else
 151     die->SetPreFilterUnlikeOnly();
 152
 153   // setup eta correction
 154   if(isESD) SetEtaCorrection();
 155
 156   return die;
 157 }
 158
 159 //______________________________________________________________________________________
 160 void SetupTrackCuts(Bool_t isESD, AliDielectron *die, Int_t cutDefinition)
 161 {
 162   //
 163   // Setup the track cuts
 164   //
 165
 166   // Quality cuts
 167   AliDielectronCutGroup* cuts = new AliDielectronCutGroup("cuts","cuts",AliDielectronCutGroup::kCompAND);
 168   die->GetTrackFilter().AddCuts(cuts);
 169
 170   //Pt cut, should make execution a bit faster
 171   AliDielectronVarCuts *pt = new AliDielectronVarCuts("Pt>.8","Pt>.8");
 172   pt->AddCut(AliDielectronVarManager::kPt,0.8,1e30);
 173   cuts->AddCut(pt);
 174
 175         // track cuts ESD and AOD
 176   AliDielectronVarCuts *varCuts = new AliDielectronVarCuts("VarCuts","VarCuts");
 177   varCuts->AddCut(AliDielectronVarManager::kImpactParXY, -1.0,   1.0);
 178   varCuts->AddCut(AliDielectronVarManager::kImpactParZ,  -3.0,   3.0);
 179   varCuts->AddCut(AliDielectronVarManager::kEta,         -0.9,   0.9);
 180   //varCuts->AddCut(AliDielectronVarManager::kTPCchi2Cl,    0.0,   4.0);
 181   varCuts->AddCut(AliDielectronVarManager::kNclsTPC,     70.0, 160.0);
 182   varCuts->AddCut(AliDielectronVarManager::kKinkIndex0,   0.0);
 183
 184    AliDielectronTrackCuts *trkCuts = new AliDielectronTrackCuts("TrkCuts","TrkCuts");
 185   if(isESD){
 186     varCuts->AddCut(AliDielectronVarManager::kTPCchi2Cl,    0.0,   4.0);
 187     switch(cutDefinition) {
 188     case kTOFTRD2: varCuts->AddCut(AliDielectronVarManager::kITSLayerFirstCls,-0.01,0.5); //ITS(0) = SPDfirst
 189       break;
 190     default:      varCuts->AddCut(AliDielectronVarManager::kITSLayerFirstCls,-0.01,1.5); //ITS(0-1) = SPDany
 191       break;
 192    }
 193   }else
 194       {
 195       switch(cutDefinition) {
 196       case kTOFTRD2: trkCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kFirst);  //ITS(0) = SPDfirst
 197       break;
 198       default:  trkCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);//ITS(0-1) = SPDany
 199       break;
 200       }
 201     }
 202
 203   cuts->AddCut(varCuts);
 204
 205   trkCuts->SetRequireITSRefit(kTRUE);
 206   trkCuts->SetRequireTPCRefit(kTRUE);
 207   cuts->AddCut(trkCuts);
 208
 209   //Do we have an MC handler?
 210   //  Bool_t hasMC=(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()!=0x0);
 211
 212   /* vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv PID CUTS vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv */
 213   AliDielectronPID *pid = new AliDielectronPID("PID","PID");
 214
 215   ////////////////////////////////// DATA
 216   if(!hasMC) {
 217     pid->AddCut(AliDielectronPID::kTPC,AliPID::kPion,-100.,3.5,0.,0.,kTRUE);
 218     pid->AddCut(AliDielectronPID::kTPC,AliPID::kProton,-100.,3.5,0.,0.,kTRUE);
 219
 220     //if(cutDefinition==kTRD || cutDefinition>=kTOFTRD || cutDefinition>=kTOFTRD2)
 221     if(cutDefinition>=kTOFTRD || cutDefinition>=kTOFTRD2)
 222       pid->AddCut(AliDielectronPID::kTRDeleEff,AliPID::kElectron,.8,1.,3.5.,6.,kFALSE,
 223                   AliDielectronPID::kIfAvailable,AliDielectronVarManager::kTRDpidQuality);
 224   }
 225
 226   ////////////////////////////////// MC
 227   if(hasMC) {
 228
 229     // electron
 230     Double_t nSigmaPi = 3.5; Double_t nSigmaP = 3.5;
 231     Double_t resolution=0.0549;
 232     Double_t BBpro[5] = {0};
 233     Double_t BBpio[5] = {0};
 234
 235     for(Int_t icent=0; icent<8; icent++) {
 236
 237       switch (icent) {
 238         case 0:  // 0-10%
 239           BBpro[0] = 0.031555;  BBpro[1] = 26.0595; BBpro[2] = 3.02422e-11;  BBpro[3] = 2.05594; BBpro[4] = 5.99848;
 240           BBpio[0] = 0.0252122; BBpio[1] = 38.8991; BBpio[2] = 4.0901e-11;   BBpio[3] = 5.27988; BBpio[4] = 4.3108;
 241           break;
 242         case 1:  // 10-20%
 243           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
 244           BBpio[0] = 0.0252127; BBpio[1] = 33.8617; BBpio[2] = 3.56866e-11;  BBpio[3] = 5.24831; BBpio[4] = 4.31093;
 245           break;
 246         case 2:  // 20-30%
 247           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
 248           BBpio[0] = 0.0263205; BBpio[1] = 37.9307; BBpio[2] = 4.29724e-11;  BBpio[3] = 5.74458; BBpio[4] = 4.32459;
 249           break;
 250         case 3:  // 30-40%
 251           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
 252           BBpio[0] = 0.026294;  BBpio[1] = 39.0346; BBpio[2] = 4.12261e-11;  BBpio[3] = 5.28808; BBpio[4] = 4.31301;
 253           break;
 254         case 4:  // 40-50%
 255           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
 256           BBpio[0] = 0.0263134; BBpio[1] = 38.2084; BBpio[2] = 3.75159e-11;  BBpio[3] = 5.78125; BBpio[4] = 4.31363;
 257           break;
 258         case 5:  // 50-60%
 259           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
 260           BBpio[0] = 0.0263134; BBpio[1] = 38.2084; BBpio[2] = 3.75159e-11;  BBpio[3] = 5.78125; BBpio[4] = 4.31363;
 261           break;
 262         case 6:  // 60-70%
 263           BBpro[0] = 0.031555;  BBpro[1] = 26.0595; BBpro[2] = 3.02422e-11;  BBpro[3] = 2.05594; BBpro[4] = 5.99848;
 264           BBpio[0] = 0.026302;  BBpio[1] = 38.6888; BBpio[2] = 3.56792e-11;  BBpio[3] = 5.2465;  BBpio[4] = 4.31094;
 265           break;
 266         case 7:  // 70-80%
 267           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
 268           BBpio[0] = 0.0263134; BBpio[1] = 38.2084; BBpio[2] = 3.75159e-11;  BBpio[3] = 5.78125; BBpio[4] = 4.31363;
 269           break;
 270         case 8:  // 80-90%
 271           BBpro[0] = 0.0313438; BBpro[1] = 25.8666; BBpro[2] = 4.5457e-11;   BBpro[3] = 2.07912; BBpro[4] = 5.99986;
 272           BBpio[0] = 0.0252127; BBpio[1] = 33.8617; BBpio[2] = 3.56866e-11;  BBpio[3] = 5.24831; BBpio[4] = 4.31093;
 273           break;
 274         case 9:  // 90-100%
 275           BBpro[0] = 0.0319126; BBpro[1] = 36.8784; BBpro[2] = 3.4274e-11;   BBpro[3] = 3.2431;  BBpro[4] = 5.93388;
 276           BBpio[0] = 0.027079;  BBpio[1] = 67.5936; BBpio[2] = 9.72548e-11;  BBpio[3] = 9.61382; BBpio[4] = 5.99372;
 277           break;
 278       }
 279
 280
 281       TF1 *ffPro=new TF1(Form("fBethe%d_c%d",AliPID::kProton,icent), Form("(%f*%f+(AliExternalTrackParam::BetheBlochAleph(x/%f,[0],[1],[2],[3],[4])-AliExternalTrackParam::BetheBlochAleph(x/%f,[0],[1],[2],[3],[4])))/%f", nSigmaP,resolution, AliPID::ParticleMass(AliPID::kProton), AliPID::ParticleMass(AliPID::kElectron), resolution), 0.05,200.);
 282
 283       TF1 *ffPio=new TF1(Form("fBethe%d_c%d",AliPID::kPion,icent), Form("(%f*%f+(AliExternalTrackParam::BetheBlochAleph(x/%f,[0],[1],[2],[3],[4])-AliExternalTrackParam::BetheBlochAleph(x/%f,[0],[1],[2],[3],[4])))/%f", nSigmaPi,resolution, AliPID::ParticleMass(AliPID::kPion), AliPID::ParticleMass(AliPID::kElectron), resolution), 0.05,200.);
 284
 285       TString list=gSystem->Getenv("LIST");
 286
 287       //LHC11a10b
 288       if (list.Contains("LHC11a10b") || list.IsNull()) {
 289         printf("LHC11a10b parameters\n");
 290         ffPro->SetParameters(BBpro[0],BBpro[1],BBpro[2],BBpro[3],BBpro[4]);
 291         ffPio->SetParameters(BBpio[0],BBpio[1],BBpio[2],BBpio[3],BBpio[4]);
 292
 293         // proton cut
 294         pid->AddCut(AliDielectronPID::kTPC,AliPID::kElectron,ffPro,10,((double)icent)*10.,((double)icent+1)*10,
 295                     kFALSE,AliDielectronPID::kRequire,AliDielectronVarManager::kCentrality);
 296         // pion cut
 297         pid->AddCut(AliDielectronPID::kTPC,AliPID::kElectron,ffPio,10,((double)icent)*10.,((double)icent+1)*10,
 298                     kFALSE,AliDielectronPID::kRequire,AliDielectronVarManager::kCentrality);
 299       }
 300     }
 301
 302     // shifts for the nSigma electrons
 303     TGraph* nSigmaCorrection = new TGraph();
 304     // LHC11a10b
 305     if (list.Contains("LHC11a10b") || list.IsNull()) {
 306       nSigmaCorrection->SetPoint(0, 137161., -0.50-(0.28));
 307       nSigmaCorrection->SetPoint(1, 139510., -0.50-(0.28));
 308       pid->SetCorrGraph(nSigmaCorrection);
 309     }
 310
 311   } //hasMC
 312
 313   ////////////////////////////////// DATA + MC
 314   // pid cuts TPC + TOF & TRD
 315   pid->AddCut(AliDielectronPID::kTPC,AliPID::kElectron,-3.,3.);
 316   if(cutDefinition==kTOFTRD || cutDefinition==kTOFTRD2)
 317     pid->AddCut(AliDielectronPID::kTOF,AliPID::kElectron,-3,3.,0.,0.,kFALSE,AliDielectronPID::kIfAvailable);
 318
 319   //if(cutDefinition!=krec)
 320   cuts->AddCut(pid);
 321   /* ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ PID CUTS ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^*/
 322
 323
 324   // exclude conversion electrons selected by the tender
 325   AliDielectronTrackCuts *noconv=new AliDielectronTrackCuts("noConv","noConv");
 326   noconv->SetV0DaughterCut(AliPID::kElectron,kTRUE);
 327   cuts->AddCut(noconv);
 328
 329 }
 330
 331 //______________________________________________________________________________________
 332 void SetupPairCuts(AliDielectron *die, Int_t cutDefinition)
 333 {
 334   //
 335   // Setup the pair cuts
 336   //
 337
 338   // conversion rejection
 339   Double_t gCut = 0.05;             // default
 340
 341   AliDielectronVarCuts *gammaCut=new AliDielectronVarCuts("gammaCut","gammaCut");
 342   gammaCut->AddCut(AliDielectronVarManager::kM,0.,gCut);
 343   die->GetPairPreFilter().AddCuts(gammaCut);
 344
 345
 346   // rapidity selection
 347   AliDielectronVarCuts *rapCut=new AliDielectronVarCuts("|Y|<.9","|Y|<.9");
 348   rapCut->AddCut(AliDielectronVarManager::kY,-0.9,0.9);
 349   die->GetPairFilter().AddCuts(rapCut);
 350
 351 }
 352
 353 //______________________________________________________________________________________
 354 void InitHistograms(AliDielectron *die, Int_t cutDefinition)
 355 {
 356   //
 357   // Initialise the histograms
 358   //
 359   //  Bool_t hasMC=(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()!=0x0);
 360
 361   //Setup histogram Manager
 362   AliDielectronHistos *histos=new AliDielectronHistos(die->GetName(),die->GetTitle());
 363
 364
 365   //add histograms to event class
 366   histos->AddClass("Event");
 367   histos->UserHistogram("Event","VtxZ","Vertex Z;z (cm)",
 368                         300,-15.,15.,
 369                         AliDielectronVarManager::kZvPrim);
 370   //histos->UserHistogram("Event","Centrality","Centrality;centrality (%)",
 371   //                      "0.,5.,10.,20.,40.,50.,60.,80.,100.",
 372   //                      AliDielectronVarManager::kCentrality);
 373   histos->UserHistogram("Event","Centrality","Centrality;centrality (%)",
 374                         20,0.,100.,AliDielectronVarManager::kCentrality);
 375
 376   histos->UserHistogram("Event","Multiplicity","Multiplicity V0;Multiplicity V0",
 377                         500,0.,25000., AliDielectronVarManager::kMultV0);
 378   histos->UserHistogram("Event","Cent_Mult","Centrality vs. Multiplicity;centrality (%);Multiplicity V0",
 379                         10,0.,100., 500,0.,25000.,AliDielectronVarManager::kCentrality,AliDielectronVarManager::kMultV0);
 380   histos->UserProfile("Event","Cent_Nacc",
 381                       "accepted tracks;centrality (%)",
 382                       AliDielectronVarManager::kNacc,
 383                       "0.,5.,10.,20.,40.,50.,60.,80.,100.",
 384                       AliDielectronVarManager::kCentrality);
 385   histos->UserProfile("Event","Cent_NVtxContrib",
 386                       "number of vertex contributors;centrality (%)",
 387                       AliDielectronVarManager::kNVtxContrib,
 388                       "0.,5.,10.,20.,40.,50.,60.,80.,100.",
 389                       AliDielectronVarManager::kCentrality);
 390
 391
 392
 393
 394   ////// MONTE CARLO //////
 395   /*
 396   if(cutDefinition == kTOFTRD && hasMC) {
 397     histos->AddClass("MCEvent");
 398     histos->UserHistogram("MCEvent","Cent_NJPsis","Centrality vs. generated incl. J/#psi per event;centrality (%);N_{J/#psi}",
 399                           10,0.,100., 21,-0.5,20.5,
 400                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kNumberOfJPsis);
 401   }
 402   */
 403
 404
 405   //Initialise histogram classes
 406   histos->SetReservedWords("Track;Pair");
 407
 408
 409   //Pair classes
 410   // to fill also mixed event histograms loop until 10
 411   for (Int_t i=0; i<3; ++i){
 412     histos->AddClass(Form("Pair_%s",AliDielectron::PairClassName(i)));
 413   }
 414
 415   if(cutDefinition <= kTOFTRD) {
 416
 417     //legs from pair
 418     for (Int_t i=0; i<3; ++i){
 419       histos->AddClass(Form("Track_Legs_%s",AliDielectron::PairClassName(i)));
 420     }
 421
 422     //Track classes
 423     //to fill also track info from 2nd event loop until 2
 424     for (Int_t i=0; i<2; ++i){
 425       histos->AddClass(Form("Track_%s",AliDielectron::TrackClassName(i)));
 426     }
 427
 428     //track rotation
 429     //   histos->AddClass(Form("Pair_%s",AliDielectron::PairClassName(AliDielectron::kEv1PMRot)));
 430     //   histos->AddClass(Form("Track_Legs_%s",AliDielectron::PairClassName(AliDielectron::kEv1PMRot)));
 431
 432     //add histograms to Track classes
 433     //histos->UserHistogram("Track","TOFbit","TOFbit;bit;#tracks",19,-9.5,9.5,AliDielectronVarManager::kTOFPIDBit);
 434     histos->UserHistogram("Track","Pt","Pt;Pt [GeV];#tracks",
 435                           400,0,20.,
 436                           AliDielectronVarManager::kPt);
 437     histos->UserHistogram("Track","TPCnCls","Number of Clusters TPC;TPC number clusteres;#tracks",
 438                           160,-0.5,159.5,
 439                           AliDielectronVarManager::kNclsTPC);
 440     histos->UserHistogram("Track","TPCsignalN","Number of Clusters TPC;TPC number clusteres;#tracks",
 441                           160,-0.5,159.5,
 442                           AliDielectronVarManager::kTPCsignalN);
 443
 444     histos->UserHistogram("Track","dXY","dXY;dXY [cm];#tracks",
 445                           500,-1.,1.,
 446                           AliDielectronVarManager::kImpactParXY);
 447     histos->UserHistogram("Track","dZ","dZ;dZ [cm];#tracks",
 448                           600,-3.,3.,
 449                           AliDielectronVarManager::kImpactParZ);
 450     histos->UserHistogram("Track","Eta_Phi","Eta Phi Map; Eta; Phi;#tracks",
 451                           200,-1,1,200,0,6.285,
 452                           AliDielectronVarManager::kEta,AliDielectronVarManager::kPhi);
 453
 454     histos->UserHistogram("Track","dEdx_P","dEdx;P [GeV];TPC signal (arb units);#tracks",
 455                           400,0.2,20.,200,0.,200.,
 456                           AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCsignal,kTRUE);
 457     histos->UserHistogram("Track","TPCnSigmaEle_P","TPC number of sigmas Electrons;P [GeV];TPC number of sigmas Electrons;#tracks",
 458                           400,0.2,20.,200,-10.,10.,
 459                           AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaEle,kTRUE);
 460
 461     histos->UserHistogram("Track","Ncl",";Number clusters TPC;Number clusters TPC",
 462                           160,-0.5,159.5,
 463                           AliDielectronVarManager::kNclsTPC);
 464     histos->UserHistogram("Track","NclFr",";Number of findable clusters (robust);Number findable clusters TPC",
 465                           160,-0.5,159.5,
 466                           AliDielectronVarManager::kNFclsTPCr);
 467     histos->UserHistogram("Track","Ncl_NclFr","Number of (findable) clusters TPC;found clusters;findable clusters",
 468                           160,-0.5,159.5,160,-0.5,159.5,
 469                           AliDielectronVarManager::kNclsTPC,AliDielectronVarManager::kNFclsTPCr);
 470     histos->UserHistogram("Track","NtrklTRD",";Number tracklets TRD for pid;Number tracklets TRD",
 471                           8,-0.5,7.5,
 472                           AliDielectronVarManager::kTRDpidQuality);
 473
 474     //add histograms to Pair classes
 475     histos->UserHistogram("Pair","InvMass","Inv.Mass;Inv. Mass [GeV];#pairs",
 476                           300,.0,300*0.04,
 477                           AliDielectronVarManager::kM); // 40MeV bins, 12GeV/c2
 478     histos->UserHistogram("Pair","Rapidity","Rapidity;Rapidity;#pairs",
 479                           100,-1.,1.,
 480                           AliDielectronVarManager::kY);
 481     histos->UserHistogram("Pair","OpeningAngle","Opening angle;angle",
 482                           100,0.,3.15,
 483                           AliDielectronVarManager::kOpeningAngle);
 484     histos->UserHistogram("Pair","Chi2NDF","#chi^{2}/NDF;#chi^{2}/NDF",
 485                           100,0.,20,
 486                           AliDielectronVarManager::kChi2NDF);
 487
 488    histos->UserHistogram("Pair","PseudoProperTime","Pseudoproper decay length; pseudoproper-decay-length[#mum];Entries/40#mum",
 489                           150,-0.3.,0.3,AliDielectronVarManager::kPseudoProperTime);
 490
 491
 492    }
 493
 494   die->SetHistogramManager(histos);
 495 }
 496
 497
 498 void InitCF(AliDielectron* die, Int_t cutDefinition)
 499 {
 500   //
 501   // Setup the CF Manager if needed
 502   //
 503   //  Bool_t hasMC=(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()!=0x0);
 504
 505   AliDielectronCF *cf=new AliDielectronCF(die->GetName(),die->GetTitle());
 506
 507   // pair variables
 508   cf->AddVariable(AliDielectronVarManager::kM,125,0.,125*.04); //40Mev Steps
 509   if(cutDefinition!=kSubRndm) cf->AddVariable(AliDielectronVarManager::kPairType,11,0,11);
 510
 511   if(cutDefinition <= kTOFTRD2) {
 512
 513     // pair and event vars
 514       cf->AddVariable(AliDielectronVarManager::kCentrality,"0.,5.,10.,20.,40.,50.,60.,80.");
 515       cf->AddVariable(AliDielectronVarManager::kPt,"0., 1., 2.5, 5., 100.0");
 516       if(!hasMC) cf->AddVariable(AliDielectronVarManager::kZvPrim,20, -10., 10.);
 517       if(hasMC) cf->AddVariable(AliDielectronVarManager::kNacc,20,0.,3000.0);
 518       if(hasMC) cf->AddVariable(AliDielectronVarManager::kNVtxContrib,20,0.,4000.);
 519       if(hasMC) cf->AddVariable(AliDielectronVarManager::kTRDpidEffPair,101,0.0,1.01);
 520       //cf->AddVariable(AliDielectronVarManager::kY,"-0.8,0.8");
 521
 522     //leg variables
 523       cf->AddVariable(AliDielectronVarManager::kPt,"0.8, 1.0, 1.1, 1.2, 1.5, 100.0",kTRUE);
 524       cf->AddVariable(AliDielectronVarManager::kTPCnSigmaEle,"-3,-2.5,-2,2,2.5,3",kTRUE);
 525       //cf->AddVariable(AliDielectronVarManager::kTPCnSigmaPio,"3.5,4.0,4.5,5.0,100",kTRUE);
 526       //    cf->AddVariable(AliDielectronVarManager::kTPCnSigmaPro,"3.5,4.0,4.5,5.0,100",kTRUE);
 527
 528       // standard vars
 529         if(hasMC) cf->AddVariable(AliDielectronVarManager::kEta,"-0.9,0.9",kTRUE);
 530         //cf->AddVariable(AliDielectronVarManager::kNclsTPC,"70, 90, 100, 120, 160",kTRUE);
 531     }
 532
 533
 534
 535   if(hasMC) {
 536     cf->AddVariable(AliDielectronVarManager::kRunNumber, GetRunNumbers() ); // LHC10h only
 537     //    cf->AddVariable(AliDielectronVarManager::kRunNumber, 170593-136831, 136831, 170593); // LHC10h -> LHC11h
 538     if(cutDefinition==kTOFTRD || cutDefinition==kTOFTRD2) cf->SetStepForMCtruth();
 539      cf->SetStepsForMCtruthOnly();
 540            // cf->SetStepsForBackground();
 541   }
 542
 543   die->SetCFManagerPair(cf);
 544 }
 545
 546 void AddMCSignals(AliDielectron *die){
 547   //Do we have an MC handler?
 548   //Bool_t hasMC=(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()!=0x0);
 549   if (!hasMC) return;
 550
 551   AliDielectronSignalMC* inclusiveJpsi = new AliDielectronSignalMC("inclusiveJpsi","Inclusive J/psi");
 552   inclusiveJpsi->SetLegPDGs(11,-11);
 553   inclusiveJpsi->SetMotherPDGs(443,443);
 554   inclusiveJpsi->SetMothersRelation(AliDielectronSignalMC::kSame);
 555   inclusiveJpsi->SetFillPureMCStep(kTRUE);
 556   inclusiveJpsi->SetCheckBothChargesLegs(kTRUE,kTRUE);
 557   inclusiveJpsi->SetCheckBothChargesMothers(kTRUE,kTRUE);
 558   die->AddSignalMC(inclusiveJpsi);
 559
 560   AliDielectronSignalMC* promptJpsi = new AliDielectronSignalMC("promptJpsi","Prompt J/psi");   // prompt J/psi (not from beauty decays)
 561   promptJpsi->SetLegPDGs(11,-11);
 562   promptJpsi->SetMotherPDGs(443,443);
 563   promptJpsi->SetGrandMotherPDGs(503,503,kTRUE,kTRUE);   // not from beauty hadrons
 564   promptJpsi->SetMothersRelation(AliDielectronSignalMC::kSame);
 565   promptJpsi->SetFillPureMCStep(kTRUE);
 566   promptJpsi->SetLegSources(AliDielectronSignalMC::kFinalState, AliDielectronSignalMC::kFinalState);
 567   promptJpsi->SetCheckBothChargesLegs(kTRUE,kTRUE);
 568   promptJpsi->SetCheckBothChargesMothers(kTRUE,kTRUE);
 569   promptJpsi->SetCheckBothChargesGrandMothers(kTRUE,kTRUE);
 570   die->AddSignalMC(promptJpsi);
 571
 572   AliDielectronSignalMC* beautyJpsi = new AliDielectronSignalMC("beautyJpsi","Beauty J/psi");
 573   beautyJpsi->SetLegPDGs(11,-11);
 574   beautyJpsi->SetMotherPDGs(443,443);
 575   beautyJpsi->SetMothersRelation(AliDielectronSignalMC::kSame);
 576   beautyJpsi->SetGrandMotherPDGs(500,500);
 577   beautyJpsi->SetFillPureMCStep(kTRUE);
 578   beautyJpsi->SetCheckBothChargesLegs(kTRUE,kTRUE);
 579   beautyJpsi->SetCheckBothChargesMothers(kTRUE,kTRUE);
 580   beautyJpsi->SetCheckBothChargesGrandMothers(kTRUE,kTRUE);
 581   die->AddSignalMC(beautyJpsi);
 582
 583   AliDielectronSignalMC* directJpsi = new AliDielectronSignalMC("directJpsi","Direct J/psi");   // embedded J/psi
 584   directJpsi->SetLegPDGs(11,-11);
 585   directJpsi->SetMotherPDGs(443,443);
 586   directJpsi->SetMothersRelation(AliDielectronSignalMC::kSame);
 587   directJpsi->SetFillPureMCStep(kTRUE);
 588   directJpsi->SetLegSources(AliDielectronSignalMC::kFinalState, AliDielectronSignalMC::kFinalState);
 589   directJpsi->SetMotherSources(AliDielectronSignalMC::kDirect, AliDielectronSignalMC::kDirect);
 590   directJpsi->SetCheckBothChargesLegs(kTRUE,kTRUE);
 591   directJpsi->SetCheckBothChargesMothers(kTRUE,kTRUE);
 592   die->AddSignalMC(directJpsi);
 593 }
 594
 595 void SetEtaCorrection()
 596 {
 597   if (AliDielectronPID::GetEtaCorrFunction()) return;
 598
 599   TString list=gSystem->Getenv("LIST");
 600
 601   TString etaMap="$TRAIN_ROOT/jpsi_JPSI/EtaCorrMaps.root";
 602   TString trainRoot=gSystem->Getenv("TRAIN_ROOT");
 603   if (trainRoot.IsNull()) etaMap="$ALICE_ROOT/PWGDQ/dielectron/files/EtaCorrMaps.root";
 604   if (gSystem->AccessPathName(gSystem->ExpandPathName(etaMap.Data()))){
 605     Error("ConfigPbPb","Eta map not found: %s",etaMap.Data());
 606     return;
 607   }
 608
 609   TFile f(etaMap.Data());
 610   if (!f.IsOpen()) return;
 611   TList *keys=f.GetListOfKeys();
 612
 613   for (Int_t i=0; i<keys->GetEntries(); ++i){
 614     TString kName=keys->At(i)->GetName();
 615     TPRegexp reg(kName);
 616     if (reg.MatchB(list)){
 617       printf("Using Eta Correction Function: %s\n",kName.Data());
 618       AliDielectronPID::SetEtaCorrFunction((TF1*)f.Get(kName.Data()));
 619     }
 620   }
 621 }
 622
 623 TVectorD *GetRunNumbers() {
 624
 625   Double_t runLHC10h[] = { // all runs
 626     136851, 136854, 136879, 137042, 137045, 137124, 137125, 137132, 137133, 137135, 137136, 137137, 137161, 137162, 137163, 137165, 137230, 137231, 137232, 137235, 137236, 137243, 137365, 137366, 137370, 137430, 137431, 137432, 137434, 137439, 137440, 137441, 137443, 137530, 137531, 137539, 137541, 137544, 137546, 137549, 137595, 137608, 137609, 137638, 137639, 137685, 137686, 137689, 137691, 137692, 137693, 137704, 137718, 137722, 137724, 137748, 137751, 137752, 137843, 137844, 137847, 137848, 138125, 138126, 138150, 138151, 138153, 138154, 138190, 138192, 138197, 138200, 138201, 138225, 138275, 138359, 138364, 138396, 138438, 138439, 138442, 138469, 138533, 138534, 138578, 138579, 138582, 138583, 138620, 138621, 138624, 138637, 138638, 138652, 138653, 138662, 138666, 138730, 138731, 138732, 138736, 138737, 138740, 138742, 138795, 138796, 138826, 138828, 138830, 138831, 138836, 138837, 138870, 138871, 138872, 138924, 138965, 138972, 138973, 138976, 138977, 138978, 138979, 138980, 138982, 138983, 139024, 139025, 139028, 139029, 139030, 139031, 139034, 139036, 139037, 139038, 139042, 139104, 139105, 139107, 139110, 139172, 139173, 139308, 139309, 139310, 139311, 139314, 139316, 139328, 139329, 139360, 139437, 139438, 139439, 139440, 139441, 139465, 139466, 139467, 139470, 139471, 139503, 139504, 139505, 139507, 139510, 139511, 139513, 139514, 139517,
 627 0
 628   };
 629
 630
 631   Int_t sizeLHC10h = (int) (sizeof(runLHC10h)/sizeof(Double_t));
 632   runLHC10h[sizeLHC10h-1] =   runLHC10h[sizeLHC10h-2] + 1.;
 633   TVectorD *vecLHC10h = new TVectorD(sizeLHC10h, runLHC10h);
 634
 635   return vecLHC10h;
 636
 637 }