]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWGCF/FLOW/macros/AddTwoParticleResonanceFlowTask.C
test task to replace phi and K* analysis
[u/mrichter/AliRoot.git] / PWGCF / FLOW / macros / AddTwoParticleResonanceFlowTask.C
1 /////////////////////////////////////////////////////////////////////////////////////////////
2 //
3 // AddTaskPhiFlow macro
4 // Author: Redmer A. Bertens, Utrecht University, 2012
5 //         rbertens@cern.ch, r.a.bertens@uu.nl
6 //         Commented where necessary
7 /////////////////////////////////////////////////////////////////////////////////////////////
8
9 class AliAnalysisDataContainer;
10 class AliFlowTrackCuts;
11 class AliFlowTrackSimpleCuts;
12 class AliFlowEventCuts;
13 class AliFlowEventSimpleCuts;
14 class AliAnalysisDataContainer;
15
16 AliAnalysisTwoParticleResonanceFlowTask* AddTwoParticleResonanceFlowTask(Bool_t SP = 1, // select flow analysis methods
17                                        Bool_t SPSUB = 0,
18                                        Bool_t QC = 0,
19                                        Bool_t EP = 0,
20                                        Bool_t EP3sub = 0,
21                                        Bool_t VZERO_SP = 0, // use vzero sp method
22                                        Float_t centrMin = 20., // centrality selection
23                                        Float_t centrMax = 30.,
24                                        Float_t ITSsigma = 0., // pid mode (see task implementation)
25                                        Float_t ITSrange = 0.,
26                                        Float_t TPCcontrol = 1.,
27                                        Float_t TPCsigma = 3.,
28                                        Float_t TPCrange = 0.,
29                                        Float_t ITScontrol = -1.,
30                                        Float_t Bpurity = 0.3,
31                                        TString suffixName = "UniqueID", // unique id for output objects
32                                        Bool_t bCentralTrigger = kTRUE, // trigger selection
33                                        Float_t EtaGap = 0., // eta gap for SPSUB
34                                        TString DCA = "pt", // dca mode (see task implementation)
35                                        Int_t harm = 2, // harmonic vn
36                                        UInt_t poi_filter = 1, // aod filterbits
37                                        UInt_t rp_filter = 1,
38                                        Float_t minMass = .71, // min mass for flow analysis
39                                        Float_t maxMaxx = 1.13, // max mass for flow analysis
40                                        Int_t speciesA = 3, // species a
41                                        Int_t speciesB = 2, // species b
42                                        Int_t chargeA = -1, // charge a
43                                        Int_t chargeB = 1, // charge b
44                                        Float_t massA = 4.93676999999999977e-01, // mass species a
45                                        Float_t massB = 1.39570e-01, // mass species b
46                                        Bool_t PhiMinusPsiMethod = 0, // use phi minus psi method
47                                        Bool_t event_mixing = kTRUE, // use event mixing        
48                                        Bool_t highPtMode = kFALSE, // use with caution !!! disables invariant mass fit method
49                                        Float_t deltaPhiMass = 0.0003, // dM in which to look for phi 
50                                        Float_t POIPtMax = 5., // max pt of daughterp particles
51                                        Bool_t shrinkSP = kTRUE, // shrink output
52                                        Bool_t debug = kTRUE) // macro debug mode, for task's debug mode see header
53 {
54    /* some defaults that are used frequently:
55     * for phi:
56     *  - minMass = 0.99
57     *  - maxMass = 1.092
58     * for kstar
59     *  - minMass = 0.71
60     *  - maxMass = 1.13
61     * kaon mass = 4.93676999999999977e-01
62     * pion mass = 1.39570e-01
63     */
64    // some defaults that have been removed as function arguments (august 30 2012)
65    Float_t POIPtMin = 0.2;  // pt of daughters
66    Float_t deltaDip = 0.;
67    Float_t deltaDipMaxPt = 0.;
68    Bool_t TPCStandAloneTracks = kFALSE; // deprecated
69    Bool_t useGlobalRPCuts = kTRUE; // deprecated
70    Float_t vertexZ = 10.;
71    Float_t POIEtaMin = -0.8;
72    Float_t POIEtaMax = 0.8;
73    // start of macro
74    Float_t PIDconfig[] = {ITSsigma, ITSrange, TPCcontrol, TPCsigma, TPCrange, ITScontrol, Bpurity};
75    // main function, create and add tasks
76    if(debug) cout << " === Adding Task PhiFlow === " << endl;
77    // set up main output container's name
78    TString fileName = AliAnalysisManager::GetCommonFileName();
79    fileName += ":PhiReconstruction";
80    suffixName += Form("%.0f", centrMin);
81    suffixName += Form("%.0f", centrMax);
82    fileName+=suffixName;
83    if(debug) cout << "    --> Reconstruction data container: " << fileName << endl;
84    // check validity of arguments
85    if((!SPSUB)&&(EtaGap > 0.)) {
86        if(debug) cout << " --> Fatal error: Eta gap is introduced but method SPSUB is not enabled !!! <-- " << endl;
87        return 0x0;
88    }
89    else if ((QC||SP||EP||EP3sub)&&(EtaGap > 0.)) {
90        if(debug) cout << " --> Fatal error: one or more of the flow analyses is not compatible with the Eta Gap !!! <--" << endl;
91        return 0x0;
92    }
93    else if(EP3sub) {
94        if(highPtMode) {
95            if(debug) cout << " --> Can't launch 3 subevent method for high pt analysis, exiting ... <--" << endl;
96            return 0x0;
97        }
98        if(debug) cout << " --> Starting 3 subevent plane method - try at your own risk !!! <-- " << endl;
99        if(!(harm!=2||harm!=3)) {
100            if(debug) cout << " --> Fatal error: can only return v2 and v3 with 3 subevent method! " << endl;
101            return 0x0;
102        }
103    }
104    if((EP3sub||PhiMinusPsiMethod)&&debug) {
105        gROOT->LoadMacro("$ALICE_ROOT/PWGCF/FLOW/macros/AddTaskVZERO.C");
106        AddTaskVZERO(0,0,0,0);
107    }
108    // get the manager and input event handler
109    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
110    if (!mgr) {
111       if(debug) cout << " Fatal error: no analysis manager found! " << endl;
112       return 0x0;
113    }
114    if (!mgr->GetInputEventHandler()) {
115       if(debug) cout << " Fatal error: no imput event handler found!" << endl;
116       return 0x0;
117    }
118    if(EP3sub&&debug) { // don't use this on train ! this is why it's only enabled in macro debug mode
119          gROOT->LoadMacro("$ALICE_ROOT/PWGCF/FLOW/macros/AddTaskVZERO.C");
120          AddTaskVZERO(0,0,0,0);
121    }
122    // create the main task
123    AliAnalysisTwoParticleResonanceFlowTask *task = new AliAnalysisTwoParticleResonanceFlowTask("TaskPhiFlow");
124    if(debug) cout << " === AliAnalysisTaskPhiFlow === " << task << endl;
125    if(!task) {
126        if(debug) cout << " --> Unexpected error occurred: NO TASK WAS CREATED! (could be a library problem!) " << endl;
127        return 0x0;
128    }
129    if(task->UsePhiMinusPsiMethod(PhiMinusPsiMethod)) {
130        cout << " --> Using phi - Psi method <-- ... " << endl;
131        Float_t dphibins[] = {0., 0.63, 1.26, 1.89, 2.52, 3.15};
132        task->SetdPhiBins(dphibins, (Int_t)(sizeof(dphibins)/sizeof(dphibins[1])-1));
133        event_mixing = kTRUE;
134    }
135    task->SetupSpeciesA(speciesA, chargeA, massA);
136    task->SetupSpeciesB(speciesB, chargeB, massB);
137
138    if(task->SetVZEROSubEvents(EP3sub)) cout << " --> Setting up VZERO subevents method ... " << endl;
139    if(event_mixing) {
140       if(debug) cout << " --> Enabeling event mixing for phi reconstruction - try at your own risk !!! <-- " << endl;
141       // set vertex and mixing bins - arrays MUST have length 20!
142       Int_t c[] = {0, 2, 4, 6, 8, 10, 20, 30, 40, 50, 60, 70, 80, 90, 101, 0, 0, 0, 0, 0};
143       Int_t v[] = {-10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10, 0, 0, 0, 0, 0, 0, 0, 0, 0};
144       if(((Int_t)(sizeof(c)/sizeof(c[1]))!=20)||((Int_t)(sizeof(v)/sizeof(v[1]))!=20)) {
145           cout << " --> Fatal error: check mixing parameters ! <-- " << endl;
146           return 0x0;
147       }
148       else task->SetMixingBins(c, v);
149    }
150    // set triggers
151    if (bCentralTrigger) task->SelectCollisionCandidates(AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral);
152    else                 task->SelectCollisionCandidates(AliVEvent::kMB);
153    if(debug) cout << "    --> Set trigger selection to ";
154    if(debug&&bCentralTrigger) cout << " kMB, kCentral, kSemiCentral " << endl;
155    if(debug&&(!bCentralTrigger)) cout << " kMB " << endl;
156    //set RP cuts for flow package analysis
157    AliFlowTrackCuts* cutsRP = new AliFlowTrackCuts("RFPcuts");
158    if(!cutsRP) {
159        if(debug) cout << " Fatal error: no RP cuts found, could be a library problem! " << endl;
160        return 0x0;
161    }
162    if(useGlobalRPCuts&&(!VZERO_SP)) {
163        AliFlowTrackCuts::trackParameterType rptype = AliFlowTrackCuts::kGlobal;
164        cutsRP->SetParamType(rptype);
165        cutsRP->SetPtRange(0.2, 5.0);
166        cutsRP->SetEtaRange(-0.8, 0.8);
167        cutsRP->SetMinNClustersTPC(70);
168        cutsRP->SetMinChi2PerClusterTPC(0.1);
169        cutsRP->SetMaxChi2PerClusterTPC(4.0);
170        cutsRP->SetRequireTPCRefit(kTRUE);
171        cutsRP->SetMaxDCAToVertexXY(0.3);
172        cutsRP->SetMaxDCAToVertexZ(0.3);
173        cutsRP->SetAcceptKinkDaughters(kFALSE);
174        cutsRP->SetMinimalTPCdedx(10.);
175        if(rp_filter < 9999 ) {
176            if(debug) cout << "  > set RP filterbit " << rp_filter << endl;     
177            cutsRP->SetAODfilterBit(rp_filter);
178        }
179        if(debug) cout << "    --> kGlobal RP's " << cutsRP << endl;
180    }
181    if(VZERO_SP) { // use vzero sub analysis
182        cutsRP = cutsRP->GetStandardVZEROOnlyTrackCuts(); // select vzero tracks
183        SP = kFALSE; // disable other methods
184        SPSUB = kTRUE; // calculate sp_qa and sp_qb
185        QC = kFALSE;
186        EP = kFALSE;
187        EP3sub = kFALSE;
188        EtaGap = 0.; // no eta gap, full tpc poi's
189        if(debug) cout << "    --> VZERO RP's " << cutsRP << endl;
190    }
191    task->SetRPCuts(cutsRP);
192    // set POI cuts for kaon selection
193    AliFlowTrackCuts* cutsPOI = new AliFlowTrackCuts("GlobalPOI");
194    if(!cutsPOI) {
195        if(debug) cout << " Fatal error: no POI cuts (could be a library problem)!" << endl;
196        return 0x0;
197    }
198    cutsPOI = cutsPOI->GetStandardGlobalTrackCuts2010();
199    cutsPOI->SetPtRange(POIPtMin, POIPtMax); // pt range of DAUGHTERS !!!
200    cutsPOI->SetMaxDCAToVertexXY(0.3); // FIXME not implemented in aod086 aod095 see PassesDCACut() in implementation
201    cutsPOI->SetMaxDCAToVertexZ(0.3);
202    if(poi_filter < 9999 ) {
203        if(debug) cout << "  > set POI filterbit " << poi_filter << endl;     
204        cutsPOI->SetAODfilterBit(poi_filter);
205    }
206    if(debug) cout << "    --> cutsPOI " << cutsPOI << endl;
207    task->SetPOICuts(cutsPOI);
208    //set POI cuts for aods XY Z - 3 distinct cases
209    Float_t POIDCA[] = {0., 0., 0., 0., 0.};
210    if(DCA == "none" ) { // 1 --- do nothing
211        if (debug) cout << " --> No DCA cut on POI's <-- " << endl;
212        for (Int_t i = 0; i < 5; i++) POIDCA[i] = 0.;
213    }
214    if(DCA == "fix" ) { // 2 --- use fixed values for xy z
215        if (debug) cout << " --> Fixed DCA cut on POI's <-- " << endl;
216        POIDCA[0] = -1.; POIDCA[1] = 0.3; POIDCA[2] = 0.3; POIDCA[3] = 0.; POIDCA[4] = 0.;
217    }
218    if(DCA == "pt" ) { // 3 --- use pt dependent cut
219        if (debug) cout << " --> Pt dependent DCA cut on POI's <-- " << endl;
220        POIDCA[0] = 1.; POIDCA[1] = 0.0105; POIDCA[2] = 0.0350; POIDCA[3] = 1.1; POIDCA[4] = 2.;
221    }
222    task->SetPOIDCAXYZ(POIDCA);
223    Float_t AddTaskMacroSummary[] = {(Float_t)0, (Float_t)SP, (Float_t)SPSUB, (Float_t)QC, (Float_t)EP, (Float_t)EP3sub, (Float_t)VZERO_SP, (Float_t)EtaGap, (Float_t)harm, (Float_t)highPtMode, (Float_t)deltaDip, (Float_t)POIPtMax}; 
224    task->SetAddTaskMacroSummary(AddTaskMacroSummary);
225
226    if(highPtMode) { // high pt loop - loop will end macro
227        Float_t _pt[] = {0., 1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 12., 15.};
228        task->SetPtBins(_pt, (Int_t)(sizeof(_pt)/sizeof(_pt[1]))-1);
229        // general approach: use kinematic filters which select kaon pairs with a certain mass window
230        AliFlowTrackSimpleCuts* HighPtSubEventFilterA = new AliFlowTrackSimpleCuts("HighPtSubEventFilterA"); 
231        HighPtSubEventFilterA->SetEtaMin(-0.8);
232        HighPtSubEventFilterA->SetEtaMax(0.0);
233        HighPtSubEventFilterA->SetMassMin(1.019445 - deltaPhiMass);
234        HighPtSubEventFilterA->SetMassMax(1.019445 + deltaPhiMass);
235        AliFlowTrackSimpleCuts* HighPtSubEventFilterB = new AliFlowTrackSimpleCuts("HighPtSubEventFilterB"); 
236        HighPtSubEventFilterB->SetEtaMin(0.0);
237        HighPtSubEventFilterB->SetEtaMax(+0.8);
238        HighPtSubEventFilterB->SetMassMin(1.019445 - deltaPhiMass);
239        HighPtSubEventFilterB->SetMassMax(1.019445 + deltaPhiMass);
240        AliFlowTrackSimpleCuts* HighPtGenericFilter = new AliFlowTrackSimpleCuts("HighPtGenericFilter");
241        HighPtGenericFilter->SetEtaMin(-0.8);
242        HighPtGenericFilter->SetEtaMax(+0.8);
243        HighPtGenericFilter->SetMassMin(1.019445 - deltaPhiMass);
244        HighPtGenericFilter->SetMassMax(1.019445 + deltaPhiMass);
245        if(debug) cout << "    --> Created poi filters " << endl;
246        // set pair and event cuts
247        if((deltaDip>0.005)&&(deltaDipMaxPt>0.005)) task->SetMaxDeltaDipAngleAndPt(deltaDip, deltaDipMaxPt);
248        else cout << " --> Disabled Delta-Dip exclusion. <-- " << endl;
249        task->SetCandidateEtaAndPt(POIEtaMin, POIEtaMax, 0., 15.);
250        task->SetCentralityParameters(centrMin, centrMax, "TRK");
251        task->SetVertexZ(vertexZ);
252        if(debug) cout << "    --> Set pair cuts and event cuts" << endl;
253        // specify the PID procedure which will be used
254        task->SetPIDConfiguration(PIDconfig);
255        AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
256        AliAnalysisDataContainer *coutHist = mgr->CreateContainer(Form("PhiV2_%s", OutputName(centrMin, centrMax,PIDconfig, suffixName, bCentralTrigger, EtaGap, POIEtaMin, POIEtaMax, 0., 15., deltaDip, deltaDipMaxPt, DCA, harm, TPCStandAloneTracks, vertexZ, debug, useGlobalRPCuts).Data()), TList::Class(), AliAnalysisManager::kOutputContainer, fileName);
257        if(debug) cout << "    --> Created IO containers " << cinput << ", " << coutHist << endl;
258        mgr->AddTask(task);
259        if(debug) cout << " === ADDING MAIN TASK == " << endl;
260        mgr->ConnectInput(task, 0, cinput);
261        mgr->ConnectOutput(task, 1, coutHist);
262        if(debug) cout << "    --> Connected IO containers " << endl;
263        if (SP || EP || QC || SPSUB) // if flow analysis should be done after reconstruction
264        {
265           Int_t mb(999);
266           if(debug) cout << " === RECEIVED REQUEST FOR FLOW ANALYSIS === " << endl;
267           AliAnalysisDataContainer *flowEvent = mgr->CreateContainer(Form("FC%s", suffixName.Data()), AliFlowEventSimple::Class(), AliAnalysisManager::kExchangeContainer);
268           mgr->ConnectOutput(task, 2, flowEvent);
269           if(debug) cout << "    --> Created IO containers " << flowEvent << endl;
270           if(debug) cout << "    --> suffixName " << suffixName << endl;
271           if (QC) {  // add qc tasks
272              AddQCmethod(Form("QCTPCMB_%d_%s", mb, suffixName.Data()), harm, flowEvent, HighPtGenericFilter, debug, 0x0, suffixName.Data());
273              if(debug) cout << "    --> Hanging QC task ... " << mb << " succes! "<< endl;
274           }
275           if (SPSUB) {  // add sp subevent tasks
276              AddSPmethod(Form("SPTPCMB_%d_%s", mb, suffixName.Data()), -0.8, -.5*EtaGap, .5*EtaGap, +0.8, "Qa", harm, flowEvent, HighPtSubEventFilterB, NULL, false, shrinkSP, debug, true, suffixName.Data());
277              if(debug) cout << "    --> Hanging SP Qa task " << mb << " succes!" << endl;
278              AddSPmethod(Form("SPTPCMB_%d_%s", mb, suffixName.Data()), -0.8, -.5*EtaGap, .5*EtaGap, +0.8, "Qb", harm, flowEvent, HighPtSubEventFilterA, NULL, false, shrinkSP, debug, true, suffixName.Data());
279              if(debug) cout << "    --> Hanging SP Qb task ..." << mb << " succes!"<< endl;
280           }
281           if (SP) { // add sp tasks
282              AddSPmethod(Form("SPTPCMB_%d_%s", mb, suffixName.Data()), -0.8, -0.0, +0.0, +0.8, "QaQb", harm, flowEvent, HighPtGenericFilter, NULL, false, shrinkSP, debug, 0x0, suffixName.Data());
283              if(debug) cout << "    --> Hanging SP task ... " << mb << " succes!" << endl;
284           }
285           if (EP) { // add ep tasks
286              AddSPmethod(Form("EPTPCMB_%d_%s", mb, suffixName.Data()), -0.8, -0.0, +0.0, +0.8, "QaQb", harm, flowEvent, HighPtGenericFilter, NULL, true, shrinkSP, debug, 0x0, suffixName.Data());
287              if(debug) cout << "    --> Hanging EP task ... " << mb << " succes!" << endl;
288           }
289        }
290        // print summary to screen
291        cout << endl << endl << "       ==== AddTaskPhiFlow launched  ===== " << endl;
292        cout << " ************ Configured PID routine ************ " << endl;
293        cout << "      0 < " << PIDconfig[1] << " p_t, ITS || TPC with s < " << PIDconfig[0] << endl;
294        if(PIDconfig[2] < 0.) cout << "    --> TPC control disabled " << endl;
295        if(PIDconfig[2] > 0.) cout << "    --> TPC control enabled " << endl;
296        cout << "    " << PIDconfig[1] << " < " << PIDconfig[4] << " p_t, TPC || ITS with s < " << PIDconfig[3] << endl;
297        if(PIDconfig[5] < 0.) cout << "    --> ITS control disabled " << endl;
298        if(PIDconfig[5] > 0.) cout << "    --> ITS control enabled " << endl;
299        cout << "      " << PIDconfig[4] << " < 7 p_t, TPC / TOF Bayesian with p < " << PIDconfig[6] << endl;
300        cout << " ************ Configured DCA setup ************** " << endl;
301        cout << "  DCA type: " << DCA;
302        if (DCA == "") cout << "default";
303        cout << endl << " ************* Task statisti:q!cs ****************** " << endl;
304        cout << "   -> Launched PHI reconstruction " << endl;
305        if(SP) cout << "   --> Launched QaQb SP filters and corresponding SP task " << endl;
306        if(EP) cout << "   --> Launched QaQb QC filters and corresponding EP task " << endl;
307        if(SPSUB) cout << "   --> Launched Qa&&Qb SP filters and corresponding SP tasks " << endl;
308        if(QC) cout << "   --> Launched QaQb QC filters and corresponding QC task " << endl;
309        if(EP3sub) cout << " --> Launched VZERO subevent analysis alongside reconstruction - USE WITH CAUTION!" << endl;
310        cout << " ************************************************ " << endl;
311        TString condit = "";
312        (task->SetQA(kFALSE)) ? condit+= " --> Enabled QA plots <-- " : condit+= " --> Disabled QA plots <-- ";
313        (task->SetIsMC(kFALSE)) ? condit+= " --> MC mode <-- " : condit+= " --> DATA mode <-- ";
314        (task->UseEventMixing(event_mixing)) ? condit+= " --> Using EVENT MIXING <--" : condit+= "--> Combinatorial background <--";
315        cout << condit << endl;
316        cout << "           --> Now go for a coffee! <-- " << endl;
317        cout << " ************************************************ " << endl;
318        return task;
319    } //  end of high pt loop - high-pt task is ready to go at this point
320    Float_t _pt[] = {0.0, 0.3, 0.6, 0.9, 1.2, 1.5, 1.8, 2.1, 2.4, 2.7, 3.0, 3.5, 4.0, 4.5, 5.0, 5.5, 6.0, 6.5, 7.0};
321    task->SetPtBins(_pt, (Int_t)(sizeof(_pt)/sizeof(_pt[1]))-1);
322    // POI filter cuts, will filter invm mass bands and subevents
323    AliFlowTrackSimpleCuts* POIfilterQC[30];
324    AliFlowTrackSimpleCuts* POIfilterSP[30][2];
325    Double_t flowBands[2][30];
326    Double_t _inc = (maxMaxx-minMass)/30.;
327    for (Int_t mb = 0; mb < 30; mb++) {
328       flowBands[0][mb] = minMass + mb * _inc;
329       flowBands[1][mb] = minMass + (mb + 1) * _inc;
330       POIfilterSP[mb][0] = new AliFlowTrackSimpleCuts(Form("FilterPOISP_MB%d_ETANEG", mb));
331       if(!POIfilterSP[mb][0]) {
332           if(debug) cout << " FAILED to initialize POI filters, could be a library problem!" << endl;
333           return 0x0;
334       }
335       POIfilterSP[mb][0]->SetEtaMin(-0.8);
336       POIfilterSP[mb][0]->SetEtaMax(0.0);
337       POIfilterSP[mb][0]->SetMassMin(flowBands[0][mb]);
338       POIfilterSP[mb][0]->SetMassMax(flowBands[1][mb]);
339       POIfilterSP[mb][1] = new AliFlowTrackSimpleCuts(Form("FilterPOISP_MB%d_ETAPOS", mb));
340       POIfilterSP[mb][1]->SetEtaMin(0.0);
341       POIfilterSP[mb][1]->SetEtaMax(+0.8);
342       POIfilterSP[mb][1]->SetMassMin(flowBands[0][mb]);
343       POIfilterSP[mb][1]->SetMassMax(flowBands[1][mb]);
344       POIfilterQC[mb] = new AliFlowTrackSimpleCuts(Form("FilterPOIQC_MB%d", mb));
345       POIfilterQC[mb]->SetEtaMin(-0.8);
346       POIfilterQC[mb]->SetEtaMax(+0.8);
347       POIfilterQC[mb]->SetMassMin(flowBands[0][mb]);
348       POIfilterQC[mb]->SetMassMax(flowBands[1][mb]);
349    }
350    if(debug) cout << "    --> Created poi filters " << endl;
351    task->SetCommonConstants(30, flowBands[0][0], flowBands[1][29]);
352    if(debug) cout << "    --> Set common constants " << endl;
353    // set pair and event cuts
354    if((deltaDip>0.005)&&(deltaDipMaxPt>0.005)) task->SetMaxDeltaDipAngleAndPt(deltaDip, deltaDipMaxPt);
355    else cout << " --> Disabled Delta-Dip exclusion. <-- " << endl;
356    task->SetCandidateEtaAndPt(POIEtaMin, POIEtaMax, 0., 10.);
357    task->SetCentralityParameters(centrMin, centrMax, "TRK");
358    task->SetVertexZ(vertexZ);
359    if(debug) cout << "    --> Set pair cuts and event cuts" << endl;
360    // set the kaon cuts, and specify the PID procedure which will be used
361    task->SetPIDConfiguration(PIDconfig);
362    if(debug) {
363        cout << " ************ Configured PID routine ************ " << endl;
364        cout << "    0 < " << PIDconfig[1] << " p_t, ITS || TPC with s < " << PIDconfig[0] << endl;
365        if(PIDconfig[2] < 0.) cout << "  --> TPC control disabled " << endl;
366        if(PIDconfig[2] > 0.) cout << "  --> TPC control enabled " << endl;
367        cout << "    " << PIDconfig[1] << " < " << PIDconfig[4] << " p_t, TPC || ITS with s < " << PIDconfig[3] << endl;
368        if(PIDconfig[5] < 0.) cout << "  --> ITS control disabled " << endl;
369        if(PIDconfig[5] > 0.) cout << "  --> ITS control enabled " << endl;
370        cout << "    " << PIDconfig[4] << " < 7 p_t, TPC / TOF Bayesian with p < " << PIDconfig[6] << endl;
371        cout << " ************************************************ " << endl;
372    }
373    if (TPCStandAloneTracks)
374    { // switch to tpc standalone tracks for POI selection
375       if(debug) cout << "    --> Switching to TPC standalone analsis " << endl;
376       task->SetRequireStrictKaonCuts();
377       task->SetRequireTPCStandAloneKaons();
378    }
379    // create and configure IO containers
380    AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
381    AliAnalysisDataContainer *coutHist = mgr->CreateContainer(Form("PhiV2_%s", OutputName(centrMin, centrMax,PIDconfig, suffixName, bCentralTrigger, EtaGap, POIEtaMin, POIEtaMax, 0., 10., deltaDip, deltaDipMaxPt, DCA, harm, TPCStandAloneTracks, vertexZ, debug, useGlobalRPCuts).Data()), TList::Class(), AliAnalysisManager::kOutputContainer, fileName);
382    if(debug) cout << "    --> Created IO containers " << cinput << ", " << coutHist << endl;
383    mgr->AddTask(task);
384    if(debug) cout << " === ADDING MAIN TASK == " << endl;
385    mgr->ConnectInput(task, 0, cinput);
386    mgr->ConnectOutput(task, 1, coutHist);
387    if(debug) cout << "    --> Connected IO containers " << endl;
388    if (SP || EP || QC || SPSUB) // if flow analysis should be done after reconstruction
389    {
390       if(debug) cout << " === RECEIVED REQUEST FOR FLOW ANALYSIS === " << endl;
391       AliAnalysisDataContainer *flowEvent = mgr->CreateContainer(Form("FC%s", suffixName.Data()), AliFlowEventSimple::Class(), AliAnalysisManager::kExchangeContainer);
392       mgr->ConnectOutput(task, 2, flowEvent);
393       if(debug) cout << "    --> Created IO containers " << flowEvent << endl;
394       if(debug) cout << "    --> suffixName " << suffixName << endl;
395       for (int mb = 0; mb != 30; ++mb) {
396          if (QC) {  // add qc tasks
397             AddQCmethod(Form("QCTPCMB_%d_%s", mb, suffixName.Data()), harm, flowEvent, POIfilterQC[mb], debug, 0x0, suffixName.Data());
398             if(debug) cout << "    --> Hanging QC task ... " << mb << " succes! "<< endl;
399          }
400          if (SPSUB) {  // add sp subevent tasks
401             AddSPmethod(Form("SPTPCMB_%d_%s", mb, suffixName.Data()), -0.8, -.5*EtaGap, .5*EtaGap, +0.8, "Qa", harm, flowEvent, POIfilterSP[mb][1], NULL, false, shrinkSP, debug, true, suffixName.Data(), VZERO_SP);
402             if(debug) cout << "    --> Hanging SP Qa task " << mb << " succes!" << endl;
403             AddSPmethod(Form("SPTPCMB_%d_%s", mb, suffixName.Data()), -0.8, -.5*EtaGap, .5*EtaGap, +0.8, "Qb", harm, flowEvent, POIfilterSP[mb][0], NULL, false, shrinkSP, debug, true, suffixName.Data(), VZERO_SP);
404             if(debug) cout << "    --> Hanging SP Qb task ..." << mb << " succes!"<< endl;
405          }
406          if (SP) { // add sp tasks
407             AddSPmethod(Form("SPTPCMB_%d_%s", mb, suffixName.Data()), -0.8, -0.0, +0.0, +0.8, "QaQb", harm, flowEvent, POIfilterQC[mb], NULL, false, shrinkSP, debug, 0x0, suffixName.Data());
408             if(debug) cout << "    --> Hanging SP task ... " << mb << " succes!" << endl;
409          }
410          if (EP) { // add ep tasks
411             AddSPmethod(Form("EPTPCMB_%d_%s", mb, suffixName.Data()), -0.8, -0.0, +0.0, +0.8, "QaQb", harm, flowEvent, POIfilterQC[mb], NULL, true, shrinkSP, debug, 0x0, suffixName.Data());
412             if(debug) cout << "    --> Hanging EP task ... " << mb << " succes!" << endl;
413          }
414       }
415    }
416    // print summary to screen
417    cout << endl << endl << "       ==== AddTaskPhiFlow launched  ===== " << endl;
418    cout << " ************ Configured PID routine ************ " << endl;
419    cout << "      0 < " << PIDconfig[1] << " p_t, ITS || TPC with s < " << PIDconfig[0] << endl;
420    if(PIDconfig[2] < 0.) cout << "    --> TPC control disabled " << endl;
421    if(PIDconfig[2] > 0.) cout << "    --> TPC control enabled " << endl;
422    cout << "    " << PIDconfig[1] << " < " << PIDconfig[4] << " p_t, TPC || ITS with s < " << PIDconfig[3] << endl;
423    if(PIDconfig[5] < 0.) cout << "    --> ITS control disabled " << endl;
424    if(PIDconfig[5] > 0.) cout << "    --> ITS control enabled " << endl;
425    cout << "      " << PIDconfig[4] << " < 7 p_t, TPC / TOF Bayesian with p < " << PIDconfig[6] << endl;
426    cout << " ************ Configured DCA setup ************** " << endl;
427    cout << "  DCA type: " << DCA;
428    if (DCA == "") cout << "default";
429    cout << endl << " ************* Task statisti:q!cs ****************** " << endl;
430    cout << "   -> Launched PHI reconstruction " << endl;
431    if(SP) cout << "   --> Launched 30 QaQb SP filters and corresponding 30 SP tasks " << endl;
432    if(EP) cout << "   --> Launched 30 QaQb QC filters and corresponding 30 EP tasks " << endl;
433    if(SPSUB) cout << "   --> Launched 30+30 Qa&&Qb SP filters and corresponding 30+30 SP tasks " << endl;
434    if(QC) cout << "   --> Launched 30 QaQb QC filters and corresponding 30 QC tasks " << endl;
435    if(EP3sub) cout << " --> Launched VZERO subevent analysis alongside reconstruction - USE WITH CAUTION!" << endl;
436    cout << " ************************************************ " << endl;
437    TString condit = "";
438    (task->SetQA(kFALSE)) ? condit+= " --> Enabled QA plots <-- " : condit+= " --> Disabled QA plots <-- ";
439    (task->SetIsMC(kFALSE)) ? condit+= " --> MC mode <-- " : condit+= " --> DATA mode <-- ";
440    (task->UseEventMixing(event_mixing)) ? condit+= " --> Using EVENT MIXING <--" : condit+= "--> Combinatorial background <--";
441    cout << condit << endl;
442    cout << "           --> Now go for a coffee! <-- " << endl;
443    cout << " ************************************************ " << endl;
444    return task;
445 }
446 //_____________________________________________________________________________
447 void AddSPmethod(char *name, double minEtaA, double maxEtaA, double minEtaB, double maxEtaB, char *Qvector, int harmonic, AliAnalysisDataContainer *flowEvent, AliFlowTrackSimpleCuts *cutsPOI = NULL, AliFlowTrackSimpleCuts *cutsRFP = NULL, bool bEP, bool shrink = false, bool debug, bool etagap = false, char *suffixName, Bool_t VZERO_SP = kFALSE)
448 {
449    // add sp task and invm filter tasks
450    if(debug) (bEP) ? cout << " ****** Reveived request for EP task ****** " << endl : cout << " ******* Switching to SP task ******* " << endl;
451    TString fileName = AliAnalysisManager::GetCommonFileName();
452    (bEP) ? fileName+=":EP" : fileName+=":SP";
453    fileName+=suffixName;
454    if(etagap) {
455        fileName+="_SUBEVENTS";
456        if(debug) cout << "    --> Setting up subevent analysis <-- " << endl;
457    }
458    if(debug) cout << "    --> fileName " << fileName << endl;
459    TString myNameSP;
460    (bEP) ? myNameSP = Form("%sEPv%d%s", name, harmonic, Qvector): myNameSP = Form("%sSPv%d%s", name, harmonic, Qvector);
461    if(debug) cout << " Task and filter name: " << myNameSP << endl;
462    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
463    AliAnalysisDataContainer *flowEvent2 = mgr->CreateContainer(Form("Filter_%s", myNameSP.Data()), AliFlowEventSimple::Class(), AliAnalysisManager::kExchangeContainer);
464    AliAnalysisTaskFilterFE *tskFilter = new AliAnalysisTaskFilterFE(Form("TaskFilter_%s", myNameSP.Data()), cutsRFP, cutsPOI);
465    tskFilter->SetSubeventEtaRange(minEtaA, maxEtaA, minEtaB, maxEtaB);
466    if(VZERO_SP) tskFilter->SetSubeventEtaRange(-10, 0, 0, 10);
467    mgr->AddTask(tskFilter);
468    mgr->ConnectInput(tskFilter, 0, flowEvent);
469    mgr->ConnectOutput(tskFilter, 1, flowEvent2);
470    AliAnalysisDataContainer *outSP = mgr->CreateContainer(myNameSP.Data(), TList::Class(), AliAnalysisManager::kOutputContainer, fileName);
471    AliAnalysisTaskScalarProduct *tskSP = new AliAnalysisTaskScalarProduct(Form("TaskScalarProduct_%s", myNameSP.Data()), kFALSE);
472    tskSP->SetApplyCorrectionForNUA(kTRUE);
473    tskSP->SetHarmonic(harmonic);
474    tskSP->SetTotalQvector(Qvector);
475    if (bEP) tskSP->SetBehaveAsEP();
476    if (shrink) tskSP->SetBookOnlyBasicCCH(kTRUE);
477    mgr->AddTask(tskSP);
478    mgr->ConnectInput(tskSP, 0, flowEvent2);
479    mgr->ConnectOutput(tskSP, 1, outSP);
480 }
481 //_____________________________________________________________________________
482 void AddQCmethod(char *name, int harmonic, AliAnalysisDataContainer *flowEvent, AliFlowTrackSimpleCuts *cutsPOI = NULL, Bool_t debug, AliFlowTrackSimpleCuts *cutsRFP = NULL, char *suffixName)
483 {
484    // add qc task and invm filter tasks
485    if(debug) cout << " ****** Received request for QC v" << harmonic << " task " << name << ", POI " << cutsPOI << ", IO ****** " << flowEvent << endl;
486    TString fileName = AliAnalysisManager::GetCommonFileName();
487    fileName+=":QC";
488    fileName+=suffixName;
489    if(debug) cout << "    --> Common filename: " << fileName << endl;
490    TString myName = Form("%s", name);
491    if(debug) cout << "    --> myName: " << myName << endl;
492    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
493    AliAnalysisDataContainer *flowEvent2 = mgr->CreateContainer(Form("Filter_%s", myName.Data()), AliFlowEventSimple::Class(),AliAnalysisManager::kExchangeContainer);
494    AliAnalysisTaskFilterFE *tskFilter = new AliAnalysisTaskFilterFE(Form("TaskFilter_%s", myName.Data()), cutsRFP, cutsPOI);
495    mgr->AddTask(tskFilter);
496    mgr->ConnectInput(tskFilter, 0, flowEvent);
497    mgr->ConnectOutput(tskFilter, 1, flowEvent2);
498    AliAnalysisDataContainer *outQC = mgr->CreateContainer(myName.Data(), TList::Class(), AliAnalysisManager::kOutputContainer, fileName);
499    AliAnalysisTaskQCumulants *tskQC = new AliAnalysisTaskQCumulants(Form("TaskQCumulants_%s", myName.Data()), kFALSE);
500    tskQC->SetApplyCorrectionForNUA(kTRUE);
501    tskQC->SetHarmonic(harmonic);
502    tskQC->SetBookOnlyBasicCCH(kTRUE);
503    mgr->AddTask(tskQC);
504    mgr->ConnectInput(tskQC, 0, flowEvent2);
505    mgr->ConnectOutput(tskQC, 1, outQC);
506 }
507 //_____________________________________________________________________________
508 TString OutputName( Float_t centrMin,
509                     Float_t centrMax,
510                     Float_t PIDconfig[7],
511                     TString suffixName,
512                     Bool_t bCentralTrigger,
513                     Float_t EtaGap,
514                     Float_t POIEtaMin,
515                     Float_t POIEtaMax,
516                     Float_t POIPtMin,
517                     Float_t POIPtMax,
518                     Float_t deltaDip,
519                     Float_t deltaDipMaxPt,
520                     TString DCA,
521                     Int_t harm,
522                     Bool_t TPCStandAloneTracks,
523                     Float_t vertexZ,
524                     Bool_t debug,
525                     Bool_t useGlobalRPCuts)
526 {
527    // generate output name
528    TString centralityName = "";
529    centralityName += suffixName;
530    centralityName += "_DCA";
531    centralityName += DCA;
532    centralityName += Form("_vZ%.f", vertexZ);
533    centralityName += "_";
534    for(Int_t i = 0; i < 7; i++) centralityName += Form("%.1f_", PIDconfig[i]);
535    centralityName += "_POIEta";
536    centralityName += Form("%.1f", POIEtaMin);
537    centralityName += Form("%.1f", POIEtaMax);
538    centralityName += "_gap";
539    centralityName += Form("%.1f", -0.5*EtaGap);
540    centralityName += Form("%.1f", 0.5*EtaGap);
541    centralityName += "_";
542    centralityName += Form("dDip%.2f", deltaDip);
543    centralityName += "-";
544    centralityName += Form("dDipPt%.2f", deltaDipMaxPt);
545    if (TPCStandAloneTracks) {
546       centralityName += "-";
547       centralityName += "TPCStandAloneTracks";
548    }
549    if (bCentralTrigger) {
550       centralityName += "-";
551       centralityName += "kMBkCkSC";
552    }
553    if (!useGlobalRPCuts) {
554        centralityName += "-";
555        centralityName += "TPCRP";
556    }
557    if(debug) cout << "    --> centralityName " << centralityName << endl;
558    return centralityName;
559 }
560 //_____________________________________________________________________________
561