Moving PbPb multiplicity in the new directory structure
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / ChargedHadrons / multPbPb / run.C
1 // TODO:
2 // 1. Check cuts for 2010 (Jochen?)
3 // 2. Run with many centrality bins at once
4 #include <string.h>
5
6 enum { kMyRunModeLocal = 0, kMyRunModeCAF, kMyRunModeGRID};
7
8 TList * listToLoad = new TList();
9
10 TChain * GetAnalysisChain(const char * incollection);
11
12 void run(Char_t* data, Long64_t nev = -1, Long64_t offset = 0, Bool_t debug = kFALSE, Int_t runMode = 0, Bool_t isMC = 0, 
13          Int_t centrBin = 0, const char * centrEstimator = "VOM", Int_t useOtherCentralityCut = 0, Int_t trackMin=0, Int_t trackMax=10000, 
14          const char* option = "",TString customSuffix = "", Int_t workers = -1, Bool_t useSingleBin=kTRUE)
15 {
16   // runMode:
17   //
18   // 0 local 
19   // 1 proof
20
21   if (nev < 0)
22     nev = 1234567890;
23
24   InitAndLoadLibs(runMode,workers,debug);
25
26   // Create the analysis manager
27   mgr = new AliAnalysisManager;
28
29   // Add ESD handler
30   AliESDInputHandler* esdH = new AliESDInputHandler;
31   // Do I need any of this? 
32   esdH->SetInactiveBranches("AliESDACORDE FMD ALIESDTZERO ALIESDZDC AliRawDataErrorLogs CaloClusters Cascades EMCALCells EMCALTrigger ESDfriend Kinks AliESDTZERO ALIESDACORDE MuonTracks TrdTracks");
33   mgr->SetInputEventHandler(esdH);
34
35   if(isMC) {
36     AliMCEventHandler* handler = new AliMCEventHandler;
37     handler->SetPreReadMode(AliMCEventHandler::kLmPreRead);
38     mgr->SetMCtruthEventHandler(handler);
39   }
40
41
42   // If we are running on grid, we need the alien handler
43   if (runMode == kMyRunModeGRID) {
44     // Create and configure the alien handler plugin
45     gROOT->LoadMacro("CreateAlienHandler.C");
46     AliAnalysisGrid *alienHandler = CreateAlienHandler(data, listToLoad, "full", isMC);  
47     if (!alienHandler) {
48       cout << "Cannot create alien handler" << endl;    
49       exit(1);
50     }
51     mgr->SetGridHandler(alienHandler);  
52   }
53
54
55
56   // physics selection
57   gROOT->ProcessLine(".L $ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C");
58   physicsSelectionTask = AddTaskPhysicsSelection(isMC);
59
60   //PID
61   gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPIDResponse.C");
62   AddTaskPIDResponse(isMC); 
63
64
65   // Centrality
66   AliCentralitySelectionTask *taskCentr = new AliCentralitySelectionTask("CentralitySelection");
67   const char * file1 = "$ALICE_ROOT/ANALYSIS/macros/test_AliCentralityBy1D.root";
68   const char * file2 = "$ALICE_ROOT/ANALYSIS/macros/test_AliCentralityByFunction.root";
69   if(isMC) taskCentr-> SetMCInput();
70   taskCentr->SetPass(2);
71
72   // const char * file1 = "$ALICE_ROOT/ANALYSIS/macros/AliCentralityBy1D_LHC10g2a_100.root";
73   // const char * file2 = "$ALICE_ROOT/ANALYSIS/macros/AliCentralityByFunction_LHC10g2a_100.root";
74   // const char * file1 = "$ALICE_ROOT/ANALYSIS/macros/AliCentralityBy1D_137161_GLAU.root";
75   // const char * file2 = "$ALICE_ROOT/ANALYSIS/macros/test_AliCentralityByFunction.root";
76   
77   // taskCentr->SetPercentileFile (file1);
78   // taskCentr->SetPercentileFile2(file2);
79   //FIXME: include back centrality estimator
80   mgr->AddTask(taskCentr);
81   mgr->ConnectInput (taskCentr,0, mgr->GetCommonInputContainer());
82
83   // Create my own centrality selector
84   AliAnalysisMultPbCentralitySelector * centrSelector = new AliAnalysisMultPbCentralitySelector();
85   centrSelector->SetIsMC(isMC);
86   centrSelector->SetCentrTaskFiles(file1,file2); // for bookkeping only
87   centrSelector->SetCentralityBin(centrBin);
88   if (!useSingleBin) centrSelector->SetCentralityBin(0); // FIXME: ok?
89   centrSelector->SetCentralityEstimator(centrEstimator);
90
91   
92   if(useOtherCentralityCut == 1){
93     cout << "Setting centrality by MULT" << endl;
94     centrSelector->SetUseMultRange();
95     centrSelector->SetMultRange(trackMin,trackMax);
96   }
97   if(useOtherCentralityCut == 2){
98     cout << "Setting centrality by V0" << endl;
99     
100     centrSelector->SetUseV0Range();
101     centrSelector->SetMultRange(trackMin,trackMax);
102   }
103   if(useOtherCentralityCut == 3){
104     cout << "Setting centrality by SPD outer" << endl;    
105     centrSelector->SetUseSPDOuterRange();
106     centrSelector->SetMultRange(trackMin,trackMax);
107   }
108
109   // Parse option strings
110   TString optionStr(option);
111   
112   // remove SAVE option if set
113   // This  is copied from a macro by Jan. The reason I kept it is that I may want to pass textual options to the new task at some point
114   Bool_t doSave = kFALSE;
115   TString optionStr(option);
116   if (optionStr.Contains("SAVE"))
117     {
118       optionStr = optionStr(0,optionStr.Index("SAVE")) + optionStr(optionStr.Index("SAVE")+4, optionStr.Length());
119       doSave = kTRUE;
120     }
121
122   AliESDtrackCuts * cuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kFALSE);  
123   TString pathsuffix = "";
124
125   if(!useSingleBin) pathsuffix += "_AllCentr";
126
127   if (optionStr.Contains("DCA")) {
128     delete cuts;
129     //    cuts = AliESDtrackCuts::GetStandardITSPureSATrackCuts2009();
130     cout << ">>>> USING DCA cut" << endl;
131     cuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kFALSE);  
132     pathsuffix+="_DCAcut";
133   }
134
135   if (optionStr.Contains("ITSsa")) {
136     delete cuts;
137     cuts = AliESDtrackCuts::GetStandardITSPureSATrackCuts2009();
138     cout << ">>>> USING ITS sa tracks" << endl;
139     pathsuffix+="_ITSsa";
140   }
141
142   if (optionStr.Contains("TPC")) {
143     delete cuts;
144     cuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
145     cout << ">>>> USING TPC only tracks" << endl;
146     pathsuffix+="_TPC";
147   }
148
149   Bool_t useMCKinematics = isMC;
150   if (optionStr.Contains("NOMCKIN")) {
151     cout << ">>>> Ignoring MC kinematics" << endl;
152     useMCKinematics=kFALSE;
153     pathsuffix+="_NOMCKIN";
154   }
155   
156   AliLog::SetClassDebugLevel("AliESDtrackCuts", AliLog::kDebug);// FIXME
157   cuts->DefineHistograms();
158   
159   // load my task
160   if (useSingleBin) {
161     gROOT->ProcessLine(".L $ALICE_ROOT/PWG0/multPbPb/AddTaskMultPbPbTracks.C");
162     AliAnalysisTaskMultPbTracks * task = AddTaskMultPbPbTracks("multPbPbtracks.root", cuts, centrSelector); 
163     task->SetIsMC(useMCKinematics);
164     task->SetOfflineTrigger(AliVEvent::kMB);
165     if(optionStr.Contains("TPC")) task->SetTPCOnly();
166     if(useMCKinematics) task->GetHistoManager()->SetSuffix("MC");
167     if(customSuffix!=""){
168       cout << "Setting custom suffix: " << customSuffix << endl;    
169       task->GetHistoManager()->SetSuffix(customSuffix);
170     }
171   } else {
172     gROOT->ProcessLine(".L $ALICE_ROOT/PWG0/multPbPb/AddTaskMultPbPbTracksAllCentrality.C");
173     centrSelector->SetUseV0Range(kTRUE);
174     Int_t ncentr = 11;
175    
176     const Float_t minCentr[] = {0 ,79 ,239,559 ,1165,2135,3555,5525,8213 ,12191,15079};
177     const Float_t maxCentr[] = {79,239,559,1165,2135,3555,5525,8213,12191,15079,21000};
178     AliAnalysisTaskMultPbTracks ** tasks = AddTaskMultPbPbTracksAllCentrality("multPbPbtracks.root", cuts, centrSelector, ncentr,minCentr,maxCentr); 
179     for(Int_t icentr = 0; icentr < ncentr; icentr++){
180       tasks[icentr]->Print();
181       cout << "MC KINEMATICS:" << useMCKinematics << endl;
182       
183       tasks[icentr]->SetIsMC(useMCKinematics);
184       tasks[icentr]->SetOfflineTrigger(AliVEvent::kMB);
185       if(optionStr.Contains("TPC")) tasks[icentr]->SetTPCOnly();
186       if(useMCKinematics) tasks[icentr]->GetHistoManager()->SetSuffix("MC");
187       if(customSuffix!=""){
188         cout << "Setting custom suffix: " << customSuffix+long(icentr) << endl;    
189         tasks[icentr]->GetHistoManager()->SetSuffix(customSuffix+long(icentr));
190       } 
191     }    
192   }
193   // Init and run the analy
194   if (!mgr->InitAnalysis()) return;
195
196   mgr->PrintStatus();
197   
198   if (runMode == kMyRunModeLocal ) {
199     // If running in local mode, create chain of ESD files
200     cout << "RUNNING LOCAL, CHAIN" << endl;    
201     TChain * chain = GetAnalysisChain(data);
202     //    chain->Print();
203     mgr->StartAnalysis("local",chain,nev);
204   } else if (runMode == kMyRunModeCAF) {
205     mgr->StartAnalysis("proof",TString(data)+"#esdTree",nev);
206   } else if (runMode == kMyRunModeGRID) {
207     mgr->StartAnalysis("grid");
208   } else {
209     cout << "ERROR: unknown run mode" << endl;        
210   }
211
212   if (!useOtherCentralityCut) {
213     pathsuffix = pathsuffix + "_" + centrEstimator + "_bin_"+long(centrBin);
214   } else if(useOtherCentralityCut==1){
215     pathsuffix = pathsuffix + "_TrackRange_" + long(trackMin) + "_" + long(trackMax);
216   } else if(useOtherCentralityCut==2){
217     pathsuffix = pathsuffix + "_V0Range_" + long(trackMin) + "_" + long(trackMax);
218   } else if(useOtherCentralityCut==3){
219     pathsuffix = pathsuffix + "_SPDOutRange_" + long(trackMin) + "_" + long(trackMax);
220   }
221   pathsuffix += customSuffix;
222
223   if (doSave) MoveOutput(data, pathsuffix.Data());
224
225   
226
227   
228 }
229
230
231 void MoveOutput(const char * data, const char * suffix = ""){
232
233   TString path("output/");
234   path = path + TString(data).Tokenize("/")->Last()->GetName() + suffix;
235   
236   TString fileName = "multPbPbtracks.root";
237   gSystem->mkdir(path, kTRUE);
238   gSystem->Rename(fileName, path + "/" + fileName);
239   for(Int_t ibin = 0; ibin < 20; ibin++){
240     TString fileBin = fileName;
241     fileBin.ReplaceAll(".root",Form("_%2.2d.root",ibin));
242     gSystem->Rename(fileBin, path + "/" + fileBin);    
243   }
244   
245   gSystem->Rename("event_stat.root", path + "/event_stat.root");      
246   Printf(">>>>> Moved files to %s", path.Data());
247 }  
248
249
250
251 TChain * GetAnalysisChain(const char * incollection){
252   // Builds a chain of esd files
253   // incollection can be
254   // - a single root file
255   // - an xml collection of files on alien
256   // - a ASCII containing a list of local root files
257   TChain* analysisChain = 0;
258   // chain
259   analysisChain = new TChain("esdTree");
260   if (TString(incollection).Contains(".root")){
261     analysisChain->Add(incollection);
262   }
263   else if (TString(incollection).Contains("xml")){
264     TGrid::Connect("alien://");
265     TAlienCollection * coll = TAlienCollection::Open (incollection);
266     while(coll->Next()){
267       analysisChain->Add(TString("alien://")+coll->GetLFN());
268     }
269   } else {
270     ifstream file_collect(incollection);
271     TString line;
272     while (line.ReadLine(file_collect) ) {
273       analysisChain->Add(line.Data());
274     }
275   }
276   analysisChain->GetListOfFiles()->Print();
277
278   return analysisChain;
279 }
280
281
282 void InitAndLoadLibs(Int_t runMode=kMyRunModeLocal, Int_t workers=0,Bool_t debug=0) {
283   // Loads libs and par files + custom task and classes
284
285   // Custom stuff to be loaded
286   listToLoad->Add(new TObjString("$ALICE_ROOT/ANALYSIS/AliCentralitySelectionTask.cxx+"));
287   listToLoad->Add(new TObjString("$ALICE_ROOT/PWGPP/background/AliHistoListWrapper.cxx+"));
288   listToLoad->Add(new TObjString("$ALICE_ROOT/PWG0/multPbPb/AliAnalysisMultPbTrackHistoManager.cxx+"));
289   listToLoad->Add(new TObjString("$ALICE_ROOT/PWG0/multPbPb/AliAnalysisMultPbCentralitySelector.cxx+"));
290   listToLoad->Add(new TObjString("$ALICE_ROOT/PWG0/multPbPb/AliAnalysisTaskMultPbTracks.cxx+"));
291
292
293   if (runMode == kMyRunModeCAF)
294     {
295       cout << "Init in CAF mode" << endl;
296     
297       gEnv->SetValue("XSec.GSI.DelegProxy", "2");
298       TProof * p = TProof::Open("alice-caf.cern.ch", workers>0 ? Form("workers=%d",workers) : "1x");
299       //      TProof * p = TProof::Open("skaf.saske.sk", workers>0 ? Form("workers=%d",workers) : "");    
300       p->Exec("TObject *o = gEnv->GetTable()->FindObject(\"Proof.UseMergers\"); gEnv->GetTable()->Remove(o);", kTRUE);
301
302       TProof::Mgr("alice-caf.cern.ch")->SetROOTVersion("VO_ALICE@ROOT::v5-28-00f");
303       //      TProof::Mgr("alice-caf.cern.ch")->SetROOTVersion("5.28/00f");
304       gProof->EnablePackage("VO_ALICE@AliRoot::v4-21-33-AN");
305       gSystem->Load("libCore.so");  
306       gSystem->Load("libTree.so");
307       gSystem->Load("libGeom.so");
308       gSystem->Load("libVMC.so");
309       gSystem->Load("libPhysics.so");
310       gSystem->Load("libSTEERBase");
311       gSystem->Load("libESD");
312       gSystem->Load("libAOD");
313       gSystem->Load("libANALYSIS");
314       gSystem->Load("libOADB");
315       gSystem->Load("libANALYSISalice");   
316
317       // Enable the needed package
318       // gProof->UploadPackage("$ALICE_ROOT/obj/STEERBase");
319       // gProof->EnablePackage("$ALICE_ROOT/obj/STEERBase");
320       // gProof->UploadPackage("$ALICE_ROOT/obj/ESD");
321       // gProof->EnablePackage("$ALICE_ROOT/obj/ESD");
322       // gProof->UploadPackage("$ALICE_ROOT/obj/AOD");
323       // gProof->EnablePackage("$ALICE_ROOT/obj/AOD");
324       // gProof->UploadPackage("$ALICE_ROOT/obj/ANALYSIS");
325       // gProof->EnablePackage("$ALICE_ROOT/obj/ANALYSIS");
326       // gProof->UploadPackage("$ALICE_ROOT/obj/OADB");
327       // gProof->EnablePackage("$ALICE_ROOT/obj/OADB");
328       // gProof->UploadPackage("$ALICE_ROOT/obj/ANALYSISalice");
329       // gProof->EnablePackage("$ALICE_ROOT/obj/ANALYSISalice");
330       // gProof->UploadPackage("$ALICE_ROOT/obj/PWG0base");
331       // gProof->EnablePackage("$ALICE_ROOT/obj/PWG0base");
332       gROOT->ProcessLine(gSystem->ExpandPathName(".include $ALICE_ROOT/PWG0/multPb"));
333       gROOT->ProcessLine(gSystem->ExpandPathName(".include $ALICE_ROOT/PWGPP/background"));
334     }
335   else
336     {
337       cout << "Init in Local or Grid mode" << endl;
338       gSystem->Load("libCore.so");  
339       gSystem->Load("libTree.so");
340       gSystem->Load("libGeom.so");
341       gSystem->Load("libVMC.so");
342       gSystem->Load("libPhysics.so");
343       gSystem->Load("libSTEERBase");
344       gSystem->Load("libESD");
345       gSystem->Load("libAOD");
346       gSystem->Load("libANALYSIS");
347       gSystem->Load("libOADB");
348       gSystem->Load("libANALYSISalice");   
349       // Use AliRoot includes to compile our task
350       gROOT->ProcessLine(".include $ALICE_ROOT/include");
351
352       // gSystem->Load("libVMC");
353       // gSystem->Load("libTree");
354       // gSystem->Load("libSTEERBase");
355       // gSystem->Load("libESD");
356       // gSystem->Load("libAOD");
357       // gSystem->Load("libANALYSIS");
358       // gSystem->Load("libANALYSISalice");
359       // gSystem->Load("libPWG0base");
360     
361       gROOT->ProcessLine(gSystem->ExpandPathName(".include $ALICE_ROOT/PWG0/multPb"));
362       gROOT->ProcessLine(gSystem->ExpandPathName(".include $ALICE_ROOT/PWGPP/background"));
363       //    gROOT->ProcessLine(gSystem->ExpandPathName(".include $ALICE_ROOT/PWGPP/background/"));
364     }
365   // Load helper classes
366   TIterator * iter = listToLoad->MakeIterator();
367   TObjString * name = 0;
368   while (name = (TObjString *)iter->Next()) {
369     gSystem->ExpandPathName(name->String());
370     cout << name->String().Data();
371     if (runMode == kMyRunModeCAF) {
372       gProof->Load(name->String()+(debug?"+g":""));   
373     } else {
374       gROOT->LoadMacro(name->String()+(debug?"+g":""));   
375     }
376   }
377
378 }