- updates for AN67
[u/mrichter/AliRoot.git] / PWGDQ / dielectron / macrosJPSI / ConfigJpsi_jb_PbPb.C
1 void InitHistograms(AliDielectron *die, Int_t cutDefinition);
2 void InitCF(AliDielectron* die, Int_t cutDefinition);
3 void InitHF(AliDielectron* die, Int_t cutDefinition);
4
5 void SetupTrackCuts(AliDielectron *die, Int_t cutDefinition);
6 void SetupPairCuts(AliDielectron *die, Int_t cutDefinition);
7
8 void AddMCSignals(AliDielectron *die);
9 void SetEtaCorrection();
10 TVectorD *GetRunNumbers();
11
12 TString names=("TPC;TOF;TRD;TOFTRD;TOFTRD2D;rec;Gam0;Gam01;Gam05;Gam10;Gam15;Gam20;EtaGap01;EtaGap02;EtaGap03;EtaGap04;EtaGap05;SubLS;SubRndm");
13 enum { kTPC=0, kTOF, kTRD, kTOFTRD, kTOFTRD2D, krec,kGam0, kGam01, kGam05, kGam10, kGam15, kGam20, kEtaGap01, kEtaGap02, kEtaGap03, kEtaGap04, kEtaGap05, kSubLS, kSubRndm };
14
15 TObjArray *arrNames=names.Tokenize(";");
16 const Int_t nDie=arrNames->GetEntries();
17
18 Bool_t  isESD = kTRUE;
19 Bool_t  hasMC = kFALSE;
20 TString list  = gSystem->Getenv("LIST");
21
22
23 AliDielectron* ConfigJpsi_jb_PbPb(Int_t cutDefinition, TString prod="")
24 {
25   //
26   // Setup the instance of AliDielectron
27   //
28
29   // find mc or not?
30   if( list.IsNull()) list=prod;
31   if( list.Contains("LHC10h")   || list.Contains("LHC11h")   ) hasMC=kFALSE;
32   if( list.Contains("LHC11a10") || list.Contains("LHC12a17") ) hasMC=kTRUE;
33
34   //ESD handler?
35   isESD=(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->IsA()==AliESDInputHandler::Class());
36
37   // switch off some configurations
38   if(hasMC) { // MONTE CARLO
39     switch(cutDefinition) {
40     case kTPC:       return 0x0;
41     case kTOF:       return 0x0;
42     case kTRD:       return 0x0;
43     case kTOFTRD:    return 0x0;
44     case kTOFTRD2D:  return 0x0;
45     case krec:       return 0x0;
46     case kGam0:      return 0x0;
47     case kGam01:     return 0x0;
48     case kGam05:     return 0x0;
49     case kGam10:     return 0x0;
50     case kGam15:     return 0x0;
51     case kGam20:     return 0x0;
52     case kEtaGap01:  return 0x0;
53     case kEtaGap02:  return 0x0;
54     case kEtaGap03:  return 0x0;
55     case kEtaGap04:  return 0x0;
56     case kEtaGap05:  return 0x0;
57     case kSubLS:     return 0x0;
58     case kSubRndm:   return 0x0;
59     case kCutStats:  return 0x0;
60     }
61   } else { // COLLISION DATA
62     switch(cutDefinition) {
63     case kTPC:       return 0x0;
64     //            case kTOF:       return 0x0;
65     case kTRD:       return 0x0;
66     //            case kTOFTRD:    return 0x0;
67     //            case kTOFTRD2D:  return 0x0;
68     case krec:       return 0x0;
69     case kGam0:      return 0x0;
70     case kGam01:     return 0x0;
71     case kGam05:     return 0x0;
72     case kGam10:     return 0x0;
73     case kGam15:     return 0x0;
74     case kGam20:     return 0x0;
75     case kEtaGap01:  return 0x0;
76     case kEtaGap02:  return 0x0;
77     case kEtaGap03:  return 0x0;
78     case kEtaGap04:  return 0x0;
79     case kEtaGap05:  return 0x0;
80     case kSubLS:     return 0x0;
81     case kSubRndm:   return 0x0;
82     }
83   }
84
85   // create the actual framework object
86   TString name=Form("%02d",cutDefinition);
87   if (cutDefinition<arrNames->GetEntriesFast()){
88     name=arrNames->At(cutDefinition)->GetName();
89   }
90   AliDielectron *die = new AliDielectron(Form("%s",name.Data()), Form("Track cuts: %s",name.Data()));
91   die->SetHasMC(hasMC);
92
93   printf(" Add %s %s config %s for %s \n",(isESD?"ESD":"AOD"),(hasMC?"MC":""),name.Data(),list.Data());
94
95   // Monte Carlo Signals and TRD efficiency tables
96   if(hasMC) {
97     AddMCSignals(die);
98     printf(" Add %d MC signals \n",die->GetMCSignals()->GetEntriesFast());
99
100     // trd tables
101     /* if (list.Contains("LHC11a") ) {
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
113   // cut setup
114   SetupTrackCuts(die,cutDefinition);
115   SetupPairCuts(die,cutDefinition);
116
117   // histogram setup
118   InitHistograms(die,cutDefinition);
119   printf(" Add %d types and %03d histos to the manager \n",die->GetHistogramList()->GetEntries(),
120          0/*die->GetHistoManager()->GetList()->GetEntries()*/);
121
122   // CF container setup, switched off
123   if(cutDefinition <  kEtaGap01 ||
124      cutDefinition == kSubRndm  ) {
125     if(0) InitCF(die,cutDefinition);
126     if(0) printf(" Add %d pair, %d leg vars, %p steps and %p bins to the container \n",
127                  die->GetCFManagerPair()->GetNvarsPair(),die->GetCFManagerPair()->GetNvarsLeg(),
128                  die->GetCFManagerPair()->GetContainer(), die->GetCFManagerPair()->GetContainer() );
129   }
130
131   // HF arrays setup
132   //  if(cutDefinition <  kEtaGap01 ) {
133     InitHF(die,cutDefinition);
134     //}
135
136
137   // bgrd estimators
138   if(!hasMC) {
139       printf(" Activate bgrd estimators: ME and ROT \n");
140
141       // rotations
142       AliDielectronTrackRotator *rot=new AliDielectronTrackRotator;
143       rot->SetIterations(10);
144       rot->SetConeAnglePhi(TMath::Pi());
145       rot->SetStartAnglePhi(TMath::Pi());
146       die->SetTrackRotator(rot);
147       // mixing
148       AliDielectronMixingHandler *mix=new AliDielectronMixingHandler;
149       mix->AddVariable(AliDielectronVarManager::kZvPrim,     30,-10.,10.);
150       //      mix->AddVariable(AliDielectronVarManager::kCentrality,  8,  0.,80.);
151       mix->SetMixType(AliDielectronMixingHandler::kAll);
152       mix->SetDepth(120);
153       die->SetMixingHandler(mix);
154
155
156       // TPC event plane configurations
157       Double_t gGap;
158       switch(cutDefinition) {
159       case kEtaGap01:   gGap=0.1;   break;
160       case kEtaGap02:   gGap=0.2;   break;
161       case kEtaGap03:   gGap=0.3;   break;
162       case kEtaGap04:   gGap=0.4;   break;
163       case kEtaGap05:   gGap=0.5;   break;
164       default: gGap=0.0;
165       }
166
167     AliDielectronVarCuts *poi = new AliDielectronVarCuts("PoI","PoI");
168     poi->AddCut(AliDielectronVarManager::kM,2.92,3.20);     // particles of interest, jpsi mass window
169     die->GetEventPlanePOIPreFilter().AddCuts(poi);
170
171     if(cutDefinition >= kEtaGap01 &&
172        cutDefinition <  kSubLS     ) {
173       AliDielectronVarCuts *etaGap = new AliDielectronVarCuts(AliDielectronVarManager::GetValueName(AliDielectronVarManager::kEta),"etaGap");
174       etaGap->AddCut(AliDielectronVarManager::kEta,-1*gGap,gGap,kTRUE);
175       die->GetEventPlanePreFilter().AddCuts(etaGap);
176     }
177
178     if(cutDefinition==kSubLS) die->SetLikeSignSubEvents();
179     die->SetPreFilterEventPlane();
180   }
181
182   // prefilter settings
183   //    die->SetPreFilterAllSigns();
184   die->SetPreFilterUnlikeOnly();
185
186   // setup eta correction
187   //  if(isESD && list.Contains("LHC10h")) SetEtaCorrection();
188
189   // VZERO calibration
190   TString trainRoot=gSystem->Getenv("TRAIN_ROOT");
191   if (!trainRoot.IsNull()) {
192     die->SetVZEROCalibrationFilename("$TRAIN_ROOT/util/dielectron/dielectron/VzeroCalibrationLHC10h.root");
193     die->SetVZERORecenteringFilename("$TRAIN_ROOT/util/dielectron/dielectron/VzeroRecenteringLHC10h.root");
194   }
195
196   return die;
197 }
198
199 //______________________________________________________________________________________
200 void SetupTrackCuts(AliDielectron *die, Int_t cutDefinition)
201 {
202   //
203   // Setup the track cuts
204   //
205
206   // Quality cuts
207   AliDielectronCutGroup* cuts = new AliDielectronCutGroup("cuts","cuts",AliDielectronCutGroup::kCompAND);
208   die->GetTrackFilter().AddCuts(cuts);
209
210   // AOD track filter (needs to be first cut)
211   AliDielectronTrackCuts *trkFilter = new AliDielectronTrackCuts("TrkFilter","TrkFilter");
212   trkFilter->SetAODFilterBit(AliDielectronTrackCuts::kTPCqualSPDanyPIDele);
213   //  if(!isESD) cuts->AddCut(trkFilter);
214
215   //Pt cut, should make execution a bit faster
216   AliDielectronVarCuts *pt = new AliDielectronVarCuts("PtCut","PtCut");
217   pt->AddCut(AliDielectronVarManager::kPt,1.1,1e30);    //0.8
218   cuts->AddCut(pt);
219
220   // track cuts ESD and AOD
221   AliDielectronVarCuts *varCuts = new AliDielectronVarCuts("VarCuts","VarCuts");
222   varCuts->AddCut(AliDielectronVarManager::kImpactParXY, -1.0,   1.0);
223   varCuts->AddCut(AliDielectronVarManager::kImpactParZ,  -3.0,   3.0);
224   varCuts->AddCut(AliDielectronVarManager::kEta,         -0.9,   0.9);
225   varCuts->AddCut(AliDielectronVarManager::kTPCchi2Cl,    0.0,   4.0);
226   varCuts->AddCut(AliDielectronVarManager::kNclsTPC,     70.0, 160.0);
227   varCuts->AddCut(AliDielectronVarManager::kKinkIndex0,   0.0);
228   if(cutDefinition==kTOF || cutDefinition==kTOFTRD || cutDefinition==kTOFTRD2D)
229     varCuts->AddCut(AliDielectronVarManager::kTOFbeta,      0.2,   0.9, kTRUE);
230   cuts->AddCut(varCuts);
231
232   AliDielectronTrackCuts *trkCuts = new AliDielectronTrackCuts("TrkCuts","TrkCuts");
233   // trkCuts->SetClusterRequirementITS(AliDielectronTrackCuts::kSPD,AliDielectronTrackCuts::kAny);
234   // trkCuts->SetClusterRequirementITS(AliDielectronTrackCuts::kSPD,AliDielectronTrackCuts::kFirst);
235   trkCuts->SetITSclusterCut(AliDielectronTrackCuts::kOneOf, 15); // ITS-4 = 1+2+4+8
236   // trkCuts->SetMaxWaivedITSNcls(1);
237   trkCuts->SetRequireITSRefit(kTRUE);
238   trkCuts->SetRequireTPCRefit(kTRUE);
239   cuts->AddCut(trkCuts);
240
241   /* vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv PID CUTS vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv */
242   AliDielectronPID *pid = new AliDielectronPID("PID","PID");
243
244   ////////////////////////////////// DATA
245   if(!hasMC) {
246     pid->AddCut(AliDielectronPID::kTPC,AliPID::kPion,-100.,4.0,0.,0.,kTRUE);
247     pid->AddCut(AliDielectronPID::kTPC,AliPID::kProton,-100.,3.5,0.,0.,kTRUE);
248
249     if(cutDefinition==kTRD || cutDefinition==kTOFTRD)
250       pid->AddCut(AliDielectronPID::kTRDeleEff,AliPID::kElectron,.8,1.,3.5.,6.,kFALSE,
251                   AliDielectronPID::kIfAvailable,AliDielectronVarManager::kTRDpidQuality);
252
253     if(cutDefinition==kTOFTRD2D)
254       pid->AddCut(AliDielectronPID::kTRDeleEff2D,AliPID::kElectron,.8,1.,3.5.,6.,kFALSE,
255                   AliDielectronPID::kIfAvailable,AliDielectronVarManager::kTRDpidQuality);
256   }
257
258   ////////////////////////////////// MC
259   if(hasMC) {
260
261     // electron
262     Double_t nSigmaPi = 3.5; Double_t nSigmaP = 3.5;
263     Double_t resolution=0.0549;
264     Double_t BBpro[5] = {0};
265     Double_t BBpio[5] = {0};
266
267     for(Int_t icent=0; icent<8; icent++) {
268
269       switch (icent) {
270         case 0:  // 0-10%
271           BBpro[0] = 0.031555;  BBpro[1] = 26.0595; BBpro[2] = 3.02422e-11;  BBpro[3] = 2.05594; BBpro[4] = 5.99848;
272           BBpio[0] = 0.0252122; BBpio[1] = 38.8991; BBpio[2] = 4.0901e-11;   BBpio[3] = 5.27988; BBpio[4] = 4.3108;
273           break;
274         case 1:  // 10-20%
275           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
276           BBpio[0] = 0.0252127; BBpio[1] = 33.8617; BBpio[2] = 3.56866e-11;  BBpio[3] = 5.24831; BBpio[4] = 4.31093;
277           break;
278         case 2:  // 20-30%
279           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
280           BBpio[0] = 0.0263205; BBpio[1] = 37.9307; BBpio[2] = 4.29724e-11;  BBpio[3] = 5.74458; BBpio[4] = 4.32459;
281           break;
282         case 3:  // 30-40%
283           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
284           BBpio[0] = 0.026294;  BBpio[1] = 39.0346; BBpio[2] = 4.12261e-11;  BBpio[3] = 5.28808; BBpio[4] = 4.31301;
285           break;
286         case 4:  // 40-50%
287           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
288           BBpio[0] = 0.0263134; BBpio[1] = 38.2084; BBpio[2] = 3.75159e-11;  BBpio[3] = 5.78125; BBpio[4] = 4.31363;
289           break;
290         case 5:  // 50-60%
291           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
292           BBpio[0] = 0.0263134; BBpio[1] = 38.2084; BBpio[2] = 3.75159e-11;  BBpio[3] = 5.78125; BBpio[4] = 4.31363;
293           break;
294         case 6:  // 60-70%
295           BBpro[0] = 0.031555;  BBpro[1] = 26.0595; BBpro[2] = 3.02422e-11;  BBpro[3] = 2.05594; BBpro[4] = 5.99848;
296           BBpio[0] = 0.026302;  BBpio[1] = 38.6888; BBpio[2] = 3.56792e-11;  BBpio[3] = 5.2465;  BBpio[4] = 4.31094;
297           break;
298         case 7:  // 70-80%
299           BBpro[0] = 0.0315171; BBpro[1] = 25.8656; BBpro[2] = 3.03896e-11;  BBpro[3] = 2.05802; BBpro[4] = 5.99999;
300           BBpio[0] = 0.0263134; BBpio[1] = 38.2084; BBpio[2] = 3.75159e-11;  BBpio[3] = 5.78125; BBpio[4] = 4.31363;
301           break;
302         case 8:  // 80-90%
303           BBpro[0] = 0.0313438; BBpro[1] = 25.8666; BBpro[2] = 4.5457e-11;   BBpro[3] = 2.07912; BBpro[4] = 5.99986;
304           BBpio[0] = 0.0252127; BBpio[1] = 33.8617; BBpio[2] = 3.56866e-11;  BBpio[3] = 5.24831; BBpio[4] = 4.31093;
305           break;
306         case 9:  // 90-100%
307           BBpro[0] = 0.0319126; BBpro[1] = 36.8784; BBpro[2] = 3.4274e-11;   BBpro[3] = 3.2431;  BBpro[4] = 5.93388;
308           BBpio[0] = 0.027079;  BBpio[1] = 67.5936; BBpio[2] = 9.72548e-11;  BBpio[3] = 9.61382; BBpio[4] = 5.99372;
309           break;
310       }
311
312
313       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.);
314
315       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.);
316
317       //LHC11a10b
318       if (list.Contains("LHC11a10b")) {
319         printf(" LHC11a10b parameters\n");
320         ffPro->SetParameters(BBpro[0],BBpro[1],BBpro[2],BBpro[3],BBpro[4]);
321         ffPio->SetParameters(BBpio[0],BBpio[1],BBpio[2],BBpio[3],BBpio[4]);
322
323         // proton cut
324         pid->AddCut(AliDielectronPID::kTPC,AliPID::kElectron,ffPro,10,((double)icent)*10.,((double)icent+1)*10,
325                     kFALSE,AliDielectronPID::kRequire,AliDielectronVarManager::kCentrality);
326         // pion cut
327         pid->AddCut(AliDielectronPID::kTPC,AliPID::kElectron,ffPio,10,((double)icent)*10.,((double)icent+1)*10,
328                     kFALSE,AliDielectronPID::kRequire,AliDielectronVarManager::kCentrality);
329       }
330     }
331
332     // shifts for the nSigma electrons
333     TGraph* nSigmaCorrection = new TGraph();
334     // LHC11a10b
335     if (list.Contains("LHC11a10b")) {
336       nSigmaCorrection->SetPoint(0, 137161., -0.50-(0.28));
337       nSigmaCorrection->SetPoint(1, 139510., -0.50-(0.28));
338       pid->SetCorrGraph(nSigmaCorrection);
339     }
340
341   } //hasMC
342
343   ////////////////////////////////// DATA + MC
344   // pid cuts TPC + TOF & TRD
345   pid->AddCut(AliDielectronPID::kTPC,AliPID::kElectron,-4.,4.);
346   
347   if(cutDefinition==kTOF || cutDefinition==kTOFTRD || cutDefinition==kTOFTRD2D)
348     pid->AddCut(AliDielectronPID::kTOF,AliPID::kElectron,-3,3.,0.,0.,kFALSE,AliDielectronPID::kIfAvailable);
349
350   if(cutDefinition!=krec) cuts->AddCut(pid);
351   /* ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ PID CUTS ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^*/
352
353
354   // exclude conversion electrons selected by the tender
355   AliDielectronTrackCuts *noconv=new AliDielectronTrackCuts("noConv","noConv");
356   noconv->SetV0DaughterCut(AliPID::kElectron,kTRUE);
357   cuts->AddCut(noconv);
358
359 }
360
361 //______________________________________________________________________________________
362 void SetupPairCuts(AliDielectron *die, Int_t cutDefinition)
363 {
364   //
365   // Setup the pair cuts
366   //
367
368   // conversion rejection
369   Double_t gCut;
370   switch(cutDefinition) {
371     case kTPC:    gCut=0.05;  break;
372     case krec:    gCut=0.05;  break;
373     case kGam10:  gCut=0.1;   break;
374     case kGam15:  gCut=0.15;  break;
375     case kGam20:  gCut=0.2;   break;
376     case kGam05:  gCut=0.05;  break;
377     case kGam01:  gCut=0.01;  break;
378     case kGam0:   gCut=0.0;   break;
379     default: gCut=0.05;             // default
380   }
381
382   AliDielectronVarCuts *gammaCuts = new AliDielectronVarCuts("GammaCuts","GammaCuts");
383 //  gammaCuts->AddCut(AliDielectronVarManager::kOpeningAngle, 0.0,   0.1,  kTRUE);
384 //  gammaCuts->AddCut(AliDielectronVarManager::kLegDist,      0.0,   0.25, kTRUE);
385 //  gammaCuts->AddCut(AliDielectronVarManager::kR,            3.0,   90.0, kTRUE);
386 //  gammaCuts->AddCut(AliDielectronVarManager::kPsiPair,      0.0,   0.05, kTRUE);
387 //  gammaCuts->AddCut(AliDielectronVarManager::kChi2NDF,      0.0,   10.0, kTRUE);
388   gammaCuts->AddCut(AliDielectronVarManager::kM,            0.0,   gCut);
389   die->GetPairPreFilter().AddCuts(gammaCuts);
390
391
392   // rapidity selection
393   //  AliDielectronVarCuts *rapCut=new AliDielectronVarCuts("|Y|<.9","|Y|<.9");
394   // rapCut->AddCut(AliDielectronVarManager::kY,-0.9,0.9);
395   // die->GetPairFilter().AddCuts(rapCut);
396
397 }
398
399 //______________________________________________________________________________________
400 void InitHistograms(AliDielectron *die, Int_t cutDefinition)
401 {
402   //
403   // Initialise the histograms
404   //
405
406   // booleans for histo selection
407   Bool_t bHistEvts = kFALSE, bHistPair = kFALSE, bHistPairME = kFALSE, bHistFlow = kFALSE, bHistFlowQA=kFALSE;
408   switch (cutDefinition) {
409   case kTPC:
410   case kTOF:      bHistEvts=kTRUE; bHistFlow=kTRUE; bHistPair=kTRUE; bHistPairME=kTRUE; break;
411   case kTRD:
412   case kTOFTRD:
413   case kTOFTRD2D: bHistFlow=kTRUE; bHistPair=kTRUE; break;
414   case krec:
415   case kITScls:
416   case kITSamy:
417   case kDCA:
418   case kChi:      break;
419   case kEtaGap01:
420   case kEtaGap02:
421   case kEtaGap03:
422   case kEtaGap04:
423   case kEtaGap05:
424   case kSubRndm:  bHistFlow=kTRUE; bHistFlowQA=kTRUE; break;
425   case kSubLS:    bHistFlow=kTRUE; bHistFlowQA=kTRUE; bHistEvts=kTRUE; break;
426   case kGam0:
427   case kGam01:
428   case kGam05:
429   case kGam10:
430   case kGam15:
431   case kGam20:    break;
432   }
433
434
435   //Setup histogram Manager
436   AliDielectronHistos *histos=new AliDielectronHistos(die->GetName(),die->GetTitle());
437
438   //add histograms to event class
439   histos->AddClass("Event");
440   histos->UserHistogram("Event","Centrality","Centrality;centrality (%);events",
441                         "0.,5.,10.,20.,40.,50.,60.,80.",
442                         AliDielectronVarManager::kCentrality);
443
444   ////// EVENT HISTOS /////
445   if(bHistEvts) {
446     histos->UserHistogram("Event","RunNumber","Events per run;run;events",
447                           GetRunNumbers(),
448                           AliDielectronVarManager::kRunNumber);
449     histos->UserHistogram("Event","VtxZ","Vertex Z;z (cm)", 300,-15.,15.,
450                           AliDielectronVarManager::kZvPrim);
451     histos->UserHistogram("Event","Multiplicity","Multiplicity V0;Multiplicity V0;events", 500,0.,25000.,
452                           AliDielectronVarManager::kMultV0);
453     histos->UserHistogram("Event","Cent_Mult","Centrality vs. Multiplicity;centrality (%);Multiplicity V0",
454                           10,0.,100., 500,0.,25000.,
455                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kMultV0);
456     histos->UserProfile("Event","Cent_Nacc", "accepted tracks;centrality (%)",
457                         AliDielectronVarManager::kNacc,
458                         "0.,5.,10.,20.,40.,50.,60.,80.,100.",
459                         AliDielectronVarManager::kCentrality);
460     histos->UserProfile("Event","Cent_NVtxContrib", "number of vertex contributors;centrality (%)",
461                         AliDielectronVarManager::kNVtxContrib,
462                       "0.,5.,10.,20.,40.,50.,60.,80.,100.",
463                         AliDielectronVarManager::kCentrality);
464   } //hist: event
465
466
467   ////// FLOW //////
468   if(bHistFlow) {
469
470     // RP angles versus centrality
471     histos->UserHistogram("Event","Cent_TPCrpH2","TPC RP;centrality (%);#Psi^{TPC} (rad.)",
472                           10,0.,100.,100,-2.,2.,
473                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kTPCrpH2);
474     histos->UserHistogram("Event","Cent_TPCsub1rpH2","TPC-1 RP;centrality (%);#Psi^{sub1} (rad.)",
475                           10,0.,100.,100,-2.,2.,
476                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kTPCsub1rpH2);
477     histos->UserHistogram("Event","Cent_TPCsub2rpH2","TPC-2 RP;centrality (%);#Psi^{sub2} (rad.)",
478                           10,0.,100.,100,-2.,2.,
479                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kTPCsub2rpH2);
480
481     histos->UserHistogram("Event","Cent_v0ACrpH2","VZERO-AC RP;centrality (%);#Psi_{2}^{V0AC} (rad.)",
482                           10,0.,100.,100,-2.0,2.0,
483                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kv0ACrpH2);
484     histos->UserHistogram("Event","Cent_v0ArpH2","VZERO-A RP;centrality (%);#Psi_{2}^{V0A} (rad.)",
485                           10,0.,100.,100,-2.0,2.0,
486                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kv0ArpH2);
487     histos->UserHistogram("Event","Cent_v0CrpH2","VZERO-C RP;centrality (%);#Psi_{2}^{V0C} (rad.)",
488                           10,0.,100.,100,-2.0,2.0,
489                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kv0CrpH2);
490   } // hist: flow
491
492   if(bHistFlowQA) {
493     // TPC event plane
494     histos->UserHistogram("Event","TPCxH2","TPC Qx component;TPCxH2",
495                           100,-1500.,1500.,
496                           AliDielectronVarManager::kTPCxH2);
497     histos->UserHistogram("Event","TPCyH2","TPC Qy component;TPCyH2",
498                           100,-1500.,1500.,
499                           AliDielectronVarManager::kTPCyH2);
500     histos->UserHistogram("Event","TPCrpH2","TPC reaction plane; #Psi^{TPC}",
501                           100,-2.,2.,
502                           AliDielectronVarManager::kTPCrpH2);
503     histos->UserHistogram("Event","TPCsub1xH2","TPC Qx component sub1;TPCsub1xH2",
504                           100,-1500.,1500.,
505                           AliDielectronVarManager::kTPCsub1xH2);
506     histos->UserHistogram("Event","TPCsub1yH2","TPC Qy component sub1;TPCsub1yH2",
507                           100,-1500.,1500.,
508                           AliDielectronVarManager::kTPCsub1yH2);
509     histos->UserHistogram("Event","TPCsub1rpH2","TPC reaction plane sub1; #Psi^{sub1}",
510                           100,-2.,2.,
511                           AliDielectronVarManager::kTPCsub1rpH2);
512     histos->UserHistogram("Event","TPCsub2xH2","TPC Qx component sub2;TPCsub2xH2",
513                           100,-1500.,1500.,
514                           AliDielectronVarManager::kTPCsub2xH2);
515     histos->UserHistogram("Event","TPCsub2yH2","TPC Qy component sub2;TPCsub2yH2",
516                           100,-1500.,1500.,
517                           AliDielectronVarManager::kTPCsub2yH2);
518     histos->UserHistogram("Event","TPCsub2rpH2","TPC reaction plane sub2; #Psi^{sub2}",
519                           100,-2.,2.,
520                           AliDielectronVarManager::kTPCsub2rpH2);
521     histos->UserHistogram("Event","TPCsub12DiffH2","TPC reaction plane diff; cos(2(#Psi^{sub1}-#Psi^{sub2}))",
522                           100,-1.,1.,
523                           AliDielectronVarManager::kTPCsub12DiffH2);
524     /* // uncorrected eventplane
525      histos->UserHistogram("Event","TPCxH2uc","TPC Qx component;TPCxH2uc",
526      100,-1500.,1500.,
527      AliDielectronVarManager::kTPCxH2uc);
528      histos->UserHistogram("Event","TPCyH2uc","TPC Qy component;TPCyH2uc",
529      100,-1500.,1500.,
530      AliDielectronVarManager::kTPCyH2uc);
531      histos->UserHistogram("Event","TPCrpH2uc","TPC reaction plane;TPCrpH2uc",
532      100,-2.,2.,
533      AliDielectronVarManager::kTPCrpH2uc);
534      histos->UserHistogram("Event","TPCsub1xH2uc","TPC Qx component sub1;TPCsub1xH2uc",
535      100,-1500.,1500.,
536      AliDielectronVarManager::kTPCsub1xH2uc);
537      histos->UserHistogram("Event","TPCsub1yH2uc","TPC Qy component sub1;TPCsub1yH2uc",
538      100,-1500.,1500.,
539      AliDielectronVarManager::kTPCsub1yH2uc);
540      histos->UserHistogram("Event","TPCsub1rpH2uc","TPC reaction plane sub1;TPCsub1rpH2uc",
541      100,-2.,2.,
542      AliDielectronVarManager::kTPCsub1rpH2uc);
543      histos->UserHistogram("Event","TPCsub2xH2uc","TPC Qx component sub2;TPCsub2xH2uc",
544      100,-1500.,1500.,
545      AliDielectronVarManager::kTPCsub2xH2uc);
546      histos->UserHistogram("Event","TPCsub2yH2uc","TPC Qy component sub2;TPCsub2yH2uc",
547      100,-1500.,1500.,
548      AliDielectronVarManager::kTPCsub2yH2uc);
549      histos->UserHistogram("Event","TPCsub2rpH2uc","TPC reaction plane sub2;TPCsub2rpH2uc",
550      100,-2.,2.,
551      AliDielectronVarManager::kTPCsub2rpH2uc);
552      histos->UserHistogram("Event","TPCsub12DiffH2uc","TPC reaction plane difference;TPCsub12DiffH2uc",
553      100,-1.,1.,
554      AliDielectronVarManager::kTPCsub12DiffH2uc);
555     */
556
557     // EP resolution calculation
558     histos->UserHistogram("Event","Cent_v0ATPCDiffH2","VZERO-A TPC diff;centrality (%);cos(2(#Psi^{V0A}-#Psi^{TPC}))",
559                           10,0.,100.,300,-1.0,1.0,
560                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kv0ATPCDiffH2);
561     histos->UserHistogram("Event","Cent_v0CTPCDiffH2","VZERO-C TPC diff;centrality (%);cos(2(#Psi^{V0C}-#Psi^{TPC}))",
562                           10,0.,100.,300,-1.0,1.0,
563                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kv0CTPCDiffH2);
564     histos->UserHistogram("Event","Cent_v0Av0CDiffH2","VZERO-A VZERO-C diff;centrality (%);cos(2(#Psi^{V0A}-#Psi^{V0C}))",
565                           10,0.,100.,300,-1.0,1.0,
566                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kv0Av0CDiffH2);
567     histos->UserHistogram("Event","Cent_TPCsub12DiffH2","TPC-sub1 TPC-sub2 diff;centrality (%);cos(2(#Psi^{sub1}-#Psi^{sub2}))",
568                           10,0.,100.,300,-1.0,1.0,
569                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kTPCsub12DiffH2);
570     // detector effects
571     histos->UserHistogram("Event","Cent_TPCsub12DiffH2Sin","TPC-sub1 TPC-sub2 diff;centrality (%);sin(2(#Psi^{sub1}-#Psi^{sub2}))",
572                           10,0.,100.,300,-1.0,1.0,
573                           AliDielectronVarManager::kCentrality,AliDielectronVarManager::kTPCsub12DiffH2Sin);
574   } //hist: flowQA
575
576   if(bHistPair) {
577     //Initialise histogram classes
578     histos->SetReservedWords("Track;Pair");
579
580     //Pair classes
581     // to fill also mixed event histograms loop until 7 or 10
582     for (Int_t i=0; i<(bHistPairME ? 8 : 3); ++i){
583       histos->AddClass(Form("Pair_%s",AliDielectron::PairClassName(i)));
584     }
585     //legs from pair (fill SE)
586     for (Int_t i=0; i<3; ++i){
587       histos->AddClass(Form("Track_Legs_%s",AliDielectron::PairClassName(i)));
588     }
589     //Track classes
590     //to fill also track info from 2nd event loop until 2
591     for (Int_t i=0; i<2; ++i){
592       histos->AddClass(Form("Track_%s",AliDielectron::TrackClassName(i)));
593     }
594     //track rotation
595     //   histos->AddClass(Form("Pair_%s",AliDielectron::PairClassName(AliDielectron::kEv1PMRot)));
596     //   histos->AddClass(Form("Track_Legs_%s",AliDielectron::PairClassName(AliDielectron::kEv1PMRot)));
597
598     // Vertex
599     histos->UserHistogram("Track","ImpactParXY",";dXY (cm);#tracks", 500,-1.,1.,
600                           AliDielectronVarManager::kImpactParXY);
601     histos->UserHistogram("Track","ImpactParZ",";dZ (cm);#tracks", 600,-3.,3.,
602                           AliDielectronVarManager::kImpactParZ);
603
604     // Kinematics
605     histos->UserHistogram("Track","Pt",";p_{T} (GeV/c);#tracks", 400,0,20.,
606                           AliDielectronVarManager::kPt);
607     histos->UserHistogram("Track","Eta_Phi",";#eta;#varphi;#tracks", 200,-1,1,200,0,6.285,
608                           AliDielectronVarManager::kEta,AliDielectronVarManager::kPhi);
609
610     // TPC
611     histos->UserHistogram("Track","TPCnCls",";N_{cls}^{TPC};#tracks", 160,-0.5,159.5,
612                           AliDielectronVarManager::kNclsTPC);
613     histos->UserHistogram("Track","TPCsignalN",";N_{cls}^{TPC};#tracks", 160,-0.5,159.5,
614                           AliDielectronVarManager::kTPCsignalN);
615     histos->UserHistogram("Track","NclFr",";N_{max.cls}^{TPC};#tracks", 160,-0.5,159.5,
616                           AliDielectronVarManager::kNFclsTPCr);
617     histos->UserHistogram("Track","Ncl_NclFr",";N_{cls}^{TPC};N_{max.cls}^{TPC};#tracks",
618                           160,-0.5,159.5,160,-0.5,159.5,
619                           AliDielectronVarManager::kNclsTPC,AliDielectronVarManager::kNFclsTPCr);
620
621     // TRD
622     histos->UserHistogram("Track","NtrklTRD",";N_{trkl}^{TRD};#tracks",
623                           8,-0.5,7.5,
624                           AliDielectronVarManager::kTRDpidQuality);
625
626     // PID
627     histos->UserHistogram("Track","dEdx_P",";p (GeV/c);TPC signal (arb units);#tracks",
628                           400,0.2,20.,200,0.,200.,
629                           AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCsignal,kTRUE);
630     histos->UserHistogram("Track","TPCnSigmaEle_P","p (GeV/c);n#sigma_{ele}^{TPC};#tracks",
631                           400,0.2,20.,200,-10.,10.,
632                           AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaEle,kTRUE);
633     histos->UserHistogram("Track","TOFbeta_P",";p (GeV/c);#beta;#tracks",
634                           250,0.0,5.0,300,0.,1.2,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTOFbeta,kTRUE);
635
636     ///// add histograms to Pair classes /////
637     histos->UserHistogram("Pair","InvMass",";m_{ee} (GeV/c^{2});#pairs",
638                           300,.0,300*0.04, AliDielectronVarManager::kM); // 40MeV bins, 12GeV/c2
639     histos->UserHistogram("Pair","Rapidity",";y;#pairs",
640                           100,-1.,1., AliDielectronVarManager::kY);
641     histos->UserHistogram("Pair","OpeningAngle","opening angle (rad.);#pairs",
642                           100,0.,3.15, AliDielectronVarManager::kOpeningAngle);
643     histos->UserHistogram("Pair","Chi2NDF",";#chi^{2}/NDF;#pairs",
644                           100,0.,20, AliDielectronVarManager::kChi2NDF);
645     histos->UserHistogram("Pair","PsiPair",";#psi;#pairs",
646                           100,0.,3.15, AliDielectronVarManager::kPsiPair);
647     histos->UserHistogram("Pair","R",";r (cm)",
648                           200,0.,100., AliDielectronVarManager::kR);
649     histos->UserHistogram("Pair","LegDist",";dca (cm)",
650                           50,0.,5., AliDielectronVarManager::kLegDist);
651     histos->UserHistogram("Pair","LegDistXY",";dca_{xy} (cm)",
652                           50,0.,5., AliDielectronVarManager::kLegDistXY);
653
654     //// FLOW results use tprofiles
655     if(bHistFlow) {
656       histos->UserProfile("Pair","M_Cent_Pt_v0ACrpH2FlowV2",
657                           "cos(2(#varphi-#Psi^{V0AC}));mass (GeV/c^{2});centrality (%);p_{T} (GeV/c)",
658                           AliDielectronVarManager::kv0ACrpH2FlowV2,
659                           125,0.,125*.04, 10, 0.,100., 200,0.,100.,
660                           AliDielectronVarManager::kM, AliDielectronVarManager::kCentrality, AliDielectronVarManager::kPt);
661       histos->UserProfile("Pair","M_Cent_Pt_v0ArpH2FlowV2",
662                           "cos(2(#varphi-#Psi^{V0A}));mass (GeV/c^{2});centrality (%);p_{T} (GeV/c)",
663                           AliDielectronVarManager::kv0ArpH2FlowV2,
664                           125,0.,125*.04, 10, 0.,100., 200,0.,100.,
665                           AliDielectronVarManager::kM, AliDielectronVarManager::kCentrality, AliDielectronVarManager::kPt);
666       histos->UserProfile("Pair","M_Cent_Pt_v0CrpH2FlowV2",
667                           "cos(2(#varphi-#Psi^{V0C}));mass (GeV/c^{2});centrality (%);p_{T} (GeV/c)",
668                           AliDielectronVarManager::kv0CrpH2FlowV2,
669                           125,0.,125*.04, 10, 0.,100., 200,0.,100.,
670                           AliDielectronVarManager::kM, AliDielectronVarManager::kCentrality, AliDielectronVarManager::kPt);
671     } //hist: flow
672   } //hist: pair
673
674
675   ////// MONTE CARLO //////
676   /*
677     if(cutDefinition == kTOFTRD && hasMC) {
678     histos->AddClass("MCEvent");
679     histos->UserHistogram("MCEvent","Cent_NJPsis","Centrality vs. generated incl. J/#psi per event;centrality (%);N_{J/#psi}",
680     10,0.,100., 21,-0.5,20.5,
681     AliDielectronVarManager::kCentrality,AliDielectronVarManager::kNumberOfJPsis);
682     }
683   */
684
685   histos->UserHistogram("Track","TRDprobEle",";P_{ele}^{TRD};#tracks", 
686                         100,0.,1.,            AliDielectronVarManager::kTRDprobEle);
687   
688   die->SetHistogramManager(histos);
689 }
690
691 void InitHF(AliDielectron* die, Int_t cutDefinition)
692 {
693   //
694   // Setup the HF arrays
695   //
696
697   AliDielectronHF *hf=new AliDielectronHF(die->GetName(),die->GetTitle());
698   //  if(hasMC) hf->SetStepForMCGenerated();
699   hf->SetPairTypes(AliDielectronHF::kAll);
700   hf->SetVariable(AliDielectronVarManager::kM, 125, 0.0, 0.04*125);
701
702   hf->AddCutVariable(AliDielectronVarManager::kCentrality,  "0.,5.,10.,20.,40.,50.,60.,80."  );
703   hf->AddCutVariable(AliDielectronVarManager::kPt,          "0.,2.5,5.,100."                 );
704   hf->AddCutVariable(AliDielectronVarManager::kDeltaPhiv0ArpH2, 8,-1.*TMath::Pi(),TMath::Pi());
705   //  hf->AddCutVariable(AliDielectronVarManager::kY,           1, -0.9, 0.9                     );
706   //  hf->AddCutVariable(AliDielectronVarManager::kPt,          "0.8, 1.0, 1.1, 1.2, 1.5, 100.0", kTRUE, AliDielectronHF::kBinToMax);
707   // hf->AddCutVariable(AliDielectronVarManager::kNclsTPC,     "70,90,100,120,160",              kTRUE, AliDielectronHF::kBinToMax);
708   //  hf->AddCutVariable(AliDielectronVarManager::kTPCnSigmaEle,"-4,-3,-2.5,-2,2,2.5,3,4",        kTRUE, AliDielectronHF::kSymBin);
709   //hf->AddCutVariable(AliDielectronVarManager::kTPCnSigmaPio,"3.,3.5,4.,100.",                 kTRUE, AliDielectronHF::kBinToMax);
710   //hf->AddCutVariable(AliDielectronVarManager::kITSLayerFirstCls,4,0.,4.,              kFALSE, kTRUE, AliDielectronHF::kBinFromMin);
711   //hf->AddCutVariable(AliDielectronVarManager::kNclsITS,         5,2.,7.,              kFALSE, kTRUE, AliDielectronHF::kBinToMax);
712   //hf->AddCutVariable(AliDielectronVarManager::kRunNumber,  GetRunNumbers());
713
714   die->SetHistogramArray(hf);
715 }
716
717 void InitCF(AliDielectron* die, Int_t cutDefinition)
718 {
719   //
720   // Setup the CF Manager if needed
721   //
722
723  AliDielectronCF *cf=new AliDielectronCF(die->GetName(),die->GetTitle());
724
725   // event variables
726   cf->AddVariable(AliDielectronVarManager::kCentrality,"0.,5.,10.,20.,40.,50.,60.,80.");
727   //    if(!hasMC) cf->AddVariable(AliDielectronVarManager::kZvPrim,20, -10., 10.);
728   if(hasMC)  cf->AddVariable(AliDielectronVarManager::kNacc,20,0.,3000.0);
729   if(hasMC)  cf->AddVariable(AliDielectronVarManager::kNVtxContrib,20,0.,4000.);
730   if(hasMC)  cf->AddVariable(AliDielectronVarManager::kRunNumber, GetRunNumbers() );
731
732   // pair variables
733   //    cf->AddVariable(AliDielectronVarManager::kY,"-0.8,0.8");
734   cf->AddVariable(AliDielectronVarManager::kM,125,0.,125*.04); //40Mev Steps
735   cf->AddVariable(AliDielectronVarManager::kPairType,11,0,11);
736   cf->AddVariable(AliDielectronVarManager::kPt,"0., 1., 2.5, 5., 100.0");
737   if(hasMC) cf->AddVariable(AliDielectronVarManager::kTRDpidEffPair,101,0.0,1.01);
738   //    if(hasMC) cf->AddVariable(AliDielectronVarManager::kThetaCS,15,-1.,1.);
739
740   // flow variables
741   cf->AddVariable(AliDielectronVarManager::kDeltaPhiv0ArpH2,4,-1.*TMath::Pi(),TMath::Pi());
742   cf->AddVariable(AliDielectronVarManager::kDeltaPhiv0CrpH2,4,-1.*TMath::Pi(),TMath::Pi());
743
744   // leg variables
745   cf->AddVariable(AliDielectronVarManager::kPt,"0.8, 1.0, 1.1, 1.2, 1.5, 100.0",kTRUE);
746   if(hasMC) cf->AddVariable(AliDielectronVarManager::kEta,"-0.9,0.9",kTRUE);
747   //    cf->AddVariable(AliDielectronVarManager::kITSLayerFirstCls,7,-1.5,5.5,kTRUE);
748   //    cf->AddVariable(AliDielectronVarManager::kNclsITS,"1,2,3,4,5,6",kTRUE);
749   cf->AddVariable(AliDielectronVarManager::kTPCnSigmaEle,"-3,-2.5,-2,2,2.5,3",kTRUE);
750   cf->AddVariable(AliDielectronVarManager::kTPCnSigmaPio,"2.5,3.0,3.5,4.0,4.5,100",kTRUE);
751   cf->AddVariable(AliDielectronVarManager::kNclsTPC,"70, 90, 100, 120, 160",kTRUE);
752   //    cf->AddVariable(AliDielectronVarManager::kTPCnSigmaPro,"3.5,4.0,4.5,5.0,100",kTRUE);
753   //    cf->AddVariable(AliDielectronVarManager::kTOFnSigmaEle,"-3,-2,2,3",kTRUE); break;
754   //    cf->AddVariable(AliDielectronVarManager::kTRDpidQuality,"3.5, 4.5, 5.5, 6.5",kTRUE);
755   //    if(!hasMC && isESD) cf->AddVariable(AliDielectronVarManager::kTRDchi2,"-1.,0.,2.,4.",kTRUE);
756
757   // mc steps
758   if(hasMC) {
759     if(cutDefinition==kTOFTRD) cf->SetStepForMCtruth();
760     cf->SetStepsForMCtruthOnly();
761     // cf->SetStepsForBackground();
762   }
763
764   die->SetCFManagerPair(cf);
765 }
766
767 void AddMCSignals(AliDielectron *die){
768   //Do we have an MC handler?
769   if (!hasMC) return;
770   
771   AliDielectronSignalMC* inclusiveJpsi = new AliDielectronSignalMC("inclusiveJpsi","Inclusive J/psi");
772   inclusiveJpsi->SetLegPDGs(11,-11);
773   inclusiveJpsi->SetMotherPDGs(443,443);
774   inclusiveJpsi->SetMothersRelation(AliDielectronSignalMC::kSame);
775   inclusiveJpsi->SetFillPureMCStep(kTRUE);
776   inclusiveJpsi->SetCheckBothChargesLegs(kTRUE,kTRUE);
777   inclusiveJpsi->SetCheckBothChargesMothers(kTRUE,kTRUE);
778   die->AddSignalMC(inclusiveJpsi);
779   
780   AliDielectronSignalMC* promptJpsi = new AliDielectronSignalMC("promptJpsi","Prompt J/psi");   // prompt J/psi (not from beauty decays)
781   promptJpsi->SetLegPDGs(11,-11);
782   promptJpsi->SetMotherPDGs(443,443);
783   promptJpsi->SetGrandMotherPDGs(503,503,kTRUE,kTRUE);   // not from beauty hadrons
784   promptJpsi->SetMothersRelation(AliDielectronSignalMC::kSame);
785   promptJpsi->SetFillPureMCStep(kTRUE);
786   promptJpsi->SetLegSources(AliDielectronSignalMC::kFinalState, AliDielectronSignalMC::kFinalState);
787   promptJpsi->SetCheckBothChargesLegs(kTRUE,kTRUE);
788   promptJpsi->SetCheckBothChargesMothers(kTRUE,kTRUE);
789   promptJpsi->SetCheckBothChargesGrandMothers(kTRUE,kTRUE);
790   die->AddSignalMC(promptJpsi);
791   
792   AliDielectronSignalMC* beautyJpsi = new AliDielectronSignalMC("beautyJpsi","Beauty J/psi");
793   beautyJpsi->SetLegPDGs(11,-11);
794   beautyJpsi->SetMotherPDGs(443,443);
795   beautyJpsi->SetMothersRelation(AliDielectronSignalMC::kSame);
796   beautyJpsi->SetGrandMotherPDGs(500,500);
797   beautyJpsi->SetFillPureMCStep(kTRUE);
798   beautyJpsi->SetCheckBothChargesLegs(kTRUE,kTRUE);
799   beautyJpsi->SetCheckBothChargesMothers(kTRUE,kTRUE);
800   beautyJpsi->SetCheckBothChargesGrandMothers(kTRUE,kTRUE);
801   die->AddSignalMC(beautyJpsi);
802   
803   AliDielectronSignalMC* directJpsi = new AliDielectronSignalMC("directJpsi","Direct J/psi");   // embedded J/psi
804   directJpsi->SetLegPDGs(11,-11);
805   directJpsi->SetMotherPDGs(443,443);
806   directJpsi->SetMothersRelation(AliDielectronSignalMC::kSame);
807   directJpsi->SetFillPureMCStep(kTRUE);
808   directJpsi->SetLegSources(AliDielectronSignalMC::kFinalState, AliDielectronSignalMC::kFinalState);
809   directJpsi->SetMotherSources(AliDielectronSignalMC::kDirect, AliDielectronSignalMC::kDirect);
810   directJpsi->SetCheckBothChargesLegs(kTRUE,kTRUE);
811   directJpsi->SetCheckBothChargesMothers(kTRUE,kTRUE);
812   die->AddSignalMC(directJpsi);
813   
814   AliDielectronSignalMC* conversionElePairs = new AliDielectronSignalMC("conversionElePairs","conversion electron pairs");  // pairs made from conversion (may be also from 2 different conversions)                                                                                                                                                                                             
815   conversionElePairs->SetLegPDGs(11,-11);
816   conversionElePairs->SetCheckBothChargesLegs(kTRUE,kTRUE);
817   conversionElePairs->SetLegSources(AliDielectronSignalMC::kSecondary, AliDielectronSignalMC::kSecondary);
818   conversionElePairs->SetMotherPDGs(22,22);
819   //   die->AddSignalMC(conversionElePairs);
820 }
821
822 void SetEtaCorrection()
823 {
824   if (AliDielectronPID::GetEtaCorrFunction()) return;
825   
826   TString etaMap="$TRAIN_ROOT/jpsi_JPSI/EtaCorrMaps.root";
827   TString trainRoot=gSystem->Getenv("TRAIN_ROOT");
828   if (trainRoot.IsNull()) etaMap="$ALICE_ROOT/PWGDQ/dielectron/files/EtaCorrMaps.root";
829   if (gSystem->AccessPathName(gSystem->ExpandPathName(etaMap.Data()))){
830     Error("ConfigPbPb","Eta map not found: %s",etaMap.Data());
831     return;
832   }
833
834   TFile f(etaMap.Data());
835   if (!f.IsOpen()) return;
836   TList *keys=f.GetListOfKeys();
837
838   for (Int_t i=0; i<keys->GetEntries(); ++i){
839     TString kName=keys->At(i)->GetName();
840     TPRegexp reg(kName);
841     if (reg.MatchB(list)){
842       printf(" Using Eta Correction Function: %s\n",kName.Data());
843       AliDielectronPID::SetEtaCorrFunction((TF1*)f.Get(kName.Data()));
844     }
845   }
846 }
847
848 TVectorD *GetRunNumbers() {
849   
850   Double_t runLHC10h[] = { // all good runs based on RCT 29.Mai
851     139510, 139507, 139505, 139503, 139465, 139438, 139437, 139360, 139329, 139328, 139314, 139310, 139309, 139173, 139107, 139105, 139038, 139037, 139036, 139029, 139028, 138872, 138871, 138870, 138837, 138732, 138730, 138666, 138662, 138653, 138652, 138638, 138624, 138621, 138583, 138582, 138579, 138578, 138534, 138469, 138442, 138439, 138438, 138396, 138364, 138275, 138225, 138201, 138197, 138192, 138190, 137848, 137844, 137752, 137751, 137724, 137722, 137718, 137704, 137693, 137692, 137691, 137686, 137685, 137639, 137638, 137608, 137595, 137549, 137546, 137544, 137541, 137539, 137531, 137530, 137443, 137441, 137440, 137439, 137434, 137432, 137431, 137430, 137366, 137243, 137236, 137235, 137232, 137231, 137230, 137162, 137161, 137135
852   };
853   
854   Double_t runLHC11h[] = { // all good runs based on RCT 29.Mai
855     170593, 170572, 170388, 170387, 170315, 170313, 170312, 170311, 170309, 170308, 170306, 170270, 170269, 170268, 170230, 170228, 170207, 170204, 170203, 170193, 170163, 170159, 170155, 170091, 170089, 170088, 170085, 170084, 170083, 170081, 170040, 170027, 169965, 169923, 169859, 169858, 169855, 169846, 169838, 169837, 169835, 169591, 169590, 169588, 169587, 169586, 169557, 169555, 169554, 169553, 169550, 169515, 169512, 169506, 169504, 169498, 169475, 169420, 169419, 169418, 169417, 169415, 169411, 169238, 169167, 169160, 169156, 169148, 169145, 169144, 169138, 169099, 169094, 169091, 169045, 169044, 169040, 169035, 168992, 168988, 168826, 168777, 168514, 168512, 168511, 168467, 168464, 168460, 168458, 168362, 168361, 168342, 168341, 168325, 168322, 168311, 168310, 168115, 168108, 168107, 168105, 168076, 168069, 167988, 167987, 167985, 167920, 167915
856   };
857   
858   // selection via environement variable (works only for gsi trains)
859
860   
861   if(list.Contains("LHC10h") || list.Contains("LHC11a10")) {
862     Int_t size = (int) (sizeof(runLHC10h)/sizeof(Double_t));
863     TVectorD *vec = new TVectorD(size+1);
864     
865     (*vec)[size] = runLHC10h[0] + 1;
866     for (int i = 0; i < size; i++) {
867       (*vec)[i] = runLHC10h[size-1-i];
868     }
869     //    vec->Print("");
870     return vec;
871   }
872
873   if( list.Contains("LHC11h") || list.Contains("LHC12a17") ) {
874     
875     Int_t size = (int) (sizeof(runLHC11h)/sizeof(Double_t));
876     TVectorD *vec = new TVectorD(size+1);
877     
878     (*vec)[size] = runLHC11h[0] + 1;
879     for (int i = 0; i < size; i++) {
880       (*vec)[i] = runLHC11h[size-1-i];
881     }
882     //   vec->Print("");
883     return vec;
884   }
885
886   TVectorD *vec = new TVectorD(2);
887   (*vec)[0] = 0;
888   (*vec)[0] = 1;
889   return vec;
890      
891 }