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