5ea4a0e0748c96b75e59351547b392869bec7150
[u/mrichter/AliRoot.git] / PWGDQ / dielectron / macrosJPSI / AddTask_jpsi_Default.C
1 void InitHistograms(AliDielectron *die, Int_t cutDefinition);
2 void InitCF(AliDielectron* die, Int_t cutDefinition);
3
4 void SetupTrackCuts(AliDielectron *die, Int_t cutDefinition);
5 void SetupPairCuts(AliDielectron *die, Int_t cutDefinition);
6
7 TString names=("default");
8 TObjArray *arrNames=names.Tokenize(";");
9
10 const Int_t nDie=arrNames->GetEntries();
11
12 //______________________________________________________________________________________
13 AliAnalysisTask* AddTask_jpsi_Default(TString prod="")
14 {
15   //get the current analysis manager
16
17   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
18   if (!mgr) {
19     Error("AddTask_jpsi_JPsi", "No analysis manager found.");
20     return 0;
21   }
22
23   //Do we have an MC handler?
24   Bool_t hasMC=(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()!=0x0);
25   //AOD input?
26   Bool_t isAOD=mgr->GetInputEventHandler()->IsA()==AliAODInputHandler::Class();
27   
28   //Get the current train configuration
29   TString trainConfig=gSystem->Getenv("CONFIG_FILE");
30   
31   TString list=gSystem->Getenv("LIST");
32   if( list.IsNull()) list=prod;
33   
34   //create task and add it to the manager
35   AliAnalysisTaskMultiDielectron *task=new AliAnalysisTaskMultiDielectron("JpsiDefault");
36   if (!hasMC) task->UsePhysicsSelection();
37   if (list.Contains("LHC11d")) task->SetTriggerMask(AliVEvent::kEMCEJE+AliVEvent::kEMC7+AliVEvent::kEMCEGA);
38   if (list.Contains("LHC11h")) {
39     task->SetTriggerMask(AliVEvent::kMB+AliVEvent::kCentral+AliVEvent::kSemiCentral);
40     printf("AddTask_jpsi_Default: LHC11h data detected. Using trigger classes: %d\n",AliVEvent::kMB+AliVEvent::kCentral+AliVEvent::kSemiCentral);
41   }
42   if (list.Contains("LHC12h")) task->SetTriggerMask(AliVEvent::kAnyINT);
43   mgr->AddTask(task);
44
45   //add dielectron analysis with different cuts to the task
46   for (Int_t i=0; i<nDie; ++i){ //nDie defined in config file
47     AliDielectron *jpsi=ConfigDefault(i);
48     if (!jpsi) continue;
49     task->AddDielectron(jpsi);
50   }
51   
52   //Add event filter
53   AliDielectronEventCuts *eventCuts=new AliDielectronEventCuts("eventCuts","Vertex Track && |vtxZ|<10 && ncontrib>0");
54   if (isAOD) eventCuts->SetVertexType(AliDielectronEventCuts::kVtxAny);
55   eventCuts->SetRequireVertex();
56   eventCuts->SetMinVtxContributors(1);
57   eventCuts->SetVertexZ(-10.,10.);
58   task->SetEventFilter(eventCuts);
59
60   //   task->SetTriggerOnV0AND();
61 //   if ( trainConfig=="pp" ) task->SetRejectPileup();
62   
63   //create output container
64   AliAnalysisDataContainer *coutput1 =
65     mgr->CreateContainer("jpsi_Default_tree",
66                          TTree::Class(),
67                          AliAnalysisManager::kExchangeContainer,
68                          "jpsi_Default_default");
69   
70   AliAnalysisDataContainer *cOutputHist1 =
71     mgr->CreateContainer("jpsi_Default_QA",
72                          TList::Class(),
73                          AliAnalysisManager::kOutputContainer,
74                          "jpsi_Default.root");
75
76   AliAnalysisDataContainer *cOutputHist2 =
77     mgr->CreateContainer("jpsi_Default_CF",
78                          TList::Class(),
79                          AliAnalysisManager::kOutputContainer,
80                          "jpsi_Default.root");
81   
82   AliAnalysisDataContainer *cOutputHist3 =
83     mgr->CreateContainer("jpsi_Default_EventStat",
84                          TH1D::Class(),
85                          AliAnalysisManager::kOutputContainer,
86                          "jpsi_Default.root");
87   
88   mgr->ConnectInput(task,  0, mgr->GetCommonInputContainer());
89   mgr->ConnectOutput(task, 0, coutput1 );
90   mgr->ConnectOutput(task, 1, cOutputHist1);
91   mgr->ConnectOutput(task, 2, cOutputHist2);
92   mgr->ConnectOutput(task, 3, cOutputHist3);
93   
94   return task;
95 }
96
97
98 //______________________________________________________________________________________
99 //______________________________________________________________________________________
100 //______________________________________________________________________________________
101 //
102 // Here the configuration part starts
103 //
104 AliDielectron* ConfigDefault(Int_t cutDefinition)
105 {
106   //
107   // Setup the instance of AliDielectron
108   //
109   
110   Bool_t hasMC=(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()!=0x0);
111   
112   // create the actual framework object
113   TString name=Form("%02d",cutDefinition);
114   if (cutDefinition<arrNames->GetEntriesFast()){
115     name=arrNames->At(cutDefinition)->GetName();
116   }
117   AliDielectron *die =
118   new AliDielectron(Form("%s",name.Data()),
119                     Form("Track cuts: %s",name.Data()));
120   
121   if (hasMC) SetupMCsignals(die);
122   // cut setup
123   SetupTrackCuts(die,cutDefinition);
124   //
125   SetupPairCuts(die,cutDefinition);
126   //
127   InitHistograms(die,cutDefinition);
128   //
129   InitCF(die,cutDefinition);
130   //
131   //   AliDielectronMixingHandler *mix=new AliDielectronMixingHandler;
132   //   mix->AddVariable(AliDielectronVarManager::kZvPrim,"-10,-8,-5,0,5,8,10");
133   //   mix->AddVariable(AliDielectronVarManager::kPhi ,"0,3.14,6.3");
134   //   mix->SetDepth(10);
135   //  die->SetMixingHandler(mix);
136   //
137   return die;
138 }
139
140 //______________________________________________________________________________________
141 void SetupTrackCuts(AliDielectron *die, Int_t cutDefinition)
142 {
143   //
144   // Setup the track cuts
145   //
146
147   AliDielectronCutGroup* cuts = new AliDielectronCutGroup("cuts","cuts",AliDielectronCutGroup::kCompAND);
148   die->GetTrackFilter().AddCuts(cuts);
149
150   //default quality cuts
151   AliDielectronTrackCuts *cut1=new AliDielectronTrackCuts("cut1","cut1");
152   cut1->SetRequireITSRefit(kTRUE);
153   cut1->SetRequireTPCRefit(kTRUE);
154   cut1->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
155   cuts->AddCut(cut1);
156
157   //pt and kink mother
158   AliDielectronVarCuts *cut2 = new AliDielectronVarCuts("cut2Cut","cut2 cut");
159   cut2->AddCut(AliDielectronVarManager::kPt,            0.8 ,     1.e30 );
160   cut2->AddCut(AliDielectronVarManager::kImpactParZ,   -3.  ,     3.    );
161   cut2->AddCut(AliDielectronVarManager::kImpactParXY,  -1.  ,     1.    );
162   cut2->AddCut(AliDielectronVarManager::kEta,          -0.9 ,     0.9   );
163   cut2->AddCut(AliDielectronVarManager::kNclsTPC,      70.  ,  160.     );
164   cut2->AddCut(AliDielectronVarManager::kTPCchi2Cl,     0.  ,    4.     );
165   //   cut2->AddCut(AliDielectronVarManager::kITSchi2Cl,     0.  ,   36.     ); // removes all tracks in AOD!! why???
166   cut2->AddCut(AliDielectronVarManager::kKinkIndex0,    0.              );
167   cut2->AddCut(AliDielectronVarManager::kTPCnSigmaEle,  -3. ,    3.     );
168   cut2->AddCut(AliDielectronVarManager::kTPCnSigmaPio,   3.5, 1000.     );
169   cut2->AddCut(AliDielectronVarManager::kTPCnSigmaPro,   3. , 1000.     );
170   cuts->AddCut(cut2);
171
172
173   //exclude conversion electrons selected by the tender
174   //   AliDielectronTrackCuts *noconv=new AliDielectronTrackCuts("noConv","noConv");
175   //   noconv->SetV0DaughterCut(AliPID::kElectron,kTRUE);
176   //   cuts->AddCut(noconv);
177 }
178
179 //______________________________________________________________________________________
180 void SetupPairCuts(AliDielectron *die, Int_t cutDefinition)
181 {
182   //
183   // Setup the pair cuts
184   //
185
186
187   // add conversion rejection
188   AliDielectronVarCuts *gammaCut=new AliDielectronVarCuts("gammaCut","gammaCut");
189   gammaCut->AddCut(AliDielectronVarManager::kM,0.,.05);
190   die->GetPairPreFilter().AddCuts(gammaCut);
191   die->SetPreFilterUnlikeOnly();
192
193 }
194
195 //______________________________________________________________________________________
196 void InitHistograms(AliDielectron *die, Int_t cutDefinition)
197 {
198   //
199   // Initialise the histograms
200   //
201
202   //Setup histogram Manager
203   AliDielectronHistos *histos=
204   new AliDielectronHistos(die->GetName(),
205                           die->GetTitle());
206
207   //Initialise histogram classes
208   histos->SetReservedWords("Track;Pair");
209
210   //Track classes
211   //to fill also track info from 2nd event loop until 2
212   for (Int_t i=0; i<2; ++i){
213     histos->AddClass(Form("Track_%s",AliDielectron::TrackClassName(i)));
214   }
215
216   //Pair classes
217   // to fill also mixed event histograms loop until 10
218   for (Int_t i=0; i<3; ++i){
219     histos->AddClass(Form("Pair_%s",AliDielectron::PairClassName(i)));
220   }
221
222   //add histograms to event class
223   if (cutDefinition==0) {
224     histos->AddClass("Event");
225     histos->UserHistogram("Event","","",300,-15.,15.,AliDielectronVarManager::kZvPrim);
226     histos->UserHistogram("Event","","",
227                           100,-2.,2.,100,-2.,2.,AliDielectronVarManager::kXvPrim,AliDielectronVarManager::kYvPrim);
228   }
229
230   //add histograms to Track classes
231   histos->UserHistogram("Track","","",200,0,20.,AliDielectronVarManager::kPt);
232   histos->UserHistogram("Track","","",
233                         400,0.2,20.,200,0.,200.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCsignal,kTRUE);
234
235   histos->UserHistogram("Track","","",
236                         100,0.2,20.,100,-10.,10.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaEle,kTRUE);
237   histos->UserHistogram("Track","","",
238                         100,0.2,20.,100,-10.,10.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaPio,kTRUE);
239   histos->UserHistogram("Track","","",
240                         100,0.2,20.,100,-10.,10.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaPro,kTRUE);
241
242   histos->UserHistogram("Track","","",
243                         100,-2,2,144,0,6.285,AliDielectronVarManager::kEta,AliDielectronVarManager::kPhi);
244   histos->UserHistogram("Track","","",
245                         160,-0.5,159.5,AliDielectronVarManager::kNclsTPC);
246   histos->UserHistogram("Track","","",
247                         100,0.,10.,AliDielectronVarManager::kTPCchi2Cl);
248   histos->UserHistogram("Track","","",
249                         150,-15,15,160,-0.5,159.5,AliDielectronVarManager::kYsignedIn,AliDielectronVarManager::kTPCsignalN);
250   histos->UserHistogram("Track","","",
251                         1000,0.,0.,AliDielectronVarManager::kKinkIndex0);
252
253   histos->UserHistogram("Track","","",500,-1.,1.,AliDielectronVarManager::kImpactParXY);
254   histos->UserHistogram("Track","","",600,-3.,3.,AliDielectronVarManager::kImpactParZ);
255
256   //add histograms to Pair classes
257   histos->UserHistogram("Pair","","",
258                         301,-.01,6.01,AliDielectronVarManager::kM);
259   histos->UserHistogram("Pair","","",
260                         100,-1.,1.,AliDielectronVarManager::kY);
261   histos->UserHistogram("Pair","","",
262                         100,0.,3.15,AliDielectronVarManager::kOpeningAngle);
263
264   die->SetHistogramManager(histos);
265
266 }
267
268 //______________________________________________________________________________________
269 void InitCF(AliDielectron* die, Int_t cutDefinition)
270 {
271   //
272   // Setupd the CF Manager if needed
273   //
274
275   AliDielectronCF *cf=new AliDielectronCF(die->GetName(),die->GetTitle());
276
277   //pair variables
278   cf->AddVariable(AliDielectronVarManager::kPt,"0.0, 1.0, 1.3, 2.0, 3.0, 5., 7.0, 10.0, 100.0");
279
280   cf->AddVariable(AliDielectronVarManager::kZvPrim,"-10,-5,0,5,10");
281   cf->AddVariable(AliDielectronVarManager::kY,"-1,-0.9,-0.8,-0.3,0,0.3,0.9,1.0");
282   cf->AddVariable(AliDielectronVarManager::kM,125,0.,125*.04); //40Mev Steps
283   cf->AddVariable(AliDielectronVarManager::kPairType,11,0,11);
284
285   //leg variables
286   cf->AddVariable(AliDielectronVarManager::kPt,"0.0, 0.8, 1.0, 1.1, 1.2, 1.3, 100.0",kTRUE);
287   cf->AddVariable(AliDielectronVarManager::kEta,"-1.,-0.9,-0.8,0,0.8,0.9,1.0",kTRUE);
288
289   //cf->AddVariable(AliDielectronVarManager::kITSLayerFirstCls,6,0,6,kTRUE);
290
291
292   Bool_t hasMC=(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()!=0x0);
293   if (hasMC){
294     cf->AddVariable(AliDielectronVarManager::kPdgCode,10000,-5000.5,4999.5,kTRUE);
295     cf->AddVariable(AliDielectronVarManager::kPdgCodeMother,10000,-5000.5,4999.5,kTRUE);
296     cf->AddVariable(AliDielectronVarManager::kPdgCodeGrandMother,10000,-5000.5,4999.5,kTRUE);
297
298     //only steps for efficiencies
299     cf->SetStepsForMCtruthOnly();
300   }
301
302   //only in this case write MC truth info
303   if (cutDefinition==0){
304     cf->SetStepForMCtruth();
305   }
306
307   cf->SetStepsForSignal();
308
309   die->SetCFManagerPair(cf);
310 }
311
312 //______________________________________________________________________________________
313 void SetupMCsignals(AliDielectron *die){
314   AliDielectronSignalMC* promptJpsi = new AliDielectronSignalMC("promptJpsi","Prompt J/psi");   // prompt J/psi (not from beauty decays)
315   promptJpsi->SetLegPDGs(11,-11);
316   promptJpsi->SetMotherPDGs(443,443);
317   promptJpsi->SetGrandMotherPDGs(503,503,kTRUE,kTRUE);   // not from beauty hadrons
318   promptJpsi->SetMothersRelation(AliDielectronSignalMC::kSame);
319   promptJpsi->SetFillPureMCStep(kTRUE);
320   promptJpsi->SetLegSources(AliDielectronSignalMC::kFinalState, AliDielectronSignalMC::kFinalState);
321   promptJpsi->SetCheckBothChargesLegs(kTRUE,kTRUE);
322   promptJpsi->SetCheckBothChargesMothers(kTRUE,kTRUE);
323   promptJpsi->SetCheckBothChargesGrandMothers(kTRUE,kTRUE);
324   die->AddSignalMC(promptJpsi);
325
326   AliDielectronSignalMC* promptJpsiNonRad = new AliDielectronSignalMC("promptJpsiNonRad","Prompt J/psi non-Radiative");   // prompt J/psi (not from beauty decays)
327   promptJpsiNonRad->SetLegPDGs(11,-11);
328   promptJpsiNonRad->SetMotherPDGs(443,443);
329   promptJpsiNonRad->SetGrandMotherPDGs(503,503,kTRUE,kTRUE);   // not from beauty hadrons
330   promptJpsiNonRad->SetMothersRelation(AliDielectronSignalMC::kSame);
331   promptJpsiNonRad->SetFillPureMCStep(kTRUE);
332   promptJpsiNonRad->SetLegSources(AliDielectronSignalMC::kFinalState, AliDielectronSignalMC::kFinalState);
333   promptJpsiNonRad->SetCheckBothChargesLegs(kTRUE,kTRUE);
334   promptJpsiNonRad->SetCheckBothChargesMothers(kTRUE,kTRUE);
335   promptJpsiNonRad->SetCheckBothChargesGrandMothers(kTRUE,kTRUE);
336   promptJpsiNonRad->SetJpsiRadiative(AliDielectronSignalMC::kIsNotRadiative);
337   die->AddSignalMC(promptJpsiNonRad);
338
339   AliDielectronSignalMC* promptJpsiRad = new AliDielectronSignalMC("promptJpsiRad","Prompt J/psi Radiative");   // prompt J/psi (not from beauty decays)
340   promptJpsiRad->SetLegPDGs(11,-11);
341   promptJpsiRad->SetMotherPDGs(443,443);
342   promptJpsiRad->SetGrandMotherPDGs(503,503,kTRUE,kTRUE);   // not from beauty hadrons
343   promptJpsiRad->SetMothersRelation(AliDielectronSignalMC::kSame);
344   promptJpsiRad->SetFillPureMCStep(kTRUE);
345   promptJpsiRad->SetLegSources(AliDielectronSignalMC::kFinalState, AliDielectronSignalMC::kFinalState);
346   promptJpsiRad->SetCheckBothChargesLegs(kTRUE,kTRUE);
347   promptJpsiRad->SetCheckBothChargesMothers(kTRUE,kTRUE);
348   promptJpsiRad->SetCheckBothChargesGrandMothers(kTRUE,kTRUE);
349   promptJpsiRad->SetJpsiRadiative(AliDielectronSignalMC::kIsRadiative);
350   die->AddSignalMC(promptJpsiRad);
351
352 }