updated for nonHFE reconstruction study
[u/mrichter/AliRoot.git] / PWGHF / hfe / macros / AddTaskFlowTPCTOFQCSP.C
1 ///////////////////////////////////////////////////////////////////\r
2 //                                                               //            \r
3 // AddTaskFlowTPCTOFQCSP macro                                     //\r
4 // Author: Andrea Dubla, Utrecht University, 2012                //\r
5 //                                                               //\r
6 ///////////////////////////////////////////////////////////////////\r
7 class AliAnalysisDataContainer;\r
8 class AliFlowTrackCuts;\r
9 class AliFlowTrackSimpleCuts;\r
10 class AliFlowEventCuts;\r
11 class AliFlowEventSimpleCuts;\r
12 class AliAnalysisDataContainer;\r
13 class AliHFEextraCuts;\r
14 \r
15 AliAnalysisTaskFlowTPCTOFQCSP* AddTaskFlowTPCTOFQCSP(\r
16                                               TString uniqueID = "",\r
17                                               Float_t centrMin ,\r
18                                               Float_t centrMax ,\r
19                                               Double_t InvmassCut,\r
20                                               Double_t pTCut,\r
21                                               Bool_t op_ang = kFALSE,\r
22                                               Double_t op_angle_cut,\r
23                                               Int_t Trigger,\r
24                                               Double_t minTPCnsigma,\r
25                                               Double_t maxTPCnsigma,\r
26                                               Double_t minTOFnSigma,\r
27                                               Double_t maxTOFnSigma,\r
28                                               Int_t minTPCCluster,\r
29                                               Int_t TPCS,\r
30                                               Int_t Vz,\r
31                                               AliHFEextraCuts::ITSPixel_t pixel,\r
32                                               Bool_t PhiCut = kFALSE,\r
33                                               Bool_t QaPidSparse = kFALSE,\r
34                                               const char *Cent = "V0M",\r
35                                               Bool_t QC = kTRUE, // use qc2 and qc4\r
36                                               Bool_t SP_TPC = kTRUE, //use tpc sp method\r
37                                               Bool_t VZERO_SP = kFALSE, // use vzero sp method\r
38                                               Int_t harmonic = 2,\r
39                                               Bool_t shrinkSP = kTRUE,\r
40                                               Bool_t debug = kFALSE,\r
41                                               Int_t RPFilterBit = 1\r
42                                               )\r
43 \r
44 {\r
45 \r
46 \r
47 \r
48     if(debug) cout << " === Adding Task ElectFlow === " << endl;    \r
49     TString fileName = AliAnalysisManager::GetCommonFileName();\r
50     fileName += ":ElectroID_";\r
51     fileName += uniqueID;\r
52     if(debug) cout << "    --> Reconstruction data container: " << fileName << endl;\r
53     AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();\r
54     if (!mgr) {\r
55         if(debug) cout << " Fatal error: no analysis manager found! " << endl;\r
56         return 0x0;\r
57     }\r
58     if (!mgr->GetInputEventHandler()) {\r
59         if(debug) cout << " Fatal error: no imput event handler found!" << endl;\r
60         return 0x0;\r
61     }\r
62 \r
63 //create a task\r
64   AliAnalysisTaskFlowTPCTOFQCSP *taskHFE = ConfigHFEStandardCuts(minTPCCluster, pixel);    //kTRUE if MC\r
65     \r
66   if(debug) cout << " === AliAnalysisTaskFlowTPCTOFQCSP === " << taskHFE << endl;\r
67   if(!taskHFE) {\r
68         if(debug) cout << " --> Unexpected error occurred: NO TASK WAS CREATED! (could be a library problem!) " << endl;\r
69         return 0x0;\r
70     }\r
71     \r
72 // Set centrality percentiles and method V0M, FMD, TRK, TKL, CL0, CL1, V0MvsFMD, TKLvsV0M, ZEMvsZDC\r
73   taskHFE->SetCentralityParameters(centrMin, centrMax, Cent);\r
74   taskHFE->SetInvariantMassCut(InvmassCut);\r
75   taskHFE->SetTrigger(Trigger);\r
76   taskHFE->SetTPCS(TPCS);\r
77   taskHFE->SetVz(Vz);\r
78   taskHFE->SetIDCuts(minTPCnsigma, maxTPCnsigma, minTOFnSigma, maxTOFnSigma);\r
79   taskHFE->SetQAPIDSparse(QaPidSparse);\r
80   taskHFE->SetPhiCut(PhiCut);\r
81   taskHFE->SetOpeningAngleflag(op_ang);\r
82   taskHFE->SetOpeningAngleCut(op_angle_cut);\r
83   taskHFE->SetpTCuttrack(pTCut);\r
84 \r
85     \r
86 //set RP cuts for flow package analysis\r
87   cutsRP = new AliFlowTrackCuts(Form("RFPcuts%s",uniqueID));\r
88     if(!cutsRP) {\r
89         if(debug) cout << " Fatal error: no RP cuts found, could be a library problem! " << endl;\r
90         return 0x0;\r
91     }\r
92 \r
93     if(!VZERO_SP) {\r
94         AliFlowTrackCuts::trackParameterType rptype = AliFlowTrackCuts::kGlobal;\r
95         cutsRP->SetParamType(rptype);\r
96         cutsRP->SetAODfilterBit(RPFilterBit);\r
97         cutsRP->SetPtRange(0.2, 5.0);\r
98         cutsRP->SetEtaRange(-0.8, 0.8);\r
99         cutsRP->SetMinNClustersTPC(70);\r
100         cutsRP->SetMinChi2PerClusterTPC(0.1);\r
101         cutsRP->SetMaxChi2PerClusterTPC(4.0);\r
102         cutsRP->SetRequireTPCRefit(kTRUE);\r
103         cutsRP->SetMaxDCAToVertexXY(0.3);\r
104         cutsRP->SetMaxDCAToVertexZ(0.3);\r
105         cutsRP->SetAcceptKinkDaughters(kFALSE);\r
106         cutsRP->SetMinimalTPCdedx(10.);\r
107         if(debug) cout << "    --> kGlobal RP's " << cutsRP << endl;\r
108     }\r
109     if(VZERO_SP) { // use vzero sub analysis\r
110         cutsRP = cutsRP->GetStandardVZEROOnlyTrackCuts(); // select vzero tracks\r
111         SP_TPC = kFALSE; // disable other methods\r
112         QC = kFALSE;\r
113         if(debug) cout << "    --> VZERO RP's " << cutsRP << endl;\r
114     }\r
115 \r
116     AliFlowTrackSimpleCuts *POIfilterLeft = new AliFlowTrackSimpleCuts();\r
117     AliFlowTrackSimpleCuts *POIfilterRight = new AliFlowTrackSimpleCuts();\r
118     if(SP_TPC){\r
119         POIfilterLeft->SetEtaMin(-0.8);\r
120         POIfilterLeft->SetEtaMax(0.0);\r
121         POIfilterLeft->SetMassMin(263731); POIfilterLeft->SetMassMax(263733);\r
122 \r
123         POIfilterRight->SetEtaMin(0.0);\r
124         POIfilterRight->SetEtaMax(0.8);\r
125         POIfilterRight->SetMassMin(263731); POIfilterRight->SetMassMax(263733);\r
126     }\r
127 \r
128     \r
129     AliFlowTrackSimpleCuts *POIfilterVZERO = new AliFlowTrackSimpleCuts();\r
130     if(VZERO_SP || QC){\r
131         POIfilterVZERO->SetEtaMin(-0.8);\r
132         POIfilterVZERO->SetEtaMax(0.8);\r
133         POIfilterVZERO->SetMassMin(263731); POIfilterVZERO->SetMassMax(263733);\r
134  \r
135     }\r
136  \r
137     \r
138     taskHFE->SetRPCuts(cutsRP);\r
139  \r
140  \r
141  AliAnalysisDataContainer *coutput3 = mgr->CreateContainer(Form("ccontainer0_%s",uniqueID.Data()),TList::Class(),AliAnalysisManager::kOutputContainer,fileName);\r
142   \r
143  mgr->ConnectInput(taskHFE,0,mgr->GetCommonInputContainer());\r
144  mgr->ConnectOutput(taskHFE,1,coutput3);\r
145     \r
146  \r
147     if(debug) cout << " === RECEIVED REQUEST FOR FLOW ANALYSIS === " << endl;\r
148     AliAnalysisDataContainer *flowEvent = mgr->CreateContainer(Form("FlowContainer_%s",uniqueID.Data()), AliFlowEventSimple::Class(), AliAnalysisManager::kExchangeContainer);\r
149     mgr->ConnectOutput(taskHFE, 2, flowEvent);\r
150     if(debug) cout << "    --> Created IO containers " << flowEvent << endl;   \r
151     \r
152     \r
153   mgr->AddTask(taskHFE);\r
154     \r
155     if (QC) {  // add qc tasks\r
156         TPCTOF::AddQCmethod(Form("QCTPCin_%s",uniqueID.Data()), harmonic, flowEvent,  debug ,uniqueID, -0.8, -0.0, 0.0, 0.8,false,POIfilterVZERO);\r
157         if(debug) cout << "    --> Hanging QC task ...succes! "<< endl;\r
158     }   \r
159     if (SP_TPC) {  // add sp subevent tasks\r
160         TPCTOF::AddSPmethod(Form("SPTPCQa_in_%s", uniqueID.Data()), -0.8, -.0, .0, +0.8, "Qa", harmonic, flowEvent, false, shrinkSP, debug,uniqueID, false, POIfilterRight);\r
161         if(debug) cout << "    --> Hanging SP Qa task ... succes!" << endl;\r
162         TPCTOF::AddSPmethod(Form("SPTPCQb_in_%s", uniqueID.Data()), -0.8, -.0, .0, +0.8, "Qb", harmonic, flowEvent,  false, shrinkSP, debug,uniqueID, false, POIfilterLeft);\r
163         if(debug) cout << "    --> Hanging SP Qb task ... succes!"<< endl;\r
164     }\r
165     if (VZERO_SP) {  // add sp subevent tasks\r
166         TPCTOF::AddSPmethod(Form("SPVZEROQa_in_%s", uniqueID.Data()), -0.8, -.0, .0, +0.8, "Qa", harmonic, flowEvent, false, shrinkSP, debug,uniqueID, true, POIfilterVZERO);\r
167         if(debug) cout << "    --> Hanging SP Qa task ... succes!" << endl;\r
168         TPCTOF::AddSPmethod(Form("SPVZEROQb_in_%s", uniqueID.Data()), -0.8, -.0, .0, +0.8, "Qb", harmonic, flowEvent,  false, shrinkSP, debug,uniqueID, true, POIfilterVZERO);\r
169         if(debug) cout << "    --> Hanging SP Qb task ... succes!"<< endl;\r
170     }\r
171 \r
172 \r
173     \r
174 return taskHFE;\r
175 \r
176 }\r
177 \r
178 //_____________________________________________________________________________\r
179 //_____________________________________________________________________________\r
180                     \r
181 AliAnalysisTaskFlowTPCTOFQCSP* ConfigHFEStandardCuts(Int_t minTPCCulster,AliHFEextraCuts::ITSPixel_t pixel){\r
182     //\r
183     // HFE standard task configuration\r
184     //\r
185     \r
186     Bool_t kAnalyseTaggedTracks = kTRUE;\r
187     \r
188     AliHFEcuts *hfecuts = new AliHFEcuts("hfeCuts","HFE Standard Cuts");  //TODO....change the cuts values to PbPb\r
189     //  hfecuts->CreateStandardCuts();\r
190     hfecuts->SetMinNClustersTPC(minTPCCulster);\r
191     hfecuts->SetMinNClustersITS(3);\r
192     hfecuts->SetMinNTrackletsTRD(0);\r
193     hfecuts->SetMinRatioTPCclusters(0.6);\r
194     \r
195     //   hfecuts->SetEtaRange(-0.9,0.9);\r
196     //   hfecuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable);\r
197     hfecuts->SetRequireITSPixel();\r
198     hfecuts->SetCutITSpixel(pixel);//kAny\r
199     hfecuts->SetMaxChi2perClusterITS(-1);\r
200     hfecuts->SetMaxChi2perClusterTPC(3.5);\r
201     hfecuts->SetCheckITSLayerStatus(kFALSE); // shud be put back\r
202     //  hfecuts->UnsetVertexRequirement();\r
203     hfecuts->SetVertexRange(10.);\r
204     hfecuts->SetRequireSigmaToVertex();\r
205     //hfecuts->SetSigmaToVertex(10);\r
206     hfecuts->SetTOFPIDStep(kFALSE);\r
207     //  hfecuts->SetQAOn();\r
208     hfecuts->SetPtRange(0, 30);\r
209     \r
210     AliAnalysisTaskFlowTPCTOFQCSP *task = new AliAnalysisTaskFlowTPCTOFQCSP("HFE_Flow_TPCTOF");\r
211     printf("*************************************************************");\r
212     printf("task -------------------------------------------- %p\n ", task);\r
213     printf("*************************************************************\n");\r
214 \r
215     \r
216     task->SetHFECuts(hfecuts);\r
217     \r
218     //   task->SetInvariantMassCut(0.05);\r
219     //  task->SetRejectKinkMother(kTRUE);\r
220     //  task->SetRemovePileUp(kTRUE);\r
221     \r
222     // Define PID\r
223  //   AliHFEpid *pid = task->GetPID();\r
224  //   if(useMC) pid->SetHasMCData(kTRUE);\r
225  //   pid->AddDetector("TPC", 0);\r
226    // pid->AddDetector("EMCAL", 1);\r
227     \r
228  //   printf("*************************************\n");\r
229  //   printf("Configuring standard Task:\n");\r
230     //  task->PrintStatus();\r
231   //  pid->PrintStatus();\r
232   //  printf("*************************************\n");\r
233     return task;\r
234     \r
235     \r
236 }\r
237 \r
238 //_____________________________________________________________________________\r
239 \r
240 namespace TPCTOF{\r
241     //_____________________________________________________________________________\r
242     void AddSPmethod(char *name, double minEtaA, double maxEtaA, double minEtaB, double maxEtaB, char *Qvector, int harmonic, AliAnalysisDataContainer *flowEvent, bool bEP, bool shrink = false, bool debug, TString uniqueID,Bool_t VZERO_SP = kFALSE,  AliFlowTrackSimpleCuts* POIfilter)\r
243     {\r
244         // add sp task and invm filter tasks\r
245         if(debug) (bEP) ? cout << " ****** Reveived request for EP task ****** " << endl : cout << " ******* Switching to SP task ******* " << endl;\r
246         TString fileName = AliAnalysisManager::GetCommonFileName();\r
247         (bEP) ? fileName+=":EP_tpctof" : fileName+=":SP_tpctof";\r
248         //          if(etagap) {\r
249         //            fileName+="_SUBEVENTS";\r
250         //          if(debug) cout << "    --> Setting up subevent analysis <-- " << endl;\r
251         //    }\r
252         if(debug) cout << "    --> fileName " << fileName << endl;\r
253         TString myFolder = fileName;\r
254         if(debug) cout << "    --> myFolder " << myFolder << endl;\r
255         TString myNameSP;\r
256         (bEP) ? myNameSP = Form("%sEPv%d%s", name, harmonic, Qvector): myNameSP = Form("%sSPv%d%s", name, harmonic, Qvector);\r
257         if(debug) cout << " myNameSP " << myNameSP << endl;\r
258         AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();\r
259         AliAnalysisDataContainer *flowEventOut = mgr->CreateContainer(Form("Filter_%s",myNameSP.Data()),AliFlowEventSimple::Class(),AliAnalysisManager::kExchangeContainer);\r
260         AliAnalysisTaskFilterFE *tskFilter = new AliAnalysisTaskFilterFE(Form("TaskFilter_%s", myNameSP.Data()), NULL, POIfilter);\r
261         tskFilter->SetSubeventEtaRange(minEtaA, maxEtaA, minEtaB, maxEtaB);\r
262         if(VZERO_SP) tskFilter->SetSubeventEtaRange(-10, 0, 0, 10);\r
263         mgr->AddTask(tskFilter);\r
264         mgr->ConnectInput(tskFilter, 0, flowEvent);\r
265         mgr->ConnectOutput(tskFilter, 1, flowEventOut);\r
266         AliAnalysisDataContainer *outSP = mgr->CreateContainer(myNameSP.Data(), TList::Class(), AliAnalysisManager::kOutputContainer, fileName);\r
267         AliAnalysisTaskScalarProduct *tskSP = new AliAnalysisTaskScalarProduct(Form("TaskScalarProduct_%s", myNameSP.Data()), kFALSE);\r
268         tskSP->SetApplyCorrectionForNUA(kTRUE);\r
269         tskSP->SetHarmonic(harmonic);\r
270         tskSP->SetTotalQvector(Qvector);\r
271         if (bEP) tskSP->SetBehaveAsEP();\r
272         if (shrink) tskSP->SetBookOnlyBasicCCH(kTRUE);\r
273         mgr->AddTask(tskSP);\r
274         mgr->ConnectInput(tskSP, 0, flowEventOut);\r
275         mgr->ConnectOutput(tskSP, 1, outSP);\r
276     }\r
277     //_____________________________________________________________________________\r
278     void AddQCmethod(char *name, int harmonic, AliAnalysisDataContainer *flowEvent, Bool_t debug, TString uniqueID,double minEtaA, double maxEtaA, double minEtaB, double maxEtaB,Bool_t VZERO_SP = kFALSE,  AliFlowTrackSimpleCuts* POIfilter)\r
279     {\r
280         // add qc task and invm filter tasks\r
281         if(debug) cout << " ****** Received request for QC v" << harmonic << " task " << name << ", IO ****** " << flowEvent << endl;\r
282         TString fileName = AliAnalysisManager::GetCommonFileName();\r
283         fileName+=":QC_tpctof";\r
284         if(debug) cout << "    --> Common filename: " << fileName << endl;\r
285         TString myFolder = Form("v%d", harmonic);\r
286         if(debug) cout << "    --> myFolder: " << myFolder << endl;\r
287         TString myName = Form("%s", name);\r
288         if(debug) cout << "    --> myName: " << myName << endl;\r
289         AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();\r
290         AliAnalysisDataContainer *flowEventOut = mgr->CreateContainer(Form("Filter_%s", myName.Data()), AliFlowEventSimple::Class(),AliAnalysisManager::kExchangeContainer);\r
291         AliAnalysisTaskFilterFE *tskFilter = new AliAnalysisTaskFilterFE(Form("TaskFilter_%s", myName.Data()), NULL, POIfilter);\r
292         tskFilter->SetSubeventEtaRange(minEtaA, maxEtaA, minEtaB, maxEtaB);\r
293         //    if(VZERO_SP) tskFilter->SetSubeventEtaRange(-10, 0, 0, 10);\r
294         mgr->AddTask(tskFilter);\r
295         mgr->ConnectInput(tskFilter, 0, flowEvent);\r
296         mgr->ConnectOutput(tskFilter, 1, flowEventOut);\r
297         \r
298         AliAnalysisDataContainer *outQC = mgr->CreateContainer(myName.Data(), TList::Class(), AliAnalysisManager::kOutputContainer, fileName);\r
299         AliAnalysisTaskQCumulants *tskQC = new AliAnalysisTaskQCumulants(Form("TaskQCumulants_%s", myName.Data()), kFALSE);\r
300         tskQC->SetApplyCorrectionForNUA(kTRUE);\r
301         tskQC->SetHarmonic(harmonic);\r
302         tskQC->SetBookOnlyBasicCCH(kTRUE);\r
303         mgr->AddTask(tskQC);\r
304         mgr->ConnectInput(tskQC, 0, flowEventOut);\r
305         mgr->ConnectOutput(tskQC, 1, outQC);\r
306     }\r
307     //_____________________________________________________________________________\r
308 }\r
309 \r
310 \r
311 \r