fe369f89cb52d5f67f2130a22763129f8681e53b
[u/mrichter/AliRoot.git] / MFT / AODtrain.C
1 // ### Settings that make sense when using the Alien plugin
2 //==============================================================================
3 Int_t       runOnData          = 0;       // Set to 1 if processing real data
4 Int_t       iCollision         = 0;       // 0=pp, 1=Pb-Pb
5 //==============================================================================
6 Bool_t      doCDBconnect        = kFALSE;
7 Bool_t      usePhysicsSelection = kTRUE; // use physics selection
8 Bool_t      useTender           = kFALSE; // use tender wagon
9 Bool_t      useCentrality       = kTRUE; // centrality
10 Bool_t      useV0tender         = kFALSE;  // use V0 correction in tender
11 Bool_t      useDBG              = kTRUE;  // activate debugging
12 Bool_t      useMC               = kTRUE;  // use MC info
13 Bool_t      useKFILTER          = kTRUE;  // use Kinematics filter
14 Bool_t      useTR               = kTRUE;  // use track references
15 Bool_t      useCORRFW           = kFALSE; // do not change
16 Bool_t      useAODTAGS          = kFALSE; // use AOD tags
17 Bool_t      useSysInfo          = kFALSE; // use sys info
18
19 // ### Analysis modules to be included. Some may not be yet fully implemented.
20 //==============================================================================
21 Int_t       iAODhandler        = 1;      // Analysis produces an AOD or dAOD's
22 Int_t       iESDMCLabelAddition= 0;
23 Int_t       iESDfilter         = 1;      // ESD to AOD filter (barrel + muon tracks)
24 Int_t       iMUONcopyAOD       = 1;      // Task that copies only muon events in a separate AOD (PWG3)
25 Int_t       iJETAN             = 0;      // Jet analysis (PWG4)
26 Int_t       iJETANdelta        = 0;      // Jet delta AODs
27 Int_t       iPWGHFvertexing     = 0;      // Vertexing HF task (PWG3)
28 Int_t       iPWGDQJPSIfilter    = 0;      // JPSI filtering (PWG3)
29 Int_t       iPWGHFd2h           = 0;      // D0->2 hadrons (PWG3)
30 Bool_t doPIDResponse  = 1;
31 Bool_t doPIDqa        = 1; //new
32
33 // ### Configuration macros used for each module
34 //==============================================================================
35  TString configPWGHFd2h = (iCollision==0)?"$ALICE_ROOT/PWGHF/vertexingHF/ConfigVertexingHF.C"
36                           :"$ALICE_ROOT/PWGHF/vertexingHF/ConfigVertexingHF_Pb_AllCent.C";
37                                                  
38 // Temporaries.
39 class AliOADBPhysicsSelection;                                                                                                                  
40 AliOADBPhysicsSelection *CreateOADBphysicsSelection();
41 void AODmerge();
42 void AddAnalysisTasks(Int_t);
43 Bool_t LoadCommonLibraries();
44 Bool_t LoadAnalysisLibraries();
45 Bool_t LoadLibrary(const char *);
46 TChain *CreateChain();
47 const char *cdbPath = "raw://";
48 Int_t run_number = 0;
49
50 //______________________________________________________________________________
51 void AODtrain(Int_t merge=0)
52 {
53 // Main analysis train macro.
54    AliLog::SetGlobalDebugLevel(3);
55   if (merge) {
56     TGrid::Connect("alien://");
57     if (!gGrid || !gGrid->IsConnected()) {
58       ::Error("QAtrain", "No grid connection");
59       return;
60     }
61   }
62   // Set temporary merging directory to current one
63   gSystem->Setenv("TMPDIR", gSystem->pwd());
64   // Set temporary compilation directory to current one
65   gSystem->SetBuildDir(gSystem->pwd(), kTRUE);
66    printf("==================================================================\n");
67    printf("===========    RUNNING FILTERING TRAIN   ==========\n");
68    printf("==================================================================\n");
69    printf("=  Configuring analysis train for:                               =\n");
70    if (usePhysicsSelection)   printf("=  Physics selection                                                =\n");
71    if (useTender)    printf("=  TENDER                                                        =\n");
72    if (iESDfilter)   printf("=  ESD filter                                                    =\n");
73    if (iMUONcopyAOD) printf("=  MUON copy AOD                                                 =\n");
74    if (iJETAN)       printf("=  Jet analysis                                                  =\n");
75    if (iJETANdelta)  printf("=     Jet delta AODs                                             =\n");
76    if (iPWGHFvertexing) printf("=  PWGHF vertexing                                                =\n");
77    if (iPWGDQJPSIfilter) printf("=  PWGDQ j/psi filter                                             =\n");
78    if (iPWGHFd2h) printf("=  PWGHF D0->2 hadrons QA                                     =\n");
79
80    // Load common libraries and set include path
81    if (!LoadCommonLibraries()) {
82       ::Error("AnalysisTrain", "Could not load common libraries");
83       return;
84    }
85     
86    // Make the analysis manager and connect event handlers
87    AliAnalysisManager *mgr  = new AliAnalysisManager("Analysis Train", "Production train");
88    if (useSysInfo) mgr->SetNSysInfo(100);
89    // Load analysis specific libraries
90    if (!LoadAnalysisLibraries()) {
91       ::Error("AnalysisTrain", "Could not load analysis libraries");
92       return;
93    }   
94
95    // Create input handler (input container created automatically)
96    // ESD input handler
97    AliESDInputHandler *esdHandler = new AliESDInputHandler();
98    esdHandler->SetNeedField();
99    mgr->SetInputEventHandler(esdHandler);       
100    // Monte Carlo handler
101    if (useMC) {
102       AliMCEventHandler* mcHandler = new AliMCEventHandler();
103       mgr->SetMCtruthEventHandler(mcHandler);
104       mcHandler->SetReadTR(useTR); 
105    }   
106    // AOD output container, created automatically when setting an AOD handler
107    if (iAODhandler) {
108       // AOD output handler
109       AliAODHandler* aodHandler   = new AliAODHandler();
110       aodHandler->SetOutputFileName("AliAOD.root");
111       mgr->SetOutputEventHandler(aodHandler);
112    }
113    // Debugging if needed
114    if (useDBG) mgr->SetDebugLevel(3);
115
116    AddAnalysisTasks(merge, cdbPath);
117    if (merge) {
118       AODmerge();
119       mgr->InitAnalysis();
120       mgr->SetGridHandler(new AliAnalysisAlien);
121       mgr->StartAnalysis("gridterminate",0);
122       return;
123    }   
124    // Run the analysis                                                                                                                     
125    //
126    TChain *chain = CreateChain();
127    if (!chain) return;
128                                                                                                                                                    
129    TStopwatch timer;
130    timer.Start();
131    mgr->SetSkipTerminate(kTRUE);
132    if (mgr->InitAnalysis()) {
133       mgr->PrintStatus();
134       mgr->StartAnalysis("local", chain);
135    }
136    timer.Print();
137 }                                                                                                                                          
138                                                                                                                                             
139 //______________________________________________________________________________                                                           
140 void AddAnalysisTasks(Int_t merge, const char *cdb_location){                                                                                                                                          
141   // Add all analysis task wagons to the train                                                                                               
142    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();                                                                     
143
144   //
145   // Tender and supplies. Needs to be called for every event.
146   //
147    AliAnalysisManager::SetCommonFileName("AODQA.root");
148    if (useTender) {
149       gROOT->LoadMacro("$ALICE_ROOT/TENDER/TenderSupplies/AddTaskTender.C");
150       // IF V0 tender needed, put kTRUE below
151       AliAnalysisTaskSE *tender = AddTaskTender(useV0tender);
152 //      tender->SetDebugLevel(2);
153    }
154
155   // CDB connection
156   //
157   if (doCDBconnect && !useTender) {
158     gROOT->LoadMacro("$ALICE_ROOT/PWGPP/PilotTrain/AddTaskCDBconnect.C");
159     AliTaskCDBconnect *taskCDB = AddTaskCDBconnect(cdb_location, run_number);
160     if (!taskCDB) return;
161     AliCDBManager *cdb = AliCDBManager::Instance();
162     cdb->SetDefaultStorage(cdb_location);
163 //    taskCDB->SetRunNumber(run_number);
164   }    
165
166    if (usePhysicsSelection) {
167    // Physics selection task
168       gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C");
169       mgr->RegisterExtraFile("event_stat.root");
170       AliPhysicsSelectionTask *physSelTask = AddTaskPhysicsSelection(kTRUE);
171 //      AliOADBPhysicsSelection * oadbDefaultPbPb = CreateOADBphysicsSelection();      
172 //      physSelTask->GetPhysicsSelection()->SetCustomOADBObjects(oadbDefaultPbPb,0,0);
173 //      if (!merge) mgr->AddStatisticsTask(AliVEvent::kAny);
174    }
175    // Centrality (only Pb-Pb)
176    if (iCollision && useCentrality) {
177       gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskCentrality.C");
178       AliCentralitySelectionTask *taskCentrality = AddTaskCentrality();
179       taskCentrality->SelectCollisionCandidates(AliVEvent::kAny);
180    }
181
182         if(iESDMCLabelAddition) {
183                 gROOT->LoadMacro("$ALICE_ROOT/PWG/muondep/AddTaskESDMCLabelAddition.C");
184                 AliAnalysisTaskESDMCLabelAddition *esdmclabel = AddTaskESDMCLabelAddition(kFALSE);
185         }
186         
187    if (iESDfilter) {
188       //  ESD filter task configuration.
189       gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/ESDfilter/macros/AddTaskESDFilter.C");     
190       if (iMUONcopyAOD) {
191          printf("Registering delta AOD file\n");
192          mgr->RegisterExtraFile("AliAOD.Muons.root");
193          mgr->RegisterExtraFile("AliAOD.Dimuons.root");
194          AliAnalysisTaskESDfilter *taskesdfilter = AddTaskESDFilter(useKFILTER, kTRUE, kFALSE, kFALSE /*usePhysicsSelection*/,kFALSE,kTRUE,kTRUE,kTRUE,1100,1);
195       } else {
196     AliAnalysisTaskESDfilter *taskesdfilter = AddTaskESDFilter(useKFILTER, kFALSE, kFALSE, kFALSE /*usePhysicsSelection*/,kFALSE,kTRUE,kTRUE,kTRUE,1100,1);
197       }   
198    }   
199
200 // ********** PWG3 wagons ******************************************************           
201    // PWGHF vertexing
202    if (iPWGHFvertexing) {
203       gROOT->LoadMacro("$ALICE_ROOT/PWGHF/vertexingHF/macros/AddTaskVertexingHF.C");
204       if (!iPWGHFd2h) TFile::Cp(gSystem->ExpandPathName(configPWGHFd2h.Data()), "file:ConfigVertexingHF.C");
205       AliAnalysisTaskSEVertexingHF *taskvertexingHF = AddTaskVertexingHF();
206       if (!taskvertexingHF) ::Warning("AnalysisTrainNew", "AliAnalysisTaskSEVertexingHF cannot run for this train conditions - EXCLUDED");
207       else mgr->RegisterExtraFile("AliAOD.VertexingHF.root");
208       taskvertexingHF->SelectCollisionCandidates(0);
209    }   
210       
211    // PWGDQ JPSI filtering (only pp)
212    if (iPWGDQJPSIfilter && (iCollision==0)) {
213       gROOT->LoadMacro("$ALICE_ROOT/PWGDQ/dielectron/macros/AddTaskJPSIFilter.C");
214       AliAnalysisTaskSE *taskJPSIfilter = AddTaskJPSIFilter();
215       if (!taskJPSIfilter) ::Warning("AnalysisTrainNew", "AliAnalysisTaskDielectronFilter cannot run for this train conditions - EXCLUDED");
216       else mgr->RegisterExtraFile("AliAOD.Dielectron.root");
217       taskJPSIfilter->SelectCollisionCandidates(0);
218    }   
219
220    // PWGHF D2h
221    if (iPWGHFd2h) {   
222      gROOT->LoadMacro("$ALICE_ROOT/PWGHF/vertexingHF/AddD2HTrain.C");
223      TFile::Cp(gSystem->ExpandPathName(configPWGHFd2h.Data()), "file:ConfigVertexingHF.C");
224      AddD2HTrain(kFALSE, 1,0,0,0,0,0,0,0,0,0,0);                                 
225    }
226    
227    // ********** PWG4 wagons ******************************************************
228    // Jet analysis
229
230    // Configurations flags, move up?
231    TString kDeltaAODJetName = "AliAOD.Jets.root"; //
232    Bool_t  kIsPbPb = (iCollision==0)?false:true; // can be more intlligent checking the name of the data set
233    TString kDefaultJetBackgroundBranch = "";
234    TString kJetSubtractBranches = "";
235    UInt_t kHighPtFilterMask = 128;// from esd filter
236    UInt_t iPhysicsSelectionFlag = AliVEvent::kMB;
237    if (iJETAN) {
238      gROOT->LoadMacro("$ALICE_ROOT/PWG4/macros/AddTaskJets.C");
239      // Default jet reconstructor running on ESD's
240      AliAnalysisTaskJets *taskjets = AddTaskJets("AOD","UA1",0.4,kHighPtFilterMask,1.,0); // no background subtraction     
241      if (!taskjets) ::Fatal("AnalysisTrainNew", "AliAnalysisTaskJets cannot run for this train conditions - EXCLUDED");
242      if(kDeltaAODJetName.Length()>0) taskjets->SetNonStdOutputFile(kDeltaAODJetName.Data());
243      if (iJETANdelta) {
244         //            AddTaskJetsDelta("AliAOD.Jets.root"); // need to modify this accordingly in the add task jets
245         mgr->RegisterExtraFile(kDeltaAODJetName.Data());
246         TString cTmp("");
247         if(kIsPbPb){
248           // UA1 intrinsic background subtraction
249           taskjets = AddTaskJets("AOD","UA1",0.4,kHighPtFilterMask,1.,2); // background subtraction
250           if(kDeltaAODJetName.Length()>0)taskjets->SetNonStdOutputFile(kDeltaAODJetName.Data());
251        }
252        // SICONE 
253        taskjets = AddTaskJets("AOD","SISCONE",0.4,kHighPtFilterMask,0.15,0); //no background subtration to be done later....                                                                                  
254        if(kDeltaAODJetName.Length()>0)taskjets->SetNonStdOutputFile(kDeltaAODJetName.Data());
255        cTmp = taskjets->GetNonStdBranch();
256        if(cTmp.Length()>0)kJetSubtractBranches += Form("%s ",cTmp.Data());
257          
258        // Add the clusters..
259        gROOT->LoadMacro("$ALICE_ROOT/PWG4/macros/AddTaskJetCluster.C");
260        AliAnalysisTaskJetCluster *taskCl = 0;
261        Float_t fCenUp = 0;
262        Float_t fCenLo = 0;
263        Float_t fTrackEtaWindow = 0.9;
264        taskCl = AddTaskJetCluster("AOD","",kHighPtFilterMask,iPhysicsSelectionFlag,"KT",0.4,0,1, kDeltaAODJetName.Data(),0.15,fTrackEtaWindow,0); // this one is for the background and random jets, random cones with no skip                                                                                 
265        taskCl->SetBackgroundCalc(kTRUE);
266        taskCl->SetNRandomCones(10);
267        taskCl->SetCentralityCut(fCenLo,fCenUp);
268        taskCl->SetGhostEtamax(fTrackEtaWindow);
269        kDefaultJetBackgroundBranch = Form("%s_%s",AliAODJetEventBackground::StdBranchName(),taskCl->GetJetOutputBranch());
270
271        taskCl = AddTaskJetCluster("AOD","",kHighPtFilterMask,iPhysicsSelectionFlag,"ANTIKT",0.4,2,1,kDeltaAODJetName.Data(),0.15);
272        taskCl->SetCentralityCut(fCenLo,fCenUp);
273        if(kIsPbPb)taskCl->SetBackgroundBranch(kDefaultJetBackgroundBranch.Data());
274        taskCl->SetNRandomCones(10);
275        kJetSubtractBranches += Form("%s ",taskCl->GetJetOutputBranch());
276
277        taskCl = AddTaskJetCluster("AOD","",kHighPtFilterMask,iPhysicsSelectionFlag,"ANTIKT",0.2,0,1,kDeltaAODJetName.Data(),0.15);
278        taskCl->SetCentralityCut(fCenLo,fCenUp);
279        if(kIsPbPb)taskCl->SetBackgroundBranch(kDefaultJetBackgroundBranch.Data());
280        kJetSubtractBranches += Form("%s ",taskCl->GetJetOutputBranch());
281          
282        // DO THE BACKGROUND SUBTRACTION
283        if(kIsPbPb&&kJetSubtractBranches.Length()){
284          gROOT->LoadMacro("$ALICE_ROOT/PWG4/macros/AddTaskJetBackgroundSubtract.C");
285          AliAnalysisTaskJetBackgroundSubtract *taskSubtract = 0;
286          taskSubtract = AddTaskJetBackgroundSubtract(kJetSubtractBranches,1,"B0","B%d");
287          taskSubtract->SetBackgroundBranch(kDefaultJetBackgroundBranch.Data());
288          if(kDeltaAODJetName.Length()>0)taskSubtract->SetNonStdOutputFile(kDeltaAODJetName.Data());
289        }
290      } 
291    }
292 }
293
294 //______________________________________________________________________________
295 Bool_t LoadCommonLibraries()
296 {
297
298   // Load common analysis libraries.
299    if (!gSystem->Getenv("ALICE_ROOT")) {
300       ::Error("AnalysisTrainNew.C::LoadCommonLibraries", "Analysis train requires that analysis libraries are compiled with a local AliRoot"); 
301       return kFALSE;
302    }   
303    Bool_t success = kTRUE;
304    // Load framework classes. Par option ignored here.
305    success &= LoadLibrary("libSTEERBase.so");
306    success &= LoadLibrary("libESD.so");
307    success &= LoadLibrary("libAOD.so");
308    success &= LoadLibrary("libANALYSIS.so");
309    success &= LoadLibrary("libOADB.so");
310    success &= LoadLibrary("libANALYSISalice.so");
311    success &= LoadLibrary("libESDfilter.so");
312    success &= LoadLibrary("libCORRFW.so");
313    gROOT->ProcessLine(".include $ALICE_ROOT/include");
314    if (success) {
315       ::Info("AnalysisTrainNew.C::LoadCommodLibraries", "Load common libraries:    SUCCESS");
316       ::Info("AnalysisTrainNew.C::LoadCommodLibraries", "Include path for Aclic compilation:\n%s",
317               gSystem->GetIncludePath());
318    } else {           
319       ::Info("AnalysisTrainNew.C::LoadCommodLibraries", "Load common libraries:    FAILED");
320    }   
321    return success;
322 }
323
324 //______________________________________________________________________________
325 Bool_t LoadAnalysisLibraries()
326 {
327
328   // Load common analysis libraries.
329   if (useTender || doCDBconnect) {
330     if (!LoadLibrary("Tender") ||
331         !LoadLibrary("TenderSupplies")) return kFALSE;
332   }       
333   // CDBconnect
334   if (doCDBconnect && !useTender) {
335     if (!LoadLibrary("PWGPP")) return kFALSE;
336   }
337   
338   if (iESDfilter || iPWGMuonTrain) {
339     if (!LoadLibrary("PWGmuon")) return kFALSE;
340   }   
341   if (iESDMCLabelAddition) {
342     if (!LoadLibrary("PWGmuondep")) return kFALSE;
343   }
344   // JETAN
345   if (iJETAN) {
346     if (!LoadLibrary("JETAN")) return kFALSE;
347   }
348   if (iJETANdelta) {
349     if (!LoadLibrary("JETAN") ||
350         !LoadLibrary("CGAL") ||
351         !LoadLibrary("fastjet") ||
352         !LoadLibrary("siscone") ||
353         !LoadLibrary("SISConePlugin") ||
354         !LoadLibrary("FASTJETAN")) return kFALSE;
355   }     
356   // PWG3 Vertexing HF
357   if (iPWGHFvertexing || iPWGHFd2h) {
358     if (!LoadLibrary("PWGflowBase") ||
359         !LoadLibrary("PWGflowTasks") ||
360         !LoadLibrary("PWGHFvertexingHF")) return kFALSE;
361   }    
362   //    if (iPWGHFvertexing || iPWG3d2h) {
363   //     if (!LoadLibrary("PWG3base") ||
364   //         !LoadLibrary("PWGHFvertexingHF")) return kFALSE;
365   //  }   
366   // PWG3 dielectron
367   if (iPWGDQJPSIfilter) {
368     if (!LoadLibrary("PWGDQdielectron")) return kFALSE;
369   }   
370   
371   ::Info("AnalysisTrainNew.C::LoadAnalysisLibraries", "Load other libraries:   SUCCESS");
372   return kTRUE;
373 }
374
375 //______________________________________________________________________________
376 Bool_t LoadLibrary(const char *module)
377 {
378 // Load a module library in a given mode. Reports success.
379    Int_t result;
380    TString mod(module);
381    if (!mod.Length()) {
382       ::Error("AnalysisTrainNew.C::LoadLibrary", "Empty module name");
383       return kFALSE;
384    }   
385    // If a library is specified, just load it
386    if (mod.EndsWith(".so")) {
387       mod.Remove(mod.Index(".so"));
388       result = gSystem->Load(mod);
389       if (result < 0) {
390          ::Error("AnalysisTrainNew.C::LoadLibrary", "Could not load library %s", module);
391          return kFALSE;
392       }
393       return kTRUE;
394    } 
395    // Check if the library is already loaded
396    if (strlen(gSystem->GetLibraries(module, "", kFALSE)) > 0) return kTRUE;    
397    result = gSystem->Load(Form("lib%s", module));
398    if (result < 0) {
399       ::Error("AnalysisTrainNew.C::LoadLibrary", "Could not load module %s", module);
400       return kFALSE;
401    }
402    return kTRUE;
403 }           
404
405
406 //______________________________________________________________________________
407 TChain *CreateChain()
408 {
409 // Create the input chain
410    chain = new TChain("esdTree");
411    if (gSystem->AccessPathName("AliESDs.root")) 
412       ::Error("AnalysisTrainNew.C::CreateChain", "File: AliESDs.root not in ./data dir");
413    else 
414        chain->Add("AliESDs.root");
415    if (chain->GetNtrees()) return chain;
416    return NULL;
417 }   
418
419 //______________________________________________________________________________
420 void AODmerge()
421 {
422 // Merging method. No staging and no terminate phase.
423   TStopwatch timer;
424   timer.Start();
425   TString outputDir = "wn.xml";
426   TString outputFiles = "EventStat_temp.root,AliAOD.root,AliAOD.Muons.root";
427   TString mergeExcludes = "";
428   TObjArray *list = outputFiles.Tokenize(",");
429   TIter *iter = new TIter(list);
430   TObjString *str;
431   TString outputFile;
432   Bool_t merged = kTRUE;
433   while((str=(TObjString*)iter->Next())) {
434     outputFile = str->GetString();
435     // Skip already merged outputs
436     if (!gSystem->AccessPathName(outputFile)) {
437        printf("Output file <%s> found. Not merging again.",outputFile.Data());
438        continue;
439     }
440     if (mergeExcludes.Contains(outputFile.Data())) continue;
441     merged = AliAnalysisAlien::MergeOutput(outputFile, outputDir, 10, 0);
442     if (!merged) {
443        printf("ERROR: Cannot merge %s\n", outputFile.Data());
444        return;
445     }
446   }
447   // all outputs merged, validate
448   ofstream out;
449   out.open("outputs_valid", ios::out);
450   out.close();
451   timer.Print();
452 }