4d68722d192acc85b999e402e1a168f2abb5b794
[u/mrichter/AliRoot.git] / PWGCF / FLOW / Documentation / examples / manual / runFlowOnDataExample.C
1 void runFlowOnDataExample() {
2     // example macro of running a flow analysis on local 
3     // data
4     //
5     // all steps are explained in detail in
6     // chapter 3.4.3 of the flow package manual
7     // $ALICE_ROOT/PWGCF/FLOW/Documentation/FlowPackageManual.pdf
8     // which is also available on the twiki page
9     // https://twiki.cern.ch/twiki/bin/viewauth/ALICE/FlowPAGFlowPackageManual
10     //
11     // author: Redmer Alexander Bertens, Utrecht University
12     // rbertens@cern.ch , rbertens@nikhef.nl , r.a.bertens@uu.nl
13
14   // load libraries
15   gSystem->Load("libCore.so");        
16   gSystem->Load("libGeom.so");
17   gSystem->Load("libVMC.so");
18   gSystem->Load("libPhysics.so");
19   gSystem->Load("libTree.so");
20   gSystem->Load("libSTEERBase.so");
21   gSystem->Load("libESD.so");
22   gSystem->Load("libAOD.so");
23   gSystem->Load("libANALYSIS.so");
24   gSystem->Load("libANALYSISalice.so");
25   gSystem->Load("libEventMixing.so");
26   gSystem->Load("libCORRFW.so");
27   gSystem->Load("libPWGTools.so");
28   gSystem->Load("libPWGCFebye.so");
29   gSystem->Load("libPWGflowBase.so");
30   gSystem->Load("libPWGflowTasks.so");
31
32     // create the analysis manager
33   AliAnalysisManager* mgr = new AliAnalysisManager("MyManager");
34
35   // create a tchain which will point to an aod tree
36   TChain* chain = new TChain("aodTree");
37   // add a few files to the chain (change this so that your local files are added)
38   chain->Add("/home/rbertens/Documents/CERN/ALICE_DATA/data/2010/LHC10h/000139510/ESDs/pass2/AOD086/0003/AliAOD.root");
39   chain->Add("/home/rbertens/Documents/CERN/ALICE_DATA/data/2010/LHC10h/000139510/ESDs/pass2/AOD086/0003/AliAOD.root");
40   chain->Add("/home/rbertens/Documents/CERN/ALICE_DATA/data/2010/LHC10h/000139510/ESDs/pass2/AOD086/0004/AliAOD.root");
41   chain->Add("/home/rbertens/Documents/CERN/ALICE_DATA/data/2010/LHC10h/000139510/ESDs/pass2/AOD086/0005/AliAOD.root");
42   chain->Add("/home/rbertens/Documents/CERN/ALICE_DATA/data/2010/LHC10h/000139510/ESDs/pass2/AOD086/0006/AliAOD.root");
43   chain->Add("/home/rbertens/Documents/CERN/ALICE_DATA/data/2010/LHC10h/000139510/ESDs/pass2/AOD086/0007/AliAOD.root");
44   chain->Add("/home/rbertens/Documents/CERN/ALICE_DATA/data/2010/LHC10h/000139510/ESDs/pass2/AOD086/0008/AliAOD.root");
45   chain->Add("/home/rbertens/Documents/CERN/ALICE_DATA/data/2010/LHC10h/000139510/ESDs/pass2/AOD086/0009/AliAOD.root");
46   chain->Add("/home/rbertens/Documents/CERN/ALICE_DATA/data/2010/LHC10h/000139510/ESDs/pass2/AOD086/0010/AliAOD.root");
47   // create an input handler
48   AliVEventHandler* inputH = new AliAODInputHandler();
49   // and connect it to the manager
50   mgr->SetInputEventHandler(inputH);
51
52    // the manager is static, so get the existing manager via the static method
53   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
54   if (!mgr) {
55       printf("No analysis manager to connect to!\n");
56       return NULL;
57   }
58         
59   // just to see if all went well, check if the input event handler has been connected
60   if (!mgr->GetInputEventHandler()) {
61       printf("This task requires an input event handler!\n");
62       return NULL;
63     }
64
65   // create instance of the class. because possible qa plots are added in a second ouptut slot,
66   // the flow analysis task must know if you want to save qa plots at the time of class construction
67   Bool_t doQA = kTRUE;
68   // craete instance of the class
69   AliAnalysisTaskFlowEvent* taskFE = new AliAnalysisTaskFlowEvent("FlowEventTask", "", doQA);
70   // add the task to the manager
71   mgr->AddTask(taskFE);
72   // set the trigger selection
73   taskFE->SelectCollisionCandidates(AliVEvent::kMB);
74
75     // define the event cuts object
76   AliFlowEventCuts* cutsEvent = new AliFlowEventCuts("EventCuts");
77   // configure some event cuts, starting with centrality
78   cutsEvent->SetCentralityPercentileRange(20., 30.);
79   // method used for centrality determination
80   cutsEvent->SetCentralityPercentileMethod(AliFlowEventCuts::kV0);
81   // vertex-z cut
82   cutsEvent->SetPrimaryVertexZrange(-10.,10.);
83   // enable the qa plots
84   cutsEvent->SetQA(doQA);
85   // explicit multiplicity outlier cut
86   cutsEvent->SetCutTPCmultiplicityOutliersAOD(kTRUE);
87   cutsEvent->SetLHC10h(kTRUE);
88   
89   
90   // and, last but not least, pass these cuts to your flow event task
91   taskFE->SetCutsEvent(cutsEvent);
92
93   //create the track cuts object using a static function of AliFlowTrackCuts
94   AliFlowTrackCuts* cutsRP = AliFlowTrackCuts::GetAODTrackCutsForFilterBit(1, "RP_cuts");
95   // specify the pt range
96   cutsRP->SetPtRange(0.2, 5.);
97   // specify eta range
98   cutsRP->SetEtaRange(-0.8, 0.8);
99   // specify track type
100   cutsRP->SetParamType(AliFlowTrackCuts::kAODFilterBit);
101   // enable saving qa histograms
102   cutsRP->SetQA(kTRUE);
103
104     //create the track cuts object using a static function of AliFlowTrackCuts
105   AliFlowTrackCuts* cutsPOI = AliFlowTrackCuts::GetAODTrackCutsForFilterBit(1, "pion selection");
106   // specify the pt range
107   cutsPOI->SetPtRange(0.2, 5.);
108   // specify eta range
109   cutsPOI->SetEtaRange(-0.8, 0.8);
110   // specify the track type
111   cutsPOI->SetParamType(AliFlowTrackCuts::kAODFilterBit);
112   // enable saving qa histograms
113   cutsPOI->SetQA(kTRUE);
114
115     // which particle do we want to identify ?
116   AliPID::EParticleType particleType=AliPID::kPion;
117   // specify the pid method that we want to use  
118   AliFlowTrackCuts::PIDsource sourcePID=AliFlowTrackCuts::kTOFbayesian;
119   // define the probability (between 0 and 1) 
120   Double_t probability = .9;
121   // pass these variables to the track cut object
122   cutsPOI->SetPID(particleType, sourcePID, probability);
123   // the bayesian pid routine uses priors tuned to an average centrality
124   cutsPOI->SetPriors(35.);
125
126     // connect the RP's to the flow event task
127   taskFE->SetCutsRP(cutsRP);
128   // connect the POI's to the flow event task
129   taskFE->SetCutsPOI(cutsPOI);
130
131     // get the default name of the output file ("AnalysisResults.root")
132   TString file = AliAnalysisManager::GetCommonFileName();
133   // get the common input container from the analysis manager
134   AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
135    // create a data container for the output of the flow event task  
136    // the output of the task is the AliFlowEventSimle class which will
137    // be passed to the flow analysis tasks. note that we use a kExchangeContainer here,
138    // which exchanges data between classes of the analysis chain, but is not
139    // written to the output file
140   AliAnalysisDataContainer *coutputFE = mgr->CreateContainer(
141       "FlowEventContainer",
142       AliFlowEventSimple::Class(),
143       AliAnalysisManager::kExchangeContainer);
144   // connect the input data to the flow event task
145   mgr->ConnectInput(taskFE,0,cinput);
146   // and connect the output to the flow event task
147   mgr->ConnectOutput(taskFE,1,coutputFE);
148   // create an additional container for the QA output of the flow event task
149   // the QA histograms will be stored in a sub-folder of the output file called 'QA'
150   TString taskFEQAname = file;
151   taskFEQAname += ":QA";
152   AliAnalysisDataContainer* coutputFEQA = mgr->CreateContainer(
153       "FlowEventContainerQA",
154        TList::Class(),
155        AliAnalysisManager::kOutputContainer,
156        taskFEQAname.Data()       
157        );
158   // and connect the qa output container to the flow event. 
159   // this container will be written to the output file
160   mgr->ConnectOutput(taskFE,2,coutputFEQA);
161
162     // declare necessary pointers
163   AliAnalysisDataContainer *coutputQC[3];
164   AliAnalysisTaskQCumulants *taskQC[3];
165
166   // the tasks will be creaated and added to the manager in a loop
167   for(Int_t i = 0; i < 3; i++) {
168       // create the flow analysis tasks
169       taskQC[i] = new AliAnalysisTaskQCumulants(Form("TaskQCumulants_n=%i", i+2));
170       // set the triggers 
171       taskQC[i]->SelectCollisionCandidates(AliVEvent::kMB);
172       // and set the correct harmonic n
173       taskQC[i]->SetHarmonic(i+2);
174
175       // connect the task to the analysis manager
176       mgr->AddTask(taskQC[i]);
177
178       // create and connect the output containers
179       TString outputQC = file;
180       // create a sub-folder in the output file for each flow analysis task's output
181       outputQC += Form(":QC_output_for_n=%i", i+2);
182       /// create the output containers
183       coutputQC[i] = mgr->CreateContainer(
184           outputQC.Data(),
185           TList::Class(),
186           AliAnalysisManager::kOutputContainer,
187           outputQC);
188       // connect the output of the flow event task to the flow analysis task
189       mgr->ConnectInput(taskQC[i], 0, coutputFE);
190       // and connect the output of the flow analysis task to the output container
191       // which will be written to the output file
192       mgr->ConnectOutput(taskQC[i], 1, coutputQC[i]);
193   }
194   
195   // check if we can initialize the manager
196   if(!mgr->InitAnalysis()) return;   
197   // print the status of the manager to screen 
198   mgr->PrintStatus();
199   // print to screen how the analysis is progressing
200   mgr->SetUseProgressBar(1, 25);
201   // start the analysis locally, reading the events from the tchain
202   mgr->StartAnalysis("local", chain);
203
204   
205 }