]>
Commit | Line | Data |
---|---|---|
8d312f00 | 1 | /************************************************************************* |
2 | * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
02914cfa | 16 | #define ALIFLOWANALYSISWITHSCALARPRODUCT_CXX |
8d312f00 | 17 | |
929098e4 | 18 | #include "TFile.h" |
d7eb18ec | 19 | #include "TList.h" |
8d312f00 | 20 | #include "TMath.h" |
02914cfa | 21 | #include "TString.h" |
8d312f00 | 22 | #include "TProfile.h" |
23 | #include "TVector2.h" | |
a554b203 | 24 | #include "TH1D.h" |
ad1b7114 | 25 | #include "TH1F.h" |
a554b203 | 26 | #include "TH2D.h" |
8d312f00 | 27 | |
929098e4 | 28 | #include "AliFlowCommonConstants.h" |
8d312f00 | 29 | #include "AliFlowEventSimple.h" |
a93de0f0 | 30 | #include "AliFlowVector.h" |
8d312f00 | 31 | #include "AliFlowTrackSimple.h" |
32 | #include "AliFlowCommonHist.h" | |
395fadba | 33 | #include "AliFlowCommonHistResults.h" |
8d312f00 | 34 | #include "AliFlowAnalysisWithScalarProduct.h" |
8d312f00 | 35 | |
12dc476e | 36 | ////////////////////////////////////////////////////////////////////////////// |
8d312f00 | 37 | // AliFlowAnalysisWithScalarProduct: |
02914cfa | 38 | // Description: Maker to analyze Flow from the Scalar Product method. |
39 | // authors: Naomi van del Kolk (kolk@nikhef.nl) | |
40 | // Ante Bilandzic (anteb@nikhef.nl) | |
41 | // mods: Carlos Perez (cperez@nikhef.nl) | |
12dc476e | 42 | ////////////////////////////////////////////////////////////////////////////// |
8d312f00 | 43 | |
44 | ClassImp(AliFlowAnalysisWithScalarProduct) | |
45 | ||
02914cfa | 46 | //----------------------------------------------------------------------- |
47 | AliFlowAnalysisWithScalarProduct::AliFlowAnalysisWithScalarProduct(): | |
48 | fDebug(0), | |
49 | fUsePhiWeights(0), | |
50 | fApplyCorrectionForNUA(0), | |
51 | fHarmonic(2), | |
52 | fNormalizationType(1), | |
53 | fTotalQvector(3), | |
54 | fWeightsList(NULL), | |
55 | fHistList(NULL), | |
56 | fHistProConfig(NULL), | |
57 | fHistProQaQbNorm(NULL), | |
58 | fHistSumOfWeights(NULL), | |
59 | fHistProNUAq(NULL), | |
60 | fHistProQNorm(NULL), | |
61 | fHistProQaQb(NULL), | |
62 | fHistProQaQbM(NULL), | |
63 | fHistMaMb(NULL), | |
64 | fHistQNormQaQbNorm(NULL), | |
65 | fHistQaNormMa(NULL), | |
66 | fHistQbNormMb(NULL), | |
67 | fResolution(NULL), | |
68 | fHistQaQb(NULL), | |
69 | fHistQaQbCos(NULL), | |
70 | fCommonHists(NULL), | |
71 | fCommonHistsuQ(NULL), | |
72 | fCommonHistsRes(NULL) | |
8d312f00 | 73 | { |
02914cfa | 74 | //ctor |
75 | for(int i=0; i!=2; ++i) { | |
76 | fPhiWeightsSub[i] = NULL; | |
77 | for(int j=0; j!=2; ++j) { | |
78 | fHistProUQ[i][j] = NULL; | |
79 | fHistProUQQaQb[i][j] = NULL; | |
80 | fHistProNUAu[i][j][0] = NULL; | |
81 | fHistProNUAu[i][j][1] = NULL; | |
82 | for(int k=0; k!=3; ++k) | |
83 | fHistSumOfWeightsu[i][j][k] = NULL; | |
29ecccee | 84 | } |
29ecccee | 85 | } |
8d312f00 | 86 | } |
1dfa3c16 | 87 | //----------------------------------------------------------------------- |
02914cfa | 88 | AliFlowAnalysisWithScalarProduct::~AliFlowAnalysisWithScalarProduct() |
1dfa3c16 | 89 | { |
02914cfa | 90 | //destructor |
91 | delete fWeightsList; | |
92 | delete fHistList; | |
1dfa3c16 | 93 | } |
8d312f00 | 94 | //----------------------------------------------------------------------- |
95 | void AliFlowAnalysisWithScalarProduct::Init() { | |
8d312f00 | 96 | //Define all histograms |
02914cfa | 97 | printf("---Analysis with the Scalar Product Method--- Init\n"); |
98 | printf("--- fNormalizationType %d ---\n", fNormalizationType); | |
489d5531 | 99 | //save old value and prevent histograms from being added to directory |
100 | //to avoid name clashes in case multiple analaysis objects are used | |
101 | //in an analysis | |
102 | Bool_t oldHistAddStatus = TH1::AddDirectoryStatus(); | |
103 | TH1::AddDirectory(kFALSE); | |
104 | ||
02914cfa | 105 | fHistList = new TList(); |
106 | fHistList->SetName("cobjSP"); | |
107 | fHistList->SetOwner(); | |
108 | ||
109 | TList *uQRelated = new TList(); | |
110 | uQRelated->SetName("uQ"); | |
111 | uQRelated->SetOwner(); | |
112 | ||
113 | TList *nuaRelated = new TList(); | |
114 | nuaRelated->SetName("NUA"); | |
115 | nuaRelated->SetOwner(); | |
116 | ||
117 | TList *errorRelated = new TList(); | |
118 | errorRelated->SetName("error"); | |
119 | errorRelated->SetOwner(); | |
120 | ||
9812922b | 121 | TList *tQARelated = new TList(); |
122 | tQARelated->SetName("QA"); | |
123 | tQARelated->SetOwner(); | |
02914cfa | 124 | |
ad1b7114 | 125 | fCommonHists = new AliFlowCommonHist("AliFlowCommonHist_SP"); |
02914cfa | 126 | (fCommonHists->GetHarmonic())->Fill(0.5,fHarmonic); // store harmonic |
127 | fHistList->Add(fCommonHists); | |
128 | ||
ad1b7114 | 129 | fCommonHistsuQ = new AliFlowCommonHist("AliFlowCommonHist_uQ"); |
02914cfa | 130 | (fCommonHistsuQ->GetHarmonic())->Fill(0.5,fHarmonic); // store harmonic |
131 | fHistList->Add(fCommonHistsuQ); | |
132 | ||
ad1b7114 | 133 | fCommonHistsRes = new AliFlowCommonHistResults("AliFlowCommonHistResults_SP","",fHarmonic); |
02914cfa | 134 | fHistList->Add(fCommonHistsRes); |
135 | ||
136 | fHistProConfig = new TProfile("FlowPro_Flags_SP","Flow_Flags_SP",4,0.5,4.5,"s"); | |
137 | fHistProConfig->GetXaxis()->SetBinLabel(1,"fApplyCorrectionForNUA"); | |
138 | fHistProConfig->GetXaxis()->SetBinLabel(2,"fNormalizationType"); | |
139 | fHistProConfig->GetXaxis()->SetBinLabel(3,"fUsePhiWeights"); | |
140 | fHistProConfig->GetXaxis()->SetBinLabel(4,"fHarmonic"); | |
141 | fHistProConfig->Fill(1,fApplyCorrectionForNUA); | |
142 | fHistProConfig->Fill(2,fNormalizationType); | |
143 | fHistProConfig->Fill(3,fUsePhiWeights); | |
144 | fHistProConfig->Fill(4,fHarmonic); | |
145 | fHistList->Add(fHistProConfig); | |
381342ec | 146 | |
a554b203 | 147 | fHistProQaQbNorm = new TProfile("FlowPro_QaQbNorm_SP","FlowPro_QaQbNorm_SP", 1, 0.5, 1.5,"s"); |
02914cfa | 148 | fHistProQaQbNorm->SetYTitle("<QaQb/Na/Nb>"); |
149 | errorRelated->Add(fHistProQaQbNorm); | |
150 | ||
ad1b7114 | 151 | fHistProNUAq = new TProfile("FlowPro_NUAq_SP","FlowPro_NUAq_SP", 6, 0.5, 6.5,"s"); |
02914cfa | 152 | fHistProNUAq->GetXaxis()->SetBinLabel( 1,"<<sin(#Phi_{a})>>"); |
153 | fHistProNUAq->GetXaxis()->SetBinLabel( 2,"<<cos(#Phi_{a})>>"); | |
154 | fHistProNUAq->GetXaxis()->SetBinLabel( 3,"<<sin(#Phi_{b})>>"); | |
155 | fHistProNUAq->GetXaxis()->SetBinLabel( 4,"<<cos(#Phi_{b})>>"); | |
156 | fHistProNUAq->GetXaxis()->SetBinLabel( 5,"<<sin(#Phi_{t})>>"); | |
157 | fHistProNUAq->GetXaxis()->SetBinLabel( 6,"<<cos(#Phi_{t})>>"); | |
158 | nuaRelated->Add(fHistProNUAq); | |
159 | ||
160 | fHistSumOfWeights = new TH1D("Flow_SumOfWeights_SP","Flow_SumOfWeights_SP",2,0.5, 2.5); | |
161 | fHistSumOfWeights->GetXaxis()->SetBinLabel( 1,"#Sigma Na*Nb"); | |
162 | fHistSumOfWeights->GetXaxis()->SetBinLabel( 2,"#Sigma (Na*Nb)^2"); | |
163 | errorRelated->Add(fHistSumOfWeights); | |
164 | ||
165 | TString sPOI[2] = {"RP","POI"}; // backward compatibility | |
166 | TString sEta[2] = {"Pt","eta"}; // backward compatibility | |
167 | TString sTitle[2] = {"p_{T} [GeV]","#eta"}; | |
168 | TString sWeights[3] = {"uQ","uQuQ","uQQaQb"}; | |
169 | Int_t iNbins[2]; | |
170 | Double_t dMin[2], dMax[2]; | |
171 | iNbins[0] = AliFlowCommonConstants::GetMaster()->GetNbinsPt(); | |
172 | iNbins[1] = AliFlowCommonConstants::GetMaster()->GetNbinsEta(); | |
173 | dMin[0] = AliFlowCommonConstants::GetMaster()->GetPtMin(); | |
174 | dMin[1] = AliFlowCommonConstants::GetMaster()->GetEtaMin(); | |
175 | dMax[0] = AliFlowCommonConstants::GetMaster()->GetPtMax(); | |
176 | dMax[1] = AliFlowCommonConstants::GetMaster()->GetEtaMax(); | |
177 | for(Int_t iPOI=0; iPOI!=2; ++iPOI) for(Int_t iSpace=0; iSpace!=2; ++iSpace) { | |
178 | // uQ | |
ad1b7114 | 179 | fHistProUQ[iPOI][iSpace] = new TProfile( Form( "FlowPro_UQ_%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ), |
02914cfa | 180 | Form( "FlowPro_UQ%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ), |
181 | iNbins[iSpace], dMin[iSpace], dMax[iSpace], "s"); | |
182 | fHistProUQ[iPOI][iSpace]->SetXTitle(sTitle[iSpace].Data()); | |
183 | fHistProUQ[iPOI][iSpace]->SetYTitle("<uQ>"); | |
184 | uQRelated->Add(fHistProUQ[iPOI][iSpace]); | |
185 | ||
186 | // NUAu | |
187 | fHistProNUAu[iPOI][iSpace][0] = new TProfile( Form("FlowPro_NUAu_%s%s_IM_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ), | |
188 | Form("FlowPro_NUAu_%s%s_IM_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ), | |
189 | iNbins[iSpace], dMin[iSpace], dMax[iSpace]); | |
190 | fHistProNUAu[iPOI][iSpace][0]->SetXTitle(sTitle[iSpace].Data()); | |
191 | nuaRelated->Add(fHistProNUAu[iPOI][iSpace][0]); | |
192 | fHistProNUAu[iPOI][iSpace][1] = new TProfile( Form("FlowPro_NUAu_%s%s_RE_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ), | |
193 | Form("FlowPro_NUAu_%s%s_RE_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ), | |
194 | iNbins[iSpace], dMin[iSpace], dMax[iSpace]); | |
195 | fHistProNUAu[iPOI][iSpace][1]->SetXTitle(sTitle[iSpace].Data()); | |
196 | nuaRelated->Add(fHistProNUAu[iPOI][iSpace][1]); | |
197 | ||
198 | // uQ QaQb | |
199 | fHistProUQQaQb[iPOI][iSpace] = new TProfile( Form("FlowPro_UQQaQb_%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ), | |
200 | Form("FlowPro_UQQaQb_%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ), | |
201 | iNbins[iSpace], dMin[iSpace], dMax[iSpace]); | |
202 | fHistProUQQaQb[iPOI][iSpace]->SetXTitle(sTitle[iSpace].Data()); | |
203 | fHistProUQQaQb[iPOI][iSpace]-> SetYTitle("<Qu QaQb>"); | |
204 | errorRelated->Add(fHistProUQQaQb[iPOI][iSpace]); | |
205 | ||
206 | // uWeights | |
207 | for(Int_t i=0; i!=3; ++i) { | |
208 | fHistSumOfWeightsu[iPOI][iSpace][i] = new TH1D( Form("Flow_SumOfWeights_%s%s_%s_SP",sWeights[i].Data(),sPOI[iPOI].Data(),sEta[iSpace].Data()), | |
209 | Form("Flow_SumOfWeights_%s%s_%s_SP",sWeights[i].Data(),sPOI[iPOI].Data(),sEta[iSpace].Data()), | |
210 | iNbins[iSpace], dMin[iSpace], dMax[iSpace]); | |
211 | fHistSumOfWeightsu[iPOI][iSpace][i]->SetYTitle(sWeights[i].Data()); | |
212 | fHistSumOfWeightsu[iPOI][iSpace][i]->SetXTitle(sTitle[iSpace].Data()); | |
213 | errorRelated->Add(fHistSumOfWeightsu[iPOI][iSpace][i]); | |
214 | } | |
12dc476e | 215 | } |
29195b69 | 216 | //weights |
217 | if(fUsePhiWeights) { | |
218 | if(!fWeightsList) { | |
02914cfa | 219 | printf( "WARNING: fWeightsList is NULL in the Scalar Product method.\n" ); |
29195b69 | 220 | exit(0); |
221 | } | |
02914cfa | 222 | fPhiWeightsSub[0] = dynamic_cast<TH1F*> |
223 | (fWeightsList->FindObject("phi_weights_sub0")); | |
224 | if(!fPhiWeightsSub[0]) { | |
225 | printf( "WARNING: phi_weights_sub0 not found in the Scalar Product method.\n" ); | |
226 | exit(0); | |
04c6b875 | 227 | } |
02914cfa | 228 | nuaRelated->Add( fPhiWeightsSub[0] ); |
229 | fPhiWeightsSub[1] = dynamic_cast<TH1F*> | |
230 | (fWeightsList->FindObject("phi_weights_sub1")); | |
231 | if(!fPhiWeightsSub[1]) { | |
232 | printf( "WARNING: phi_weights_sub1 not found in the Scalar Product method.\n" ); | |
233 | exit(0); | |
234 | } | |
235 | nuaRelated->Add( fPhiWeightsSub[1] ); | |
236 | } | |
04c6b875 | 237 | |
e2d51347 | 238 | |
02914cfa | 239 | fHistProQNorm = new TProfile("FlowPro_QNorm_SP","FlowPro_QNorm_SP", 1,0.5,1.5,"s"); |
240 | fHistProQNorm->SetYTitle("<|Qa+Qb|>"); | |
9812922b | 241 | tQARelated->Add(fHistProQNorm); |
8d312f00 | 242 | |
02914cfa | 243 | fHistProQaQb = new TProfile("FlowPro_QaQb_SP","FlowPro_QaQb_SP", 1,0.5,1.5,"s"); |
244 | fHistProQaQb->SetYTitle("<QaQb>"); | |
9812922b | 245 | tQARelated->Add(fHistProQaQb); |
02914cfa | 246 | |
247 | fHistProQaQbM = new TProfile("FlowPro_QaQbvsM_SP","FlowPro_QaQbvsM_SP",1000,0.0,10000); | |
248 | fHistProQaQbM->SetYTitle("<QaQb>"); | |
249 | fHistProQaQbM->SetXTitle("M"); | |
250 | fHistProQaQbM->Sumw2(); | |
9812922b | 251 | tQARelated->Add(fHistProQaQbM); |
02914cfa | 252 | |
253 | fHistMaMb = new TH2D("Flow_MavsMb_SP","Flow_MavsMb_SP",100,0.,100.,100,0.,100.); | |
254 | fHistMaMb->SetYTitle("Ma"); | |
255 | fHistMaMb->SetXTitle("Mb"); | |
9812922b | 256 | tQARelated->Add(fHistMaMb); |
02914cfa | 257 | |
258 | fHistQNormQaQbNorm = new TH2D("Flow_QNormvsQaQbNorm_SP","Flow_QNormvsQaQbNorm_SP",88,-1.1,1.1,22,0.,1.1); | |
259 | fHistQNormQaQbNorm->SetYTitle("|Q/Mq|"); | |
260 | fHistQNormQaQbNorm->SetXTitle("QaQb/MaMb"); | |
9812922b | 261 | tQARelated->Add(fHistQNormQaQbNorm); |
02914cfa | 262 | |
263 | fHistQaNormMa = new TH2D("Flow_QaNormvsMa_SP","Flow_QaNormvsMa_SP",100,0.,100.,22,0.,1.1); | |
264 | fHistQaNormMa->SetYTitle("|Qa/Ma|"); | |
265 | fHistQaNormMa->SetXTitle("Ma"); | |
9812922b | 266 | tQARelated->Add(fHistQaNormMa); |
02914cfa | 267 | |
268 | fHistQbNormMb = new TH2D("Flow_QbNormvsMb_SP","Flow_QbNormvsMb_SP",100,0.,100.,22,0.,1.1); | |
269 | fHistQbNormMb->SetYTitle("|Qb/Mb|"); | |
270 | fHistQbNormMb->SetXTitle("Mb"); | |
9812922b | 271 | tQARelated->Add(fHistQbNormMb); |
02914cfa | 272 | |
273 | fResolution = new TH1D("Flow_resolution_SP","Flow_resolution_SP",100,-1.0,1.0); | |
274 | fResolution->SetYTitle("dN/d(Cos2(#phi_a - #phi_b))"); | |
275 | fResolution->SetXTitle("Cos2(#phi_a - #phi_b)"); | |
9812922b | 276 | tQARelated->Add(fResolution); |
1b302085 | 277 | |
02914cfa | 278 | fHistQaQb = new TH1D("Flow_QaQb_SP","Flow_QaQb_SP",200,-100.,100.); |
279 | fHistQaQb->SetYTitle("dN/dQaQb"); | |
280 | fHistQaQb->SetXTitle("dQaQb"); | |
9812922b | 281 | tQARelated->Add(fHistQaQb); |
1b302085 | 282 | |
02914cfa | 283 | fHistQaQbCos = new TH1D("Flow_QaQbCos_SP","Flow_QaQbCos_SP",63,0.,TMath::Pi()); |
284 | fHistQaQbCos->SetYTitle("dN/d#phi"); | |
285 | fHistQaQbCos->SetXTitle("#phi"); | |
9812922b | 286 | tQARelated->Add(fHistQaQbCos); |
1b302085 | 287 | |
02914cfa | 288 | fHistList->Add(uQRelated); |
289 | fHistList->Add(nuaRelated); | |
290 | fHistList->Add(errorRelated); | |
9812922b | 291 | fHistList->Add(tQARelated); |
1b302085 | 292 | |
02914cfa | 293 | TH1::AddDirectory(oldHistAddStatus); |
1b302085 | 294 | } |
295 | ||
296 | //----------------------------------------------------------------------- | |
02914cfa | 297 | void AliFlowAnalysisWithScalarProduct::Make(AliFlowEventSimple* anEvent) { |
298 | // Scalar Product method | |
299 | if (!anEvent) return; // for coverity | |
300 | ||
301 | // Get Q vectors for the subevents | |
302 | AliFlowVector* vQarray = new AliFlowVector[2]; | |
303 | if (fUsePhiWeights) | |
304 | anEvent->Get2Qsub(vQarray,fHarmonic,fWeightsList,kTRUE); | |
305 | else | |
306 | anEvent->Get2Qsub(vQarray,fHarmonic); | |
307 | // Subevent a | |
308 | AliFlowVector vQa = vQarray[0]; | |
309 | // Subevent b | |
310 | AliFlowVector vQb = vQarray[1]; | |
311 | delete [] vQarray; | |
312 | ||
313 | Double_t dMa = vQa.GetMult(); | |
314 | if( dMa < 2 ) return; | |
315 | Double_t dMb = vQb.GetMult(); | |
316 | if( dMb < 2 ) return; | |
317 | //fill control histograms | |
318 | if (fUsePhiWeights) { | |
319 | fCommonHists->FillControlHistograms(anEvent,fWeightsList,kTRUE); | |
320 | } else { | |
321 | fCommonHists->FillControlHistograms(anEvent); | |
322 | } | |
1b302085 | 323 | |
02914cfa | 324 | //Normalizing: weight the Q vectors for the subevents |
325 | Double_t dNa = fNormalizationType ? dMa: vQa.Mod(); // SP corresponds to true | |
326 | Double_t dNb = fNormalizationType ? dMb: vQb.Mod(); // SP corresponds to true | |
327 | Double_t dWa = fNormalizationType ? dMa: 1; // SP corresponds to true | |
328 | Double_t dWb = fNormalizationType ? dMb: 1; // SP corresponds to true | |
12dc476e | 329 | |
02914cfa | 330 | //Scalar product of the two subevents vectors |
331 | Double_t dQaQb = (vQa*vQb); | |
8d312f00 | 332 | |
02914cfa | 333 | //printf("==============\n"); |
334 | //printf("vQa: { %f, %f }\n",vQa.X(),vQa.Y()); | |
335 | //printf("QaQb/|Qa||Qb|: %f\n",dQaQb/vQa.Mod()/vQb.Mod()); | |
336 | //printf("QaQb/|Ma||Mb|: %f\n",dQaQb/dMa/dMb); | |
337 | ||
338 | // 01 10 11 <=== fTotalQVector | |
339 | // Q ?= Qa or Qb or QaQb | |
340 | AliFlowVector vQm; | |
341 | vQm.Set(0.0,0.0); | |
342 | Double_t dNq=0; | |
343 | if( (fTotalQvector%2)>0 ) { // 01 or 11 | |
344 | vQm += vQa; | |
345 | dNq += dMa; | |
346 | } | |
347 | if( fTotalQvector>1 ) { // 10 or 11 | |
348 | vQm += vQb; | |
349 | dNq += dMb; | |
350 | } | |
351 | Double_t dWq = fNormalizationType ? dNq: 1; // SP corresponds to true | |
352 | dNq = fNormalizationType ? dNq: vQm.Mod(); // SP corresponds to true | |
353 | ||
354 | //fill some EP control histograms | |
355 | fHistProQaQbNorm->Fill(1.,dQaQb/dNa/dNb,dWa*dWb); //Fill (QaQb/NaNb) with weight (WaWb). | |
356 | //needed for the error calculation: | |
357 | fHistSumOfWeights -> Fill(1.,dWa*dWb); | |
358 | fHistSumOfWeights -> Fill(2.,pow(dWa*dWb,2.)); | |
359 | //needed for correcting non-uniform acceptance: | |
360 | fHistProNUAq->Fill(1.,vQa.Y()/dNa,dWa); // to get <<sin(phi_a)>> | |
361 | fHistProNUAq->Fill(2.,vQa.X()/dNa,dWa); // to get <<cos(phi_a)>> | |
362 | fHistProNUAq->Fill(3.,vQb.Y()/dNb,dWb); // to get <<sin(phi_b)>> | |
363 | fHistProNUAq->Fill(4.,vQb.X()/dNb,dWb); // to get <<cos(phi_b)>> | |
364 | fHistProNUAq->Fill(5.,vQm.Y()/dNq,dWq); | |
365 | fHistProNUAq->Fill(6.,vQm.X()/dNq,dWq); | |
366 | ||
367 | //loop over the tracks of the event | |
368 | AliFlowTrackSimple* pTrack = NULL; | |
369 | Int_t iNumberOfTracks = anEvent->NumberOfTracks(); | |
370 | for (Int_t i=0;i<iNumberOfTracks;i++) { | |
371 | pTrack = anEvent->GetTrack(i) ; | |
372 | if (!pTrack) continue; | |
373 | Double_t dPhi = pTrack->Phi(); | |
374 | Double_t dPt = pTrack->Pt(); | |
375 | Double_t dEta = pTrack->Eta(); | |
376 | ||
377 | //calculate vU | |
378 | TVector2 vU; | |
379 | Double_t dUX = TMath::Cos(fHarmonic*dPhi); | |
380 | Double_t dUY = TMath::Sin(fHarmonic*dPhi); | |
381 | vU.Set(dUX,dUY); | |
382 | ||
383 | // 01 10 11 <=== fTotalQVector | |
384 | // Q ?= Qa or Qb or QaQb | |
385 | vQm.Set(0.0,0.0); //start the loop fresh | |
386 | Double_t dMq=0; | |
387 | if( (fTotalQvector%2)>0 ) { // 01 or 11 | |
388 | vQm += vQa; | |
389 | dMq += dMa; | |
29195b69 | 390 | } |
02914cfa | 391 | if( fTotalQvector>1 ) { // 10 or 11 |
392 | vQm += vQb; | |
393 | dMq += dMb; | |
1b302085 | 394 | } |
1b302085 | 395 | |
02914cfa | 396 | //remove track if in subevent |
397 | for(Int_t inSubEvent=0; inSubEvent!=2; ++inSubEvent) { | |
398 | if( !pTrack->InSubevent( inSubEvent ) ) | |
399 | continue; | |
400 | if(inSubEvent==0) | |
401 | if( (fTotalQvector%2)!=1 ) | |
402 | continue; | |
403 | if(inSubEvent==1) | |
404 | if( fTotalQvector>0 ) | |
405 | continue; | |
406 | Double_t dW=1; | |
407 | Double_t dPhiCenter = dPhi; | |
408 | //if phi weights are used | |
409 | if(fUsePhiWeights && fPhiWeightsSub[inSubEvent]) { | |
410 | Int_t iNBinsPhiSub = fPhiWeightsSub[inSubEvent]->GetNbinsX(); | |
411 | Int_t phiBin = 1+(Int_t)(TMath::Floor(dPhi*iNBinsPhiSub/TMath::TwoPi())); | |
412 | dW = fPhiWeightsSub[inSubEvent]->GetBinContent(phiBin); | |
413 | dPhiCenter = fPhiWeightsSub[inSubEvent]->GetBinCenter(phiBin); | |
1b302085 | 414 | } |
02914cfa | 415 | Double_t dQmX = vQm.X() - dW*(pTrack->Weight())* TMath::Cos(fHarmonic*dPhiCenter); |
416 | Double_t dQmY = vQm.Y() - dW*(pTrack->Weight())* TMath::Sin(fHarmonic*dPhiCenter); | |
417 | vQm.Set(dQmX,dQmY); | |
418 | dMq = dMq-dW*pTrack->Weight(); | |
419 | } | |
420 | dNq = fNormalizationType ? dMq : vQm.Mod(); | |
421 | dWq = fNormalizationType ? dMq : 1; | |
422 | ||
423 | //Filling QA (for compatibility with previous version) | |
424 | fHistProQaQb->Fill(1,vQa*vQb,1); | |
425 | fHistProQaQbM->Fill(anEvent->GetNumberOfRPs()+0.5,(vQa*vQb)/dMa/dMb,dMa*dMb); | |
426 | fHistQaQbCos->Fill(TMath::ACos((vQa*vQb)/vQa.Mod()/vQb.Mod())); | |
427 | fResolution->Fill( TMath::Cos( vQa.Phi()-vQb.Phi() ) ); | |
428 | fHistQaQb->Fill(vQa*vQb); | |
429 | fHistMaMb->Fill(dMb,dMa); | |
430 | fHistProQNorm->Fill(1,vQm.Mod()/dMq,dMq); | |
431 | fHistQNormQaQbNorm->Fill((vQa*vQb)/dMa/dMb,vQm.Mod()/dMq); | |
432 | fHistQaNormMa->Fill(dMa,vQa.Mod()/dMa); | |
433 | fHistQbNormMb->Fill(dMb,vQb.Mod()/dMb); | |
434 | ||
435 | Double_t dUQ = vU*vQm; | |
436 | ||
437 | //fill the profile histograms | |
438 | for(Int_t iPOI=0; iPOI!=2; ++iPOI) { | |
439 | if( (iPOI==0)&&(!pTrack->InRPSelection()) ) | |
440 | continue; | |
441 | if( (iPOI==1)&&(!pTrack->InPOISelection()) ) | |
442 | continue; | |
443 | fHistProUQ[iPOI][0]->Fill(dPt ,dUQ/dNq,dWq); //Fill (uQ/Nq') with weight (Nq') | |
444 | fHistProUQ[iPOI][1]->Fill(dEta,dUQ/dNq,dWq); //Fill (uQ/Nq') with weight (Nq') | |
445 | //needed for the error calculation: | |
446 | fHistProUQQaQb[iPOI][0]-> Fill(dPt ,dUQ/dNq*dQaQb/dNa/dNb,dWq*dWa*dWb); //Fill [Qu/Nq']*[QaQb/NaNb] with weight (Nq')NaNb | |
447 | fHistProUQQaQb[iPOI][1]-> Fill(dEta,dUQ/dNq*dQaQb/dNa/dNb,dWq*dWa*dWb); //Fill [Qu/Nq']*[QaQb/NaNb] with weight (Nq')NaNb | |
448 | fHistSumOfWeightsu[iPOI][0][0]->Fill(dPt ,dWq); // sum of Nq' | |
449 | fHistSumOfWeightsu[iPOI][0][1]->Fill(dPt ,pow(dWq,2.));// sum of Nq'^2 | |
450 | fHistSumOfWeightsu[iPOI][0][2]->Fill(dPt ,dWq*dWa*dWb);// sum of Nq'*Na*Nb | |
451 | fHistSumOfWeightsu[iPOI][1][0]->Fill(dEta,dWq); // sum of Nq' | |
452 | fHistSumOfWeightsu[iPOI][1][1]->Fill(dEta,pow(dWq,2.));// sum of Nq'^2 | |
453 | fHistSumOfWeightsu[iPOI][1][2]->Fill(dEta,dNq*dWa*dWb);// sum of Nq'*Na*Nb | |
454 | //NUA: | |
455 | fHistProNUAu[iPOI][0][0]->Fill(dPt,dUY,1.); //sin u | |
456 | fHistProNUAu[iPOI][0][1]->Fill(dPt,dUX,1.); //cos u | |
457 | fHistProNUAu[iPOI][1][0]->Fill(dEta,dUY,1.); //sin u | |
458 | fHistProNUAu[iPOI][1][1]->Fill(dEta,dUX,1.); //cos u | |
459 | } | |
460 | }//loop over tracks | |
1b302085 | 461 | |
02914cfa | 462 | } |
8d312f00 | 463 | |
12dc476e | 464 | //-------------------------------------------------------------------- |
fd46c3dd | 465 | void AliFlowAnalysisWithScalarProduct::GetOutputHistograms(TList *outputListHistos){ |
fd46c3dd | 466 | //get pointers to all output histograms (called before Finish()) |
02914cfa | 467 | fHistList = outputListHistos; |
468 | ||
469 | fCommonHists = (AliFlowCommonHist*) fHistList->FindObject("AliFlowCommonHist_SP"); | |
470 | fCommonHistsuQ = (AliFlowCommonHist*) fHistList->FindObject("AliFlowCommonHist_uQ"); | |
471 | fCommonHistsRes = (AliFlowCommonHistResults*) fHistList->FindObject("AliFlowCommonHistResults_SP"); | |
472 | fHistProConfig = (TProfile*) fHistList->FindObject("FlowPro_Flags_SP"); | |
473 | if(!fHistProConfig) printf("Error loading fHistProConfig\n"); | |
474 | TList *uQ = (TList*) fHistList->FindObject("uQ"); | |
475 | TList *nua = (TList*) fHistList->FindObject("NUA"); | |
476 | TList *error = (TList*) fHistList->FindObject("error"); | |
477 | ||
478 | fHistProQaQbNorm = (TProfile*) error->FindObject("FlowPro_QaQbNorm_SP"); | |
479 | if(!fHistProQaQbNorm) printf("Error loading fHistProQaQbNorm\n"); | |
ad1b7114 | 480 | fHistProNUAq = (TProfile*) nua->FindObject("FlowPro_NUAq_SP"); |
02914cfa | 481 | if(!fHistProNUAq) printf("Error loading fHistProNUAq\n"); |
482 | fHistSumOfWeights = (TH1D*) error->FindObject("Flow_SumOfWeights_SP"); | |
483 | if(!fHistSumOfWeights) printf("Error loading fHistSumOfWeights\n"); | |
484 | ||
485 | TString sPOI[2] = {"RP","POI"}; | |
486 | TString sEta[2] = {"Pt","eta"}; | |
487 | TString sWeights[3] = {"uQ","uQuQ","uQQaQb"}; | |
488 | Int_t iNbins[2]; | |
489 | Double_t dMin[2], dMax[2]; | |
490 | ||
491 | iNbins[0] = AliFlowCommonConstants::GetMaster()->GetNbinsPt(); | |
492 | iNbins[1] = AliFlowCommonConstants::GetMaster()->GetNbinsEta(); | |
493 | dMin[0] = AliFlowCommonConstants::GetMaster()->GetPtMin(); | |
494 | dMin[1] = AliFlowCommonConstants::GetMaster()->GetEtaMin(); | |
495 | dMax[0] = AliFlowCommonConstants::GetMaster()->GetPtMax(); | |
496 | dMax[1] = AliFlowCommonConstants::GetMaster()->GetEtaMax(); | |
497 | for(Int_t iPOI=0; iPOI!=2; ++iPOI) for(Int_t iSpace=0; iSpace!=2; ++iSpace) { | |
498 | fHistProUQ[iPOI][iSpace] = (TProfile*) uQ->FindObject( Form( "FlowPro_UQ_%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ) ); | |
499 | if(!fHistProUQ[iPOI][iSpace]) printf("Error loading fHistProUQ[%d][%d]\n",iPOI,iSpace); | |
500 | fHistProNUAu[iPOI][iSpace][0] = (TProfile*) nua->FindObject( Form("FlowPro_NUAu_%s%s_IM_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ) ); | |
501 | if(!fHistProNUAu[iPOI][iSpace][0]) printf("Error loading fHistProNUAu[%d][%d][0]\n",iPOI,iSpace); | |
502 | fHistProNUAu[iPOI][iSpace][1] = (TProfile*) nua->FindObject( Form("FlowPro_NUAu_%s%s_RE_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ) ); | |
503 | if(!fHistProNUAu[iPOI][iSpace][1]) printf("Error loading fHistProNUAu[%d][%d][1]\n",iPOI,iSpace); | |
504 | fHistProUQQaQb[iPOI][iSpace] = (TProfile*) error->FindObject( Form("FlowPro_UQQaQb_%s%s_SP", sEta[iSpace].Data(), sPOI[iPOI].Data() ) ); | |
505 | for(Int_t i=0; i!=3; ++i){ | |
506 | fHistSumOfWeightsu[iPOI][iSpace][i] = (TH1D*) error->FindObject( Form("Flow_SumOfWeights_%s%s_%s_SP",sWeights[i].Data(),sPOI[iPOI].Data(),sEta[iSpace].Data()) ); | |
507 | if(!fHistSumOfWeightsu[iPOI][iSpace][i]) printf("Error loading fHistSumOfWeightsu[%d][%d][%d]\n",iPOI,iSpace,i); | |
508 | } | |
509 | } | |
0006c329 | 510 | if(fHistProConfig) { |
511 | fApplyCorrectionForNUA = (Int_t) fHistProConfig->GetBinContent(1); | |
512 | fNormalizationType = (Int_t) fHistProConfig->GetBinContent(2); | |
513 | fUsePhiWeights = (Int_t) fHistProConfig->GetBinContent(3); | |
514 | fHarmonic = (Int_t) fHistProConfig->GetBinContent(4); | |
515 | } | |
fd46c3dd | 516 | } |
517 | ||
518 | //-------------------------------------------------------------------- | |
8d312f00 | 519 | void AliFlowAnalysisWithScalarProduct::Finish() { |
489d5531 | 520 | //calculate flow and fill the AliFlowCommonHistResults |
02914cfa | 521 | printf("AliFlowAnalysisWithScalarProduct::Finish()\n"); |
29ecccee | 522 | |
50fe33aa | 523 | // access harmonic: |
9812922b | 524 | fApplyCorrectionForNUA = fHistProConfig->GetBinContent(1); |
525 | fNormalizationType = fHistProConfig->GetBinContent(2); | |
526 | fHarmonic = fHistProConfig->GetBinContent(4); | |
50fe33aa | 527 | |
02914cfa | 528 | printf("*************************************\n"); |
529 | printf("*************************************\n"); | |
530 | printf(" Integrated flow from \n"); | |
ad1b7114 | 531 | printf(" Scalar Product \n\n"); |
532 | if(!fNormalizationType) | |
533 | printf(" (BehaveAsEP) \n\n"); | |
e4cecfa0 | 534 | |
02914cfa | 535 | //Calculate reference flow |
12dc476e | 536 | //---------------------------------- |
02914cfa | 537 | //weighted average over (QaQb/NaNb) with weight (NaNb) |
538 | Double_t dEntriesQaQb = fHistProQaQbNorm->GetEntries(); | |
539 | if( dEntriesQaQb < 1 ) | |
540 | return; | |
a554b203 | 541 | Double_t dQaQb = fHistProQaQbNorm->GetBinContent(1); |
02914cfa | 542 | if( dQaQb < 0 ) |
543 | return; | |
a93de0f0 | 544 | Double_t dSpreadQaQb = fHistProQaQbNorm->GetBinError(1); |
02914cfa | 545 | |
546 | //NUA qq: | |
547 | Double_t dImQa = fHistProNUAq->GetBinContent(1); | |
548 | Double_t dReQa = fHistProNUAq->GetBinContent(2); | |
549 | Double_t dImQb = fHistProNUAq->GetBinContent(3); | |
550 | Double_t dReQb = fHistProNUAq->GetBinContent(4); | |
29ecccee | 551 | if(fApplyCorrectionForNUA) |
02914cfa | 552 | dQaQb = dQaQb - dImQa*dImQb - dReQa*dReQb; |
553 | printf("QaQb = %f +- %f\n", dQaQb, (dSpreadQaQb/TMath::Sqrt(dEntriesQaQb)) ); | |
554 | Double_t dV = TMath::Sqrt(dQaQb); | |
555 | ||
556 | printf("ResSub = %f\n", dV ); | |
557 | printf("fTotalQvector %d \n",fTotalQvector); | |
558 | if(!fNormalizationType) { | |
ad1b7114 | 559 | if(fTotalQvector>2) { |
9812922b | 560 | dV = ComputeResolution( TMath::Sqrt2()*FindXi(dV,1e-6) ); |
ad1b7114 | 561 | printf("An estimate of the event plane resolution is: %f\n", dV ); |
562 | } | |
524798bf | 563 | } |
02914cfa | 564 | printf("ResTot = %f\n", dV ); |
12dc476e | 565 | //statistical error of dQaQb: |
566 | // statistical error = term1 * spread * term2: | |
567 | // term1 = sqrt{sum_{i=1}^{N} w^2}/(sum_{i=1}^{N} w) | |
568 | // term2 = 1/sqrt(1-term1^2) | |
02914cfa | 569 | Double_t dSumOfLinearWeights = fHistSumOfWeights->GetBinContent(1); |
570 | Double_t dSumOfQuadraticWeights = fHistSumOfWeights->GetBinContent(2); | |
12dc476e | 571 | Double_t dTerm1 = 0.; |
572 | Double_t dTerm2 = 0.; | |
02914cfa | 573 | if(dSumOfLinearWeights) |
12dc476e | 574 | dTerm1 = pow(dSumOfQuadraticWeights,0.5)/dSumOfLinearWeights; |
02914cfa | 575 | if(1.-pow(dTerm1,2.) > 0.) |
12dc476e | 576 | dTerm2 = 1./pow(1-pow(dTerm1,2.),0.5); |
12dc476e | 577 | Double_t dStatErrorQaQb = dTerm1 * dSpreadQaQb * dTerm2; |
12dc476e | 578 | Double_t dVerr = 0.; |
02914cfa | 579 | if(dQaQb > 0.) |
12dc476e | 580 | dVerr = (1./(2.*pow(dQaQb,0.5)))*dStatErrorQaQb; |
02914cfa | 581 | fCommonHistsRes->FillIntegratedFlow(dV,dVerr); |
582 | printf("v%d(subevents) = %f +- %f\n",fHarmonic,dV,dVerr); | |
583 | ||
584 | Int_t iNbins[2]; | |
585 | iNbins[0] = AliFlowCommonConstants::GetMaster()->GetNbinsPt(); | |
586 | iNbins[1] = AliFlowCommonConstants::GetMaster()->GetNbinsEta(); | |
587 | ||
12dc476e | 588 | //Calculate differential flow and integrated flow (RP, POI) |
589 | //--------------------------------------------------------- | |
12dc476e | 590 | //v as a function of eta for RP selection |
02914cfa | 591 | for(Int_t iRFPorPOI=0; iRFPorPOI != 2; ++iRFPorPOI) |
592 | for(Int_t iPTorETA=0; iPTorETA != 2; ++iPTorETA) | |
593 | for(Int_t b=1; b != iNbins[iPTorETA]+1; ++b) { | |
594 | Double_t duQpro = fHistProUQ[iRFPorPOI][iPTorETA]->GetBinContent(b); | |
595 | if(fApplyCorrectionForNUA) | |
a93de0f0 | 596 | duQpro = duQpro |
02914cfa | 597 | - fHistProNUAu[iRFPorPOI][iPTorETA][1]->GetBinContent(b)*fHistProNUAq->GetBinContent(6) |
598 | - fHistProNUAu[iRFPorPOI][iPTorETA][0]->GetBinContent(b)*fHistProNUAq->GetBinContent(5); | |
a554b203 | 599 | Double_t dv2pro = -999.; |
02914cfa | 600 | if( TMath::Abs(dV!=0.) ) { dv2pro = duQpro/dV; } |
12dc476e | 601 | //calculate the statistical error |
02914cfa | 602 | Double_t dv2ProErr = CalculateStatisticalError(iRFPorPOI, iPTorETA, b, dStatErrorQaQb); |
603 | if( (iRFPorPOI==0)&&(iPTorETA==0) ) | |
604 | fCommonHistsRes->FillDifferentialFlowPtRP( b, dv2pro, dv2ProErr); | |
605 | if( (iRFPorPOI==0)&&(iPTorETA==1) ) | |
606 | fCommonHistsRes->FillDifferentialFlowEtaRP( b, dv2pro, dv2ProErr); | |
607 | if( (iRFPorPOI==1)&&(iPTorETA==0) ) | |
608 | fCommonHistsRes->FillDifferentialFlowPtPOI( b, dv2pro, dv2ProErr); | |
609 | if( (iRFPorPOI==1)&&(iPTorETA==1) ) | |
610 | fCommonHistsRes->FillDifferentialFlowEtaPOI(b, dv2pro, dv2ProErr); | |
611 | //printf("POI %d | PT %d >>> %f +- %f\n",iRFPorPOI,iPTorETA,dv2pro,dv2ProErr); | |
80f72ed6 | 612 | } |
12dc476e | 613 | |
02914cfa | 614 | printf("\n"); |
615 | printf("*************************************\n"); | |
616 | printf("*************************************\n"); | |
12dc476e | 617 | } |
8232a5ec | 618 | |
02914cfa | 619 | //----------------------------------------------------------------------- |
9812922b | 620 | void AliFlowAnalysisWithScalarProduct::WriteHistograms(TDirectoryFile *outputFileName) const { |
02914cfa | 621 | //store the final results in output .root file |
622 | outputFileName->Add(fHistList); | |
623 | outputFileName->Write(outputFileName->GetName(), TObject::kSingleKey); | |
624 | } | |
8232a5ec | 625 | |
12dc476e | 626 | //-------------------------------------------------------------------- |
9812922b | 627 | Double_t AliFlowAnalysisWithScalarProduct::CalculateStatisticalError(Int_t iRFPorPOI, Int_t iPTorETA, Int_t b, Double_t aStatErrorQaQb) const { |
12dc476e | 628 | //calculate the statistical error for differential flow for bin b |
02914cfa | 629 | Double_t duQproSpread = fHistProUQ[iRFPorPOI][iPTorETA]->GetBinError(b); |
630 | Double_t sumOfMq = fHistSumOfWeightsu[iRFPorPOI][iPTorETA][0]->GetBinContent(b); | |
631 | Double_t sumOfMqSquared = fHistSumOfWeightsu[iRFPorPOI][iPTorETA][1]->GetBinContent(b); | |
632 | Double_t dQaQb = fHistProQaQbNorm->GetBinContent(1); | |
633 | ||
ccad6081 | 634 | //non-isotropic terms: |
02914cfa | 635 | if(fApplyCorrectionForNUA) { |
636 | Double_t dImQa = fHistProNUAq->GetBinContent(1); // <<sin(phi_a)>> | |
637 | Double_t dReQa = fHistProNUAq->GetBinContent(2); // <<cos(phi_a)>> | |
638 | Double_t dImQb = fHistProNUAq->GetBinContent(3); // <<sin(phi_b)>> | |
639 | Double_t dReQb = fHistProNUAq->GetBinContent(4); // <<cos(phi_b)>> | |
640 | dQaQb = dQaQb - dImQa*dImQb - dReQa*dReQb; | |
641 | } | |
642 | ||
12dc476e | 643 | Double_t dTerm1 = 0.; |
644 | Double_t dTerm2 = 0.; | |
645 | if(sumOfMq) { | |
646 | dTerm1 = (pow(sumOfMqSquared,0.5)/sumOfMq); | |
647 | } | |
648 | if(1.-pow(dTerm1,2.)>0.) { | |
649 | dTerm2 = 1./pow(1.-pow(dTerm1,2.),0.5); | |
650 | } | |
651 | Double_t duQproErr = dTerm1*duQproSpread*dTerm2; | |
12dc476e | 652 | // covariances: |
02914cfa | 653 | Double_t dTerm1Cov = fHistSumOfWeightsu[iRFPorPOI][iPTorETA][2]->GetBinContent(b); |
654 | Double_t dTerm2Cov = fHistSumOfWeights->GetBinContent(1); | |
12dc476e | 655 | Double_t dTerm3Cov = sumOfMq; |
656 | Double_t dWeightedCovariance = 0.; | |
657 | if(dTerm2Cov*dTerm3Cov>0.) { | |
658 | Double_t dDenominator = 1.-dTerm1Cov/(dTerm2Cov*dTerm3Cov); | |
659 | Double_t dPrefactor = dTerm1Cov/(dTerm2Cov*dTerm3Cov); | |
660 | if(dDenominator!=0) { | |
02914cfa | 661 | Double_t dCovariance = ( fHistProUQQaQb[iRFPorPOI][iPTorETA]->GetBinContent(b)-dQaQb*fHistProUQ[iRFPorPOI][iPTorETA]->GetBinContent(b))/dDenominator; |
12dc476e | 662 | dWeightedCovariance = dCovariance*dPrefactor; |
663 | } | |
664 | } | |
665 | Double_t dv2ProErr = 0.; // final statitical error: | |
a554b203 | 666 | if(dQaQb>0.) { |
667 | Double_t dv2ProErrorSquared = (1./4.)*pow(dQaQb,-3.)* | |
02914cfa | 668 | (pow(fHistProUQ[iRFPorPOI][iPTorETA]->GetBinContent(b),2.)*pow(aStatErrorQaQb,2.) |
a554b203 | 669 | + 4.*pow(dQaQb,2.)*pow(duQproErr,2.) |
02914cfa | 670 | - 4.*dQaQb*fHistProUQ[iRFPorPOI][iPTorETA]->GetBinContent(b)*dWeightedCovariance); |
12dc476e | 671 | if(dv2ProErrorSquared>0.) dv2ProErr = pow(dv2ProErrorSquared,0.5); |
02914cfa | 672 | } |
12dc476e | 673 | return dv2ProErr; |
674 | } | |
29ecccee | 675 | |
9812922b | 676 | Double_t AliFlowAnalysisWithScalarProduct::ComputeResolution( Double_t x ) const { |
677 | // Computes resolution for Event Plane method | |
02914cfa | 678 | if(x > 51.3) { |
679 | printf("Warning: Estimation of total resolution might be WRONG. Please check!"); | |
680 | return 0.99981; | |
681 | } | |
682 | Double_t a = x*x/4; | |
683 | Double_t b = TMath::Exp(-a)*TMath::BesselI0(a)+TMath::Exp(-a)*TMath::BesselI1(a); | |
684 | return TMath::Sqrt(TMath::PiOver2())/2*x*b; | |
685 | } | |
29ecccee | 686 | |
9812922b | 687 | Double_t AliFlowAnalysisWithScalarProduct::FindXi( Double_t res, Double_t prec ) const { |
688 | // Computes x(res) for Event Plane method | |
02914cfa | 689 | if(res > 0.99981) { |
690 | printf("Warning: Resolution for subEvent is high. You reached the precision limit."); | |
691 | return 51.3; | |
692 | } | |
693 | int nSteps =0; | |
694 | Double_t xtmp=0, xmin=0, xmax=51.3, rtmp=0, delta=2*prec; | |
695 | while( delta > prec ) { | |
696 | xtmp = 0.5*(xmin+xmax); | |
9812922b | 697 | rtmp = ComputeResolution(xtmp); |
02914cfa | 698 | delta = TMath::Abs( res-rtmp ); |
699 | if(rtmp>res) xmax = xtmp; | |
700 | if(rtmp<res) xmin = xtmp; | |
701 | nSteps++; | |
702 | } | |
703 | return xtmp; | |
704 | } | |
29ecccee | 705 |