d17b6774eb64336c09b5dbdfc269c1620cea326e
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / macros / AddTask_GammaConvNeutralMesonPiPlPiMiPiZero_MixedMode_pPb.C
1 void AddTask_GammaConvNeutralMesonPiPlPiMiPiZero_MixedMode_pPb(    
2                                                                                 Int_t trainConfig = 1,
3                                                                                 Bool_t isMC       = kFALSE, //run MC 
4                                                                                 Bool_t enableQAMesonTask = kTRUE, //enable QA in AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero
5                                                                                 TString fileNameInputForWeighting = "MCSpectraInput.root", // path to file for weigting input
6                                                                                 Bool_t doWeighting = kFALSE,  //enable Weighting
7                                                                                 TString generatorName = "HIJING",                               
8                                                                                 TString cutnumberAODBranch = "0000000060084001001500000"
9                                                                                 ) {
10
11         // ================= Load Librariers =================================
12         gSystem->Load("libCore.so");  
13         gSystem->Load("libTree.so");
14         gSystem->Load("libGeom.so");
15         gSystem->Load("libVMC.so");
16         gSystem->Load("libPhysics.so");
17         gSystem->Load("libMinuit");
18         gSystem->Load("libSTEERBase");
19         gSystem->Load("libESD");
20         gSystem->Load("libAOD");
21         gSystem->Load("libANALYSIS");
22         gSystem->Load("libANALYSISalice");  
23         gSystem->Load("libPWGGAGammaConv.so");
24         gSystem->Load("libCDB.so");
25         gSystem->Load("libSTEER.so");
26         gSystem->Load("libSTEERBase.so");
27         gSystem->Load("libTENDER.so");
28         gSystem->Load("libTENDERSupplies.so");
29
30         Int_t isHeavyIon = 2;
31         Int_t neutralPionMode = 1;
32         
33         // ================== GetAnalysisManager ===============================
34         AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
35         if (!mgr) {
36                 Error(Form("AddTask_GammaConvNeutralMesonPiPlPiMiPiZero_pPb_%i",trainConfig), "No analysis manager found.");
37                 return ;
38         }
39
40         // ================== GetInputEventHandler =============================
41         AliVEventHandler *inputHandler=mgr->GetInputEventHandler();
42         
43         //========= Add PID Reponse to ANALYSIS manager ====
44         if(!(AliPIDResponse*)mgr->GetTask("PIDResponseTask")){
45                 gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPIDResponse.C");
46                 AddTaskPIDResponse(isMC);
47         }
48         
49         //=========  Set Cutnumber for V0Reader ================================
50         TString cutnumberPhoton = "060084001001500000000";
51         TString cutnumberEvent = "8000000";
52         TString PionCuts      = "000000200";            //Electron Cuts
53         
54
55         
56         Bool_t doEtaShift = kFALSE;
57
58         AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
59         
60         //========= Add V0 Reader to  ANALYSIS manager if not yet existent =====
61         if( !(AliV0ReaderV1*)mgr->GetTask("V0ReaderV1") ){
62                 AliV0ReaderV1 *fV0ReaderV1 = new AliV0ReaderV1("V0ReaderV1");
63                 
64                 fV0ReaderV1->SetUseOwnXYZCalculation(kTRUE);
65                 fV0ReaderV1->SetCreateAODs(kFALSE);// AOD Output
66                 fV0ReaderV1->SetUseAODConversionPhoton(kTRUE);
67                 
68                 if (!mgr) {
69                         Error("AddTask_V0ReaderV1", "No analysis manager found.");
70                         return;
71                 }
72
73                 AliConvEventCuts *fEventCuts=NULL;
74                 if(cutnumberEvent!=""){
75                         fEventCuts= new AliConvEventCuts(cutnumberEvent.Data(),cutnumberEvent.Data());
76                         fEventCuts->SetPreSelectionCutFlag(kTRUE);
77                         if(fEventCuts->InitializeCutsFromCutString(cutnumberEvent.Data())){
78                                 fEventCuts->DoEtaShift(doEtaShift);
79                                 fV0ReaderV1->SetEventCuts(fEventCuts);
80                                 fEventCuts->SetFillCutHistograms("",kTRUE);
81                         }
82                 }
83
84                 // Set AnalysisCut Number
85                 AliConversionPhotonCuts *fCuts=NULL;
86                 if(cutnumberPhoton!=""){
87                         fCuts= new AliConversionPhotonCuts(cutnumberPhoton.Data(),cutnumberPhoton.Data());
88                         fCuts->SetPreSelectionCutFlag(kTRUE);
89                         fCuts->SetIsHeavyIon(isHeavyIon);
90                         if(fCuts->InitializeCutsFromCutString(cutnumberPhoton.Data())){
91                                 fV0ReaderV1->SetConversionCuts(fCuts);
92                                 fCuts->SetFillCutHistograms("",kTRUE);
93                         }
94                 }
95
96                 if(inputHandler->IsA()==AliAODInputHandler::Class()){
97                 // AOD mode
98                         fV0ReaderV1->SetDeltaAODBranchName(Form("GammaConv_%s_gamma",cutnumberAODBranch.Data()));
99                 }
100                 fV0ReaderV1->Init();
101
102                 AliLog::SetGlobalLogLevel(AliLog::kInfo);
103
104                 //connect input V0Reader
105                 mgr->AddTask(fV0ReaderV1);
106                 mgr->ConnectInput(fV0ReaderV1,0,cinput);
107         }
108
109         //================================================
110         //========= Add Electron Selector ================
111
112
113         if( !(AliPrimaryPionSelector*)mgr->GetTask("PionSelector") ){
114
115                 AliPrimaryPionSelector *fPionSelector = new AliPrimaryPionSelector("PionSelector");
116                 // Set AnalysisCut Number
117
118                 AliPrimaryPionCuts *fPionCuts=0;
119                 if( PionCuts!=""){
120                         fPionCuts= new AliPrimaryPionCuts(PionCuts.Data(),PionCuts.Data());
121                         if(fPionCuts->InitializeCutsFromCutString(PionCuts.Data())){
122                                 fPionSelector->SetPrimaryPionCuts(fPionCuts);
123                                 fPionCuts->SetFillCutHistograms("",kTRUE);
124
125                         }
126                 }
127
128                 fPionSelector->Init();
129                 mgr->AddTask(fPionSelector);
130                 
131                 AliAnalysisDataContainer *cinput1  = mgr->GetCommonInputContainer();
132
133                 //connect input V0Reader
134                 mgr->ConnectInput (fPionSelector,0,cinput1);
135
136         }
137
138         
139         
140         AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero *task=NULL;
141
142         task= new AliAnalysisTaskNeutralMesonToPiPlPiMiPiZero(Form("GammaConvNeutralMesonPiPlPiMiPiZero_%i_%i",neutralPionMode, trainConfig));
143
144         task->SetIsHeavyIon(2);
145         task->SetIsMC(isMC);
146
147         // Cut Numbers to use in Analysis
148         Int_t numberOfCuts = 1;
149
150         TString *eventCutArray                  = new TString[numberOfCuts];
151         TString *ClusterCutarray                = new TString[numberOfCuts];
152         TString *ConvCutarray                   = new TString[numberOfCuts];
153         TString *PionCutarray                   = new TString[numberOfCuts];
154         TString *NeutralPionCutarray    = new TString[numberOfCuts];
155         TString *MesonCutarray                  = new TString[numberOfCuts];
156         
157         Bool_t doEtaShiftIndCuts = kFALSE;
158         TString stringShift = "";
159
160         // Shifting in pPb direction
161
162         doEtaShiftIndCuts = kTRUE;
163         stringShift = "pPb";
164
165         // EMCAL mode
166         if( trainConfig == 1 ) {
167                 // everything open
168                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "10000040022030000"; PionCutarray[0] = "000010400"; NeutralPionCutarray[0] = "01035030000000"; MesonCutarray[0] = "01035030000000"; 
169         } else if( trainConfig == 2 ) {
170                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, min pt charged pi = 100 MeV
171                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "10000040022030000"; PionCutarray[0] = "002010700"; NeutralPionCutarray[0] = "01035030000000"; MesonCutarray[0] = "01035030000000"; 
172         } else if( trainConfig == 3 ) {
173                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 5 sigma, min pt charged pi = 100 MeV
174                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "10000040022030000"; PionCutarray[0] = "002013700"; NeutralPionCutarray[0] = "01035030000000"; MesonCutarray[0] = "01035030000000"; 
175         } else if( trainConfig == 4 ) {
176                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 4 sigma, min pt charged pi = 100 MeV
177                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "10000040022030000"; PionCutarray[0] = "002016700"; NeutralPionCutarray[0] = "01035030000000"; MesonCutarray[0] = "01035030000000"; 
178         } else if( trainConfig == 5 ) {
179                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 4 sigma, min pt charged pi = 100 MeV
180                 // closing neural pion cuts, 0.1 < M_gamma,gamma < 0.145
181                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "10000040022030000"; PionCutarray[0] = "002016700"; NeutralPionCutarray[0] = "01035031000000"; MesonCutarray[0] = "01035030000000";      
182         } else if( trainConfig == 6 ) {
183                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 4 sigma, min pt charged pi = 100 MeV
184                 // closing neural pion cuts, 0.11 < M_gamma,gamma < 0.145
185                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "10000040022030000"; PionCutarray[0] = "002016700"; NeutralPionCutarray[0] = "01035032000000"; MesonCutarray[0] = "01035030000000";      
186         } else if( trainConfig == 7 ) {
187                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 4 sigma, min pt charged pi = 100 MeV
188                 // closing neural pion cuts, 0.12 < M_gamma,gamma < 0.145
189                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "10000040022030000"; PionCutarray[0] = "002016700"; NeutralPionCutarray[0] = "01035033000000"; MesonCutarray[0] = "01035030000000";      
190         } else if( trainConfig == 8 ) {
191                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 5 sigma, min pt charged pi = 100 MeV
192                 // closing neural pion cuts, 0.12 < M_gamma,gamma < 0.145
193                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "10000040022030000"; PionCutarray[0] = "002013700"; NeutralPionCutarray[0] = "01035033000000"; MesonCutarray[0] = "01035030000000";                      
194         } else if( trainConfig == 9 ) {
195                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, pi+pi- mass Cut at 0.75, min pt charged pi = 100 MeV
196                 // closing neural pion cuts, 0.1 < M_gamma,gamma < 0.145
197                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "10000040022030000"; PionCutarray[0] = "002010702"; NeutralPionCutarray[0] = "01035031000000"; MesonCutarray[0] = "01035030000000"; 
198         }
199
200         // PHOS mode
201         if( trainConfig == 31 ) {
202                 // everything open
203                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "20000030022000000"; PionCutarray[0] = "000010400"; NeutralPionCutarray[0] = "01035030000000"; MesonCutarray[0] = "01035030000000"; 
204         } else if( trainConfig == 32 ) {
205                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, min pt charged pi = 100 MeV
206                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "20000030022000000"; PionCutarray[0] = "002010700"; NeutralPionCutarray[0] = "01035030000000"; MesonCutarray[0] = "01035030000000"; 
207         } else if( trainConfig == 33 ) {
208                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 5 sigma, min pt charged pi = 100 MeV
209                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "20000030022000000"; PionCutarray[0] = "002013700"; NeutralPionCutarray[0] = "01035030000000"; MesonCutarray[0] = "01035030000000"; 
210         } else if( trainConfig == 34 ) {
211                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 4 sigma, min pt charged pi = 100 MeV
212                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "20000030022000000"; PionCutarray[0] = "002016700"; NeutralPionCutarray[0] = "01035030000000"; MesonCutarray[0] = "01035030000000"; 
213         } else if( trainConfig == 35 ) {
214                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 4 sigma, min pt charged pi = 100 MeV
215                 // closing neural pion cuts, 0.1 < M_gamma,gamma < 0.145
216                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "20000030022000000"; PionCutarray[0] = "002016700"; NeutralPionCutarray[0] = "01035031000000"; MesonCutarray[0] = "01035030000000";      
217         } else if( trainConfig == 36 ) {
218                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 4 sigma, min pt charged pi = 100 MeV
219                 // closing neural pion cuts, 0.11 < M_gamma,gamma < 0.145
220                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "20000030022000000"; PionCutarray[0] = "002016700"; NeutralPionCutarray[0] = "01035032000000"; MesonCutarray[0] = "01035030000000";      
221         } else if( trainConfig == 37 ) {
222                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 4 sigma, min pt charged pi = 100 MeV
223                 // closing neural pion cuts, 0.12 < M_gamma,gamma < 0.145
224                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "20000030022000000"; PionCutarray[0] = "002016700"; NeutralPionCutarray[0] = "01035033000000"; MesonCutarray[0] = "01035030000000";      
225         } else if( trainConfig == 38 ) {
226                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, ITS dEdx = \pm 5 sigma, min pt charged pi = 100 MeV
227                 // closing neural pion cuts, 0.12 < M_gamma,gamma < 0.145
228                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "20000030022000000"; PionCutarray[0] = "002013700"; NeutralPionCutarray[0] = "01035033000000"; MesonCutarray[0] = "01035030000000";                      
229         } else if( trainConfig == 39 ) {
230                 // closing charged pion cuts, minimum TPC cluster = 80, TPC dEdx pi = \pm 3 sigma, pi+pi- mass Cut at 0.75, min pt charged pi = 100 MeV
231                 // closing neural pion cuts, 0.1 < M_gamma,gamma < 0.145
232                 eventCutArray[ 0] = "8000011"; ConvCutarray[0] = "002091170008260400000"; ClusterCutarray[0] = "20000030022000000"; PionCutarray[0] = "002010702"; NeutralPionCutarray[0] = "01035031000000"; MesonCutarray[0] = "01035030000000"; 
233         }
234
235         TList *EventCutList = new TList();
236         TList *ConvCutList  = new TList();
237         TList *ClusterCutList  = new TList();
238         TList *NeutralPionCutList = new TList();
239         TList *MesonCutList = new TList();
240         TList *PionCutList  = new TList();
241
242         TList *HeaderList = new TList();
243         TObjString *Header1 = new TObjString("pi0_1");
244         HeaderList->Add(Header1);
245         TObjString *Header3 = new TObjString("eta_2");
246         HeaderList->Add(Header3);
247         
248         EventCutList->SetOwner(kTRUE);
249         AliConvEventCuts **analysisEventCuts = new AliConvEventCuts*[numberOfCuts];
250         ConvCutList->SetOwner(kTRUE);
251         AliConversionPhotonCuts **analysisCuts = new AliConversionPhotonCuts*[numberOfCuts];
252         ClusterCutList->SetOwner(kTRUE);
253         AliCaloPhotonCuts **analysisClusterCuts = new AliCaloPhotonCuts*[numberOfCuts];
254         NeutralPionCutList->SetOwner(kTRUE);
255         AliConversionMesonCuts **analysisNeutralPionCuts   = new AliConversionMesonCuts*[numberOfCuts];
256         MesonCutList->SetOwner(kTRUE);
257         AliConversionMesonCuts **analysisMesonCuts   = new AliConversionMesonCuts*[numberOfCuts];
258         PionCutList->SetOwner(kTRUE);
259         AliPrimaryPionCuts **analysisPionCuts     = new AliPrimaryPionCuts*[numberOfCuts];
260
261         for(Int_t i = 0; i<numberOfCuts; i++){
262                 analysisEventCuts[i] = new AliConvEventCuts();   
263                 analysisEventCuts[i]->InitializeCutsFromCutString(eventCutArray[i].Data());
264                 EventCutList->Add(analysisEventCuts[i]);
265                 analysisEventCuts[i]->SetFillCutHistograms("",kFALSE);
266
267                 analysisCuts[i] = new AliConversionPhotonCuts();
268                 if( ! analysisCuts[i]->InitializeCutsFromCutString(ConvCutarray[i].Data()) ) {
269                                 cout<<"ERROR: analysisCuts [" <<i<<"]"<<endl;
270                                 return 0;
271                 } else {                                
272                         ConvCutList->Add(analysisCuts[i]);
273                         analysisCuts[i]->SetFillCutHistograms("",kFALSE);       
274                 }
275
276                 analysisClusterCuts[i] = new AliCaloPhotonCuts();
277                 if( ! analysisClusterCuts[i]->InitializeCutsFromCutString(ClusterCutarray[i].Data()) ) {
278                                 cout<<"ERROR: analysisClusterCuts [" <<i<<"]"<<endl;
279                                 return 0;
280                 } else {                                
281                         ClusterCutList->Add(analysisClusterCuts[i]);
282                         analysisClusterCuts[i]->SetFillCutHistograms("");                       
283                 }
284
285                 analysisNeutralPionCuts[i] = new AliConversionMesonCuts();
286                 if( ! analysisNeutralPionCuts[i]->InitializeCutsFromCutString(NeutralPionCutarray[i].Data()) ) {
287                         cout<<"ERROR: analysisMesonCuts [ " <<i<<" ] "<<endl;
288                         return 0;
289                 } else {
290                         NeutralPionCutList->Add(analysisNeutralPionCuts[i]);
291                         analysisNeutralPionCuts[i]->SetFillCutHistograms("");
292                 }
293         
294                 analysisMesonCuts[i] = new AliConversionMesonCuts();
295                 if( ! analysisMesonCuts[i]->InitializeCutsFromCutString(MesonCutarray[i].Data()) ) {
296                         cout<<"ERROR: analysisMesonCuts [ " <<i<<" ] "<<endl;
297                         return 0;
298                 } else {
299                         MesonCutList->Add(analysisMesonCuts[i]);
300                         analysisMesonCuts[i]->SetFillCutHistograms("");
301                 }
302                 analysisEventCuts[i]->SetAcceptedHeader(HeaderList);
303                 
304                 TString cutName( Form("%s_%s_%s_%s_%s_%s",eventCutArray[i].Data(), ConvCutarray[i].Data(), ClusterCutarray[i].Data(),PionCutarray[i].Data(),NeutralPionCutarray[i].Data(), MesonCutarray[i].Data() ) );
305                 analysisPionCuts[i] = new AliPrimaryPionCuts();
306                 if( !analysisPionCuts[i]->InitializeCutsFromCutString(PionCutarray[i].Data())) {
307                         cout<< "ERROR:  analysisPionCuts [ " <<i<<" ] "<<endl;
308                         return 0;
309                 } else { 
310                         PionCutList->Add(analysisPionCuts[i]);
311                         analysisPionCuts[i]->SetFillCutHistograms("",kFALSE,cutName); 
312                 }
313         }
314
315         task->SetNeutralPionMode(neutralPionMode);
316         task->SetEventCutList(numberOfCuts,EventCutList);
317         task->SetConversionCutList(ConvCutList);
318         task->SetClusterCutList(ClusterCutList);
319         task->SetNeutralPionCutList(NeutralPionCutList);
320         task->SetMesonCutList(MesonCutList);
321         task->SetPionCutList(PionCutList);
322
323         task->SetMoveParticleAccordingToVertex(kTRUE);
324
325         if(enableQAMesonTask) task->SetDoMesonQA(kTRUE);
326
327         //connect containers
328         AliAnalysisDataContainer *coutput =
329         mgr->CreateContainer(Form("GammaConvNeutralMesonPiPlPiMiPiZero_%i_%i",neutralPionMode, trainConfig), TList::Class(),
330                                                         AliAnalysisManager::kOutputContainer,Form("GammaConvNeutralMesonPiPlPiMiPiZero_%i_%i.root",neutralPionMode, trainConfig));
331
332         mgr->AddTask(task);
333         mgr->ConnectInput(task,0,cinput);
334         mgr->ConnectOutput(task,1,coutput);
335
336         return;
337
338 }