]>
Commit | Line | Data |
---|---|---|
1 | class AliAnalysisDataContainer; | |
2 | class AliFlowEventCuts; | |
3 | class AliFlowEventTrackCuts; | |
4 | class AliFlowEventTrackSimpleCuts; | |
5 | class AliRDHFCutsD0toKpi; | |
6 | class AliRDHFCutsDStartoKpipi; | |
7 | ||
8 | void AddTaskFlowD2H(TString fileNameCuts, TString folderName, Int_t nDmeson, Int_t myHarmonic, | |
9 | Bool_t bDoQC, Bool_t bDoSPTPC, Bool_t bDoSPVZERO, Bool_t bDoEPTPC, Bool_t bDoEPVZERO, | |
10 | Int_t ptBinWidth, Double_t gapTPC, Double_t etaVZERO1, Double_t etaVZERO2, Double_t etaVZERO3, Double_t etaVZERO4, | |
11 | Bool_t bOldApproach=kFALSE, Bool_t shrinkSP=kFALSE, Bool_t swapAssumption=kFALSE ) { | |
12 | TFile *filecuts = TFile::Open( fileNameCuts.Data() ); | |
13 | if( (!filecuts) || ( filecuts && !filecuts->IsOpen()) ){ | |
14 | AliFatal("Could not open cuts file."); | |
15 | } | |
16 | ||
17 | TString fileName = AliAnalysisManager::GetCommonFileName(); | |
18 | fileName.ReplaceAll(".root",""); | |
19 | ||
20 | AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); | |
21 | AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer(); | |
22 | ||
23 | TString thecuts = folderName; | |
24 | TString sgn = DMesonName( nDmeson ) + Form("w%d",ptBinWidth); | |
25 | ||
26 | //********************************************************************** | |
27 | // FLOW TRACK CUTS | |
28 | AliFlowTrackCuts* cutsRFPTPC = new AliFlowTrackCuts( "GlobalRFPTPC" ); | |
29 | cutsRFPTPC->SetParamType(AliFlowTrackCuts::kGlobal); | |
30 | cutsRFPTPC->SetPtRange(0.2,5.); | |
31 | cutsRFPTPC->SetEtaRange(-0.8,0.8); | |
32 | cutsRFPTPC->SetMinNClustersTPC(70); | |
33 | cutsRFPTPC->SetMinChi2PerClusterTPC(0.2); | |
34 | cutsRFPTPC->SetMaxChi2PerClusterTPC(4.0); | |
35 | cutsRFPTPC->SetAcceptKinkDaughters(kFALSE); | |
36 | cutsRFPTPC->SetMinimalTPCdedx(10.); | |
37 | cutsRFPTPC->SetAODfilterBit(1); | |
38 | // cutsRFPTPC->SetQA(kTRUE); | |
39 | AliFlowTrackCuts* cutsRFPVZE = new AliFlowTrackCuts( "GlobalRFPVZE" ); | |
40 | cutsRFPVZE->SetParamType(AliFlowTrackCuts::kV0); | |
41 | cutsRFPVZE->SetEtaRange(-10,+10); | |
42 | cutsRFPVZE->SetPhiMin(0); | |
43 | cutsRFPVZE->SetPhiMax(TMath::TwoPi()); | |
44 | // cutsRFPVZE->SetQA(kTRUE); | |
45 | AliFlowTrackCuts* cutsPOI = new AliFlowTrackCuts( "GlobalPOI" ); | |
46 | cutsPOI->SetParamType(AliFlowTrackCuts::kGlobal); | |
47 | cutsPOI->SetPtRange(0.2,-5.); | |
48 | cutsPOI->SetEtaRange(0.8,-0.8); | |
49 | ||
50 | //********************************************************************** | |
51 | // POI FILTER CUTS | |
52 | AliFlowTrackSimpleCuts *filterPOIQC[50]; // MASS BANDS | |
53 | AliFlowTrackSimpleCuts *filterPOISP[50][2]; // MASS BANDS || ETA || SUBEVENT GAP | |
54 | int myMassBands = MassBands( nDmeson, bOldApproach ); | |
55 | for(int mb=0; mb!=myMassBands; ++mb) { | |
56 | filterPOISP[mb][0] = new AliFlowTrackSimpleCuts( Form("FilterPOISP_MB%d_ETANEG",mb) ); | |
57 | filterPOISP[mb][0]->SetEtaMin( -0.8 ); | |
58 | filterPOISP[mb][0]->SetEtaMax( 0. ); | |
59 | filterPOISP[mb][0]->SetMassMin( MassBandLowEdge(nDmeson,mb,bOldApproach) ); | |
60 | filterPOISP[mb][0]->SetMassMax( MassBandLowEdge(nDmeson,mb+1,bOldApproach) ); | |
61 | filterPOISP[mb][1] = new AliFlowTrackSimpleCuts( Form("FilterPOISP_MB%d_ETAPOS",mb) ); | |
62 | filterPOISP[mb][1]->SetEtaMin( 0. ); | |
63 | filterPOISP[mb][1]->SetEtaMax( +0.8 ); | |
64 | filterPOISP[mb][1]->SetMassMin( MassBandLowEdge(nDmeson,mb,bOldApproach) ); | |
65 | filterPOISP[mb][1]->SetMassMax( MassBandLowEdge(nDmeson,mb+1,bOldApproach) ); | |
66 | filterPOIQC[mb] = new AliFlowTrackSimpleCuts( Form("FilterPOIQC_MB%d",mb) ); | |
67 | filterPOIQC[mb]->SetEtaMin( -0.8 ); | |
68 | filterPOIQC[mb]->SetEtaMax( +0.8 ); | |
69 | filterPOIQC[mb]->SetMassMin( MassBandLowEdge(nDmeson,mb,bOldApproach) ); | |
70 | filterPOIQC[mb]->SetMassMax( MassBandLowEdge(nDmeson,mb+1,bOldApproach) ); | |
71 | } | |
72 | ||
73 | // * DMESON SELECTOR ************************************************* | |
74 | AliAnalysisTaskFlowD2H *taskSel; | |
75 | switch (nDmeson) { | |
76 | case ( AliRDHFCuts::kD0toKpiCuts ): | |
77 | AliRDHFCutsD0toKpi *myCutsD0 = (AliRDHFCutsD0toKpi*)filecuts->Get("D0toKpiCuts"); | |
78 | if(!myCutsD0) { | |
79 | AliFatal("Problems reaching D0toKpiCuts"); | |
80 | } | |
81 | taskSel = new AliAnalysisTaskFlowD2H( Form("TaskD0Selector_%s",thecuts.Data()), | |
82 | cutsRFPTPC, cutsRFPVZE, myCutsD0, nDmeson ); | |
83 | break; | |
84 | case ( AliRDHFCuts::kDstarCuts ): | |
85 | AliRDHFCutsDStartoKpipi *myCutsDStar = (AliRDHFCutsDStartoKpipi*)filecuts->Get("DStartoKpipiCuts"); | |
86 | if(!myCutsDStar) { | |
87 | AliFatal("Problems reaching DStarToKpipiCuts"); | |
88 | } | |
89 | taskSel = new AliAnalysisTaskFlowD2H( Form("TaskDStarSelector_%s",thecuts.Data()), | |
90 | cutsRFPTPC, cutsRFPVZE, myCutsDStar, nDmeson); | |
91 | break; | |
92 | case (AliRDHFCuts::kDplusCuts): | |
93 | AliRDHFCutsDplustoKpipi *myCutsDplus = (AliRDHFCutsDplustoKpipi*)filecuts->Get("AnalysisCuts"); | |
94 | if(!myCutsDplus) { | |
95 | AliFatal("Problems reaching AnalysisCuts"); | |
96 | } | |
97 | taskSel = new AliAnalysisTaskFlowD2H( Form("TaskDplusSelector_%s",thecuts.Data()), | |
98 | cutsRFPTPC, cutsRFPVZE, myCutsDplus, nDmeson ); | |
99 | break; | |
100 | } | |
101 | taskSel->SetCommonConstants( MassBins(nDmeson), MinMass(nDmeson), MaxMass(nDmeson), ptBinWidth ); | |
102 | if (swapAssumption) | |
103 | taskSel->SetSwapAsumption(); | |
104 | ||
105 | //taskSel->SelectCollisionCandidates(trigger); | |
106 | //taskSel->SetDebug(); | |
107 | AliAnalysisDataContainer *coutHist = mgr->CreateContainer( Form("%sSelector_%s",sgn.Data(),thecuts.Data()), | |
108 | TList::Class(),AliAnalysisManager::kOutputContainer, | |
109 | Form("%s.root:FlowD2H_%s",fileName.Data(),thecuts.Data()) ); | |
110 | AliAnalysisDataContainer *exc_TPC = mgr->CreateContainer( Form("TPCEventWithCandidates_%s_%s",sgn.Data(),thecuts.Data()), | |
111 | AliFlowEventSimple::Class(), | |
112 | AliAnalysisManager::kExchangeContainer ); | |
113 | AliAnalysisDataContainer *exc_VZE = mgr->CreateContainer( Form("VZEEventWithCandidates_%s_%s",sgn.Data(),thecuts.Data()), | |
114 | AliFlowEventSimple::Class(), | |
115 | AliAnalysisManager::kExchangeContainer ); | |
116 | mgr->AddTask(taskSel); | |
117 | mgr->ConnectOutput(taskSel,1,coutHist); | |
118 | mgr->ConnectOutput(taskSel,2,exc_TPC); | |
119 | mgr->ConnectOutput(taskSel,3,exc_VZE); | |
120 | mgr->ConnectInput (taskSel,0,cinput1); | |
121 | ||
122 | // * HANGING ANALYSIS TASKS ******************************************** | |
123 | int harm = myHarmonic; | |
124 | for(int mb=0; mb!=myMassBands; ++mb) { | |
125 | if(bDoQC) { | |
126 | AddQCmethod( Form("%sQCTPCMB%d",sgn.Data(),mb), fileName.Data(), thecuts.Data(), harm, exc_TPC, filterPOIQC[mb]); | |
127 | } | |
128 | if(bDoSPVZERO) { | |
129 | AddSPmethod( Form("%sSPVZEMB%d",sgn.Data(),mb), fileName.Data(), thecuts.Data(), etaVZERO1, etaVZERO2, etaVZERO3, etaVZERO4, | |
130 | "Qa", harm, exc_VZE, 0, filterPOIQC[mb], NULL, false, shrinkSP ); | |
131 | AddSPmethod( Form("%sSPVZEMB%d",sgn.Data(),mb), fileName.Data(), thecuts.Data(), etaVZERO1, etaVZERO2, etaVZERO3, etaVZERO4, | |
132 | "Qb", harm, exc_VZE, 0, filterPOIQC[mb], NULL, false, shrinkSP ); | |
133 | } | |
134 | if(bDoEPVZERO) { | |
135 | AddSPmethod( Form("%sEPVZEMB%d",sgn.Data(),mb), fileName.Data(), thecuts.Data(), etaVZERO1, etaVZERO2, etaVZERO3, etaVZERO4, | |
136 | "QaQb", harm, exc_VZE, 0, filterPOIQC[mb], NULL, true, shrinkSP ); | |
137 | } | |
138 | if(bDoSPTPC) { | |
139 | for(int eg=0; eg!=2; ++eg) { | |
140 | AddSPmethod( Form("%sSPTPCMB%d",sgn.Data(),mb), fileName.Data(), thecuts.Data(), -0.8, -gapTPC*eg, +gapTPC*eg, +0.8, | |
141 | "Qa", harm, exc_TPC, eg, filterPOISP[mb][1], NULL, false, shrinkSP ); | |
142 | AddSPmethod( Form("%sSPTPCMB%d",sgn.Data(),mb), fileName.Data(), thecuts.Data(), -0.8, -gapTPC*eg, +gapTPC*eg, +0.8, | |
143 | "Qb", harm, exc_TPC, eg, filterPOISP[mb][0], NULL, false, shrinkSP ); | |
144 | } | |
145 | } | |
146 | if(bDoEPTPC) { | |
147 | AddSPmethod( Form("%sEPTPCMB%d",sgn.Data(),mb), fileName.Data(), thecuts.Data(), -0.8, -0.0, +0.0, +0.8, | |
148 | "QaQb", harm, exc_TPC, 0, filterPOIQC[mb], NULL, true, shrinkSP ); | |
149 | } | |
150 | } | |
151 | } | |
152 | ||
153 | void AddSPmethod(char *name, char *fileName, char *thecuts, | |
154 | double minEtaA, double maxEtaA, double minEtaB, double maxEtaB, | |
155 | char *Qvector, int harmonic, AliAnalysisDataContainer *flowEvent, int eg, | |
156 | AliFlowTrackSimpleCuts *cutsPOI=NULL, AliFlowTrackSimpleCuts *cutsRFP=NULL, | |
157 | bool bEP, bool shrink=false ) { | |
158 | TString myFolder = Form("v%d_%s",harmonic,thecuts); | |
159 | TString myNameSP = Form("%sSPv%d%sGAP%d_%s",name,harmonic,Qvector,eg,thecuts); | |
160 | AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); | |
161 | AliAnalysisDataContainer *flowEvent2 = mgr->CreateContainer( Form("Filter_%s", myNameSP.Data()), | |
162 | AliFlowEventSimple::Class(), | |
163 | AliAnalysisManager::kExchangeContainer ); | |
164 | AliAnalysisTaskFilterFE *tskFilter = new AliAnalysisTaskFilterFE( Form("TaskFilter_%s",myNameSP.Data()), | |
165 | cutsRFP, cutsPOI); | |
166 | tskFilter->SetSubeventEtaRange( minEtaA, maxEtaA, minEtaB, maxEtaB ); | |
167 | mgr->AddTask(tskFilter); | |
168 | mgr->ConnectInput( tskFilter,0,flowEvent); | |
169 | mgr->ConnectOutput(tskFilter,1,flowEvent2); | |
170 | AliAnalysisDataContainer *outSP = mgr->CreateContainer( myNameSP.Data(), | |
171 | TList::Class(),AliAnalysisManager::kOutputContainer, | |
172 | Form("%s.root:FlowD2H_SP_%s",fileName,myFolder.Data()) ); | |
173 | AliAnalysisTaskScalarProduct *tskSP = new AliAnalysisTaskScalarProduct( Form("TaskScalarProduct_%s", | |
174 | myNameSP.Data()),kFALSE); | |
175 | tskSP->SetApplyCorrectionForNUA(kTRUE); | |
176 | tskSP->SetHarmonic(harmonic); | |
177 | tskSP->SetTotalQvector(Qvector); | |
178 | if(bEP) tskSP->SetBehaveAsEP(); | |
179 | if(shrink) tskSP->SetBookOnlyBasicCCH(kTRUE); | |
180 | mgr->AddTask(tskSP); | |
181 | mgr->ConnectInput( tskSP,0,flowEvent2); | |
182 | mgr->ConnectOutput(tskSP,1,outSP); | |
183 | } | |
184 | ||
185 | void AddQCmethod(char *name, char *fileName, char *thecuts, | |
186 | int harmonic, AliAnalysisDataContainer *flowEvent, | |
187 | AliFlowTrackSimpleCuts *cutsPOI=NULL, AliFlowTrackSimpleCuts *cutsRFP=NULL) { | |
188 | TString myFolder = Form("v%d_%s",harmonic,thecuts); | |
189 | TString myName = Form("%sv%d_%s",name,harmonic,thecuts); | |
190 | AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); | |
191 | AliAnalysisDataContainer *flowEvent2 = mgr->CreateContainer( Form("Filter_%s", myName.Data()), | |
192 | AliFlowEventSimple::Class(), | |
193 | AliAnalysisManager::kExchangeContainer ); | |
194 | AliAnalysisTaskFilterFE *tskFilter = new AliAnalysisTaskFilterFE( Form("TaskFilter_%s",myName.Data()), | |
195 | cutsRFP, cutsPOI); | |
196 | mgr->AddTask(tskFilter); | |
197 | mgr->ConnectInput( tskFilter,0,flowEvent); | |
198 | mgr->ConnectOutput(tskFilter,1,flowEvent2); | |
199 | AliAnalysisDataContainer *outQC = mgr->CreateContainer( myName.Data(), | |
200 | TList::Class(),AliAnalysisManager::kOutputContainer, | |
201 | Form("%s.root:FlowD2H_QC_%s",fileName,myFolder.Data()) ); | |
202 | AliAnalysisTaskQCumulants *tskQC = new AliAnalysisTaskQCumulants( Form("TaskQCumulants_%s", | |
203 | myName.Data()),kFALSE); | |
204 | tskQC->SetApplyCorrectionForNUA(kTRUE); | |
205 | tskQC->SetHarmonic(harmonic); | |
206 | tskQC->SetBookOnlyBasicCCH(kTRUE); | |
207 | mgr->AddTask(tskQC); | |
208 | mgr->ConnectInput( tskQC,0,flowEvent2); | |
209 | mgr->ConnectOutput(tskQC,1,outQC); | |
210 | } | |
211 | ||
212 | int MassBands( int nDmeson, bool bOldApproach=false ) { | |
213 | switch (nDmeson) { | |
214 | case ( AliRDHFCuts::kDplusCuts ): | |
215 | case ( AliRDHFCuts::kD0toKpiCuts ): | |
216 | if(bOldApproach) return 5; | |
217 | else return 26; | |
218 | case ( AliRDHFCuts::kDstarCuts ): | |
219 | return 13; | |
220 | } | |
221 | } | |
222 | ||
223 | double MassBandLowEdge( int nDmeson, int mb, bool bOldApproach=false ) { | |
224 | switch (nDmeson) { | |
225 | case ( AliRDHFCuts::kDplusCuts ): | |
226 | case ( AliRDHFCuts::kD0toKpiCuts ): // 2 + 20 + 4 | |
227 | double lowEdgeMinimal[5+1] = {1.75,1.80,1.83,1.90,1.93,2.03}; | |
228 | double lowEdge[26+1] = { 1.66, 1.71, 1.76, 1.77, 1.78, 1.79, 1.80, 1.81, 1.82, 1.83, | |
229 | 1.84, 1.85, 1.86, 1.87, 1.88, 1.89, 1.90, 1.91, 1.92, 1.93, | |
230 | 1.94, 1.95, 1.96, 2.01, 2.06, 2.11, 2.16 }; | |
231 | if(bOldApproach) return lowEdgeMinimal[mb]; | |
232 | else return lowEdge[mb]; | |
233 | case ( AliRDHFCuts::kDstarCuts ): // 2 + 10 + 3 | |
234 | double lowEdge[13+1]={0.138, 0.139, 0.141, 0.143, 0.144, 0.145, 0.146, 0.147, 0.148, 0.150, 0.152, 0.154, 0.156, 0.158}; | |
235 | return lowEdge[mb]; | |
236 | } | |
237 | } | |
238 | ||
239 | double MinMass( int nDmeson ) { | |
240 | switch (nDmeson) { | |
241 | case ( AliRDHFCuts::kDplusCuts ): | |
242 | case ( AliRDHFCuts::kD0toKpiCuts ): | |
243 | return 1.66; | |
244 | case ( AliRDHFCuts::kDstarCuts ): | |
245 | return 0.138; | |
246 | } | |
247 | } | |
248 | ||
249 | double MaxMass( int nDmeson ) { | |
250 | switch (nDmeson) { | |
251 | case ( AliRDHFCuts::kDplusCuts ): | |
252 | case ( AliRDHFCuts::kD0toKpiCuts ): | |
253 | return 2.16; | |
254 | case ( AliRDHFCuts::kDstarCuts ): | |
255 | return 0.158; | |
256 | } | |
257 | } | |
258 | ||
259 | int MassBins( int nDmeson ) { | |
260 | switch (nDmeson) { | |
261 | case ( AliRDHFCuts::kDplusCuts ): | |
262 | case ( AliRDHFCuts::kD0toKpiCuts ): | |
263 | return 50; | |
264 | case ( AliRDHFCuts::kDstarCuts ): | |
265 | return 40; | |
266 | } | |
267 | } | |
268 | ||
269 | TString DMesonName( int nDmeson ) { | |
270 | TString toReturn; | |
271 | switch (nDmeson) { | |
272 | case ( AliRDHFCuts::kDplusCuts ): | |
273 | toReturn = "DPlus"; break; | |
274 | case ( AliRDHFCuts::kD0toKpiCuts ): | |
275 | toReturn = "D0"; break; | |
276 | case ( AliRDHFCuts::kDstarCuts ): | |
277 | toReturn = "DStar"; break; | |
278 | } | |
279 | return toReturn; | |
280 | } |