hooks for PMD flow analysis
[u/mrichter/AliRoot.git] / PWG2 / FLOW / macros / runFlowTaskCentralityTrain4Candidates.C
1 enum anaModes {mLocal,mLocalPAR,mPROOF,mGrid,mGridPAR};
2 //mLocal: Analyze locally files in your computer using aliroot
3 //mLocalPAR: Analyze locally files in your computer using root + PAR files
4 //mPROOF: Analyze CAF files with PROOF
5 //mGrid: Analyze files on Grid via AliEn plug-in and using precompiled FLOW libraries
6 //       (Remark: When using this mode set also Bool_t bUseParFiles = kFALSE; in CreateAlienHandler.C)
7 //mGridPAR: Analyze files on Grid via AliEn plug-in and using par files for FLOW package
8 //          (Remark: when using this mode set also Bool_t bUseParFiles = kTRUE; in CreateAlienHandler.C)
9  
10 // CENTRALITY DEFINITION
11 //Int_t binfirst = 4;  //where do we start numbering bins
12 //Int_t binlast = 6;  //where do we stop numbering bins
13 //const Int_t numberOfCentralityBins = 9;
14 Int_t binfirst = 0;  //where do we start numbering bins
15 Int_t binlast = 8;  //where do we stop numbering bins
16 const Int_t numberOfCentralityBins = 9;
17 Float_t centralityArray[numberOfCentralityBins+1] = {0.,5.,10.,20.,30.,40.,50.,60.,70.,80.}; // in centrality percentile
18 //Int_t centralityArray[numberOfCentralityBins+1] = {41,80,146,245,384,576,835,1203,1471,10000}; // in terms of TPC only reference multiplicity
19
20 TString commonOutputFileName = "outputCentrality"; // e.g.: result for centrality bin 0 will be in the file "outputCentrality0.root", etc
21
22
23 //void runFlowTaskCentralityTrain(Int_t mode=mLocal, Int_t nRuns = 10, 
24 //Bool_t DATA = kFALSE, const Char_t* dataDir="/Users/snelling/alice_data/Therminator_midcentral", Int_t offset = 0)
25
26 void runFlowTaskCentralityTrain4Candidates(Int_t mode = mGridPAR, Int_t nRuns = 50000000, 
27                  Bool_t DATA = kTRUE, const Char_t* dataDir="/alice/data/LHC10h_000137161_p1_plusplusplus", Int_t offset=0) 
28 //void runFlowTaskCentralityTrain4Candidates(Int_t mode = mLocal, Int_t nRuns = 50000000, 
29 //                              Bool_t DATA = kTRUE, const Char_t* dataDir="./data/", Int_t offset=0) 
30 //void runFlowTaskCentralityTrain(Int_t mode = mGridPAR, Bool_t DATA = kTRUE)
31 {
32   // Time:
33   TStopwatch timer;
34   timer.Start();
35   // Cross-check user settings before starting:
36   //  CrossCheckUserSettings(DATA);
37   // Load needed libraries:
38   LoadLibraries(mode);
39   // Create and configure the AliEn plug-in:
40   if(mode == mGrid || mode == mGridPAR) 
41     {    
42       gROOT->LoadMacro("CreateAlienHandler.C");
43       AliAnalysisGrid *alienHandler = CreateAlienHandler();  
44       if(!alienHandler) return;
45     }
46   // Chains: 
47   if(mode == mLocal || mode == mLocalPAR) {
48 //    TChain* chain = CreateESDChain(dataDir, nRuns, offset);
49     TChain* chain = new TChain();
50     chain->Add("~/alice/datasets/Pb/AliESDs.root/esdTree");
51     //TChain* chain = CreateAODChain(dataDir, nRuns, offset);
52   }
53   
54   // Create analysis manager:
55   AliAnalysisManager *mgr = new AliAnalysisManager("FlowAnalysisManager"); 
56   // Connect plug-in to the analysis manager:
57   if(mode == mGrid || mode == mGridPAR) 
58     { 
59       mgr->SetGridHandler(alienHandler);
60     }
61   
62   // Event handlers:
63   AliVEventHandler* esdH = new AliESDInputHandler;
64   mgr->SetInputEventHandler(esdH);
65   if (!DATA) {
66     AliMCEventHandler *mc = new AliMCEventHandler();
67     mgr->SetMCtruthEventHandler(mc); 
68   }
69
70   // Task to check the offline trigger:
71   gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C"); 
72   AddTaskPhysicsSelection(!DATA);
73
74   //Add the centrality determination task
75   gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskCentrality.C");
76   AddTaskCentrality();
77
78   // Setup analysis per centrality bin:
79   gROOT->LoadMacro("AddTaskFlowCentrality4Candidates.C");
80   for (Int_t i=binfirst; i<binlast+1; i++)
81   {
82     Float_t lowCentralityBinEdge = centralityArray[i];
83     Float_t highCentralityBinEdge = centralityArray[i+1];
84     Printf("\nWagon for centrality bin %i: %.0f-%.0f",i,lowCentralityBinEdge,highCentralityBinEdge);
85     AddTaskFlowCentrality4Candidates( lowCentralityBinEdge,
86                                       highCentralityBinEdge,
87                                       commonOutputFileName,
88                                       0.474, 0.490 );
89     AddTaskFlowCentrality4Candidates( lowCentralityBinEdge,
90                                       highCentralityBinEdge,
91                                       commonOutputFileName,
92                                       0.490, 0.506 );
93     AddTaskFlowCentrality4Candidates( lowCentralityBinEdge,
94                                       highCentralityBinEdge,
95                                       commonOutputFileName,
96                                       0.506, 0.522 );
97   } // end of for (Int_t i=0; i<numberOfCentralityBins; i++)
98
99   // Enable debug printouts:
100   mgr->SetDebugLevel(2);
101   // Run the analysis:
102   if(!mgr->InitAnalysis()) return;  
103   mgr->PrintStatus();
104   if(mode == mLocal || mode == mLocalPAR) {
105     mgr->StartAnalysis("local",chain);
106   } else if(mode == mPROOF) {
107     mgr->StartAnalysis("proof",dataDir,nRuns,offset);
108   } else if(mode == mGrid || mode == mGridPAR) { 
109     mgr->StartAnalysis("grid");
110   }
111
112   // Print real and CPU time used for analysis:
113   timer.Stop();
114   timer.Print();  
115   
116 } // end of void runFlowTaskCentralityTrain(...)
117
118 //===============================================================================================
119 /*
120 void CrossCheckUserSettings(Bool_t bData) 
121 {
122  // Check in this method if the user settings make sense. 
123  if(LYZ1SUM && LYZ2SUM) {cout<<" WARNING: You cannot run LYZ1 and LYZ2 at the same time! LYZ2 needs the output from LYZ1 !!!!"<<endl; exit(0); }
124  if(LYZ1PROD && LYZ2PROD) {cout<<" WARNING: You cannot run LYZ1 and LYZ2 at the same time! LYZ2 needs the output from LYZ1 !!!!"<<endl; exit(0); }
125  if(LYZ2SUM && LYZEP) {cout<<" WARNING: You cannot run LYZ2 and LYZEP at the same time! LYZEP needs the output from LYZ2 !!!!"<<endl; exit(0); }
126  if(LYZ1SUM && LYZEP) {cout<<" WARNING: You cannot run LYZ1 and LYZEP at the same time! LYZEP needs the output from LYZ2 !!!!"<<endl; exit(0); }
127 } // end of void CrossCheckUserSettings()
128 */
129 //===============================================================================================
130
131 void LoadLibraries(const anaModes mode) 
132 {
133   //--------------------------------------
134   // Load the needed libraries most of them already loaded by aliroot
135   //--------------------------------------
136
137   gSystem->Load("libCore");  
138   gSystem->Load("libTree");
139   gSystem->Load("libGeom");
140   gSystem->Load("libVMC");
141   gSystem->Load("libXMLIO");
142   gSystem->Load("libPhysics");
143   gSystem->Load("libXMLParser");
144   gSystem->Load("libProof");
145   
146   if (mode==mLocal || mode==mGrid || mode == mGridPAR || mode == mLocalPAR )
147   {
148     gSystem->Load("libSTEERBase");
149     gSystem->Load("libCDB");
150     gSystem->Load("libRAWDatabase");
151     gSystem->Load("libRAWDatarec");
152     gSystem->Load("libESD");
153     gSystem->Load("libAOD");
154     gSystem->Load("libSTEER");
155     gSystem->Load("libANALYSIS");
156     gSystem->Load("libANALYSISalice");   
157     gSystem->Load("libTOFbase"); 
158     gSystem->Load("libTOFrec"); 
159
160     if (mode == mLocal || mode == mGrid)
161     {
162       gSystem->Load("libPWG2flowCommon"); 
163       gSystem->Load("libPWG2flowTasks"); 
164     }
165     if (mode == mLocalPAR || mode == mGridPAR )
166     {
167       AliAnalysisAlien::SetupPar("PWG2flowCommon");
168       AliAnalysisAlien::SetupPar("PWG2flowTasks");
169     }
170   }
171   
172   //---------------------------------------------------------
173   // <<<<<<<<<< PROOF mode >>>>>>>>>>>>
174   //---------------------------------------------------------
175   else if (mode==mPROOF) {
176     //  set to debug root versus if needed
177     //TProof::Mgr("alicecaf")->SetROOTVersion("v5-24-00a_dbg");
178     //TProof::Mgr("alicecaf")->SetROOTVersion("v5-24-00a");
179     //TProof::Reset("proof://snelling@alicecaf.cern.ch");     
180     // Connect to proof
181     printf("*** Connect to PROOF ***\n");
182     gEnv->SetValue("XSec.GSI.DelegProxy","2");
183     TProof::Open("mkrzewic@alice-caf.cern.ch");
184     //TProof::Open("mkrzewic@skaf.saske.sk");
185      // list the data available
186     //gProof->ShowDataSets("/*/*"); 
187     //gProof->ShowDataSets("/alice/sim/"); //for MC Data
188     //gProof->ShowDataSets("/alice/data/"); //for REAL Data 
189  
190     // Clear the Packages
191     /*    
192     gProof->ClearPackage("STEERBase.par");
193     gProof->ClearPackage("ESD.par");
194     gProof->ClearPackage("AOD.par");
195     */
196     //gProof->ClearPackage("ANALYSIS.par");
197     //gProof->ClearPackage("ANALYSISalice.par");
198     //gProof->ClearPackage("CORRFW.par");
199     
200     gProof->ClearPackage("PWG2flowCommon");
201     gProof->ClearPackage("PWG2flowTasks");
202     
203     // Upload the Packages
204     //gProof->UploadPackage("STEERBase.par");
205     //gProof->UploadPackage("ESD.par");    
206     //gProof->UploadPackage("AOD.par");
207        
208     //gProof->UploadPackage("ANALYSIS.par"); 
209     //gProof->UploadPackage("ANALYSISalice.par");
210     gProof->UploadPackage("CORRFW.par");
211     gProof->UploadPackage("PWG2flowCommon.par");
212     gProof->UploadPackage("PWG2flowTasks.par");
213     gProof->UploadPackage("ALIRECO.par");
214
215     // Enable the Packages 
216     // The global package
217     TList* list = new TList();
218     list->Add(new TNamed("ALIROOT_EXTRA_INCLUDES","RAW:OCDB:STEER:TOF"));
219     gProof->EnablePackage("VO_ALICE@AliRoot::v4-21-07-AN",list);
220     gProof->EnablePackage("ALIRECO");
221     //gProof->EnablePackage("ANALYSIS");
222     //gProof->EnablePackage("ANALYSISalice");
223     //gProof->EnablePackage("CORRFW");
224     gProof->EnablePackage("PWG2flowCommon");
225     gProof->EnablePackage("PWG2flowTasks");
226
227     // Show enables Packages
228     gProof->ShowEnabledPackages();
229   }  
230   
231 } // end of void LoadLibraries(const anaModes mode) 
232
233 // Helper macros for creating chains
234 // from: CreateESDChain.C,v 1.10 jgrosseo Exp
235
236 TChain* CreateESDChain(const char* aDataDir, Int_t aRuns, Int_t offset)
237 {
238   // creates chain of files in a given directory or file containing a list.
239   // In case of directory the structure is expected as:
240   // <aDataDir>/<dir0>/AliESDs.root
241   // <aDataDir>/<dir1>/AliESDs.root
242   // ...
243   
244   if (!aDataDir)
245     return 0;
246   
247   Long_t id, size, flags, modtime;
248   if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
249     {
250       printf("%s not found.\n", aDataDir);
251       return 0;
252     }
253   
254   TChain* chain = new TChain("esdTree");
255   TChain* chaingAlice = 0;
256   
257   if (flags & 2)
258     {
259       TString execDir(gSystem->pwd());
260       TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
261       TList* dirList            = baseDir->GetListOfFiles();
262       Int_t nDirs               = dirList->GetEntries();
263       gSystem->cd(execDir);
264       
265       Int_t count = 0;
266       
267       for (Int_t iDir=0; iDir<nDirs; ++iDir)
268         {
269           TSystemFile* presentDir = (TSystemFile*) dirList->At(iDir);
270           if (!presentDir || !presentDir->IsDirectory() || strcmp(presentDir->GetName(), ".") == 0 || strcmp(presentDir->GetName(), "..") == 0)
271             continue;
272           
273           if (offset > 0)
274             {
275               --offset;
276               continue;
277             }
278           
279           if (count++ == aRuns)
280             break;
281           
282           TString presentDirName(aDataDir);
283           presentDirName += "/";
284           presentDirName += presentDir->GetName();        
285           chain->Add(presentDirName + "/AliESDs.root/esdTree");
286           //  cerr<<presentDirName<<endl;
287         }
288       
289     }
290   else
291     {
292       // Open the input stream
293       ifstream in;
294       in.open(aDataDir);
295       
296       Int_t count = 0;
297       
298       // Read the input list of files and add them to the chain
299       TString esdfile;
300       while(in.good()) {
301         in >> esdfile;
302         if (!esdfile.Contains("root")) continue; // protection
303         
304         if (offset > 0)
305           {
306             --offset;
307             continue;
308           }
309         
310         if (count++ == aRuns)
311           break;
312         
313         // add esd file
314         chain->Add(esdfile);
315       }
316       
317       in.close();
318     }
319   
320   return chain;
321
322 } // end of TChain* CreateESDChain(const char* aDataDir, Int_t aRuns, Int_t offset)
323
324 //===============================================================================================
325
326 TChain* CreateAODChain(const char* aDataDir, Int_t aRuns, Int_t offset)
327 {
328   // creates chain of files in a given directory or file containing a list.
329   // In case of directory the structure is expected as:
330   // <aDataDir>/<dir0>/AliAOD.root
331   // <aDataDir>/<dir1>/AliAOD.root
332   // ...
333   
334   if (!aDataDir)
335     return 0;
336   
337   Long_t id, size, flags, modtime;
338   if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
339     {
340       printf("%s not found.\n", aDataDir);
341       return 0;
342     }
343   
344   TChain* chain = new TChain("aodTree");
345   TChain* chaingAlice = 0;
346   
347   if (flags & 2)
348     {
349       TString execDir(gSystem->pwd());
350       TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
351       TList* dirList            = baseDir->GetListOfFiles();
352       Int_t nDirs               = dirList->GetEntries();
353       gSystem->cd(execDir);
354       
355       Int_t count = 0;
356       
357       for (Int_t iDir=0; iDir<nDirs; ++iDir)
358         {
359           TSystemFile* presentDir = (TSystemFile*) dirList->At(iDir);
360           if (!presentDir || !presentDir->IsDirectory() || strcmp(presentDir->GetName(), ".") == 0 || strcmp(presentDir->GetName(), "..") == 0)
361             continue;
362           
363           if (offset > 0)
364             {
365               --offset;
366               continue;
367             }
368           
369           if (count++ == aRuns)
370             break;
371           
372           TString presentDirName(aDataDir);
373           presentDirName += "/";
374           presentDirName += presentDir->GetName();        
375           chain->Add(presentDirName + "/AliAOD.root/aodTree");
376           // cerr<<presentDirName<<endl;
377         }
378       
379     }
380   else
381     {
382       // Open the input stream
383       ifstream in;
384       in.open(aDataDir);
385       
386       Int_t count = 0;
387       
388       // Read the input list of files and add them to the chain
389       TString aodfile;
390       while(in.good()) {
391         in >> aodfile;
392         if (!aodfile.Contains("root")) continue; // protection
393         
394         if (offset > 0)
395           {
396             --offset;
397             continue;
398           }
399         
400         if (count++ == aRuns)
401           break;
402         
403         // add aod file
404         chain->Add(aodfile);
405       }
406       
407       in.close();
408     }
409   
410   return chain;
411
412 } // end of TChain* CreateAODChain(const char* aDataDir, Int_t aRuns, Int_t offset)
413