- added the ITS analysis task and updated the running macro
[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'("phos")' 2>&1 | tee task.log
12  *   aliroot -b -q -l compare_HLT_offline_local.C'("phos tpc")' 2>&1 | tee task.log
13  *   aliroot -b -q -l compare_HLT_offline_local.C 2>&1 | tee task.log
14  * </pre>
15  *
16  * If no argument is specified, ALL detector tasks are run.
17  *
18  * @ingroup alihlt_tpc
19  * @author zbyin@mail.ccnu.edu.cn, Kalliopi.Kanaki@ift.uib.no
20  */
21
22 void compare_HLT_offline_local(const char* detectorTask="all"){
23  
24   TStopwatch timer;
25   timer.Start();
26
27   gSystem->Load("libTree.so");
28   gSystem->Load("libGeom.so");
29   gSystem->Load("libVMC.so");
30   gSystem->Load("libPhysics.so");
31  
32   //----------- Loading the required libraries ---------//
33
34   gSystem->Load("libSTEERBase.so");
35   gSystem->Load("libESD.so");
36   gSystem->Load("libAOD.so");
37   gSystem->Load("libANALYSIS.so");
38   gSystem->Load("libANALYSISalice.so");
39   gSystem->Load("libHLTbase.so");
40   gROOT->ProcessLine(".include $ALICE_ROOT/include");
41
42   
43   Bool_t bAll=kFALSE, bTPC=kFALSE, bPHOS=kFALSE, bITS=kFALSE;
44  
45   TString allArgs = detectorTask;
46   TString argument;
47  
48   TObjArray *pTokens = allArgs.Tokenize(" ");
49   if(pTokens){
50      for(int i=0; i<pTokens->GetEntries(); i++){
51          argument=((TObjString*)pTokens->At(i))->GetString();
52          if(argument.IsNull()) continue;
53
54          if(argument.CompareTo("tpc", TString::kIgnoreCase)==0){
55             bTPC = kTRUE;
56             continue;
57          } 
58        
59          if(argument.CompareTo("phos", TString::kIgnoreCase)==0){
60             bPHOS = kTRUE;
61             continue;
62          }
63          
64          if(argument.CompareTo("its", TString::kIgnoreCase)==0){
65             bITS = kTRUE;
66             continue;
67          }
68         
69          if(argument.CompareTo("all",TString::kIgnoreCase)==0){
70             bTPC  = kTRUE;
71             bPHOS = kTRUE;
72             bITS  = kTRUE;
73             bAll  = kTRUE;
74             continue;
75          }
76          else break;
77     }
78   }
79     
80   
81   //-------------- Compile the analysis tasks ---------- //
82   if(bTPC)  gROOT->LoadMacro("AliAnalysisTaskHLTTPC.cxx+"); 
83   if(bPHOS) gROOT->LoadMacro("AliAnalysisTaskHLTPHOS.cxx+"); 
84   if(bITS)  gROOT->LoadMacro("AliAnalysisTaskHLTITS.cxx+");
85
86   
87   AliTagAnalysis *TagAna = new AliTagAnalysis("ESD"); 
88   //TagAna->ChainLocalTags("../Tags");
89
90 //   AliRunTagCuts *runCuts = new AliRunTagCuts();
91 //   AliLHCTagCuts *lhcCuts = new AliLHCTagCuts();
92 //   AliDetectorTagCuts *detCuts = new AliDetectorTagCuts();
93 //   AliEventTagCuts *evCuts = new AliEventTagCuts();
94 //   evCuts->SetMultiplicityRange(11,12);  
95   
96
97   TChain *chain = 0x0;
98   //chain = TagAna->QueryTags(runCuts,lhcCuts,detCuts,evCuts);
99   chain = new TChain("esdTree");
100  
101   //gROOT->LoadMacro("$ALICE_ROOT/PWG0/CreateESDChain.C");
102   //chain = CreateESDChain("esd_run84254.txt", 2);
103   
104   //chain->Add("/afs/.alihlt.cern.ch/public/rec/79876/AliESDs.root");
105   //chain->Add("/opt/HLT-public/rec/83683/rec/09000083683000.10/AliESDs.root");
106   
107   //chain->Add("/opt/HLT-public/rec/82762/ESDs/pass1/09000082762002.10/AliESDs.root");
108   //chain->Add("/opt/HLT-public/rec/82762/ESDs/pass1/09000082762004.10/AliESDs.root");
109   //chain->Add("/opt/HLT-public/rec/82762/ESDs/pass1/09000082762005.10/AliESDs.root");
110   //chain->Add("/opt/HLT-public/rec/82762/ESDs/pass1/09000082762006.10/AliESDs.root");
111   
112   //chain->Add("/opt/HLT-public/rec/84254/alien/09000084254009.200/AliESDs.root");
113   //chain->Add("/opt/HLT-public/rec/84254/alien/09000084254009.40/AliESDs.root");
114   //chain->Add("/opt/HLT-public/rec/84254/alien/09000084254009.10/AliESDs.root");
115
116   chain->Add("~/7TeV/115322/10000115322040.110/AliESDs.root");
117   //chain->SetBranchStatus("*Calo*",0);
118
119  
120   //-------- Make the analysis manager ---------------//
121  
122   AliAnalysisManager *mgr  = new AliAnalysisManager("TestManager");
123   AliESDInputHandler *esdH = new AliESDInputHandler;
124   esdH->SetReadHLT();
125   mgr->SetInputEventHandler(esdH);  
126   mgr->SetNSysInfo(1000);
127  
128   //-------------- define the tasks ------------//
129   
130   if(bTPC){ 
131      AliAnalysisTaskHLTTPC *taskTPC = new AliAnalysisTaskHLTTPC("offhlt_comparison_TPC");
132      mgr->AddTask(taskTPC);
133      AliAnalysisDataContainer *coutput1 =  mgr->CreateContainer("tpc_histograms", TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-TPC-comparison.root");  
134      mgr->ConnectInput(taskTPC,0,mgr->GetCommonInputContainer());
135      //mgr->ConnectOutput (taskTPC, 0, mgr->GetCommonOutputContainer());
136      mgr->ConnectOutput(taskTPC,1,coutput1);
137   }
138
139   if(bPHOS){
140      AliAnalysisTaskHLTPHOS *taskPHOS = new AliAnalysisTaskHLTPHOS("offhlt_comparison_PHOS");
141      mgr->AddTask(taskPHOS);
142      AliAnalysisDataContainer *coutput2 =  mgr->CreateContainer("phos_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-PHOS-comparison.root");  
143      mgr->ConnectInput(taskPHOS,0,mgr->GetCommonInputContainer());
144      mgr->ConnectOutput(taskPHOS,1,coutput2);
145   }
146   
147   if(bITS){
148      AliAnalysisTaskHLTITS *taskITS = new AliAnalysisTaskHLTITS("offhlt_comparison_ITS");
149      mgr->AddTask(taskITS);
150      AliAnalysisDataContainer *coutput3 =  mgr->CreateContainer("its_histograms",TList::Class(), AliAnalysisManager::kOutputContainer, "HLT-OFFLINE-ITS-comparison.root");  
151      mgr->ConnectInput(taskITS,0,mgr->GetCommonInputContainer());
152      mgr->ConnectOutput(taskITS,1,coutput3);
153   }
154   
155   if (!mgr->InitAnalysis()) return;
156   mgr->PrintStatus();
157   mgr->StartAnalysis("local",chain);
158
159   timer.Stop();
160   timer.Print();
161 }