correcting a typo in the include path
[u/mrichter/AliRoot.git] / HLT / QA / tasks / macros / compare-HLT-offline-local.C
1 // $Id$
2 /*
3  * Example macro to run locally an analysis task for comparing the offline
4  * with the HLT esd tree.
5  *
6  * The output is a root file containing the histograms defined in the
7  * analysis task. There is one output file per detector.
8  *
9  * Run without arguments to get a few examples how to use the macro.
10  *
11  * Usage:
12  * <pre>
13  *   aliroot -b -l -q compare-HLT-offline-local.C'("/home/blabla/AliESDs.root","global","./",kTRUE,10)' 2>&1 | tee task.log
14  *   aliroot -b -l -q compare-HLT-offline-local.C'("/home/blabla/AliESDs.root","phos global cb","./",kTRUE,100)' 2>&1 | tee task.log
15  *   aliroot -b -l -q compare-HLT-offline-local.C'("alien:///alice/data/2010/LHC10b/000115322/ESDs/pass1/10000115322040.20/AliESDs.root","global")' 2>&1 | tee log
16  * </pre>
17  * 
18  * If alien:// is contained in the name of the file, then the macro connects to the grid to access the file.
19  * 
20  * In case you want to run over many ESD files, then prepare a list of them in a .txt file and they will be chained for the analysis.
21  * The .txt file takes the place of the first argument in that case.
22  *
23  * @ingroup alihlt_qa
24  * @author Kalliopi.Kanaki@ift.uib.no, Hege.Erdal@student.uib.no
25  */
26
27 void compare_HLT_offline_local(TString file, 
28                                const char* detectorTask="global",
29                                TString taskFolder="$ALICE_ROOT/HLT/QA/tasks/", 
30                                bool fUseHLTTrigger=kFALSE, 
31                                Long64_t nEvents=1234567890
32                               )
33 {
34
35   TStopwatch timer;
36   timer.Start();
37
38   gSystem->Load("libTree.so");
39   gSystem->Load("libGeom.so");
40   gSystem->Load("libVMC.so");
41   gSystem->Load("libPhysics.so");
42  
43   //----------- Loading the required libraries ---------//
44
45   gSystem->Load("libSTEERBase.so");
46   gSystem->Load("libESD.so");
47   gSystem->Load("libAOD.so");
48   gSystem->Load("libANALYSIS.so");
49   gSystem->Load("libANALYSISalice.so");
50   gSystem->Load("libHLTbase.so");
51  
52   gSystem->AddIncludePath("-I$ALICE_ROOT/HLT/BASE -I$ALICE_ROOT/PWG1/TPC -I.");
53   
54   gSystem->Load("libTPCcalib.so");
55   gSystem->Load("libTRDbase.so");
56   gSystem->Load("libTRDrec.so");
57   gSystem->Load("libITSbase.so");
58   gSystem->Load("libITSrec.so");
59   gSystem->Load("libTENDER.so");
60   gSystem->Load("libPWG1.so");
61  
62   gROOT->ProcessLine(".include $ALICE_ROOT/include");
63   //gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C");
64   
65   Bool_t bPHOS = kFALSE, bGLOBAL = kFALSE, bEMCAL = kFALSE, bPWG1 = kFALSE, bD0 = kFALSE, bCB = kFALSE;
66  
67   TString allArgs = detectorTask;
68   TString argument;
69  
70   TObjArray *pTokens = allArgs.Tokenize(" ");
71   if(pTokens){
72     for(int i=0; i<pTokens->GetEntries(); i++){
73       argument=((TObjString*)pTokens->At(i))->GetString();
74       if(argument.IsNull()) continue;
75
76       if(argument.CompareTo("phos", TString::kIgnoreCase)==0){
77         bPHOS = kTRUE;
78         continue;
79       }         
80       else if(argument.CompareTo("emcal", TString::kIgnoreCase)==0){
81         bEMCAL = kTRUE;
82         continue;
83       }         
84       if(argument.CompareTo("global", TString::kIgnoreCase)==0){
85         bGLOBAL = kTRUE;
86         continue;
87       }      
88       if(argument.CompareTo("pwg1", TString::kIgnoreCase)==0){
89         bPWG1 = kTRUE;
90         continue;
91       }
92       if(argument.CompareTo("D0", TString::kIgnoreCase)==0){
93         bD0 = kTRUE;
94         continue;
95       }
96       if(argument.CompareTo("cb", TString::kIgnoreCase)==0){
97         bCB = kTRUE;
98         continue;
99       } 
100       if(argument.CompareTo("all",TString::kIgnoreCase)==0){
101         bPHOS   = kTRUE;
102         bEMCAL  = kTRUE;
103         bGLOBAL = kTRUE; 
104         bD0     = kTRUE;  
105         bCB     = kTRUE; 
106         continue;
107       }
108       else break;
109     }
110   }
111       
112   //-------------- Compile the analysis tasks ---------- //
113   
114   if(bPHOS){
115     //gSystem->Load("libHLTbase");
116     gSystem->Load("libAliHLTUtil");
117     gSystem->Load("libAliHLTGlobal");
118     TString strTask1("AliAnalysisTaskHLTCalo.cxx+");
119     TString strTask2("AliAnalysisTaskHLTPHOS.cxx+");    
120     gROOT->LoadMacro(taskFolder+strTask1); 
121     gROOT->LoadMacro(taskFolder+strTask2); 
122     cout << "\n========= You are loading the following tasks --> "<< (taskFolder+strTask1).Chop()  << " and " <<  (taskFolder+strTask2).Chop() << endl;
123   }
124   
125   if(bEMCAL){
126     //gSystem->Load("libHLTbase");
127     gSystem->Load("libAliHLTUtil");
128     gSystem->Load("libAliHLTGlobal");
129     TString strTask1("AliAnalysisTaskHLTCalo.cxx+");
130     TString strTask2("AliAnalysisTaskHLTEMCAL.cxx+");
131     gROOT->LoadMacro(taskFolder+strTask1); 
132     gROOT->LoadMacro(taskFolder+strTask2); 
133     cout << "\n========= You are loading the following tasks --> "<< (taskFolder+strTask1).Chop()  << " and " <<  (taskFolder+strTask2).Chop() << endl;
134   }  
135   
136   if(bGLOBAL){
137      TString strTask("AliAnalysisTaskHLT.cxx+");
138      gROOT->LoadMacro(taskFolder+strTask);
139      cout << "\n========= You are loading the following task --> "<< (taskFolder+strTask).Chop()  << endl;
140   }
141   if(bD0){
142      TString strTask("AliAnalysisTaskD0Trigger.cxx+");
143      gROOT->LoadMacro(taskFolder+strTask); 
144      cout << "\n========= You are loading the following task --> "<< (taskFolder+strTask).Chop()  << endl;
145   }
146  
147   if(bCB){
148     TString strTask("AliAnalysisTaskHLTCentralBarrel.cxx+");
149     gROOT->LoadMacro(taskFolder+strTask);
150     cout << "\n========= You are loading the following task --> "<< (taskFolder+strTask).Chop()  << endl;
151   }
152   
153   if(bPWG1) gROOT->LoadMacro("$ALICE_ROOT/HLT/QA/tasks/macros/AddTaskPerformance.C");
154    
155   if(file.Contains("alien")) TGrid::Connect("alien://");
156     
157   if(file.Contains("AliESDs.root")){
158     TChain *chain = new TChain("esdTree"); 
159     chain->Add(file);
160   }
161   
162   //Constructs chain from filenames in *.txt
163   //on the form $DIR/AliESDs.root
164   else if(file.Contains(".txt")){
165     gROOT->LoadMacro("$ALICE_ROOT/PWG0/CreateESDChain.C");
166     chain=CreateESDChain(file.Data());
167   }
168   
169   else if(!file){
170     printf("File %s does not exist or is corrupted.\n",file.Data());
171     return;  
172   }
173
174   if(!chain){
175     Printf("Chain is empty");
176     return;
177   }
178    
179   //-------- Make the analysis manager ---------------//
180  
181   AliAnalysisManager *mgr  = new AliAnalysisManager("TestManager");
182   AliESDInputHandler *esdH = new AliESDInputHandler;
183
184   //For the PWG1 task, setting HLT is handled inside AliPerformanceTask.C
185   if(!bPWG1)  esdH->SetReadHLT();
186   esdH->SetReadFriends(kFALSE);
187   mgr->SetInputEventHandler(esdH);  
188   mgr->SetNSysInfo(1000);
189
190   //AliPhysicsSelectionTask *physSelTask = AddTaskPhysicsSelection(kFALSE,kTRUE);
191  
192   //-------------- define the tasks ------------//
193   
194   if(bPHOS){
195     AliAnalysisTaskHLTPHOS *taskPHOS = new AliAnalysisTaskHLTPHOS("offhlt_comparison_PHOS");
196     taskPHOS->SetUseHLTTriggerDecision(fUseHLTTrigger);
197     if(fUseHLTTrigger==kTRUE) printf("\n\nOnly HLT triggered events will be used to fill the distributions for task %s.\n\n", taskPHOS->GetName());
198     mgr->AddTask(taskPHOS);
199     if(fUseHLTTrigger==kFALSE)
200        AliAnalysisDataContainer *coutputPHOS =  mgr->CreateContainer("phos_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-PHOS-comparison.root");  
201     else 
202        AliAnalysisDataContainer *coutputPHOS =  mgr->CreateContainer("phos_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-PHOS-comparison_triggered.root");      
203     mgr->ConnectInput(taskPHOS,0,mgr->GetCommonInputContainer());
204     mgr->ConnectOutput(taskPHOS,1,coutputPHOS);
205   }
206
207   if(bEMCAL){
208     AliAnalysisTaskHLTEMCAL *taskEMCAL = new AliAnalysisTaskHLTEMCAL("offhlt_comparison_EMCAL");
209     taskEMCAL->SetUseHLTTriggerDecision(fUseHLTTrigger);
210     if(fUseHLTTrigger==kTRUE) printf("\n\nOnly HLT triggered events will be used to fill the distributions for task %s.\n\n", taskEMCAL->GetName());
211     mgr->AddTask(taskEMCAL);
212     if(fUseHLTTrigger==kFALSE)
213        AliAnalysisDataContainer *coutputEMCAL =  mgr->CreateContainer("emcal_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-EMCAL-comparison.root");  
214     else
215       AliAnalysisDataContainer *coutputEMCAL =  mgr->CreateContainer("emcal_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-EMCAL-comparison_triggered.root");      
216     mgr->ConnectInput(taskEMCAL,0,mgr->GetCommonInputContainer());
217     mgr->ConnectOutput(taskEMCAL,1,coutputEMCAL);
218   }
219   
220   if(bGLOBAL){
221     AliAnalysisTaskHLT *taskGLOBAL = new AliAnalysisTaskHLT("offhlt_comparison_GLOBAL");
222     taskGLOBAL->SetUseHLTTriggerDecision(fUseHLTTrigger);
223     if(fUseHLTTrigger==kTRUE) printf("\n\nOnly HLT triggered events will be used to fill the distributions for task %s.\n\n", taskGLOBAL->GetName());
224     //taskGLOBAL->SelectCollisionCandidates();
225     mgr->AddTask(taskGLOBAL);
226     if(fUseHLTTrigger==kFALSE)
227       AliAnalysisDataContainer *coutputGLOBAL =  mgr->CreateContainer("global_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-GLOBAL-comparison.root");  
228     else 
229       AliAnalysisDataContainer *coutputGLOBAL =  mgr->CreateContainer("global_histograms",TList::Class(), AliAnalysisManager::kOutputContainer,"HLT-OFFLINE-GLOBAL-comparison_triggered.root");  
230     mgr->ConnectInput(taskGLOBAL,0,mgr->GetCommonInputContainer());
231     mgr->ConnectOutput(taskGLOBAL,1,coutputGLOBAL);
232   }
233
234   if(bPWG1){
235     Bool_t hasMC=kFALSE;  
236     // -- Add Task for HLT and Offline
237     AliPerformanceTask *HLTtpcQA = AddTaskPerformance(hasMC,kFALSE,kTRUE);
238     AliPerformanceTask *tpcQA = AddTaskPerformance(hasMC,kFALSE); 
239     if(!HLTtpcQA || !tpcQA) {
240       Error("RunPerformanceTrain","AliPerformanceTask not created!");
241       return;
242     }
243   }
244   if(bD0){
245     float cuts[7]={0.5,0.04,0.7,0.8,0.05,-0.00025,0.7};
246     AliAnalysisTaskD0Trigger *taskD0 = new AliAnalysisTaskD0Trigger("offhlt_comparison_D0",cuts);
247     mgr->AddTask(taskD0);
248     AliAnalysisDataContainer *coutputD0 =  mgr->CreateContainer("D0_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-D0-comparison.root");  
249     mgr->ConnectInput(taskD0,0,mgr->GetCommonInputContainer());
250     mgr->ConnectOutput(taskD0,1,coutputD0);
251   }  
252   
253   if(bCB){
254      AliAnalysisTaskHLTCentralBarrel *taskCB = new AliAnalysisTaskHLTCentralBarrel("offhlt_comparison_CB");    
255      AliAnalysisDataContainer *coutputCB =  mgr->CreateContainer("esd_thnsparse", TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-CentralBarrel-comparison.root");  
256      mgr->ConnectInput(taskCB,0,mgr->GetCommonInputContainer());
257      mgr->ConnectOutput(taskCB,1,coutputCB);
258   }
259   
260   if (!mgr->InitAnalysis()) return;
261   mgr->PrintStatus();
262   mgr->StartAnalysis("local",chain, nEvents);
263
264   timer.Stop();
265   timer.Print();
266 }
267
268 void compare_HLT_offline_local(){
269   cout << " " << endl;
270   cout << " Usage examples:" << endl;
271   cout << "    compare-HLT-offline-local.C'(file, taskOption, taskFolder, fUseHLTTrigger, nEvents)' 2>&1 | tee log" << endl;
272   cout << "    compare-HLT-offline-local.C'(\"AliESDs.root\",\"global\")' 2>&1 | tee log" << endl;
273   cout << "    compare-HLT-offline-local.C'(\"AliESDs.root\",\"global\",\"./\",kFALSE,nEvents)' 2>&1 | tee log" << endl;
274   cout << "    compare-HLT-offline-local.C'(\"AliESDs.root\",\"global phos cb D0\", \"./\", kTRUE, nEvents)' 2>&1 | tee log" << endl;
275   cout << "    compare-HLT-offline-local.C'(\"files.txt\",\"all\")' 2>&1 | tee log" << endl;
276   cout << "    compare-HLT-offline-local.C'(\"alien:///alice/data/2010/LHC10b/000115322/ESDs/pass1/10000115322040.20/AliESDs.root\",\"global\")' 2>&1 | tee log" << endl;
277   cout << " " << endl;
278 }