- removed the option to run over the first 500 events of the file, (committed by...
[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  * Usage:
10  * <pre>
11  *   aliroot -b -q -l compare_HLT_offline_local.C'("/home/blabla/AliESDs.root","phos")' 2>&1 | tee task.log
12  *   aliroot -b -q -l compare_HLT_offline_local.C'("/home/blabla/AliESDs.root","phos tpc global")' 2>&1 | tee task.log
13  *   aliroot -q compare-HLT-offline-local.C'("alien:///alice/data/2010/LHC10b/000115322/ESDs/pass1/10000115322040.20/AliESDs.root","global")' 2>&1 | tee log
14  * </pre>
15  * 
16  * If alien:// is contained in the name of the file, then the macro connects to the grid to access the file.
17  * 
18  * 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.
19  * The .txt file takes the place of the first argument in that case.
20  *
21  * @ingroup alihlt_qa
22  * @author zbyin@mail.ccnu.edu.cn, Kalliopi.Kanaki@ift.uib.no
23  */
24
25 void compare_HLT_offline_local(TString file, const char* detectorTask="global"){
26
27   TStopwatch timer;
28   timer.Start();
29
30   gSystem->Load("libTree.so");
31   gSystem->Load("libGeom.so");
32   gSystem->Load("libVMC.so");
33   gSystem->Load("libPhysics.so");
34  
35   //----------- Loading the required libraries ---------//
36
37   gSystem->Load("libSTEERBase.so");
38   gSystem->Load("libESD.so");
39   gSystem->Load("libAOD.so");
40   gSystem->Load("libANALYSIS.so");
41   gSystem->Load("libANALYSISalice.so");
42   gSystem->Load("libHLTbase.so");
43  
44   gROOT->ProcessLine(".include $ALICE_ROOT/include");
45
46   
47   Bool_t bTPC=kFALSE, bPHOS=kFALSE, bITS=kFALSE, bGLOBAL=kFALSE, bEMCAL = kFALSE;
48  
49   TString allArgs = detectorTask;
50   TString argument;
51  
52   TObjArray *pTokens = allArgs.Tokenize(" ");
53   if(pTokens){
54      for(int i=0; i<pTokens->GetEntries(); i++){
55          argument=((TObjString*)pTokens->At(i))->GetString();
56          if(argument.IsNull()) continue;
57
58          if(argument.CompareTo("tpc", TString::kIgnoreCase)==0){
59             bTPC = kTRUE;
60             continue;
61          }        
62          if(argument.CompareTo("phos", TString::kIgnoreCase)==0){
63             bPHOS = kTRUE;
64             continue;
65          }         
66          else if(argument.CompareTo("emcal", TString::kIgnoreCase)==0){
67            bEMCAL = kTRUE;
68            continue;
69          }         
70
71           if(argument.CompareTo("its", TString::kIgnoreCase)==0){
72             bITS = kTRUE;
73             continue;
74          }      
75           if(argument.CompareTo("global", TString::kIgnoreCase)==0){
76             bGLOBAL = kTRUE;
77             continue;
78          }        
79          if(argument.CompareTo("all",TString::kIgnoreCase)==0){
80             bTPC    = kTRUE;
81             bPHOS   = kTRUE;
82             bEMCAL  = kTRUE;
83             bITS    = kTRUE;
84             bGLOBAL = kTRUE;    
85             continue;
86          }
87          else break;
88     }
89   }
90     
91   
92   //-------------- Compile the analysis tasks ---------- //
93   if(bTPC) gROOT->LoadMacro("AliAnalysisTaskHLTTPC.cxx+"); 
94   
95   if(bPHOS) {
96     AliHLTSystem * pHLT = AliHLTPluginBase::GetInstance();
97     pHLT->LoadComponentLibraries("libHLTbase");
98     pHLT->LoadComponentLibraries("libAliHLTUtil");
99     pHLT->LoadComponentLibraries("libAliHLTGlobal");
100     gROOT->LoadMacro("AliAnalysisTaskHLTCalo.cxx+"); 
101     gROOT->LoadMacro("AliAnalysisTaskHLTPHOS.cxx+"); 
102   }
103   
104   if(bEMCAL) {
105     AliHLTSystem * pHLT = AliHLTPluginBase::GetInstance();
106     pHLT->LoadComponentLibraries("libHLTbase");
107     pHLT->LoadComponentLibraries("libAliHLTUtil");
108     pHLT->LoadComponentLibraries("libAliHLTGlobal");
109     gROOT->LoadMacro("AliAnalysisTaskHLTCalo.cxx+"); 
110     gROOT->LoadMacro("AliAnalysisTaskHLTEMCAL.cxx+"); 
111   }  
112   
113   if(bITS)    gROOT->LoadMacro("AliAnalysisTaskHLTITS.cxx+");
114   if(bGLOBAL) gROOT->LoadMacro("AliAnalysisTaskHLT.cxx+");
115   
116   //if(!AliAnalysisGrid::CreateToken()) return NULL;
117   
118   if(file.Contains("alien")) TGrid::Connect("alien://");
119     
120   if(file.Contains("AliESDs.root")){
121     TChain *chain = new TChain("esdTree"); 
122     chain->Add(file);
123   }
124   
125   //Constructs chain from filenames in *.txt
126   //on the form $DIR/AliESDs.root
127   else if(file.Contains(".txt")){
128     gROOT->LoadMacro("$ALICE_ROOT/PWG0/CreateESDChain.C");
129     chain=CreateESDChain(file.Data());
130   }
131
132   if(!chain){
133     Printf("Chain is empty");
134     return;
135   }
136
137
138    
139   //-------- Make the analysis manager ---------------//
140  
141   AliAnalysisManager *mgr  = new AliAnalysisManager("TestManager");
142   AliESDInputHandler *esdH = new AliESDInputHandler;
143   esdH->SetReadHLT();
144   esdH->SetReadFriends(kFALSE);
145   mgr->SetInputEventHandler(esdH);  
146   mgr->SetNSysInfo(1000);
147  
148   //-------------- define the tasks ------------//
149   
150   if(bTPC){ 
151      AliAnalysisTaskHLTTPC *taskTPC = new AliAnalysisTaskHLTTPC("offhlt_comparison_TPC");
152      mgr->AddTask(taskTPC);
153      AliAnalysisDataContainer *coutput1 =  mgr->CreateContainer("tpc_histograms", TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-TPC-comparison.root");  
154      mgr->ConnectInput(taskTPC,0,mgr->GetCommonInputContainer());
155      mgr->ConnectOutput(taskTPC,1,coutput1);
156   }
157
158   if(bPHOS){
159      AliAnalysisTaskHLTPHOS *taskPHOS = new AliAnalysisTaskHLTPHOS("offhlt_comparison_PHOS");
160      mgr->AddTask(taskPHOS);
161      AliAnalysisDataContainer *coutput2 =  mgr->CreateContainer("phos_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-PHOS-comparison.root");  
162      mgr->ConnectInput(taskPHOS,0,mgr->GetCommonInputContainer());
163      mgr->ConnectOutput(taskPHOS,1,coutput2);
164   }
165
166   if(bEMCAL){
167      AliAnalysisTaskHLTEMCAL *taskEMCAL = new AliAnalysisTaskHLTEMCAL("offhlt_comparison_EMCAL");
168      mgr->AddTask(taskEMCAL);
169      AliAnalysisDataContainer *coutput5 =  mgr->CreateContainer("emcal_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-EMCAL-comparison.root");  
170      mgr->ConnectInput(taskEMCAL,0,mgr->GetCommonInputContainer());
171      mgr->ConnectOutput(taskEMCAL,1,coutput5);
172   }
173   
174   if(bITS){
175      AliAnalysisTaskHLTITS *taskITS = new AliAnalysisTaskHLTITS("offhlt_comparison_ITS");
176      mgr->AddTask(taskITS);
177      AliAnalysisDataContainer *coutput3 =  mgr->CreateContainer("its_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-ITS-comparison.root");  
178      mgr->ConnectInput(taskITS,0,mgr->GetCommonInputContainer());
179      mgr->ConnectOutput(taskITS,1,coutput3);
180   }
181  
182   if(bGLOBAL){
183      AliAnalysisTaskHLT *taskGLOBAL = new AliAnalysisTaskHLT("offhlt_comparison_GLOBAL");
184      mgr->AddTask(taskGLOBAL);
185      AliAnalysisDataContainer *coutput4 =  mgr->CreateContainer("global_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-GLOBAL-comparison.root");  
186      mgr->ConnectInput(taskGLOBAL,0,mgr->GetCommonInputContainer());
187      mgr->ConnectOutput(taskGLOBAL,1,coutput4);
188   }
189   
190   if (!mgr->InitAnalysis()) return;
191   mgr->PrintStatus();
192   mgr->StartAnalysis("local",chain);
193
194   timer.Stop();
195   timer.Print();
196 }