]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/FLOW/AliFlowCommon/AliFlowAnalysisWithCumulants.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
cleanup of GFC cumulants (some well deserved attention)
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowCommon / AliFlowAnalysisWithCumulants.cxx
CommitLineData
2188af53 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
d6130938 16/* $Id$ */
17
18/*************************************************
19 * Flow analysis with cumulants. In this class *
20 * cumulants are calculated by making use of the *
21 * formalism of generating functions proposed by *
22 * Ollitrault et al. *
23 * *
24 * Author: Ante Bilandzic *
25 * (abilandzic@gmail.com) *
26 *************************************************/
2188af53 27
f1d945a1 28#define AliFlowAnalysisWithCumulants_cxx
29
30#include "Riostream.h"
929098e4 31#include "TMath.h"
f1d945a1 32#include "TFile.h"
929098e4 33#include "TList.h"
aaebd73d 34#include "TProfile.h"
35#include "TProfile2D.h"
36#include "TProfile3D.h"
929098e4 37#include "TH1F.h"
38#include "TH1D.h"
d6130938 39#include "TMatrixD.h"
40#include "TVectorD.h"
929098e4 41
42#include "AliFlowCommonConstants.h"
43#include "AliFlowCommonHist.h"
44#include "AliFlowCommonHistResults.h"
f1d945a1 45#include "AliFlowEventSimple.h"
46#include "AliFlowTrackSimple.h"
47#include "AliFlowAnalysisWithCumulants.h"
929098e4 48#include "AliFlowVector.h"
f1d945a1 49
2188af53 50//================================================================================================================
f1d945a1 51
52ClassImp(AliFlowAnalysisWithCumulants)
53
d6130938 54AliFlowAnalysisWithCumulants::AliFlowAnalysisWithCumulants():
2188af53 55 fHistList(NULL),
d6130938 56 fHistListName(NULL),
57 fAnalysisSettings(NULL),
58 fCommonHists(NULL),
52021ae2 59 fCommonHistsResults2nd(NULL),
60 fCommonHistsResults4th(NULL),
61 fCommonHistsResults6th(NULL),
62 fCommonHistsResults8th(NULL),
d6130938 63 fnBinsPhi(0),
64 fPhiMin(0),
65 fPhiMax(0),
66 fPhiBinWidth(0),
67 fnBinsPt(0),
68 fPtMin(0),
69 fPtMax(0),
70 fPtBinWidth(0),
71 fnBinsEta(0),
72 fEtaMin(0),
73 fEtaMax(0),
74 fEtaBinWidth(0),
75 fHarmonic(2),
76 fMultiple(1),
77 fR0(2.2),
78 fWeightsList(NULL),
e5e75b58 79 fUsePhiWeights(kFALSE),
80 fUsePtWeights(kFALSE),
9dd53ff2 81 fUseEtaWeights(kFALSE),
d6130938 82 fUseParticleWeights(NULL),
83 fPhiWeights(NULL),
84 fPtWeights(NULL),
85 fEtaWeights(NULL),
86 fMultiplicityWeight(NULL),
87 fReferenceFlowList(NULL),
88 fReferenceFlowProfiles(NULL),
89 fReferenceFlowResults(NULL),
90 fReferenceFlowFlags(NULL),
91 fnBinsMult(10000),
92 fMinMult(0.),
93 fMaxMult(10000.),
94 fGEBE(NULL),
95 fReferenceFlowGenFun(NULL),
96 fQvectorComponents(NULL),
97 fAverageOfSquaredWeight(NULL),
98 fAvM(0.),
99 fnEvts(0),
100 fReferenceFlowCumulants(NULL),
101 fReferenceFlow(NULL),
102 fChi(NULL),
103 fDiffFlowList(NULL),
104 fDiffFlowProfiles(NULL),
105 fDiffFlowResults(NULL),
106 fDiffFlowFlags(NULL),
107 fTuningList(NULL),
108 fTuningProfiles(NULL),
109 fTuningResults(NULL),
110 fTuningFlags(NULL),
111 fTuneParameters(kFALSE)
52021ae2 112 {
d6130938 113 // Constructor.
114
115 // Base list to hold all output objects:
116 fHistList = new TList();
117 fHistListName = new TString("cobjGFC");
118 fHistList->SetName(fHistListName->Data());
119 fHistList->SetOwner(kTRUE);
120
121 // List to hold histograms with phi, pt and eta weights:
122 fWeightsList = new TList();
123 fWeightsList->SetName("Weights");
124 fWeightsList->SetOwner(kTRUE);
125 fHistList->Add(fWeightsList);
126
127 // Multiplicity weight:
128 fMultiplicityWeight = new TString("unit");
129
130 // Initialize all arrays:
131 this->InitializeArrays();
813a4157 132
d6130938 133} // end of AliFlowAnalysisWithCumulants::AliFlowAnalysisWithCumulants()
134
135//================================================================================================================
f1d945a1 136
2188af53 137AliFlowAnalysisWithCumulants::~AliFlowAnalysisWithCumulants()
138{
d6130938 139 // Desctructor.
140
2188af53 141 delete fHistList;
e5e75b58 142 delete fWeightsList;
d6130938 143
144} // end of AliFlowAnalysisWithCumulants::~AliFlowAnalysisWithCumulants()
729ec982 145
2188af53 146//================================================================================================================
147
e5e75b58 148void AliFlowAnalysisWithCumulants::Init()
2188af53 149{
d6130938 150 // Initialize and book all objects.
151
152 // a) Cross check if the user settings make sense before starting;
153 // b) Access all common constants;
154 // c) Book and fill weights histograms;
155 // d) Book and nest all lists in the base list fHistList;
156 // e) Book and fill profile holding analysis settings;
157 // f) Book common control and results histograms;
158 // g) Store flags for reference flow;
159 // h) Store flags for differential flow;
160 // i) Book all objects needed for tuning.
aaebd73d 161
489d5531 162 //save old value and prevent histograms from being added to directory
163 //to avoid name clashes in case multiple analaysis objects are used
164 //in an analysis
165 Bool_t oldHistAddStatus = TH1::AddDirectoryStatus();
166 TH1::AddDirectory(kFALSE);
167
d6130938 168 this->CrossCheckSettings();
169 this->AccessConstants();
170 this->BookAndFillWeightsHistograms();
171 this->BookAndNestAllLists();
172 this->BookProfileHoldingSettings();
173 this->BookCommonHistograms();
174 this->BookEverythingForReferenceFlow();
175 this->BookEverythingForDiffFlow();
176 this->StoreReferenceFlowFlags();
177 this->StoreDiffFlowFlags();
178 if(fTuneParameters){this->BookEverythingForTuning();}
52021ae2 179
d6130938 180 (fCommonHists->GetHarmonic())->Fill(0.5,fHarmonic); // to be improved (moved somewhere else?)
52021ae2 181
489d5531 182 TH1::AddDirectory(oldHistAddStatus);
d6130938 183
184} // end of void AliFlowAnalysisWithCumulants::Init()
2188af53 185
186//================================================================================================================
187
188void AliFlowAnalysisWithCumulants::Make(AliFlowEventSimple* anEvent)
189{
d6130938 190 // Running over data only in this method.
cb308e83 191
d6130938 192 // a) Check all pointers used in this method;
193 // b) If tuning enabled, fill generating functions for different values of tuning parameters;
194 // c) For default values of tuning parameters (r0 = 2.2 and cutoff at 10th order):
195 // c1) Fill common control histograms;
196 // c2) Fill generating function for reference flow;
197 // c3) Fill profile holding average values of various Q-vector components;
198 // c4) Fill generating function for differential flow.
e5e75b58 199
d6130938 200 this->CheckPointersUsedInMake();
201 if(fTuneParameters) {this->FillGeneratingFunctionsForDifferentTuningParameters(anEvent);}
202 Int_t nRP = anEvent->GetEventNSelTracksRP(); // number of RPs (i.e. number of particles used to determine the reaction plane)
203 if(nRP<10) {return;} // generating function formalism make sense only for nRPs >= 10 for default settings
204 fCommonHists->FillControlHistograms(anEvent);
205 this->FillGeneratingFunctionForReferenceFlow(anEvent);
206 this->FillQvectorComponents(anEvent);
207 this->FillGeneratingFunctionForDiffFlow(anEvent);
208
209} // end of void AliFlowAnalysisWithCumulants::Make()
210
211//================================================================================================================
212
213void AliFlowAnalysisWithCumulants::Finish()
214{
215 // Calculate the final results.
9dd53ff2 216
d6130938 217 // a) Check all pointers used in this method;
218 // b) Access settings for analysis with Generating Function Cumulants;
219 // c) From relevant common control histogram get average multiplicity of RPs and number of events;
220 // d) Calculate cumulants for reference flow;
221 // e) Calculate from isotropic cumulants reference flow;
222 // f) Calculate error for reference flow estimates;
223 // g) Store the final results for reference flow in common hist results;
224 // h) Print on the screen the final results for reference flow;
225 // i) Calculate cumulants for differential flow;
226 // j) Calculate differential flow for RPs/POIs vs pt/eta from cumulants;
227 // k) Calculate integrated flow of RPs and POIs;
228 // l) Print on the screen the final results for integrated flow of RPs and POIs;
229 // m) If tuning enabled, calculate results for different tuning parameters.
e5e75b58 230
d6130938 231 this->CheckPointersUsedInFinish();
232 this->AccessSettings();
233 this->GetAvMultAndNoOfEvts();
234 this->CalculateCumulantsForReferenceFlow();
235 this->CalculateReferenceFlow();
236 this->CalculateReferenceFlowError();
237 this->FillCommonHistResultsForReferenceFlow();
238 if(fPrintFinalResults[0]){this->PrintFinalResults("RF");}
239 this->CalculateCumulantsForDiffFlow("RP","pt");
240 this->CalculateCumulantsForDiffFlow("RP","eta");
241 this->CalculateCumulantsForDiffFlow("POI","pt");
242 this->CalculateCumulantsForDiffFlow("POI","eta");
243 this->CalculateDifferentialFlow("RP","pt");
244 this->CalculateDifferentialFlow("RP","eta");
245 this->CalculateDifferentialFlow("POI","pt");
246 this->CalculateDifferentialFlow("POI","eta");
247 this->CalculateDifferentialFlowErrors("RP","pt");
248 this->CalculateDifferentialFlowErrors("RP","eta");
249 this->CalculateDifferentialFlowErrors("POI","pt");
250 this->CalculateDifferentialFlowErrors("POI","eta");
251 this->FillCommonHistResultsForDifferentialFlow("RP");
252 this->FillCommonHistResultsForDifferentialFlow("POI");
253 this->CalculateIntegratedFlow("RP");
254 this->CalculateIntegratedFlow("POI");
255 if(fPrintFinalResults[1]){this->PrintFinalResults("RP");}
256 if(fPrintFinalResults[2]){this->PrintFinalResults("POI");}
257 if(fTuneParameters){this->FinalizeTuning();}
258
259} // end of void AliFlowAnalysisWithCumulants::Finish()
260
261//================================================================================================================
262
263void AliFlowAnalysisWithCumulants::FinalizeTuning()
264{
265 // Finalize results with tuned inerpolating parameters.
266
267 for(Int_t r=0;r<10;r++)
268 {
269 if(TMath::Abs(fTuningR0[r])<1.e-10) continue; // protection against division by r0 bellow
270 for(Int_t pq=0;pq<5;pq++)
271 {
272 Int_t pMax = fTuningGenFun[r][pq]->GetXaxis()->GetNbins();
273 Int_t qMax = fTuningGenFun[r][pq]->GetYaxis()->GetNbins();
274 // <G[p][q]>
275 TMatrixD dAvG(pMax,qMax);
276 dAvG.Zero();
277 Bool_t someAvGEntryIsNegative = kFALSE;
278 for(Int_t p=0;p<pMax;p++)
279 {
280 for(Int_t q=0;q<qMax;q++)
281 {
282 dAvG(p,q) = fTuningGenFun[r][pq]->GetBinContent(fTuningGenFun[r][pq]->GetBin(p+1,q+1));
283 if(dAvG(p,q)<0.)
284 {
285 someAvGEntryIsNegative = kTRUE;
286 cout<<endl;
287 cout<<" WARNING: "<<Form("<G[%d][%d]> is negative !!!! GFC results are meaningless for r0 = %f, pq = %i.",p,q,fTuningR0[r],pq)<<endl;
288 cout<<endl;
289 }
290 }
291 }
292 // C[p][q] (generating function for the cumulants)
293 TMatrixD dC(pMax,qMax);
294 dC.Zero();
295 if(fAvM>0. && !someAvGEntryIsNegative)
296 {
297 for(Int_t p=0;p<pMax;p++)
298 {
299 for(Int_t q=0;q<qMax;q++)
300 {
301 dC(p,q) = fAvM*(pow(dAvG(p,q),(1./fAvM))-1.);
302 }
303 }
304 }
305 // Averaging the generating function for cumulants over azimuth
306 // in order to eliminate detector effects.
307 // <C[p][q]> (Remark: here <> stands for average over azimuth):
308 TVectorD dAvC(pMax);
309 dAvC.Zero();
310 for(Int_t p=0;p<pMax;p++)
311 {
312 Double_t temp = 0.;
313 for(Int_t q=0;q<qMax;q++)
314 {
315 temp += 1.*dC(p,q);
316 }
317 dAvC[p] = temp/qMax;
318 }
319 // Finally, the isotropic cumulants for reference flow and reference flow itself:
320 TVectorD cumulant(pMax);
321 cumulant.Zero();
322 TVectorD flow(pMax);
323 flow.Zero();
324 if(pMax==2)
325 {
326 cumulant[0]=(1./(fTuningR0[r]*fTuningR0[r]))*(2.*dAvC[0]-(1./2.)*dAvC[1]);
327 cumulant[1]=(2./pow(fTuningR0[r],4.))*((-2.)*dAvC[0]+1.*dAvC[1]);
328 if(cumulant[0]>=0.) {flow[0] = pow(cumulant[0],1./2.);}
329 if(cumulant[1]<=0.) {flow[1] = pow(-1.*cumulant[1],1./4.);}
330 }
331 else if(pMax==3)
332 {
333 cumulant[0] = (1./(fTuningR0[r]*fTuningR0[r]))*(3.*dAvC[0]-(3./2.)*dAvC[1]+(1./3.)*dAvC[2]);
334 cumulant[1] = (2./pow(fTuningR0[r],4.))*((-5.)*dAvC[0]+4.*dAvC[1]-1.*dAvC[2]);
335 cumulant[2] = (6./pow(fTuningR0[r],6.))*(3.*dAvC[0]-3.*dAvC[1]+1.*dAvC[2]);
336 if(cumulant[0]>=0.) {flow[0] = pow(cumulant[0],1./2.);}
337 if(cumulant[1]<=0.) {flow[1] = pow(-1.*cumulant[1],1./4.);}
338 if(cumulant[2]>=0.) {flow[2] = pow((1./4.)*cumulant[2],1./6.);}
339 }
340 else if(pMax==4)
341 {
342 cumulant[0] = (1./(fTuningR0[r]*fTuningR0[r]))*(4.*dAvC[0]-3.*dAvC[1]+(4./3.)*dAvC[2]-(1./4.)*dAvC[3]);
343 cumulant[1] = (1./pow(fTuningR0[r],4.))*((-52./3.)*dAvC[0]+19.*dAvC[1]-(28./3.)*dAvC[2]+(11./6.)*dAvC[3]);
344 cumulant[2] = (3./pow(fTuningR0[r],6.))*(18.*dAvC[0]-24.*dAvC[1]+14.*dAvC[2]-3.*dAvC[3]);
345 cumulant[3] = (24./pow(fTuningR0[r],8.))*((-4.)*dAvC[0]+6.*dAvC[1]-4.*dAvC[2]+1.*dAvC[3]);
346 if(cumulant[0]>=0.) {flow[0] = pow(cumulant[0],1./2.);}
347 if(cumulant[1]<=0.) {flow[1] = pow(-1.*cumulant[1],1./4.);}
348 if(cumulant[2]>=0.) {flow[2] = pow((1./4.)*cumulant[2],1./6.);}
349 if(cumulant[3]<=0.) {flow[3] = pow((-1./33.)*cumulant[3],1./8.);}
350 }
351 else if(pMax==5)
352 {
353 cumulant[0] = (-1./(60*fTuningR0[r]*fTuningR0[r]))*((-300.)*dAvC[0]+300.*dAvC[1]-200.*dAvC[2]+75.*dAvC[3]-12.*dAvC[4]);
354 cumulant[1] = (-1./(6.*pow(fTuningR0[r],4.)))*(154.*dAvC[0]-214.*dAvC[1]+156.*dAvC[2]-61.*dAvC[3]+10.*dAvC[4]);
355 cumulant[2] = (3./(2.*pow(fTuningR0[r],6.)))*(71.*dAvC[0]-118.*dAvC[1]+98.*dAvC[2]-41.*dAvC[3]+7.*dAvC[4]);
356 cumulant[3] = (-24./pow(fTuningR0[r],8.))*(14.*dAvC[0]-26.*dAvC[1]+24.*dAvC[2]-11.*dAvC[3]+2.*dAvC[4]);
357 cumulant[4] = (120./pow(fTuningR0[r],10.))*(5.*dAvC[0]-10.*dAvC[1]+10.*dAvC[2]-5.*dAvC[3]+1.*dAvC[4]);
358 if(cumulant[0]>=0.) {flow[0] = pow(cumulant[0],1./2.);}
359 if(cumulant[1]<=0.) {flow[1] = pow(-1.*cumulant[1],1./4.);}
360 if(cumulant[2]>=0.) {flow[2] = pow((1./4.)*cumulant[2],1./6.);}
361 if(cumulant[3]<=0.) {flow[3] = pow((-1./33.)*cumulant[3],1./8.);}
362 if(cumulant[4]>=0.) {flow[4] = pow((1./456.)*cumulant[4],1./10.);}
363 }
364 else if(pMax==8)
365 {
366 cumulant[0] = (1./(fTuningR0[r]*fTuningR0[r]))*(8.*dAvC[0]-14.*dAvC[1]+(56./3.)*dAvC[2]-(35./2.)*dAvC[3]
367 + (56./5.)*dAvC[4]-(14./3.)*dAvC[5]+(8./7.)*dAvC[6]-(1./8.)*dAvC[7]);
368 cumulant[1] = (1./pow(fTuningR0[r],4.))*((-1924./35.)*dAvC[0]+(621./5.)*dAvC[1]-(8012./45.)*dAvC[2]
369 + (691./4.)*dAvC[3]-(564./5.)*dAvC[4]+(2143./45.)*dAvC[5]-(412./35.)*dAvC[6]+(363./280.)*dAvC[7]);
370 cumulant[2] = (1./pow(fTuningR0[r],6.))*(349.*dAvC[0]-(18353./20.)*dAvC[1]+(7173./5.)*dAvC[2]
371 - 1457.*dAvC[3]+(4891./5.)*dAvC[4]-(1683./4.)*dAvC[5]+(527./5.)*dAvC[6]-(469./40.)*dAvC[7]);
372 cumulant[3] = (1./pow(fTuningR0[r],8.))*((-10528./5.)*dAvC[0]+(30578./5.)*dAvC[1]-(51456./5.)*dAvC[2]
373 + 10993.*dAvC[3]-(38176./5.)*dAvC[4]+(16818./5.)*dAvC[5]-(4288./5.)*dAvC[6]+(967./10.)*dAvC[7]);
374 cumulant[4] = (1./pow(fTuningR0[r],10.))*(11500.*dAvC[0]-35800.*dAvC[1]+63900.*dAvC[2]-71600.*dAvC[3]
375 + 51620.*dAvC[4]-23400.*dAvC[5]+6100.*dAvC[6]-700.*dAvC[7]);
376 cumulant[5] = (1./pow(fTuningR0[r],12.))*(-52560.*dAvC[0]+172080.*dAvC[1]-321840.*dAvC[2]+376200.*dAvC[3]
377 - 281520.*dAvC[4]+131760.*dAvC[5]-35280.*dAvC[6]+4140.*dAvC[7]);
378 cumulant[6] = (1./pow(fTuningR0[r],14.))*(176400.*dAvC[0]-599760.*dAvC[1]+1164240.*dAvC[2]-1411200.*dAvC[3]
379 + 1093680.*dAvC[4]-529200.*dAvC[5]+146160.*dAvC[6]-17640.*dAvC[7]);
380 cumulant[7] = (1./pow(fTuningR0[r],16.))*(-322560*dAvC[0]+1128960.*dAvC[1]-2257920.*dAvC[2]+2822400.*dAvC[3]
381 - 2257920.*dAvC[4]+1128960.*dAvC[5]-322560.*dAvC[6]+40320.*dAvC[7]);
382 if(cumulant[0]>=0.) {flow[0] = pow(cumulant[0],1./2.);}
383 if(cumulant[1]<=0.) {flow[1] = pow(-1.*cumulant[1],1./4.);}
384 if(cumulant[2]>=0.) {flow[2] = pow((1./4.)*cumulant[2],1./6.);}
385 if(cumulant[3]<=0.) {flow[3] = pow((-1./33.)*cumulant[3],1./8.);}
386 if(cumulant[4]>=0.) {flow[4] = pow((1./456.)*cumulant[4],1./10.);}
387 if(cumulant[5]<=0.) {flow[5] = pow((-1./9460.)*cumulant[5],1./12.);}
388 if(cumulant[6]>=0.) {flow[6] = pow((1./274800.)*cumulant[6],1./14.);}
389 if(cumulant[7]<=0.) {flow[7] = pow((-1./10643745.)*cumulant[7],1./16.);}
390 }
391 // Store cumulants and reference flow:
392 for(Int_t co=0;co<pMax;co++) // cumulant order
393 {
394 fTuningCumulants[r][pq]->SetBinContent(co+1,cumulant[co]);
395 fTuningFlow[r][pq]->SetBinContent(co+1,flow[co]);
396 }
397 } // end of for(Int_t pq=0;pq<5;pq++)
398 } // end of for(Int_t r=0;r<10;r++)
399
400} // end of void AliFlowAnalysisWithCumulants::FinalizeTuning()
e5e75b58 401
d6130938 402//================================================================================================================
403
404void AliFlowAnalysisWithCumulants::FillGeneratingFunctionsForDifferentTuningParameters(AliFlowEventSimple *anEvent)
405{
406 // Fill generating function for reference flow evaluated for different tuning parameters.
407
408 Int_t pMax[5] = {2,3,4,5,8};
409 Int_t qMax[5] = {5,7,9,11,17};
410
411 // Particle variables and weights:
412 Double_t dPhi = 0.; // azimuthal angle in the laboratory frame
413 Double_t dPt = 0.; // transverse momentum
414 Double_t dEta = 0.; // pseudorapidity
2f5e6dea 415 Double_t wPhi = 1.; // phi weight
416 Double_t wPt = 1.; // pt weight
417 Double_t wEta = 1.; // eta weight
d6130938 418
419 Int_t nPrim = anEvent->NumberOfTracks(); // nPrim = total number of primary tracks, i.e. nPrim = nRP + nPOI + rest, where:
420 // nRP = # of particles used to determine the reaction plane;
421 // nPOI = # of particles of interest for a detailed flow analysis;
422 // rest = # of particles which are not niether RPs nor POIs.
423
424 Int_t nRP = anEvent->GetEventNSelTracksRP(); // nRP = # of particles used to determine the reaction plane;
425
426 Double_t tuningGenFunEBE[10][5][8][17] = {{{{0.}}}};
427 for(Int_t r=0;r<10;r++)
e5e75b58 428 {
d6130938 429 for(Int_t pq=0;pq<5;pq++)
e5e75b58 430 {
d6130938 431 for(Int_t p=0;p<pMax[pq];p++)
e5e75b58 432 {
d6130938 433 for(Int_t q=0;q<qMax[pq];q++)
434 {
435 tuningGenFunEBE[r][pq][p][q] = 1.;
436 }
437 }
e5e75b58 438 }
d6130938 439 }
440
441 // Looping over tracks:
442 for(Int_t i=0;i<nPrim;i++)
443 {
444 AliFlowTrackSimple *aftsTrack = anEvent->GetTrack(i);
445 if(aftsTrack && aftsTrack->InRPSelection())
e5e75b58 446 {
d6130938 447 // Access particle variables and weights:
448 dPhi = aftsTrack->Phi();
449 dPt = aftsTrack->Pt();
450 dEta = aftsTrack->Eta();
451 if(fUsePhiWeights && fnBinsPhi) // determine phi weight for this particle:
e5e75b58 452 {
d6130938 453 wPhi = fPhiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*fnBinsPhi/TMath::TwoPi())));
e5e75b58 454 }
d6130938 455 if(fUsePtWeights && fnBinsPt) // determine pt weight for this particle:
456 {
457 wPt = fPtWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-fPtMin)/fPtBinWidth)));
458 }
459 if(fUseEtaWeights && fEtaBinWidth) // determine eta weight for this particle:
460 {
461 wEta = fEtaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-fEtaMin)/fEtaBinWidth)));
462 }
463 // Fill the generating functions:
464 for(Int_t r=0;r<10;r++) // 10 different values for interpolating parameter r0
465 {
466 if(TMath::Abs(fTuningR0[r])<1.e-10) continue;
467 for(Int_t pq=0;pq<5;pq++) // 5 different values for set (pMax,qMax)
468 {
469 for(Int_t p=0;p<pMax[pq];p++)
470 {
471 for(Int_t q=0;q<qMax[pq];q++)
472 {
473 tuningGenFunEBE[r][pq][p][q] *= (1.+wPhi*wPt*wEta*(2.*fTuningR0[r]*sqrt(p+1.)/nRP)*cos(fHarmonic*dPhi-2.*q*TMath::Pi()/qMax[pq]));
474 } // end of for(Int_t q=0;q<qMax[pq];q++)
475 } // end of for(Int_t p=0;p<pMax[pq];p++)
476 } // end for(Int_t pq=0;pq<5;pq++) // 5 different values for set (pMax,qMax)
477 } // end of for(Int_t r=0;r<10;r++) // 10 different values for interpolating parameter r0
478 } // end of if(aftsTrack && aftsTrack->InRPSelection())
479 } // end of for(Int_t i=0;i<nPrim;i++)
aaebd73d 480
d6130938 481 // Store G[p][q]:
482 for(Int_t r=0;r<10;r++)
483 {
484 for(Int_t pq=0;pq<5;pq++)
485 {
486 for(Int_t p=0;p<pMax[pq];p++)
487 {
488 for(Int_t q=0;q<qMax[pq];q++)
489 {
490 if(fTuningGenFun[r][pq]) {fTuningGenFun[r][pq]->Fill((Double_t)p,(Double_t)q,tuningGenFunEBE[r][pq][p][q],1.);}
491 }
492 }
493 }
494 }
f1d945a1 495
d6130938 496} // end of void AliFlowAnalysisWithCumulants::FillGeneratingFunctionsForDifferentTuningParameters(AliFlowEventSimple *anEvent)
497
498//================================================================================================================
499
500void AliFlowAnalysisWithCumulants::FillGeneratingFunctionForReferenceFlow(AliFlowEventSimple *anEvent)
501{
502 // Fill generating function for reference flow for current event.
503
504 // Particle variables and weights:
505 Double_t dPhi = 0.; // azimuthal angle in the laboratory frame
506 Double_t dPt = 0.; // transverse momentum
507 Double_t dEta = 0.; // pseudorapidity
508 Double_t wPhi = 1.; // phi weight
509 Double_t wPt = 1.; // pt weight
510 Double_t wEta = 1.; // eta weight
511
512 Int_t nPrim = anEvent->NumberOfTracks(); // nPrim = total number of primary tracks, i.e. nPrim = nRP + nPOI + rest, where:
513 // nRP = # of particles used to determine the reaction plane;
514 // nPOI = # of particles of interest for a detailed flow analysis;
515 // rest = # of particles which are not niether RPs nor POIs.
b6cd16a9 516
d6130938 517 Int_t nRP = anEvent->GetEventNSelTracksRP(); // nRP = # of particles used to determine the reaction plane;
518
519 // Initializing the generating function G[p][q] for reference flow for current event:
520 Int_t pMax = fGEBE->GetNrows();
521 Int_t qMax = fGEBE->GetNcols();
522 for(Int_t p=0;p<pMax;p++)
17100608 523 {
d6130938 524 for(Int_t q=0;q<qMax;q++)
b6cd16a9 525 {
d6130938 526 (*fGEBE)(p,q) = 1.;
17100608 527 }
528 }
d6130938 529
530 // Cross-checking the number of RPs in current event:
531 Int_t crossCheckRP = 0;
b6cd16a9 532
d6130938 533 // Looping over tracks:
b6cd16a9 534 for(Int_t i=0;i<nPrim;i++)
535 {
d6130938 536 AliFlowTrackSimple *aftsTrack = anEvent->GetTrack(i);
537 if(aftsTrack && aftsTrack->InRPSelection())
b6cd16a9 538 {
d6130938 539 crossCheckRP++;
540 // Access particle variables and weights:
541 dPhi = aftsTrack->Phi();
542 dPt = aftsTrack->Pt();
543 dEta = aftsTrack->Eta();
544 if(fUsePhiWeights && fnBinsPhi) // determine phi weight for this particle:
545 {
546 wPhi = fPhiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*fnBinsPhi/TMath::TwoPi())));
2f5e6dea 547 }
d6130938 548 if(fUsePtWeights && fnBinsPt) // determine pt weight for this particle:
549 {
550 wPt = fPtWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-fPtMin)/fPtBinWidth)));
551 }
552 if(fUseEtaWeights && fEtaBinWidth) // determine eta weight for this particle:
553 {
554 wEta = fEtaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-fEtaMin)/fEtaBinWidth)));
555 }
556 // Fill the generating function:
557 for(Int_t p=0;p<pMax;p++)
b6cd16a9 558 {
d6130938 559 for(Int_t q=0;q<qMax;q++)
b6cd16a9 560 {
d6130938 561 (*fGEBE)(p,q) *= (1.+wPhi*wPt*wEta*(2.*fR0*sqrt(p+1.)/nRP)*cos(fHarmonic*dPhi-2.*q*TMath::Pi()/qMax));
f1d945a1 562 }
9dd53ff2 563 }
d6130938 564 // Fill the profile to calculate <<w^2>>:
565 fAverageOfSquaredWeight->Fill(0.5,pow(wPhi*wPt*wEta,2.),1.);
566 } // end of if(aftsTrack && aftsTrack->InRPSelection())
2f5e6dea 567 } // end of for(Int_t i=0;i<nPrim;i++)
b6cd16a9 568
d6130938 569 // Cross check # of RPs:
570 if(anEvent && (crossCheckRP != anEvent->GetEventNSelTracksRP()))
571 {
572 cout<<endl;
573 cout<<"WARNING (GFC): crossCheckRP != nRP in GFC::Make(). Something is wrong with RP flagging !!!!"<<endl;
574 cout<<endl;
575 exit(0);
576 }
577
578 // Storing the value of G[p][q] in 2D profile in order to get eventually the avarage <G[p][q]>:
579 // Determine first the event weight for G[p][q]:
580 // (to be improved - this can be implemented much better, this shall be executed only once out of Make(), eventWeight should be a data member)
581 Double_t eventWeight = 0.;
582 if(!strcmp(fMultiplicityWeight->Data(),"unit"))
583 {
584 eventWeight = 1.;
585 } else if(!strcmp(fMultiplicityWeight->Data(),"multiplicity"))
586 {
587 eventWeight = anEvent->GetEventNSelTracksRP();
588 }
589 // Store G[p][q] weighted appropriately:
590 for(Int_t p=0;p<pMax;p++)
b6cd16a9 591 {
d6130938 592 for(Int_t q=0;q<qMax;q++)
b6cd16a9 593 {
d6130938 594 fReferenceFlowGenFun->Fill((Double_t)p,(Double_t)q,(*fGEBE)(p,q),eventWeight);
b6cd16a9 595 }
596 }
2188af53 597
d6130938 598} // end of void AliFlowAnalysisWithCumulants::FillGeneratingFunctionForReferenceFlow(AliFlowEventSimple* anEvent)
599
600//================================================================================================================
601
602void AliFlowAnalysisWithCumulants::FillQvectorComponents(AliFlowEventSimple* anEvent)
603{
604 // Fill components of Q-vector for current event (needed for error calculation).
605
606 // Remark: Components are stored in profile fQvectorComponents whose binning is organized as follows:
607 // 1st bin: Q_x
608 // 2nd bin: Q_y
609 // 3rd bin: (Q_x)^2
610 // 4th bin: (Q_y)^2
611
612 AliFlowVector afv;
613 afv.Set(0.,0.);
614 afv.SetMult(0);
615
616 Int_t n = 2; // to be removed
617
618 if(anEvent)
2188af53 619 {
d6130938 620 afv = anEvent->GetQ(1*n,fWeightsList,fUsePhiWeights,fUsePtWeights,fUseEtaWeights); // get the Q-vector for this event
621 fQvectorComponents->Fill(0.5,afv.X(),1.); // in the 1st bin fill Q_x
622 fQvectorComponents->Fill(1.5,afv.Y(),1.); // in the 2nd bin fill Q_y
623 fQvectorComponents->Fill(2.5,pow(afv.X(),2.),1.); // in the 3rd bin fill (Q_x)^2
624 fQvectorComponents->Fill(3.5,pow(afv.Y(),2.),1.); // in the 4th bin fill (Q_y)^2
625 }
626
627} // end of void AliFlowAnalysisWithCumulants::FillQvectorComponents(AliFlowEventSimple* anEvent)
628
629//================================================================================================================
630
631void AliFlowAnalysisWithCumulants::FillGeneratingFunctionForDiffFlow(AliFlowEventSimple* anEvent)
632{
633 // Fill generating function for differential flow for the current event.
634
635 // Remark 0: Generating function D[b][p][q] is a complex number => real and imaginary part are calculated separately
636 // (b denotes pt or eta bin);
637 // Remark 1: Note that bellow G[p][q] is needed, the value of generating function for reference flow for the CURRENT event.
638 // This values is obtained in method FillGeneratingFunctionForReferenceFlow() as TMatrixD fGEBE;
639 // Remark 2: Results for D[b][p][q] are stored in 3D profiles fDiffFlowGenFun[0=Re,1=Im][0=RP,1=POI][0=pt,1=eta] in order to
640 // automatically get <Re(D[b][p][q])> and <Im(D[b][p][q])> at the end of the day.
641
642 // Particle variables and weights:
643 Double_t dPhi = 0.; // azimuthal angle in the laboratory frame
644 Double_t dPt = 0.; // transverse momentum
645 Double_t dEta = 0.; // pseudorapidity
646 Double_t wPhi = 1.; // phi weight
647 Double_t wPt = 1.; // pt weight
648 Double_t wEta = 1.; // eta weight
649
650 // pMax and qMax:
651 Int_t pMax = fGEBE->GetNrows();
652 Int_t qMax = fGEBE->GetNcols();
653
654 Int_t nPrim = anEvent->NumberOfTracks(); // nPrim = total number of primary tracks, i.e. nPrim = nRP + nPOI + rest, where:
655 // nRP = # of particles used to determine the reaction plane;
656 // nPOI = # of particles of interest for a detailed flow analysis;
657 // rest = # of particles which are not niether RPs nor POIs.
2188af53 658
d6130938 659 Int_t nRP = anEvent->GetEventNSelTracksRP(); // nRP = # of particles used to determine the reaction plane
660
661 // Start the second loop over event in order to evaluate the generating function D[b][p][q] for differential flow:
b6cd16a9 662 for(Int_t i=0;i<nPrim;i++)
663 {
d6130938 664 AliFlowTrackSimple *aftsTrack = anEvent->GetTrack(i);
665 if(aftsTrack)
b6cd16a9 666 {
d6130938 667 if(!(aftsTrack->InRPSelection() || aftsTrack->InPOISelection())) continue;
668 // Differential flow of POIs:
669 if(aftsTrack->InPOISelection())
b6cd16a9 670 {
d6130938 671 // Get azimuthal angle, momentum and pseudorapidity of a particle:
672 dPhi = aftsTrack->Phi();
673 dPt = aftsTrack->Pt();
674 dEta = aftsTrack->Eta();
675 Double_t ptEta[2] = {dPt,dEta};
b9639e20 676
d6130938 677 // Count number of POIs in pt/eta bin:
678 for(Int_t pe=0;pe<2;pe++)
679 {
680 fNoOfParticlesInBin[1][pe]->Fill(ptEta[pe],ptEta[pe],1.);
f1d945a1 681 }
d6130938 682
683 if(!(aftsTrack->InRPSelection())) // particle was flagged only as POI
b6cd16a9 684 {
d6130938 685 // Fill generating function:
686 for(Int_t p=0;p<pMax;p++)
b6cd16a9 687 {
d6130938 688 for(Int_t q=0;q<qMax;q++)
689 {
690 for(Int_t ri=0;ri<2;ri++)
691 {
692 for(Int_t pe=0;pe<2;pe++)
693 {
694 if(ri==0) // Real part (to be improved - this can be implemented better)
695 {
696 fDiffFlowGenFun[ri][1][pe]->Fill(ptEta[pe],(Double_t)p,(Double_t)q, // to be improved - hardwired weight 1. in the line bellow
697 (*fGEBE)(p,q)*cos(fMultiple*fHarmonic*dPhi),1.);
698 }
699 else if(ri==1) // Imaginary part (to be improved - this can be implemented better)
700 {
701 fDiffFlowGenFun[ri][1][pe]->Fill(ptEta[pe],(Double_t)p,(Double_t)q, // to be improved - hardwired weight 1. in the line bellow
702 (*fGEBE)(p,q)*sin(fMultiple*fHarmonic*dPhi),1.);
703 }
704 } // end of for(Int_t pe=0;pe<2;pe++)
705 } // end of for(Int_t ri=0;ri<2;ri++)
706 } // end of for(Int_t q=0;q<qMax;q++)
707 } // end of for(Int_t p=0;p<pMax;p++)
708 } // end of if(!(aftsTrack->InRPSelection())) // particle was flagged only as POI
709 else if(aftsTrack->InRPSelection()) // particle was flagged both as RP and POI
b9639e20 710 {
d6130938 711 // If particle weights were used, get them:
712 if(fUsePhiWeights && fnBinsPhi) // determine phi weight for this particle:
b9639e20 713 {
d6130938 714 wPhi = fPhiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*fnBinsPhi/TMath::TwoPi())));
b9639e20 715 }
d6130938 716 if(fUsePtWeights && fnBinsPt) // determine pt weight for this particle:
b9639e20 717 {
d6130938 718 wPt = fPtWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-fPtMin)/fPtBinWidth)));
719 }
720 if(fUseEtaWeights && fEtaBinWidth) // determine eta weight for this particle:
b9639e20 721 {
d6130938 722 wEta = fEtaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-fEtaMin)/fEtaBinWidth)));
723 }
724 // Fill generating function:
725 for(Int_t p=0;p<pMax;p++)
726 {
727 for(Int_t q=0;q<qMax;q++)
728 {
729 for(Int_t ri=0;ri<2;ri++)
730 {
731 for(Int_t pe=0;pe<2;pe++)
732 {
733 if(ri==0) // Real part (to be improved - this can be implemented better)
734 {
735 fDiffFlowGenFun[ri][1][pe]->Fill(ptEta[pe],(Double_t)p,(Double_t)q, // to be improved - hardwired weight 1. in the line bellow
736 (*fGEBE)(p,q)*cos(fMultiple*fHarmonic*dPhi)/(1.+wPhi*wPt*wEta*(2.*fR0*sqrt(p+1.)/nRP)*cos(fHarmonic*dPhi-2.*q*TMath::Pi()/qMax)),1.);
737 }
738 else if(ri==1) // Imaginary part (to be improved - this can be implemented better)
739 {
740 fDiffFlowGenFun[ri][1][pe]->Fill(ptEta[pe],(Double_t)p,(Double_t)q, // to be improved - hardwired weight 1. in the line bellow
741 (*fGEBE)(p,q)*sin(fMultiple*fHarmonic*dPhi)/(1.+wPhi*wPt*wEta*(2.*fR0*sqrt(p+1.)/nRP)*cos(fHarmonic*dPhi-2.*q*TMath::Pi()/qMax)),1.);
742 }
743 } // end of for(Int_t pe=0;pe<2;pe++)
744 } // end of for(Int_t ri=0;ri<2;ri++)
745 } // end of for(Int_t q=0;q<qMax;q++)
746 } // end of for(Int_t p=0;p<pMax;p++)
747 } // end of else if (aftsTrack->InRPSelection()) // particle was flagged both as RP and POI
748 } // end of if(aftsTrack->InPOISelection())
749 // Differential flow of RPs:
750 if(aftsTrack->InRPSelection())
751 {
752 // Get azimuthal angle, momentum and pseudorapidity of a particle:
753 dPhi = aftsTrack->Phi();
754 dPt = aftsTrack->Pt();
755 dEta = aftsTrack->Eta();
756 Double_t ptEta[2] = {dPt,dEta};
757
758 // Count number of RPs in pt/eta bin:
759 for(Int_t pe=0;pe<2;pe++)
760 {
761 fNoOfParticlesInBin[0][pe]->Fill(ptEta[pe],ptEta[pe],1.);
b9639e20 762 }
763
d6130938 764 // If particle weights were used, get them:
765 if(fUsePhiWeights && fnBinsPhi) // determine phi weight for this particle:
813a4157 766 {
d6130938 767 wPhi = fPhiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*fnBinsPhi/TMath::TwoPi())));
813a4157 768 }
d6130938 769 if(fUsePtWeights && fnBinsPt) // determine pt weight for this particle:
b6cd16a9 770 {
d6130938 771 wPt = fPtWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-fPtMin)/fPtBinWidth)));
772 }
773 if(fUseEtaWeights && fEtaBinWidth) // determine eta weight for this particle:
774 {
775 wEta = fEtaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-fEtaMin)/fEtaBinWidth)));
776 }
777 // Fill generating function:
778 for(Int_t p=0;p<pMax;p++)
b6cd16a9 779 {
d6130938 780 for(Int_t q=0;q<qMax;q++)
b6cd16a9 781 {
d6130938 782 for(Int_t ri=0;ri<2;ri++)
cb308e83 783 {
d6130938 784 for(Int_t pe=0;pe<2;pe++)
b6cd16a9 785 {
d6130938 786 if(ri==0) // Real part (to be improved - this can be implemented better)
cb308e83 787 {
d6130938 788 fDiffFlowGenFun[ri][0][pe]->Fill(ptEta[pe],(Double_t)p,(Double_t)q, // to be improved - hardwired weight 1. in the line bellow
789 (*fGEBE)(p,q)*cos(fMultiple*fHarmonic*dPhi)/(1.+wPhi*wPt*wEta*(2.*fR0*sqrt(p+1.)/nRP)*cos(fHarmonic*dPhi-2.*q*TMath::Pi()/qMax)),1.);
790 }
791 else if(ri==1) // Imaginary part (to be improved - this can be implemented better)
b6cd16a9 792 {
d6130938 793 fDiffFlowGenFun[ri][0][pe]->Fill(ptEta[pe],(Double_t)p,(Double_t)q, // to be improved - hardwired weight 1. in the line bellow
794 (*fGEBE)(p,q)*sin(fMultiple*fHarmonic*dPhi)/(1.+wPhi*wPt*wEta*(2.*fR0*sqrt(p+1.)/nRP)*cos(fHarmonic*dPhi-2.*q*TMath::Pi()/qMax)),1.);
b6cd16a9 795 }
d6130938 796 } // end of for(Int_t pe=0;pe<2;pe++)
797 } // end of for(Int_t ri=0;ri<2;ri++)
798 } // end of for(Int_t q=0;q<qMax;q++)
799 } // end of for(Int_t p=0;p<pMax;p++)
800 } // end of if(aftsTrack->InRPSelection())
801 } // end of if(aftsTrack)
802 } // end of for(Int_t i=0;i<nPrim;i++)
803
804} // end of void AliFlowAnalysisWithCumulants::FillGeneratingFunctionForDiffFlow(AliFlowEventSimple* anEvent)
805
806//================================================================================================================
807
808void AliFlowAnalysisWithCumulants::GetOutputHistograms(TList *outputListHistos)
809{
810 // Get pointers to all objects saved in the output file.
811
812 if(outputListHistos)
813 {
814 this->SetHistList(outputListHistos);
815 if(!fHistList)
816 {
817 cout<<endl;
818 cout<<" WARNING (GFC): fHistList is NULL in AFAWGFC::GOH() !!!!"<<endl;
819 cout<<endl;
820 exit(0);
821 }
822 this->GetPointersForBaseHistograms();
823 this->GetPointersForCommonControlHistograms();
824 this->GetPointersForCommonResultsHistograms();
825 this->GetPointersForParticleWeightsHistograms();
826 this->GetPointersForReferenceFlowObjects();
827 this->GetPointersForDiffFlowObjects();
828 } else
829 {
830 cout<<endl;
831 cout<<" WARNING (GFC): outputListHistos is NULL in AFAWGFC::GOH() !!!!"<<endl;
832 cout<<endl;
833 exit(0);
834 }
835
836} // end of void AliFlowAnalysisWithCumulants::GetOutputHistograms(TList *outputListHistos)
837
838//================================================================================================================
839
840void AliFlowAnalysisWithCumulants::GetPointersForBaseHistograms()
841{
842 // Get pointers to base histograms.
843
844 TString analysisSettingsName = "fAnalysisSettings";
845 TProfile *analysisSettings = dynamic_cast<TProfile*>(fHistList->FindObject(analysisSettingsName.Data()));
846 if(analysisSettings)
847 {
848 this->SetAnalysisSettings(analysisSettings);
849 } else
850 {
851 cout<<endl;
852 cout<<" WARNING (GFC): analysisSettings is NULL in AFAWGFC::GPFBH() !!!!"<<endl;
853 cout<<endl;
854 exit(0);
855 }
856
857} // end of void AliFlowAnalysisWithCumulants::GetPointersForBaseHistograms()
858
859//================================================================================================================
860
861void AliFlowAnalysisWithCumulants::GetPointersForCommonControlHistograms()
862{
863 // Get pointers for common control histograms.
864
865 TString commonHistsName = "AliFlowCommonHistGFC";
866 AliFlowCommonHist *commonHist = dynamic_cast<AliFlowCommonHist*>(fHistList->FindObject(commonHistsName.Data()));
867 if(commonHist)
868 {
869 this->SetCommonHists(commonHist);
870 } else
871 {
872 cout<<endl;
873 cout<<" WARNING (GFC): commonHist is NULL in AFAWGFC::GPFCH() !!!!"<<endl;
874 cout<<endl;
875 exit(0);
876 }
877
878} // end of void AliFlowAnalysisWithCumulants::GetPointersForCommonControlHistograms()
879
880//================================================================================================================
881
882void AliFlowAnalysisWithCumulants::GetPointersForCommonResultsHistograms()
883{
884 // Get pointers for common results histograms.
885
886 TString commonHistResults2ndOrderName = "AliFlowCommonHistResults2ndOrderGFC";
887 AliFlowCommonHistResults *commonHistRes2nd = dynamic_cast<AliFlowCommonHistResults*>
888 (fHistList->FindObject(commonHistResults2ndOrderName.Data()));
889 if(commonHistRes2nd)
890 {
891 this->SetCommonHistsResults2nd(commonHistRes2nd);
892 } else
893 {
894 cout<<endl;
895 cout<<" WARNING (GFC): commonHistRes2nd is NULL in AFAWGFC::GPFCRH() !!!!"<<endl;
896 cout<<endl;
897 exit(0);
898 }
899 TString commonHistResults4thOrderName = "AliFlowCommonHistResults4thOrderGFC";
900 AliFlowCommonHistResults *commonHistRes4th = dynamic_cast<AliFlowCommonHistResults*>
901 (fHistList->FindObject(commonHistResults4thOrderName.Data()));
902 if(commonHistRes4th)
903 {
904 this->SetCommonHistsResults4th(commonHistRes4th);
905 } else
906 {
907 cout<<endl;
908 cout<<" WARNING (GFC): commonHistRes4th is NULL in AFAWGFC::GPFCRH() !!!!"<<endl;
909 cout<<endl;
910 exit(0);
911 }
912 TString commonHistResults6thOrderName = "AliFlowCommonHistResults6thOrderGFC";
913 AliFlowCommonHistResults *commonHistRes6th = dynamic_cast<AliFlowCommonHistResults*>
914 (fHistList->FindObject(commonHistResults6thOrderName.Data()));
915 if(commonHistRes6th)
916 {
917 this->SetCommonHistsResults6th(commonHistRes6th);
918 } else
919 {
920 cout<<endl;
921 cout<<" WARNING (GFC): commonHistRes6th is NULL in AFAWGFC::GPFCRH() !!!!"<<endl;
922 cout<<endl;
923 exit(0);
924 }
925 TString commonHistResults8thOrderName = "AliFlowCommonHistResults8thOrderGFC";
926 AliFlowCommonHistResults *commonHistRes8th = dynamic_cast<AliFlowCommonHistResults*>
927 (fHistList->FindObject(commonHistResults8thOrderName.Data()));
928 if(commonHistRes8th)
929 {
930 this->SetCommonHistsResults8th(commonHistRes8th);
931 } else
932 {
933 cout<<endl;
934 cout<<" WARNING (GFC): commonHistRes8th is NULL in AFAWGFC::GPFCRH() !!!!"<<endl;
935 cout<<endl;
936 exit(0);
937 }
938
939} // end of void AliFlowAnalysisWithCumulants::GetPointersForCommonResultsHistograms()
940
941//================================================================================================================
942
943void AliFlowAnalysisWithCumulants::GetPointersForParticleWeightsHistograms()
944{
945 // Get pointers for histograms holding particle weights.
946
947 TList *weightsList = dynamic_cast<TList*>(fHistList->FindObject("Weights"));
948 if(weightsList) this->SetWeightsList(weightsList);
949 TString fUseParticleWeightsName = "fUseParticleWeights";
950 TProfile *useParticleWeights = dynamic_cast<TProfile*>(weightsList->FindObject(fUseParticleWeightsName.Data()));
951 if(useParticleWeights)
952 {
953 this->SetUseParticleWeights(useParticleWeights);
954 fUsePhiWeights = (Int_t)fUseParticleWeights->GetBinContent(1);
955 fUsePtWeights = (Int_t)fUseParticleWeights->GetBinContent(2);
956 fUseEtaWeights = (Int_t)fUseParticleWeights->GetBinContent(3);
957 }
958
959} // end of void AliFlowAnalysisWithCumulants::GetPointersForParticleWeightsHistograms();
960
961//================================================================================================================
962
963void AliFlowAnalysisWithCumulants::GetPointersForReferenceFlowObjects()
964{
965 // Get pointers for all objects relevant for calculation of reference flow.
966
967 // a) Get pointers to all lists relevant for reference flow;
968 // b) Get pointer to profile holding flags;
969 // c) Get pointers to all objects in the list fReferenceFlowProfiles;
970 // d) Get pointers to all objects in the list fReferenceFlowResults.
971
972 // a) Get pointers to all lists relevant for reference flow:
973 TList *referenceFlowList = dynamic_cast<TList*>(fHistList->FindObject("Reference Flow"));
974 if(!referenceFlowList)
975 {
976 cout<<endl;
977 cout<<"WARNING (GFC): referenceFlowList is NULL in AFAWGFC::GPFRFO() !!!!"<<endl;
978 cout<<endl;
979 exit(0);
980 }
981 TList *referenceFlowProfiles = dynamic_cast<TList*>(referenceFlowList->FindObject("Profiles"));
982 if(!referenceFlowProfiles)
983 {
984 cout<<endl;
985 cout<<"WARNING (GFC): referenceFlowProfiles is NULL in AFAWGFC::GPFRFO() !!!!"<<endl;
986 cout<<endl;
987 exit(0);
988 }
989 TList *referenceFlowResults = dynamic_cast<TList*>(referenceFlowList->FindObject("Results"));
990 if(!referenceFlowResults)
991 {
992 cout<<endl;
993 cout<<"WARNING (GFC): referenceFlowResults is NULL in AFAWGFC::GPFRFO() !!!!"<<endl;
994 cout<<endl;
995 exit(0);
996 }
997
998 // b) Get pointer to profile holding flags:
999 TString referenceFlowFlagsName = "fReferenceFlowFlags";
1000 TProfile *referenceFlowFlags = dynamic_cast<TProfile*>(referenceFlowList->FindObject(referenceFlowFlagsName.Data()));
1001 if(referenceFlowFlags)
1002 {
1003 this->SetReferenceFlowFlags(referenceFlowFlags);
1004 } else
1005 {
1006 cout<<endl;
1007 cout<<"WARNING (GFC): referenceFlowFlags is NULL in AFAWGFC::GPFRFO() !!!!"<<endl;
1008 cout<<endl;
1009 exit(0);
1010 }
1011
1012 // c) Get pointers to all objects in the list fReferenceFlowProfiles:
1013 TString referenceFlowGenFunName = "fReferenceFlowGenFun";
1014 TProfile2D *referenceFlowGenFun = dynamic_cast<TProfile2D*>(referenceFlowProfiles->FindObject(referenceFlowGenFunName.Data()));
1015 if(referenceFlowGenFun)
1016 {
1017 this->SetReferenceFlowGenFun(referenceFlowGenFun);
1018 } else
1019 {
1020 cout<<endl;
1021 cout<<"WARNING (GFC): referenceFlowGenFun is NULL in AFAWGFC::GPFRFO() !!!!"<<endl;
1022 cout<<endl;
1023 exit(0);
1024 }
1025 TString qvectorComponentsName = "fQvectorComponents";
1026 TProfile *qvectorComponents = dynamic_cast<TProfile*>(referenceFlowProfiles->FindObject(qvectorComponentsName.Data()));
1027 if(qvectorComponents)
1028 {
1029 this->SetQvectorComponents(qvectorComponents);
1030 } else
1031 {
1032 cout<<endl;
1033 cout<<"WARNING (GFC): qvectorComponents is NULL in AFAWGFC::GPFRFO() !!!!"<<endl;
1034 cout<<endl;
1035 exit(0);
1036 }
1037 // <<w^2>>, where w = wPhi*wPt*wEta:
1038 TString averageOfSquaredWeightName = "fAverageOfSquaredWeight";
1039 TProfile *averageOfSquaredWeight = dynamic_cast<TProfile*>(referenceFlowProfiles->FindObject(averageOfSquaredWeightName.Data()));
1040 if(averageOfSquaredWeight)
1041 {
1042 this->SetAverageOfSquaredWeight(averageOfSquaredWeight);
1043 } else
1044 {
1045 cout<<endl;
1046 cout<<"WARNING (GFC): averageOfSquaredWeight is NULL in AFAWGFC::GPFRFO() !!!!"<<endl;
1047 cout<<endl;
1048 exit(0);
1049 }
1050
1051 // d) Get pointers to all objects in the list fReferenceFlowResults:
1052 // Final results for isotropic cumulants for reference flow:
1053 TString referenceFlowCumulantsName = "fReferenceFlowCumulants";
1054 TH1D *referenceFlowCumulants = dynamic_cast<TH1D*>(referenceFlowResults->FindObject(referenceFlowCumulantsName.Data()));
1055 if(referenceFlowCumulants)
1056 {
1057 this->SetReferenceFlowCumulants(referenceFlowCumulants);
1058 } else
1059 {
1060 cout<<endl;
1061 cout<<"WARNING (GFC): referenceFlowCumulants is NULL in AFAWGFC::GPFRFO() !!!!"<<endl;
1062 cout<<endl;
1063 exit(0);
1064 }
1065 // Final results for reference flow:
1066 TString referenceFlowName = "fReferenceFlow";
1067 TH1D *referenceFlow = dynamic_cast<TH1D*>(referenceFlowResults->FindObject(referenceFlowName.Data()));
1068 if(referenceFlow)
1069 {
1070 this->SetReferenceFlow(referenceFlow);
1071 } else
1072 {
1073 cout<<endl;
1074 cout<<"WARNING (GFC): referenceFlow is NULL in AFAWGFC::GPFRFO() !!!!"<<endl;
1075 cout<<endl;
1076 exit(0);
1077 }
1078 // Final results for resolution:
1079 TString chiName = "fChi";
1080 TH1D *chi = dynamic_cast<TH1D*>(referenceFlowResults->FindObject(chiName.Data()));
1081 if(chi)
1082 {
1083 this->SetChi(chi);
1084 } else
1085 {
1086 cout<<endl;
1087 cout<<"WARNING (GFC): chi is NULL in AFAWGFC::GPFRFO() !!!!"<<endl;
1088 cout<<endl;
1089 exit(0);
1090 }
1091
1092} // end of void AliFlowAnalysisWithCumulants::GetPointersForReferenceFlowObjects()
1093
1094//================================================================================================================
1095
1096void AliFlowAnalysisWithCumulants::GetPointersForDiffFlowObjects()
1097{
1098 // Get pointers to all objects relevant for differential flow.
1099
1100 // a) Define flags locally (to be improved: should I promote flags to data members?);
1101 // b) Get pointer to base list for differential flow fDiffFlowList and nested lists fDiffFlowListProfiles and fDiffFlowListResults;
1102 // c) Get pointer to profile fDiffFlowFlags holding all flags for differential flow;
1103 // d) Get pointers to all profiles in the list fDiffFlowProfiles;
1104 // e) Get pointers to all profiles in the list fDiffFlowResults.
1105
1106 // a) Define flags locally (to be improved: should I promote flags to data members?):
1107 TString reIm[2] = {"Re","Im"};
1108 TString rpPoi[2] = {"RP","POI"};
1109 TString ptEta[2] = {"p_{t}","#eta"};
1110 TString order[4] = {"2nd order","4th order","6th order","8th order"};
1111 //TString differentialFlowIndex[4] = {"v'{2}","v'{4}","v'{6}","v'{8}"};
1112
1113 // b) Get pointer to base list for differential flow fDiffFlowList and nested lists fDiffFlowListProfiles and fDiffFlowListResults:
1114 TList *diffFlowList = dynamic_cast<TList*>(fHistList->FindObject("Differential Flow")); // to be improved (hardwired name)
1115 if(!diffFlowList)
1116 {
1117 cout<<endl;
1118 cout<<"WARNING: diffFlowList is NULL in AFAWQC::GPFDFO() !!!!"<<endl;
1119 cout<<endl;
1120 exit(0);
1121 }
1122 TList *diffFlowProfiles = dynamic_cast<TList*>(diffFlowList->FindObject("Profiles")); // to be improved (hardwired name)
1123 if(!diffFlowProfiles)
1124 {
1125 cout<<endl;
1126 cout<<"WARNING: diffFlowProfiles is NULL in AFAWQC::GPFDFO() !!!!"<<endl;
1127 cout<<endl;
1128 exit(0);
1129 }
1130 TList *diffFlowResults = dynamic_cast<TList*>(diffFlowList->FindObject("Results")); // to be improved (hardwired name)
1131 if(!diffFlowResults)
1132 {
1133 cout<<endl;
1134 cout<<"WARNING: diffFlowResults is NULL in AFAWQC::GPFDFO() !!!!"<<endl;
1135 cout<<endl;
1136 exit(0);
1137 }
1138
1139 // c) Get pointer to profile holding flags:
1140 TString diffFlowFlagsName = "fDiffFlowFlags";
1141 TProfile *diffFlowFlags = dynamic_cast<TProfile*>(diffFlowList->FindObject(diffFlowFlagsName.Data()));
1142 if(diffFlowFlags)
1143 {
1144 this->SetDiffFlowFlags(diffFlowFlags);
1145 } else
1146 {
1147 cout<<endl;
1148 cout<<"WARNING (GFC): diffFlowFlags is NULL in AFAWGFC::GPFDFO() !!!!"<<endl;
1149 cout<<endl;
1150 exit(0);
1151 }
1152
1153 // d) Get pointers to all profiles in the list fDiffFlowListProfiles:
1154 // Generating functions for differential flow:
1155 TProfile3D *diffFlowGenFun[2][2][2] = {{{NULL}}};
1156 for(Int_t ri=0;ri<2;ri++)
1157 {
1158 for(Int_t rp=0;rp<2;rp++)
1159 {
1160 for(Int_t pe=0;pe<2;pe++)
1161 {
1162 diffFlowGenFun[ri][rp][pe] = dynamic_cast<TProfile3D*> // to be improved - harwired name fDiffFlowGenFun in the line bellow
1163 (diffFlowProfiles->FindObject(Form("fDiffFlowGenFun (%s, %s, %s)",reIm[ri].Data(),rpPoi[rp].Data(),ptEta[pe].Data())));
1164 if(diffFlowGenFun[ri][rp][pe])
1165 {
1166 this->SetDiffFlowGenFun(diffFlowGenFun[ri][rp][pe],ri,rp,pe);
1167 } else
1168 {
1169 cout<<endl;
1170 cout<<"WARNING (GFC): "<<Form("diffFlowGenFun[%d][%d][%d]",ri,rp,pe)<<" is NULL in AFAWGFC::GPFDFO() !!!!"<<endl;
1171 cout<<endl;
1172 exit(0);
1173 }
1174 }
1175 }
1176 }
1177 // Number of particles in pt/eta bin for RPs/POIs:
1178 TProfile *noOfParticlesInBin[2][2] = {{NULL}};
1179 for(Int_t rp=0;rp<2;rp++)
1180 {
1181 for(Int_t pe=0;pe<2;pe++)
1182 {
1183 noOfParticlesInBin[rp][pe] = dynamic_cast<TProfile*> // to be improved - harwired name fNoOfParticlesInBin in the line bellow
1184 (diffFlowProfiles->FindObject(Form("fNoOfParticlesInBin (%s, %s)",rpPoi[rp].Data(),ptEta[pe].Data())));
1185 if(noOfParticlesInBin[rp][pe])
1186 {
1187 this->SetNoOfParticlesInBin(noOfParticlesInBin[rp][pe],rp,pe);
1188 } else
1189 {
1190 cout<<endl;
1191 cout<<"WARNING (GFC): "<<Form("noOfParticlesInBin[%d][%d]",rp,pe)<<" is NULL in AFAWGFC::GPFDFO() !!!!"<<endl;
1192 cout<<endl;
1193 exit(0);
1194 }
1195 }
1196 }
1197 // Differential cumulants per pt/eta bin for RPs/POIs:
1198 TH1D *diffFlowCumulants[2][2][4] = {{{NULL}}};
1199 for(Int_t rp=0;rp<2;rp++)
1200 {
1201 for(Int_t pe=0;pe<2;pe++)
1202 {
1203 for(Int_t co=0;co<4;co++)
1204 {
1205 diffFlowCumulants[rp][pe][co] = dynamic_cast<TH1D*> // to be improved - hardwired name fDiffFlowCumulants in the line bellow
1206 (diffFlowResults->FindObject(Form("fDiffFlowCumulants (%s, %s, %s)",rpPoi[rp].Data(),ptEta[pe].Data(),order[co].Data())));
1207 if(diffFlowCumulants[rp][pe][co])
1208 {
1209 this->SetDiffFlowCumulants(diffFlowCumulants[rp][pe][co],rp,pe,co);
1210 } else
1211 {
1212 cout<<endl;
1213 cout<<"WARNING (GFC): "<<Form("diffFlowCumulants[%d][%d][%d]",rp,pe,co)<<" is NULL in AFAWGFC::GPFDFO() !!!!"<<endl;
1214 cout<<endl;
1215 exit(0);
b6cd16a9 1216 }
d6130938 1217 }
1218 }
1219 }
1220 // Differential flow per pt/eta bin for RPs/POIs:
1221 TH1D *diffFlow[2][2][4] = {{{NULL}}};
1222 for(Int_t rp=0;rp<2;rp++)
1223 {
1224 for(Int_t pe=0;pe<2;pe++)
1225 {
1226 for(Int_t co=0;co<4;co++)
1227 {
1228 diffFlow[rp][pe][co] = dynamic_cast<TH1D*> // to be improved - hardwired name fDiffFlow in the line bellow
1229 (diffFlowResults->FindObject(Form("fDiffFlow (%s, %s, %s)",rpPoi[rp].Data(),ptEta[pe].Data(),order[co].Data())));
1230 if(diffFlow[rp][pe][co])
1231 {
1232 this->SetDiffFlow(diffFlow[rp][pe][co],rp,pe,co);
1233 } else
1234 {
1235 cout<<endl;
1236 cout<<"WARNING (GFC): "<<Form("diffFlow[%d][%d][%d]",rp,pe,co)<<" is NULL in AFAWGFC::GPFDFO() !!!!"<<endl;
1237 cout<<endl;
1238 exit(0);
1239 }
1240 }
1241 }
1242 }
1243
1244} // end of void AliFlowAnalysisWithCumulants::GetPointersForDiffFlowObjects()
1245
1246//================================================================================================================
1247
1248void AliFlowAnalysisWithCumulants::CalculateIntegratedFlow(TString rpPoi)
1249{
1250 // Calculate final results for integrated flow of RPs and POIs.
1251 // (to be improved - this method can be implemented much better)
1252
1253 Int_t rp = -1;
1254
1255 if(rpPoi == "RP")
1256 {
1257 rp = 0;
1258 } else if(rpPoi == "POI")
1259 {
1260 rp = 1;
1261 }
1262
1263 // pt yield:
1264 TH1F *yieldPt = NULL;
1265
1266 if(rpPoi == "POI")
1267 {
1268 yieldPt = (TH1F*)(fCommonHists->GetHistPtPOI())->Clone();
1269 } else if(rpPoi == "RP")
1270 {
1271 yieldPt = (TH1F*)(fCommonHists->GetHistPtRP())->Clone();
1272 }
1273
1274 TH1D *flow2ndPt = (TH1D*)fDiffFlow[rp][0][0]->Clone();
1275 TH1D *flow4thPt = (TH1D*)fDiffFlow[rp][0][1]->Clone();
1276 TH1D *flow6thPt = (TH1D*)fDiffFlow[rp][0][2]->Clone();
1277 TH1D *flow8thPt = (TH1D*)fDiffFlow[rp][0][3]->Clone();
1278 Double_t dvn2nd = 0., dvn4th = 0., dvn6th = 0., dvn8th = 0.; // differential flow
1279 Double_t dErrvn2nd = 0., dErrvn4th = 0., dErrvn6th = 0., dErrvn8th = 0.; // error on differential flow
1280 Double_t dVn2nd = 0., dVn4th = 0., dVn6th = 0., dVn8th = 0.; // integrated flow
1281 Double_t dErrVn2nd = 0., dErrVn4th = 0., dErrVn6th = 0., dErrVn8th = 0.; // error on integrated flow
1282 Double_t dYield = 0.; // pt yield
1283 Double_t dSum2nd = 0., dSum4th = 0., dSum6th = 0., dSum8th = 0.; // needed for normalizing integrated flow
1284 fnBinsPt = flow2ndPt->GetXaxis()->GetNbins();
1285 // looping over pt bins:
1286 for(Int_t p=1;p<=fnBinsPt;p++)
1287 {
1288 dvn2nd = flow2ndPt->GetBinContent(p);
1289 dvn4th = flow4thPt->GetBinContent(p);
1290 dvn6th = flow6thPt->GetBinContent(p);
1291 dvn8th = flow8thPt->GetBinContent(p);
1292
1293 dErrvn2nd = flow2ndPt->GetBinError(p);
1294 dErrvn4th = flow4thPt->GetBinError(p);
1295 dErrvn6th = flow6thPt->GetBinError(p);
1296 dErrvn8th = flow8thPt->GetBinError(p);
1297
1298 dYield = yieldPt->GetBinContent(p);
1299
1300 dVn2nd += dvn2nd*dYield;
1301 dVn4th += dvn4th*dYield;
1302 dVn6th += dvn6th*dYield;
1303 dVn8th += dvn8th*dYield;
1304
1305 dSum2nd += dYield;
1306 dSum4th += dYield;
1307 dSum6th += dYield;
1308 dSum8th += dYield;
1309
1310 dErrVn2nd += dYield*dYield*dErrvn2nd*dErrvn2nd;
1311 dErrVn4th += dYield*dYield*dErrvn4th*dErrvn4th;
1312 dErrVn6th += dYield*dYield*dErrvn6th*dErrvn6th;
1313 dErrVn8th += dYield*dYield*dErrvn8th*dErrvn8th;
1314 } // end of for(Int_t p=1;p<=fnBinsPt;p++)
1315
1316 // normalizing the results for integrated flow:
1317 if(dSum2nd)
1318 {
1319 dVn2nd /= dSum2nd;
1320 dErrVn2nd /= (dSum2nd*dSum2nd);
1321 dErrVn2nd = TMath::Sqrt(dErrVn2nd);
1322 }
1323 if(dSum4th)
1324 {
1325 dVn4th /= dSum4th;
1326 dErrVn4th /= (dSum4th*dSum4th);
1327 dErrVn4th = TMath::Sqrt(dErrVn4th);
1328 }
1329 if(dSum6th)
1330 {
1331 dVn6th /= dSum6th;
1332 dErrVn6th /= (dSum6th*dSum6th);
1333 dErrVn6th = TMath::Sqrt(dErrVn6th);
1334 }
1335 if(dSum8th)
1336 {
1337 dVn8th /= dSum8th;
1338 dErrVn8th /= (dSum8th*dSum8th);
1339 dErrVn8th = TMath::Sqrt(dErrVn8th);
1340 }
1341
1342 // storing the results for integrated flow in common hist results:
1343 if(rpPoi == "POI")
1344 {
1345 fCommonHistsResults2nd->FillIntegratedFlowPOI(dVn2nd,dErrVn2nd);
1346 fCommonHistsResults4th->FillIntegratedFlowPOI(dVn4th,dErrVn4th);
1347 fCommonHistsResults6th->FillIntegratedFlowPOI(dVn6th,dErrVn6th);
1348 fCommonHistsResults8th->FillIntegratedFlowPOI(dVn8th,dErrVn8th);
1349 }
1350 else if(rpPoi == "RP")
1351 {
1352 fCommonHistsResults2nd->FillIntegratedFlowRP(dVn2nd,dErrVn2nd);
1353 fCommonHistsResults4th->FillIntegratedFlowRP(dVn4th,dErrVn4th);
1354 fCommonHistsResults6th->FillIntegratedFlowRP(dVn6th,dErrVn6th);
1355 fCommonHistsResults8th->FillIntegratedFlowRP(dVn8th,dErrVn8th);
1356 }
1357
1358 delete flow2ndPt;
1359 delete flow4thPt;
1360 delete flow6thPt;
1361 delete flow8thPt;
1362 delete yieldPt;
1363
1364} // end of void AliFlowAnalysisWithCumulants::CalculateIntegratedFlow(TString rpPoi)
1365
1366//================================================================================================================
1367
1368void AliFlowAnalysisWithCumulants::FillCommonHistResultsForDifferentialFlow(TString rpPoi)
1369{
1370 // Fill common result histograms for differential flow.
1371 // (to be improved - this method can be implemented much better)
1372
1373 Int_t rp = -1;
1374
1375 if(rpPoi == "RP")
1376 {
1377 rp = 0;
1378 } else if(rpPoi == "POI")
1379 {
1380 rp = 1;
1381 }
1382
1383 // pt:
1384 for(Int_t p=1;p<=fnBinsPt;p++)
1385 {
1386 // Result:
1387 Double_t v2 = fDiffFlow[rp][0][0]->GetBinContent(p);
1388 Double_t v4 = fDiffFlow[rp][0][1]->GetBinContent(p);
1389 Double_t v6 = fDiffFlow[rp][0][2]->GetBinContent(p);
1390 Double_t v8 = fDiffFlow[rp][0][3]->GetBinContent(p);
1391 // Error:
1392 Double_t v2Error = fDiffFlow[rp][0][0]->GetBinError(p);
1393 Double_t v4Error = fDiffFlow[rp][0][1]->GetBinError(p);
1394 Double_t v6Error = fDiffFlow[rp][0][2]->GetBinError(p);
1395 Double_t v8Error = fDiffFlow[rp][0][3]->GetBinError(p);
1396 // Fill common hist results:
1397 if(rpPoi == "RP")
1398 {
1399 fCommonHistsResults2nd->FillDifferentialFlowPtRP(p,v2,v2Error);
1400 fCommonHistsResults4th->FillDifferentialFlowPtRP(p,v4,v4Error);
1401 fCommonHistsResults6th->FillDifferentialFlowPtRP(p,v6,v6Error);
1402 fCommonHistsResults8th->FillDifferentialFlowPtRP(p,v8,v8Error);
1403 } else if(rpPoi == "POI")
1404 {
1405 fCommonHistsResults2nd->FillDifferentialFlowPtPOI(p,v2,v2Error);
1406 fCommonHistsResults4th->FillDifferentialFlowPtPOI(p,v4,v4Error);
1407 fCommonHistsResults6th->FillDifferentialFlowPtPOI(p,v6,v6Error);
1408 fCommonHistsResults8th->FillDifferentialFlowPtPOI(p,v8,v8Error);
1409 }
1410 } // end of for(Int_t p=1;p<=fnBinsPt;p++)
1411
1412 // eta:
1413 for(Int_t e=1;e<=fnBinsEta;e++)
1414 {
1415 // Results:
1416 Double_t v2 = fDiffFlow[rp][1][0]->GetBinContent(e);
1417 Double_t v4 = fDiffFlow[rp][1][1]->GetBinContent(e);
1418 Double_t v6 = fDiffFlow[rp][1][2]->GetBinContent(e);
1419 Double_t v8 = fDiffFlow[rp][1][3]->GetBinContent(e);
1420 // Errors:
1421 Double_t v2Error = fDiffFlow[rp][1][0]->GetBinError(e);
1422 Double_t v4Error = fDiffFlow[rp][1][1]->GetBinError(e);
1423 Double_t v6Error = fDiffFlow[rp][1][2]->GetBinError(e);
1424 Double_t v8Error = fDiffFlow[rp][1][3]->GetBinError(e);
1425 // Fill common hist results:
1426 if(rpPoi == "RP")
1427 {
1428 fCommonHistsResults2nd->FillDifferentialFlowEtaRP(e,v2,v2Error);
1429 fCommonHistsResults4th->FillDifferentialFlowEtaRP(e,v4,v4Error);
1430 fCommonHistsResults6th->FillDifferentialFlowEtaRP(e,v6,v6Error);
1431 fCommonHistsResults8th->FillDifferentialFlowEtaRP(e,v8,v8Error);
1432 } else if(rpPoi == "POI")
1433 {
1434 fCommonHistsResults2nd->FillDifferentialFlowEtaPOI(e,v2,v2Error);
1435 fCommonHistsResults4th->FillDifferentialFlowEtaPOI(e,v4,v4Error);
1436 fCommonHistsResults6th->FillDifferentialFlowEtaPOI(e,v6,v6Error);
1437 fCommonHistsResults8th->FillDifferentialFlowEtaPOI(e,v8,v8Error);
1438 }
1439 } // end of for(Int_t e=1;e<=fnBinsEta;e++)
1440
1441} // end of void AliFlowAnalysisWithCumulants::FillCommonHistResultsForDifferentialFlow(TString rpPoi)
1442
1443//================================================================================================================
1444
1445void AliFlowAnalysisWithCumulants::CalculateDifferentialFlow(TString rpPoi, TString ptEta)
1446{
1447 // Calculate differential flow for RPs/POIs vs pt/eta from cumulants.
1448
1449 Int_t rp = -1; // RP or POI
1450 Int_t pe = -1; // pt or eta
1451
1452 if(rpPoi == "RP")
1453 {
1454 rp = 0;
1455 } else if(rpPoi == "POI")
1456 {
1457 rp = 1;
1458 }
1459 if(ptEta == "pt")
1460 {
1461 pe = 0;
1462 } else if(ptEta == "eta")
1463 {
1464 pe = 1;
1465 }
1466
1467 // Reference cumulants:
1468 Double_t gfc2 = fReferenceFlowCumulants->GetBinContent(1); // reference 2nd order cumulant
1469 Double_t gfc4 = fReferenceFlowCumulants->GetBinContent(2); // reference 4th order cumulant
1470 Double_t gfc6 = fReferenceFlowCumulants->GetBinContent(3); // reference 6th order cumulant
1471 Double_t gfc8 = fReferenceFlowCumulants->GetBinContent(4); // reference 8th order cumulant
1472
1473 Int_t nBins = fDiffFlowCumulants[rp][pe][0]->GetXaxis()->GetNbins();
1474
1475 for(Int_t b=1;b<=nBins;b++)
1476 {
1477 // Differential cumulants:
1478 Double_t gfd2 = fDiffFlowCumulants[rp][pe][0]->GetBinContent(b); // differential 2nd order cumulant
1479 Double_t gfd4 = fDiffFlowCumulants[rp][pe][1]->GetBinContent(b); // differential 4th order cumulant
1480 Double_t gfd6 = fDiffFlowCumulants[rp][pe][2]->GetBinContent(b); // differential 6th order cumulant
1481 Double_t gfd8 = fDiffFlowCumulants[rp][pe][3]->GetBinContent(b); // differential 8th order cumulant
1482 // Differential flow:
1483 Double_t v2 = 0.; // v'{2,GFC}
1484 Double_t v4 = 0.; // v'{4,GFC}
1485 Double_t v6 = 0.; // v'{6,GFC}
1486 Double_t v8 = 0.; // v'{8,GFC}
1487 // 2nd order:
1488 if(gfc2>0.)
1489 {
1490 v2 = gfd2/pow(gfc2,0.5);
1491 fDiffFlow[rp][pe][0]->SetBinContent(b,v2);
1492 }
1493 // 4th order:
1494 if(gfc4<0.)
1495 {
1496 v4 = -gfd4/pow(-gfc4,.75);
1497 fDiffFlow[rp][pe][1]->SetBinContent(b,v4);
1498 }
1499 // 6th order:
1500 if(gfc6>0.)
1501 {
1502 v6 = gfd6/(4.*pow((1./4.)*gfc6,(5./6.)));
1503 fDiffFlow[rp][pe][2]->SetBinContent(b,v6);
1504 }
1505 // 8th order:
1506 if(gfc8<0.)
1507 {
1508 v8 = -gfd8/(33.*pow(-(1./33.)*gfc8,(7./8.)));
1509 fDiffFlow[rp][pe][3]->SetBinContent(b,v8);
1510 }
1511 } // end of for(Int_t b=1;b<=nBins;b++)
1512
1513} // end of void AliFlowAnalysisWithCumulants::CalculateDifferentialFlow(TString rpPoi,TString ptEta)
1514
1515//================================================================================================================
1516
1517void AliFlowAnalysisWithCumulants::CalculateDifferentialFlowErrors(TString rpPoi,TString ptEta)
1518{
1519 // Calculate errors of differential flow.
1520
1521 Int_t rp = -1; // RP or POI
1522 Int_t pe = -1; // pt or eta
1523
1524 if(rpPoi == "RP")
1525 {
1526 rp = 0;
1527 } else if(rpPoi == "POI")
1528 {
1529 rp = 1;
1530 }
1531 if(ptEta == "pt")
1532 {
1533 pe = 0;
1534 } else if(ptEta == "eta")
1535 {
1536 pe = 1;
1537 }
1538
1539 // Resolution chi:
1540 Double_t chi2 = fChi->GetBinContent(1);
1541 Double_t chi4 = fChi->GetBinContent(2);
1542 //Double_t chi6 = fChi->GetBinContent(3);
1543 //Double_t chi8 = fChi->GetBinContent(4);
1544
1545 Int_t nBins = fNoOfParticlesInBin[rp][pe]->GetXaxis()->GetNbins();
1546 for(Int_t b=1;b<=nBins;b++)
1547 {
1548 Int_t nParticles = (Int_t)fNoOfParticlesInBin[rp][pe]->GetBinEntries(b);
1549 // Error of 2nd order estimate:
1550 if(chi2>0. && nParticles>0.)
1551 {
1552 Double_t v2Error = pow((1./(2.*nParticles))*((1.+pow(chi2,2.))/pow(chi2,2.)),0.5);
1553 fDiffFlow[rp][pe][0]->SetBinError(b,v2Error);
1554 }
1555 // Error of 4th order estimate:
1556 if(chi4>0. && nParticles>0.)
1557 {
1558 Double_t v4Error = pow((1./(2.*nParticles))*((2.+6.*pow(chi4,2.)+pow(chi4,4.)+pow(chi4,6.))/pow(chi4,6.)),0.5);
1559 fDiffFlow[rp][pe][1]->SetBinError(b,v4Error);
1560 }
1561 // Error of 6th order estimate:
1562 //if(chi6>0. && nParticles>0.)
1563 //{
1564 // Double_t v6Error = ... // to be improved - yet to be calculated
1565 fDiffFlow[rp][pe][2]->SetBinError(b,0.);
1566 //}
1567 // Error of 8th order estimate:
1568 //if(chi8>0. && nParticles>0.)
1569 //{
1570 // Double_t v8Error = ... // to be improved - yet to be calculated
1571 fDiffFlow[rp][pe][3]->SetBinError(b,0.);
1572 //}
1573 } // end of for(Int_t b=1;b<=nBins;b++)
1574
1575} // end of void AliFlowAnalysisWithCumulants::CalculateDifferentialFlowErrors(TString rpPoi,TString ptEta)
1576
1577//================================================================================================================
1578
1579void AliFlowAnalysisWithCumulants::CalculateCumulantsForDiffFlow(TString rpPoi,TString ptEta)
1580{
1581 // Calculate cumulants for differential flow.
1582
1583 Int_t rp = -1; // RP or POI
1584 Int_t pe = -1; // pt or eta
1585
1586 if(rpPoi == "RP")
1587 {
1588 rp = 0;
1589 } else if(rpPoi == "POI")
1590 {
1591 rp = 1;
1592 }
1593 if(ptEta == "pt")
1594 {
1595 pe = 0;
1596 } else if(ptEta == "eta")
1597 {
1598 pe = 1;
1599 }
1600
1601 // [nBins][pMax][qMax]:
1602 Int_t nBins = fDiffFlowGenFun[0][rp][pe]->GetXaxis()->GetNbins();
1603 Int_t pMax = fDiffFlowGenFun[0][rp][pe]->GetYaxis()->GetNbins();
1604 Int_t qMax = fDiffFlowGenFun[0][rp][pe]->GetZaxis()->GetNbins();
1605 // <G[p][q]>
1606 TMatrixD dAvG(pMax,qMax);
1607 dAvG.Zero();
1608 for(Int_t p=0;p<pMax;p++)
1609 {
1610 for(Int_t q=0;q<qMax;q++)
1611 {
1612 dAvG(p,q) = fReferenceFlowGenFun->GetBinContent(fReferenceFlowGenFun->GetBin(p+1,q+1));
1613 }
1614 }
1615 // Loop over pt/eta bins and calculate differential cumulants:
1616 for(Int_t b=0;b<nBins;b++)
1617 {
1618 Double_t gfc[5] = {0.}; // to be improved (hardwired 5)
1619 Double_t D[5] = {0.}; // D_{p} in Eq. (11) in Practical guide // to be improved (hardwired 5)
1620 // ptBinRPNoOfParticles[b]=fPtBinRPNoOfParticles->GetBinEntries(b+1);
1621 for(Int_t p=0;p<pMax;p++)
1622 {
1623 Double_t tempSum = 0.;
1624 for(Int_t q=0;q<qMax;q++)
1625 {
1626 if(TMath::Abs(dAvG(p,q))>1.e-44)
1627 {
1628 Double_t X = fDiffFlowGenFun[0][rp][pe]->GetBinContent(fDiffFlowGenFun[0][rp][pe]->GetBin(b+1,p+1,q+1))/dAvG(p,q); // see Ollitrault's Practical guide (Eq. 11)
1629 Double_t Y = fDiffFlowGenFun[1][rp][pe]->GetBinContent(fDiffFlowGenFun[0][rp][pe]->GetBin(b+1,p+1,q+1))/dAvG(p,q); // see Ollitrault's Practical guide (Eq. 11)
1630 tempSum += cos(fMultiple*2.*q*TMath::Pi()/qMax)*X
1631 + sin(fMultiple*2.*q*TMath::Pi()/qMax)*Y;
1632 }
1633 }
1634 D[p] = (pow(fR0*pow(p+1.0,0.5),fMultiple)/qMax)*tempSum;
1635 }
1636 gfc[0] = (1./(fR0*fR0))*(5.*D[0]-5.*D[1]+(10./3.)*D[2]-(5./4.)*D[3]+(1./5.)*D[4]);
1637 gfc[1] = (1./pow(fR0,4.))*((-77./6.)*D[0]+(107./6.)*D[1]-(13./1.)*D[2]+(61./12.)*D[3]-(5./6.)*D[4]);
1638 gfc[2] = (1./pow(fR0,6.))*((71./2.)*D[0]-59.*D[1]+49.*D[2]-(41./2.)*D[3]+(7./2.)*D[4]);
1639 gfc[3] = (1./pow(fR0,8.))*(-84.*D[0]+156.*D[1]-144.*D[2]+66.*D[3]-12.*D[4]);
1640 // gfc[4] = (1./pow(fR0,10.))*(120.*D[0]-240.*D[1]+240.*D[2]-120.*D[3]+24.*D[4]); // 10th order cumulant (to be improved - where to store it?)
1641 // Store cumulants:
1642 for(Int_t co=0;co<4;co++)
1643 {
1644 fDiffFlowCumulants[rp][pe][co]->SetBinContent(b+1,gfc[co]);
1645 }
1646 }
1647
1648} // end of void AliFlowAnalysisWithCumulants::CalculateCumulantsForDiffFlow(TString rpPoi, TString ptEta)
1649
1650//================================================================================================================
1651
1652void AliFlowAnalysisWithCumulants::PrintFinalResults(TString type)
1653{
1654 // Printing on the screen the final results for reference flow and for integrated flow of RPs and POIs.
1655
1656 Int_t n = fHarmonic;
1657
1658 Double_t dVn[4] = {0.}; // array to hold Vn{2}, Vn{4}, Vn{6} and Vn{8}
1659 Double_t dVnErr[4] = {0.}; // array to hold errors of Vn{2}, Vn{4}, Vn{6} and Vn{8}
1660
1661 if(type == "RF")
1662 {
1663 dVn[0] = (fCommonHistsResults2nd->GetHistIntFlow())->GetBinContent(1);
1664 dVnErr[0] = (fCommonHistsResults2nd->GetHistIntFlow())->GetBinError(1);
1665 dVn[1] = (fCommonHistsResults4th->GetHistIntFlow())->GetBinContent(1);
1666 dVnErr[1] = (fCommonHistsResults4th->GetHistIntFlow())->GetBinError(1);
1667 dVn[2] = (fCommonHistsResults6th->GetHistIntFlow())->GetBinContent(1);
1668 dVnErr[2] = (fCommonHistsResults6th->GetHistIntFlow())->GetBinError(1);
1669 dVn[3] = (fCommonHistsResults8th->GetHistIntFlow())->GetBinContent(1);
1670 dVnErr[3] = (fCommonHistsResults8th->GetHistIntFlow())->GetBinError(1);
1671 } else if(type == "RP")
1672 {
1673 dVn[0] = (fCommonHistsResults2nd->GetHistIntFlowRP())->GetBinContent(1);
1674 dVnErr[0] = (fCommonHistsResults2nd->GetHistIntFlowRP())->GetBinError(1);
1675 dVn[1] = (fCommonHistsResults4th->GetHistIntFlowRP())->GetBinContent(1);
1676 dVnErr[1] = (fCommonHistsResults4th->GetHistIntFlowRP())->GetBinError(1);
1677 dVn[2] = (fCommonHistsResults6th->GetHistIntFlowRP())->GetBinContent(1);
1678 dVnErr[2] = (fCommonHistsResults6th->GetHistIntFlowRP())->GetBinError(1);
1679 dVn[3] = (fCommonHistsResults8th->GetHistIntFlowRP())->GetBinContent(1);
1680 dVnErr[3] = (fCommonHistsResults8th->GetHistIntFlowRP())->GetBinError(1);
1681 } else if(type == "POI")
1682 {
1683 dVn[0] = (fCommonHistsResults2nd->GetHistIntFlowPOI())->GetBinContent(1);
1684 dVnErr[0] = (fCommonHistsResults2nd->GetHistIntFlowPOI())->GetBinError(1);
1685 dVn[1] = (fCommonHistsResults4th->GetHistIntFlowPOI())->GetBinContent(1);
1686 dVnErr[1] = (fCommonHistsResults4th->GetHistIntFlowPOI())->GetBinError(1);
1687 dVn[2] = (fCommonHistsResults6th->GetHistIntFlowPOI())->GetBinContent(1);
1688 dVnErr[2] = (fCommonHistsResults6th->GetHistIntFlowPOI())->GetBinError(1);
1689 dVn[3] = (fCommonHistsResults8th->GetHistIntFlowPOI())->GetBinContent(1);
1690 dVnErr[3] = (fCommonHistsResults8th->GetHistIntFlowPOI())->GetBinError(1);
1691 } else
1692 {
1693 cout<<endl;
1694 cout<<" WARNING: Impossible type (can be RF, RP or POI) !!!!"<<endl;
1695 cout<<" Results will not be printed on the screen."<<endl;
1696 cout<<endl;
1697 exit(0);
1698 }
1699
1700 TString title = " flow estimates from GF-cumulants";
1701 TString subtitle = " (";
1702
1703 if(!(fUsePhiWeights||fUsePtWeights||fUseEtaWeights))
1704 {
1705 subtitle.Append(type);
1706 subtitle.Append(", without weights)");
1707 } else
1708 {
1709 subtitle.Append(type);
1710 subtitle.Append(", with weights)");
1711 }
1712
1713 cout<<endl;
1714 cout<<"*************************************"<<endl;
1715 cout<<"*************************************"<<endl;
1716 cout<<title.Data()<<endl;
1717 cout<<subtitle.Data()<<endl;
1718 cout<<endl;
1719
1720 for(Int_t i=0;i<4;i++)
1721 {
1722 cout<<" v_"<<n<<"{"<<2*(i+1)<<"} = "<<dVn[i]<<" +/- "<<dVnErr[i]<<endl;
1723 }
1724
1725 cout<<endl;
1726 if(type == "RF")
1727 {
1728 cout<<" nEvts = "<<(Int_t)fCommonHists->GetHistMultRP()->GetEntries()<<", <M> = "<<(Double_t)fCommonHists->GetHistMultRP()->GetMean()<<endl;
1729 }
1730 else if (type == "RP")
1731 {
1732 cout<<" nEvts = "<<(Int_t)fCommonHists->GetHistMultRP()->GetEntries()<<", <M> = "<<(Double_t)fCommonHists->GetHistMultRP()->GetMean()<<endl;
1733 }
1734 else if (type == "POI")
1735 {
1736 cout<<" nEvts = "<<(Int_t)fCommonHists->GetHistMultPOI()->GetEntries()<<", <M> = "<<(Double_t)fCommonHists->GetHistMultPOI()->GetMean()<<endl;
1737 }
1738 cout<<"*************************************"<<endl;
1739 cout<<"*************************************"<<endl;
1740 cout<<endl;
1741
1742} // end of AliFlowAnalysisWithCumulants::PrintFinalResults(TString type);
1743
1744//================================================================================================================
1745
1746void AliFlowAnalysisWithCumulants::FillCommonHistResultsForReferenceFlow()
1747{
1748 // Fill in AliFlowCommonHistResults dedicated histograms for reference flow.
1749
1750 // Results:
1751 Double_t v2 = fReferenceFlow->GetBinContent(1);
1752 Double_t v4 = fReferenceFlow->GetBinContent(2);
1753 Double_t v6 = fReferenceFlow->GetBinContent(3);
1754 Double_t v8 = fReferenceFlow->GetBinContent(4);
1755 // Errors:
1756 Double_t v2Error = fReferenceFlow->GetBinError(1);
1757 Double_t v4Error = fReferenceFlow->GetBinError(2);
1758 Double_t v6Error = fReferenceFlow->GetBinError(3);
1759 Double_t v8Error = fReferenceFlow->GetBinError(4);
1760 // Fill results end errors in common hist results:
1761 fCommonHistsResults2nd->FillIntegratedFlow(v2,v2Error);
1762 fCommonHistsResults4th->FillIntegratedFlow(v4,v4Error);
1763 fCommonHistsResults6th->FillIntegratedFlow(v6,v6Error);
1764 fCommonHistsResults8th->FillIntegratedFlow(v8,v8Error);
1765 // Chi:
1766 Double_t chi2 = fChi->GetBinContent(1);
1767 Double_t chi4 = fChi->GetBinContent(2);
1768 Double_t chi6 = fChi->GetBinContent(3);
1769 Double_t chi8 = fChi->GetBinContent(4);
1770 // Fill resolution chi in common hist results:
1771 fCommonHistsResults2nd->FillChi(chi2);
1772 fCommonHistsResults4th->FillChi(chi4);
1773 fCommonHistsResults6th->FillChi(chi6);
1774 fCommonHistsResults8th->FillChi(chi8);
1775
1776} // end of AliFlowAnalysisWithCumulants::FillCommonHistResultsForReferenceFlow()
1777
1778//================================================================================================================
1779
1780void AliFlowAnalysisWithCumulants::CalculateReferenceFlowError()
1781{
1782 // Calculate error of reference flow harmonics.
1783
1784 // Generating Function Cumulants:
1785 Double_t gfc2 = fReferenceFlowCumulants->GetBinContent(1); // GFC{2}
1786 Double_t gfc4 = fReferenceFlowCumulants->GetBinContent(2); // GFC{4}
1787 Double_t gfc6 = fReferenceFlowCumulants->GetBinContent(3); // GFC{6}
1788 Double_t gfc8 = fReferenceFlowCumulants->GetBinContent(4); // GFC{8}
1789 // Reference flow estimates:
1790 Double_t v2 = fReferenceFlow->GetBinContent(1); // v{2,GFC}
1791 Double_t v4 = fReferenceFlow->GetBinContent(2); // v{4,GFC}
1792 Double_t v6 = fReferenceFlow->GetBinContent(3); // v{6,GFC}
1793 Double_t v8 = fReferenceFlow->GetBinContent(4); // v{8,GFC}
1794 // Statistical errors of reference flow estimates:
1795 Double_t v2Error = 0.; // statistical error of v{2,GFC}
1796 Double_t v4Error = 0.; // statistical error of v{4,GFC}
1797 Double_t v6Error = 0.; // statistical error of v{6,GFC}
1798 Double_t v8Error = 0.; // statistical error of v{8,GFC}
1799 // Chi:
1800 Double_t chi2 = 0.;
1801 Double_t chi4 = 0.;
1802 Double_t chi6 = 0.;
1803 Double_t chi8 = 0.;
1804 // <Q-vector stuff>:
1805 Double_t dAvQx = fQvectorComponents->GetBinContent(1); // <Q_x>
1806 Double_t dAvQy = fQvectorComponents->GetBinContent(2); // <Q_y>
1807 Double_t dAvQ2x = fQvectorComponents->GetBinContent(3); // <(Q_x)^2>
1808 Double_t dAvQ2y = fQvectorComponents->GetBinContent(4); // <(Q_y)^2>
1809 // <w^2>:
1810 Double_t dAvw2 = 1.;
1811 if(fnEvts>0)
1812 {
1813 dAvw2 = fAverageOfSquaredWeight->GetBinContent(1);
1814 if(TMath::Abs(dAvw2)<1.e-44)
1815 {
1816 cout<<endl;
1817 cout<<" WARNING (GFC): Average of squared weight is 0 in GFC. Most probably one of the histograms"<<endl;
1818 cout<<" in the file \"weights.root\" was empty. Nothing will be calculated !!!!"<<endl;
1819 cout<<endl;
1820 }
1821 }
1822 // Calculating statistical error of v{2,GFC}:
1823 if(fnEvts>0. && fAvM>0. && dAvw2>0. && gfc2>=0.)
1824 {
1825 if((dAvQ2x+dAvQ2y-pow(dAvQx,2.)-pow(dAvQy,2.)-pow(pow(gfc2,(1./2.))*(fAvM/dAvw2),2.)>0.))
1826 {
1827 chi2 = (fAvM*v2)/pow(dAvQ2x+dAvQ2y-pow(dAvQx,2.)-pow(dAvQy,2.)-pow(v2*fAvM/dAvw2,2.),0.5);
1828 }
1829 if(TMath::Abs(chi2)>1.e-44)
1830 {
1831 v2Error = pow(((1./(2.*fAvM*fnEvts))*((1.+2.*pow(chi2,2))/(2.*pow(chi2,2)))),0.5);
1832 }
1833 }
1834 // Calculating statistical error of v{4,GFC}:
1835 if(fnEvts>0 && fAvM>0 && dAvw2>0 && gfc4<=0.)
1836 {
1837 if((dAvQ2x+dAvQ2y-pow(dAvQx,2.)-pow(dAvQy,2.)-pow(pow(-gfc4,(1./4.))*(fAvM/dAvw2),2.)>0.))
1838 {
1839 chi4 = (fAvM*v4)/pow(dAvQ2x+dAvQ2y-pow(dAvQx,2.)-pow(dAvQy,2.)-pow(v4*fAvM/dAvw2,2.),0.5);
1840 }
1841 if(TMath::Abs(chi4)>1.e-44)
1842 {
1843 v4Error = (1./(pow(2.*fAvM*fnEvts,0.5)))*pow((1.+4.*pow(chi4,2)+1.*pow(chi4,4.)+2.*pow(chi4,6.))/(2.*pow(chi4,6.)),0.5);
1844 }
1845 }
1846 // Calculating statistical error of v{6,GFC}:
1847 if(fnEvts>0 && fAvM>0 && dAvw2>0 && gfc6>=0.)
1848 {
1849 if((dAvQ2x+dAvQ2y-pow(dAvQx,2.)-pow(dAvQy,2.)-pow(pow((1./4.)*gfc6,(1./6.))*(fAvM/dAvw2),2.)>0.))
1850 {
1851 chi6 = (fAvM*v6)/pow(dAvQ2x+dAvQ2y-pow(dAvQx,2.)-pow(dAvQy,2.)-pow(v6*fAvM/dAvw2,2.),0.5);
1852 }
1853 if(TMath::Abs(chi6)>1.e-44)
1854 {
1855 v6Error = (1./(pow(2.*fAvM*fnEvts,0.5)))*pow((3.+18.*pow(chi6,2)+9.*pow(chi6,4.)+28.*pow(chi6,6.)
1856 +12.*pow(chi6,8.)+24.*pow(chi6,10.))/(24.*pow(chi6,10.)),0.5);
1857 }
1858 }
1859 // Calculating statistical error of v{8,GFC}:
1860 if(fnEvts>0 && fAvM>0 && dAvw2>0 && gfc8<=0.)
1861 {
1862 if((dAvQ2x+dAvQ2y-pow(dAvQx,2.)-pow(dAvQy,2.)-pow(pow(-(1./33.)*gfc8,(1./8.))*(fAvM/dAvw2),2.)>0.))
1863 {
1864 chi8=(fAvM*v8)/pow(dAvQ2x+dAvQ2y-pow(dAvQx,2.)-pow(dAvQy,2.)-pow(v8*fAvM/dAvw2,2.),0.5);
1865 }
1866 if(TMath::Abs(chi8)>1.e-44)
1867 {
1868 v8Error = (1./(pow(2.*fAvM*fnEvts,0.5)))*pow((12.+96.*pow(chi8,2.)+72.*pow(chi8,4.)+304.*pow(chi8,6.)
1869 +257.*pow(chi8,8.)+804.*pow(chi8,10.)+363.*pow(chi8,12.)+726.*pow(chi8,14.))/(726.*pow(chi8,14.)),0.5);
1870 }
1871 }
1872
1873 // Store errors for reference flow:
1874 fReferenceFlow->SetBinError(1,v2Error);
1875 fReferenceFlow->SetBinError(2,v4Error);
1876 fReferenceFlow->SetBinError(3,v6Error);
1877 fReferenceFlow->SetBinError(4,v8Error);
1878 // Store resolution chi:
1879 fChi->SetBinContent(1,chi2);
1880 fChi->SetBinContent(2,chi4);
1881 fChi->SetBinContent(3,chi6);
1882 fChi->SetBinContent(4,chi8);
1883
1884} // end of void AliFlowAnalysisWithCumulants::CalculateReferenceFlowError()
1885
1886//================================================================================================================
1887
1888void AliFlowAnalysisWithCumulants::CalculateReferenceFlow()
1889{
1890 // Calculate from isotropic cumulants reference flow.
1891
1892 // Generating Function Cumulants:
1893 Double_t gfc2 = fReferenceFlowCumulants->GetBinContent(1); // GFC{2}
1894 Double_t gfc4 = fReferenceFlowCumulants->GetBinContent(2); // GFC{4}
1895 Double_t gfc6 = fReferenceFlowCumulants->GetBinContent(3); // GFC{6}
1896 Double_t gfc8 = fReferenceFlowCumulants->GetBinContent(4); // GFC{8}
1897 // Double_t gfc10 = fReferenceFlowCumulants->GetBinContent(5); // GFC{10} overflow
1898 // Reference flow estimates:
1899 Double_t v2 = 0.; // v{2,GFC}
1900 Double_t v4 = 0.; // v{4,GFC}
1901 Double_t v6 = 0.; // v{6,GFC}
1902 Double_t v8 = 0.; // v{8,GFC}
1903 // Double_t v10 = 0.; // v{10,GFC} overflow
1904 // Calculate reference flow estimates from Q-cumulants:
1905 if(gfc2>=0.) v2 = pow(gfc2,1./2.);
1906 if(gfc4<=0.) v4 = pow(-1.*gfc4,1./4.);
1907 if(gfc6>=0.) v6 = pow((1./4.)*gfc6,1./6.);
1908 if(gfc8<=0.) v8 = pow((-1./33.)*gfc8,1./8.);
1909 // if(gfc10>=0.) v10 = pow((1./456.)*gfc10,1./10.);
1910 // Store results for reference flow:
1911 fReferenceFlow->SetBinContent(1,v2);
1912 fReferenceFlow->SetBinContent(2,v4);
1913 fReferenceFlow->SetBinContent(3,v6);
1914 fReferenceFlow->SetBinContent(4,v8);
1915 // fReferenceFlow->SetBinContent(5,v10); // overflow
1916
1917} // end of void AliFlowAnalysisWithCumulants::CalculateReferenceFlow()
1918
1919//================================================================================================================
1920
1921void AliFlowAnalysisWithCumulants::CalculateCumulantsForReferenceFlow()
1922{
1923 // Calculate cumulants for reference flow.
1924
1925 Int_t pMax = fReferenceFlowGenFun->GetXaxis()->GetNbins();
1926 Int_t qMax = fReferenceFlowGenFun->GetYaxis()->GetNbins();
1927
1928 // <G[p][q]>
1929 TMatrixD dAvG(pMax,qMax);
1930 dAvG.Zero();
1931 Bool_t someAvGEntryIsNegative = kFALSE;
1932 for(Int_t p=0;p<pMax;p++)
1933 {
1934 for(Int_t q=0;q<qMax;q++)
1935 {
1936 dAvG(p,q) = fReferenceFlowGenFun->GetBinContent(fReferenceFlowGenFun->GetBin(p+1,q+1));
1937 if(dAvG(p,q)<0.)
1938 {
1939 someAvGEntryIsNegative = kTRUE;
1940 cout<<endl;
1941 cout<<" WARNING: "<<Form("<G[%d][%d]> is negative !!!! GFC results are meaningless.",p,q)<<endl;
1942 cout<<endl;
1943 }
1944 }
1945 }
1946
1947 // C[p][q] (generating function for the cumulants)
1948 TMatrixD dC(pMax,qMax);
1949 dC.Zero();
1950 if(fAvM>0. && !someAvGEntryIsNegative)
1951 {
1952 for(Int_t p=0;p<pMax;p++)
1953 {
1954 for(Int_t q=0;q<qMax;q++)
1955 {
1956 dC(p,q) = fAvM*(pow(dAvG(p,q),(1./fAvM))-1.);
1957 }
1958 }
1959 }
1960
1961 // Averaging the generating function for cumulants over azimuth
1962 // in order to eliminate detector effects.
1963 // <C[p][q]> (Remark: here <> stands for average over azimuth):
1964 TVectorD dAvC(pMax);
1965 dAvC.Zero();
1966 for(Int_t p=0;p<pMax;p++)
1967 {
1968 Double_t temp = 0.;
1969 for(Int_t q=0;q<qMax;q++)
1970 {
1971 temp += 1.*dC(p,q);
1972 }
1973 dAvC[p] = temp/qMax;
1974 }
1975
1976 // Finally, the isotropic cumulants for reference flow:
1977 TVectorD cumulant(pMax);
1978 cumulant.Zero();
1979 cumulant[0] = (-1./(60*fR0*fR0))*((-300.)*dAvC[0]+300.*dAvC[1]-200.*dAvC[2]+75.*dAvC[3]-12.*dAvC[4]);
1980 cumulant[1] = (-1./(6.*pow(fR0,4.)))*(154.*dAvC[0]-214.*dAvC[1]+156.*dAvC[2]-61.*dAvC[3]+10.*dAvC[4]);
1981 cumulant[2] = (3./(2.*pow(fR0,6.)))*(71.*dAvC[0]-118.*dAvC[1]+98.*dAvC[2]-41.*dAvC[3]+7.*dAvC[4]);
1982 cumulant[3] = (-24./pow(fR0,8.))*(14.*dAvC[0]-26.*dAvC[1]+24.*dAvC[2]-11.*dAvC[3]+2.*dAvC[4]);
1983 cumulant[4] = (120./pow(fR0,10.))*(5.*dAvC[0]-10.*dAvC[1]+10.*dAvC[2]-5.*dAvC[3]+1.*dAvC[4]);
1984
1985 // Store cumulants:
1986 // Remark: the highest order cumulant is on purpose in the overflow.
1987 for(Int_t co=0;co<pMax;co++) // cumulant order
1988 {
1989 fReferenceFlowCumulants->SetBinContent(co+1,cumulant[co]);
1990 }
1991
1992} // end of void AliFlowAnalysisWithCumulants::CalculateCumulantsForReferenceFlow()
1993
1994//================================================================================================================
1995
1996void AliFlowAnalysisWithCumulants::GetAvMultAndNoOfEvts()
1997{
1998 // From relevant common control histogram get average multiplicity of RPs and number of events.
1999
2000 fAvM = (Double_t)fCommonHists->GetHistMultRP()->GetMean();
2001 fnEvts = (Int_t)fCommonHists->GetHistMultRP()->GetEntries();
2002
2003} // end of void AliFlowAnalysisWithCumulants::GetAvMultAndNoOfEvts()
2004
2005//================================================================================================================
2006
2007void AliFlowAnalysisWithCumulants::InitializeArrays()
2008{
2009 // Initialize all arrays.
2010
2011 for(Int_t ri=0;ri<2;ri++)
2012 {
2013 for(Int_t rp=0;rp<2;rp++)
2014 {
2015 for(Int_t pe=0;pe<2;pe++)
2016 {
2017 fDiffFlowGenFun[ri][rp][pe] = NULL;
2018 }
2019 }
2020 }
2021 for(Int_t rp=0;rp<2;rp++)
2022 {
2023 for(Int_t pe=0;pe<2;pe++)
2024 {
2025 fNoOfParticlesInBin[rp][pe] = NULL;
2026 }
2027 }
2028 for(Int_t rp=0;rp<2;rp++)
2029 {
2030 for(Int_t pe=0;pe<2;pe++)
2031 {
2032 for(Int_t co=0;co<4;co++)
2033 {
2034 fDiffFlowCumulants[rp][pe][co] = NULL;
2035 fDiffFlow[rp][pe][co] = NULL;
2036 }
2037 }
2038 }
2039 for(Int_t i=0;i<3;i++)
2040 {
2041 fPrintFinalResults[i] = kTRUE;
2042 }
2043 for(Int_t r=0;r<10;r++)
2044 {
2045 fTuningR0[r] = 0.;
2046 for(Int_t pq=0;pq<5;pq++)
2047 {
2048 fTuningGenFun[r][pq] = NULL;
2049 fTuningCumulants[r][pq] = NULL;
2050 fTuningFlow[r][pq] = NULL;
2051 }
2052 }
2053
2054} // end of void AliFlowAnalysisWithCumulants::InitializeArrays()
2055
2056//================================================================================================================
2057
2058void AliFlowAnalysisWithCumulants::CrossCheckSettings()
2059{
2060 // Cross-check the user settings before starting.
2061
2062 // a) Cross check if the choice for multiplicity weight make sense.
2063
2064 // a) Cross check if the choice for multiplicity weight make sense:
2065 if(strcmp(fMultiplicityWeight->Data(),"unit") &&
2066 strcmp(fMultiplicityWeight->Data(),"multiplicity"))
2067 {
2068 cout<<endl;
2069 cout<<"WARNING (GFC): Multiplicity weight can be either \"unit\" or \"multiplicity\"."<<endl;
2070 cout<<" Certainly not \""<<fMultiplicityWeight->Data()<<"\"."<<endl;
2071 cout<<endl;
2072 exit(0);
2073 }
2074
2075
2076
2077} // end of void AliFlowAnalysisWithCumulants::CrossCheckSettings()
2078
2079//================================================================================================================
2080
2081void AliFlowAnalysisWithCumulants::AccessConstants()
2082{
2083 // Access needed common constants from AliFlowCommonConstants.
b6cd16a9 2084
d6130938 2085 fnBinsPhi = AliFlowCommonConstants::GetMaster()->GetNbinsPhi();
2086 fPhiMin = AliFlowCommonConstants::GetMaster()->GetPhiMin();
2087 fPhiMax = AliFlowCommonConstants::GetMaster()->GetPhiMax();
2088 if(fnBinsPhi) fPhiBinWidth = (fPhiMax-fPhiMin)/fnBinsPhi;
2089 fnBinsPt = AliFlowCommonConstants::GetMaster()->GetNbinsPt();
2090 fPtMin = AliFlowCommonConstants::GetMaster()->GetPtMin();
2091 fPtMax = AliFlowCommonConstants::GetMaster()->GetPtMax();
2092 if(fnBinsPt) fPtBinWidth = (fPtMax-fPtMin)/fnBinsPt;
2093 fnBinsEta = AliFlowCommonConstants::GetMaster()->GetNbinsEta();
2094 fEtaMin = AliFlowCommonConstants::GetMaster()->GetEtaMin();
2095 fEtaMax = AliFlowCommonConstants::GetMaster()->GetEtaMax();
2096 if(fnBinsEta) fEtaBinWidth = (fEtaMax-fEtaMin)/fnBinsEta;
2097
2098} // end of void AliFlowAnalysisWithCumulants::AccessConstants()
2099
2100//================================================================================================================
2101
2102void AliFlowAnalysisWithCumulants::BookAndFillWeightsHistograms()
2103{
2104 // book and fill histograms which hold phi, pt and eta weights
2105
2106 if(!fWeightsList)
2107 {
2108 cout<<"WARNING (GFC): fWeightsList is NULL in AFAWGFC::BAFWH() !!!!"<<endl;
2109 exit(0);
2110 }
2111
2112 TString fUseParticleWeightsName = "fUseParticleWeights";
2113 fUseParticleWeights = new TProfile(fUseParticleWeightsName.Data(),"0 = particle weight not used, 1 = particle weight used ",3,0,3);
2114 fUseParticleWeights->SetLabelSize(0.06);
2115 (fUseParticleWeights->GetXaxis())->SetBinLabel(1,"w_{#phi}");
2116 (fUseParticleWeights->GetXaxis())->SetBinLabel(2,"w_{p_{T}}");
2117 (fUseParticleWeights->GetXaxis())->SetBinLabel(3,"w_{#eta}");
2118 fUseParticleWeights->Fill(0.5,(Int_t)fUsePhiWeights);
2119 fUseParticleWeights->Fill(1.5,(Int_t)fUsePtWeights);
2120 fUseParticleWeights->Fill(2.5,(Int_t)fUseEtaWeights);
2121 fWeightsList->Add(fUseParticleWeights);
2122 if(fUsePhiWeights)
b6cd16a9 2123 {
d6130938 2124 if(fWeightsList->FindObject("phi_weights"))
b6cd16a9 2125 {
d6130938 2126 fPhiWeights = dynamic_cast<TH1F*>(fWeightsList->FindObject("phi_weights"));
2127 if(TMath::Abs(fPhiWeights->GetBinWidth(1)-fPhiBinWidth)>pow(10.,-6.))
b6cd16a9 2128 {
d6130938 2129 cout<<endl;
2130 cout<<"WARNING (GFC): Inconsistent binning in histograms for phi-weights throughout the code."<<endl;
2131 cout<<endl;
2132 exit(0);
2133 }
2134 } else
b6cd16a9 2135 {
d6130938 2136 cout<<endl;
2137 cout<<"WARNING (GFC): fWeightsList->FindObject(\"phi_weights\") is NULL in AFAWGFC::BAFWH() !!!!"<<endl;
2138 cout<<endl;
2139 exit(0);
b6cd16a9 2140 }
d6130938 2141 } // end of if(fUsePhiWeights)
2142
2143 if(fUsePtWeights)
2144 {
2145 if(fWeightsList->FindObject("pt_weights"))
2146 {
2147 fPtWeights = dynamic_cast<TH1D*>(fWeightsList->FindObject("pt_weights"));
2148 if(TMath::Abs(fPtWeights->GetBinWidth(1)-fPtBinWidth)>pow(10.,-6.))
2149 {
2150 cout<<endl;
2151 cout<<"WARNING (GFC): Inconsistent binning in histograms for pt-weights throughout the code."<<endl;
2152 cout<<endl;
2153 exit(0);
2154 }
2155 } else
2156 {
2157 cout<<endl;
2158 cout<<"WARNING (GFC): fWeightsList->FindObject(\"pt_weights\") is NULL in AFAWGFC::BAFWH() !!!!"<<endl;
2159 cout<<endl;
2160 exit(0);
2161 }
2162 } // end of if(fUsePtWeights)
2163
2164 if(fUseEtaWeights)
2165 {
2166 if(fWeightsList->FindObject("eta_weights"))
2167 {
2168 fEtaWeights = dynamic_cast<TH1D*>(fWeightsList->FindObject("eta_weights"));
2169 if(TMath::Abs(fEtaWeights->GetBinWidth(1)-fEtaBinWidth)>pow(10.,-6.))
2170 {
2171 cout<<endl;
2172 cout<<"WARNING (GFC): Inconsistent binning in histograms for eta-weights throughout the code."<<endl;
2173 cout<<endl;
2174 exit(0);
2175 }
2176 } else
2177 {
2178 cout<<endl;
2179 cout<<"WARNING (GFC): fUseEtaWeights && fWeightsList->FindObject(\"eta_weights\") is NULL in AFAWGFC::BAFWH() !!!!"<<endl;
2180 cout<<endl;
2181 exit(0);
2182 }
2183 } // end of if(fUseEtaWeights)
2184
2185} // end of AliFlowAnalysisWithCumulants::BookAndFillWeightsHistograms()
2186
2187//================================================================================================================
2188
2189void AliFlowAnalysisWithCumulants::BookEverythingForReferenceFlow()
2190{
2191 // Book all objects relevant for calculation of reference flow.
2192
2193 // a) Define static constants for array's boundaries;
2194 // b) Book profile to hold all flags for reference flow;
2195 // c) Book all event-by-event quantities;
2196 // d) Book all profiles;
2197 // e) Book all histograms.
2198
2199 // a) Define static constants for array's boundaries:
2200 static const Int_t pMax = 5;
2201 static const Int_t qMax = 11;
2202
2203 // b) Book profile to hold all flags for reference flow:
2204 TString referenceFlowFlagsName = "fReferenceFlowFlags";
2205 fReferenceFlowFlags = new TProfile(referenceFlowFlagsName.Data(),"Flags for Reference Flow",2,0,2);
2206 fReferenceFlowFlags->SetTickLength(-0.01,"Y");
2207 fReferenceFlowFlags->SetMarkerStyle(25);
2208 fReferenceFlowFlags->SetLabelSize(0.05);
2209 fReferenceFlowFlags->SetLabelOffset(0.02,"Y");
2210 fReferenceFlowFlags->GetXaxis()->SetBinLabel(1,"Particle weights");
2211 fReferenceFlowFlags->GetXaxis()->SetBinLabel(2,"Event weights");
2212 fReferenceFlowList->Add(fReferenceFlowFlags);
2213
2214 // c) Book all event-by-event quantities:
2215 fGEBE = new TMatrixD(pMax,qMax);
2216
2217 // d) Book all profiles:
2218 // Average of the generating function for reference flow <G[p][q]>:
2219 fReferenceFlowGenFun = new TProfile2D("fReferenceFlowGenFun","#LTG[p][q]#GT",pMax,0.,(Double_t)pMax,qMax,0.,(Double_t)qMax);
2220 fReferenceFlowGenFun->SetXTitle("p");
2221 fReferenceFlowGenFun->SetYTitle("q");
2222 fReferenceFlowProfiles->Add(fReferenceFlowGenFun);
2223 // Averages of Q-vector components:
2224 fQvectorComponents = new TProfile("fQvectorComponents","Averages of Q-vector components",4,0.,4.);
2225 fQvectorComponents->SetLabelSize(0.06);
2226 fQvectorComponents->SetMarkerStyle(25);
2227 fQvectorComponents->GetXaxis()->SetBinLabel(1,"#LTQ_{x}#GT"); // Q_{x}
2228 fQvectorComponents->GetXaxis()->SetBinLabel(2,"#LTQ_{y}#GT"); // Q_{y}
2229 fQvectorComponents->GetXaxis()->SetBinLabel(3,"#LTQ_{x}^{2}#GT"); // Q_{x}^{2}
2230 fQvectorComponents->GetXaxis()->SetBinLabel(4,"#LTQ_{y}^{2}#GT"); // Q_{y}^{2}
2231 fReferenceFlowProfiles->Add(fQvectorComponents);
2232 // <<w^2>>, where w = wPhi*wPt*wEta:
2233 fAverageOfSquaredWeight = new TProfile("fAverageOfSquaredWeight","#LT#LTw^{2}#GT#GT",1,0,1);
2234 fAverageOfSquaredWeight->SetLabelSize(0.06);
2235 fAverageOfSquaredWeight->SetMarkerStyle(25);
2236 fAverageOfSquaredWeight->SetLabelOffset(0.01);
2237 fAverageOfSquaredWeight->GetXaxis()->SetBinLabel(1,"#LT#LTw^{2}#GT#GT");
2238 fReferenceFlowProfiles->Add(fAverageOfSquaredWeight);
2239
2240 // e) Book all histograms:
2241 // Final results for isotropic cumulants for reference flow:
2242 TString referenceFlowCumulantsName = "fReferenceFlowCumulants";
2243 fReferenceFlowCumulants = new TH1D(referenceFlowCumulantsName.Data(),"Isotropic Generating Function Cumulants for reference flow",4,0,4); // to be improved (hw 4)
2244 fReferenceFlowCumulants->SetLabelSize(0.05);
2245 fReferenceFlowCumulants->SetMarkerStyle(25);
2246 fReferenceFlowCumulants->GetXaxis()->SetBinLabel(1,"GFC{2}");
2247 fReferenceFlowCumulants->GetXaxis()->SetBinLabel(2,"GFC{4}");
2248 fReferenceFlowCumulants->GetXaxis()->SetBinLabel(3,"GFC{6}");
2249 fReferenceFlowCumulants->GetXaxis()->SetBinLabel(4,"GFC{8}");
2250 fReferenceFlowResults->Add(fReferenceFlowCumulants);
2251 // Final results for reference flow:
2252 fReferenceFlow = new TH1D("fReferenceFlow","Reference flow",4,0,4); // to be improved (hardwired 4)
2253 fReferenceFlow->SetLabelSize(0.06);
2254 fReferenceFlow->SetMarkerStyle(25);
2255 fReferenceFlow->GetXaxis()->SetBinLabel(1,"v_{n}{2,GFC}");
2256 fReferenceFlow->GetXaxis()->SetBinLabel(2,"v_{n}{4,GFC}");
2257 fReferenceFlow->GetXaxis()->SetBinLabel(3,"v_{n}{6,GFC}");
2258 fReferenceFlow->GetXaxis()->SetBinLabel(4,"v_{n}{8,GFC}");
2259 fReferenceFlowResults->Add(fReferenceFlow);
2260 // Final results for resolution:
2261 fChi = new TH1D("fChi","Resolution",4,0,4); // to be improved (hardwired 4)
2262 fChi->SetLabelSize(0.06);
2263 fChi->SetMarkerStyle(25);
2264 fChi->GetXaxis()->SetBinLabel(1,"#chi_{2}");
2265 fChi->GetXaxis()->SetBinLabel(2,"#chi_{4}");
2266 fChi->GetXaxis()->SetBinLabel(3,"#chi_{6}");
2267 fChi->GetXaxis()->SetBinLabel(4,"#chi_{8}");
2268 fReferenceFlowResults->Add(fChi);
2269
2270} // end of void AliFlowAnalysisWithCumulants::BookEverythingForReferenceFlow()
2271
2272//================================================================================================================
2273
2274void AliFlowAnalysisWithCumulants::BookEverythingForTuning()
2275{
2276 // Book all objects relevant for tuning.
2277
2278 // a) Define pMax's and qMax's:
2279 // b) Book profile to hold all tuning parameters and flags;
2280 // c) Book all profiles;
2281 // d) Book all histograms.
2282
2283 // a) Define pMax's and qMax's:
2284 Int_t pMax[5] = {2,3,4,5,8};
2285 Int_t qMax[5] = {5,7,9,11,17};
2286
2287 // b) Book profile to hold all tuning parameters and flags:
2288 TString tuningFlagsName = "fTuningFlags";
2289 fTuningFlags = new TProfile(tuningFlagsName.Data(),"Tuning parameters",10,0,10);
2290 // fTuningFlags->SetTickLength(-0.01,"Y");
2291 fTuningFlags->SetMarkerStyle(25);
2292 fTuningFlags->SetLabelSize(0.05);
2293 fTuningFlags->SetLabelOffset(0.02,"X");
2294 for(Int_t r=1;r<=10;r++)
2295 {
2296 fTuningFlags->GetXaxis()->SetBinLabel(r,Form("r_{0,%d}",r-1));
2297 fTuningFlags->Fill(r-0.5,fTuningR0[r-1],1.);
2298 }
2299 fTuningList->Add(fTuningFlags);
2300
2301 // c) Book all profiles:
2302 // Average of the generating function for reference flow <G[p][q]> for different tuning parameters:
2303 for(Int_t r=0;r<10;r++)
2304 {
2305 for(Int_t pq=0;pq<5;pq++)
2306 {
2307 fTuningGenFun[r][pq] = new TProfile2D(Form("fTuningGenFun (r_{0,%i}, pq set %i)",r,pq),
2308 Form("#LTG[p][q]#GT for r_{0} = %f, p_{max} = %i, q_{max} = %i",fTuningR0[r],pMax[pq],qMax[pq]),
2309 pMax[pq],0.,(Double_t)pMax[pq],qMax[pq],0.,(Double_t)qMax[pq]);
2310 fTuningGenFun[r][pq]->SetXTitle("p");
2311 fTuningGenFun[r][pq]->SetYTitle("q");
2312 fTuningProfiles->Add(fTuningGenFun[r][pq]);
b6cd16a9 2313 }
2314 }
d6130938 2315
2316 // d) Book all histograms:
2317 // Final results for isotropic cumulants for reference flow for different tuning parameters:
2318 for(Int_t r=0;r<10;r++)
2319 {
2320 for(Int_t pq=0;pq<5;pq++)
2321 {
2322 fTuningCumulants[r][pq] = new TH1D(Form("fTuningCumulants (r_{0,%i}, pq set %i)",r,pq),
2323 Form("GFC for r_{0} = %f, p_{max} = %i, q_{max} = %i",fTuningR0[r],pMax[pq],qMax[pq]),
2324 pMax[pq],0,pMax[pq]);
2325 // fTuningCumulants[r][pq]->SetLabelSize(0.05);
2326 fTuningCumulants[r][pq]->SetMarkerStyle(25);
2327 for(Int_t b=1;b<=pMax[pq];b++)
2328 {
2329 fTuningCumulants[r][pq]->GetXaxis()->SetBinLabel(b,Form("GFC{%i}",2*b));
2330 }
2331 fTuningResults->Add(fTuningCumulants[r][pq]);
2332 }
2333 }
2334 // Final results for reference flow for different tuning parameters:
2335 for(Int_t r=0;r<10;r++)
2336 {
2337 for(Int_t pq=0;pq<5;pq++)
2338 {
2339 fTuningFlow[r][pq] = new TH1D(Form("fTuningFlow (r_{0,%i}, pq set %i)",r,pq),
2340 Form("Reference flow for r_{0} = %f, p_{max} = %i, q_{max} = %i",fTuningR0[r],pMax[pq],qMax[pq]),
2341 pMax[pq],0,pMax[pq]);
2342 // fTuningFlow[r][pq]->SetLabelSize(0.06);
2343 fTuningFlow[r][pq]->SetMarkerStyle(25);
2344 for(Int_t b=1;b<=pMax[pq];b++)
2345 {
2346 fTuningFlow[r][pq]->GetXaxis()->SetBinLabel(b,Form("v{%i,GFC}",2*b));
2347 }
2348 fTuningResults->Add(fTuningFlow[r][pq]);
2349 }
2350 }
2351
2352} // end of void AliFlowAnalysisWithCumulants::BookEverythingForTuning()
f1d945a1 2353
2188af53 2354//================================================================================================================
f1d945a1 2355
d6130938 2356void AliFlowAnalysisWithCumulants::BookEverythingForDiffFlow()
fd46c3dd 2357{
d6130938 2358 // Book all objects relevant for calculation of differential flow.
2359
2360 // a) Define static constants for array's boundaries;
2361 // b) Define local variables and local flags for booking;
2362 // c) Book profile to hold all flags for differential flow;
2363 // d) Book all event-by-event quantities;
2364 // e) Book all profiles;
2365 // f) Book all histograms.
2366
2367 // a) Define static constants for array's boundaries:
2368 static const Int_t pMax = 5;
2369 static const Int_t qMax = 11;
2370
2371 // b) Define local variables and local flags for booking:
2372 Int_t nBinsPtEta[2] = {fnBinsPt,fnBinsEta};
2373 Double_t minPtEta[2] = {fPtMin,fEtaMin};
2374 Double_t maxPtEta[2] = {fPtMax,fEtaMax};
2375 TString reIm[2] = {"Re","Im"};
2376 TString rpPoi[2] = {"RP","POI"};
2377 TString ptEta[2] = {"p_{t}","#eta"};
2378 TString order[4] = {"2nd order","4th order","6th order","8th order"};
2379
2380 // c) Book profile to hold all flags for differential flow:
2381 TString diffFlowFlagsName = "fDiffFlowFlags";
2382 fDiffFlowFlags = new TProfile(diffFlowFlagsName.Data(),"Flags for Differential Flow",1,0,1);
2383 fDiffFlowFlags->SetTickLength(-0.01,"Y");
2384 fDiffFlowFlags->SetMarkerStyle(25);
2385 fDiffFlowFlags->SetLabelSize(0.05);
2386 fDiffFlowFlags->SetLabelOffset(0.02,"Y");
2387 fDiffFlowFlags->GetXaxis()->SetBinLabel(1,"...");
2388 fDiffFlowList->Add(fDiffFlowFlags);
2389
2390 // d) Book all event-by-event quantities:
2391 // ... (to be improved - perhaps not needed)
2392
2393 // e) Book all profiles:
2394 // Generating functions for differential flow:
2395 for(Int_t ri=0;ri<2;ri++)
2396 {
2397 for(Int_t rp=0;rp<2;rp++)
2398 {
2399 for(Int_t pe=0;pe<2;pe++)
2400 {
2401 fDiffFlowGenFun[ri][rp][pe] = new TProfile3D(Form("fDiffFlowGenFun (%s, %s, %s)",reIm[ri].Data(),rpPoi[rp].Data(),ptEta[pe].Data()),
2402 Form("#LT%s[D[%s-bin][p][q]]#GT for %ss",reIm[ri].Data(),ptEta[pe].Data(),rpPoi[rp].Data()),
2403 nBinsPtEta[pe],minPtEta[pe],maxPtEta[pe],pMax,0.,(Double_t)pMax,qMax,0.,(Double_t)qMax);
2404 fDiffFlowGenFun[ri][rp][pe]->SetXTitle(ptEta[pe].Data());
2405 fDiffFlowGenFun[ri][rp][pe]->SetYTitle("p");
2406 fDiffFlowGenFun[ri][rp][pe]->SetZTitle("q");
2407 fDiffFlowGenFun[ri][rp][pe]->SetTitleOffset(1.44,"X");
2408 fDiffFlowGenFun[ri][rp][pe]->SetTitleOffset(1.44,"Y");
2409 fDiffFlowProfiles->Add(fDiffFlowGenFun[ri][rp][pe]);
2410 // to be improved - alternative // nBinsPtEta[pe],(Double_t)(fPtMin/fPtBinWidth),(Double_t)(fPtMax/fPtBinWidth),pMax,0.,(Double_t)pMax,qMax,0.,(Double_t)qMax);
2411 }
2412 }
2413 }
2414 // Number of particles in pt/eta bin for RPs/POIs:
2415 for(Int_t rp=0;rp<2;rp++)
2416 {
2417 for(Int_t pe=0;pe<2;pe++)
2418 {
2419 fNoOfParticlesInBin[rp][pe] = new TProfile(Form("fNoOfParticlesInBin (%s, %s)",rpPoi[rp].Data(),ptEta[pe].Data()),
2420 Form("Number of %ss per %s bin",rpPoi[rp].Data(),ptEta[pe].Data()),
2421 nBinsPtEta[pe],minPtEta[pe],maxPtEta[pe]);
2422 fNoOfParticlesInBin[rp][pe]->SetXTitle(ptEta[pe].Data());
2423 fDiffFlowProfiles->Add(fNoOfParticlesInBin[rp][pe]);
2424 }
2425 }
2426 // Differential cumulants per pt/eta bin for RPs/POIs:
2427 for(Int_t rp=0;rp<2;rp++)
2428 {
2429 for(Int_t pe=0;pe<2;pe++)
2430 {
2431 for(Int_t co=0;co<4;co++)
2432 {
2433 fDiffFlowCumulants[rp][pe][co] = new TH1D(Form("fDiffFlowCumulants (%s, %s, %s)",rpPoi[rp].Data(),ptEta[pe].Data(),order[co].Data()),
2434 Form("Differential %s cumulant for %ss vs %s",order[co].Data(),rpPoi[rp].Data(),ptEta[pe].Data()),
2435 nBinsPtEta[pe],minPtEta[pe],maxPtEta[pe]);
2436 fDiffFlowCumulants[rp][pe][co]->SetXTitle(ptEta[pe].Data());
2437 fDiffFlowResults->Add(fDiffFlowCumulants[rp][pe][co]);
2438 }
2439 }
2440 }
2441 // Differential flow per pt/eta bin for RPs/POIs:
2442 for(Int_t rp=0;rp<2;rp++)
2443 {
2444 for(Int_t pe=0;pe<2;pe++)
2445 {
2446 for(Int_t co=0;co<4;co++)
2447 {
2448 fDiffFlow[rp][pe][co] = new TH1D(Form("fDiffFlow (%s, %s, %s)",rpPoi[rp].Data(),ptEta[pe].Data(),order[co].Data()),
2449 Form("Differential flow from %s cumulant for %ss vs %s",order[co].Data(),rpPoi[rp].Data(),ptEta[pe].Data()),
2450 nBinsPtEta[pe],minPtEta[pe],maxPtEta[pe]);
2451 fDiffFlow[rp][pe][co]->SetXTitle(ptEta[pe].Data());
2452 fDiffFlowResults->Add(fDiffFlow[rp][pe][co]);
2453 }
2454 }
2455 }
2456
2457}// end of void AliFlowAnalysisWithCumulants::BookEverythingForDiffFlow()
2458
2459//================================================================================================================
2460
2461void AliFlowAnalysisWithCumulants::StoreReferenceFlowFlags()
2462{
2463 // Store all flags for reference flow in profile fReferenceFlowFlags.
2464
2465 if(!fReferenceFlowFlags)
2466 {
2467 cout<<endl;
2468 cout<<"WARNING: !fReferenceFlowFlags is NULL in AFAWC::SRFF() !!!!"<<endl;
2469 cout<<endl;
2470 exit(0);
2471 }
2472
2473 // Particle weights used or not:
2474 fReferenceFlowFlags->Fill(0.5,(Double_t)fUsePhiWeights||fUsePtWeights||fUseEtaWeights);
2475 // Which event weight was used to weight generating function event-by-event:
2476 if(strcmp(fMultiplicityWeight->Data(),"unit"))
2477 {
2478 fReferenceFlowFlags->Fill(1.5,0.); // 0 = "unit" (default)
2479 } else if(strcmp(fMultiplicityWeight->Data(),"multiplicity"))
2480 {
2481 fReferenceFlowFlags->Fill(1.5,1.); // 1 = "multiplicity"
2482 }
2483
2484} // end of void AliFlowAnalysisWithCumulants::StoreReferenceFlowFlags()
2485
2486//================================================================================================================
2487
2488void AliFlowAnalysisWithCumulants::StoreDiffFlowFlags()
2489{
2490 // Store all flags for differential flow in profile fDiffFlowFlags.
2491
2492 if(!fDiffFlowFlags)
2493 {
2494 cout<<endl;
2495 cout<<"WARNING: !fDiffFlowFlags is NULL in AFAWC::SRFF() !!!!"<<endl;
2496 cout<<endl;
2497 exit(0);
2498 }
2499
2500 // fDiffFlags->Fill(0.5,(Double_t) ... );
2501
2502} // end of void AliFlowAnalysisWithCumulants::StoreDiffFlowFlags()
2503
2504//================================================================================================================
2505
2506void AliFlowAnalysisWithCumulants::BookAndNestAllLists()
2507{
2508 // Book and nest all list in base list fHistList.
2509
2510 // a) Book and nest lists for reference flow;
2511 // b) Book and nest lists for differential flow;
2512 // c) Book and nest lists for tuning.
fd46c3dd 2513
d6130938 2514 // a) Book and nest all lists for reference flow:
2515 fReferenceFlowList = new TList();
2516 fReferenceFlowList->SetName("Reference Flow");
2517 fReferenceFlowList->SetOwner(kTRUE);
2518 fHistList->Add(fReferenceFlowList);
2519 fReferenceFlowProfiles = new TList();
2520 fReferenceFlowProfiles->SetName("Profiles");
2521 fReferenceFlowProfiles->SetOwner(kTRUE);
2522 fReferenceFlowList->Add(fReferenceFlowProfiles);
2523 fReferenceFlowResults = new TList();
2524 fReferenceFlowResults->SetName("Results");
2525 fReferenceFlowResults->SetOwner(kTRUE);
2526 fReferenceFlowList->Add(fReferenceFlowResults);
2527 // b) Book and nest lists for differential flow:
2528 fDiffFlowList = new TList();
2529 fDiffFlowList->SetName("Differential Flow");
2530 fDiffFlowList->SetOwner(kTRUE);
2531 fHistList->Add(fDiffFlowList);
2532 fDiffFlowProfiles = new TList();
2533 fDiffFlowProfiles->SetName("Profiles");
2534 fDiffFlowProfiles->SetOwner(kTRUE);
2535 fDiffFlowList->Add(fDiffFlowProfiles);
2536 fDiffFlowResults = new TList();
2537 fDiffFlowResults->SetName("Results");
2538 fDiffFlowResults->SetOwner(kTRUE);
2539 fDiffFlowList->Add(fDiffFlowResults);
2540 // c) Book and nest lists for tuning:
2541 if(fTuneParameters)
2542 {
2543 fTuningList = new TList();
2544 fTuningList->SetName("Tuning");
2545 fTuningList->SetOwner(kTRUE);
2546 fHistList->Add(fTuningList);
2547 fTuningProfiles = new TList();
2548 fTuningProfiles->SetName("Profiles");
2549 fTuningProfiles->SetOwner(kTRUE);
2550 fTuningList->Add(fTuningProfiles);
2551 fTuningResults = new TList();
2552 fTuningResults->SetName("Results");
2553 fTuningResults->SetOwner(kTRUE);
2554 fTuningList->Add(fTuningResults);
2555 }
fd46c3dd 2556
d6130938 2557} // end of void AliFlowAnalysisWithCumulants::BookAndNestAllLists()
2558
2559//================================================================================================================
2560
2561void AliFlowAnalysisWithCumulants::BookProfileHoldingSettings()
2562{
2563 // Book profile to hold all analysis settings.
2564
2565 TString analysisSettingsName = "fAnalysisSettings";
2566 fAnalysisSettings = new TProfile(analysisSettingsName.Data(),"Settings for analysis with Generating Function Cumulants",6,0.,6.);
2567 fAnalysisSettings->GetXaxis()->SetLabelSize(0.035);
2568 fAnalysisSettings->GetXaxis()->SetBinLabel(1,"Harmonic");
2569 fAnalysisSettings->Fill(0.5,fHarmonic);
2570 fAnalysisSettings->GetXaxis()->SetBinLabel(2,"Multiple");
2571 fAnalysisSettings->Fill(1.5,fMultiple);
2572 fAnalysisSettings->GetXaxis()->SetBinLabel(3,"r_{0}");
2573 fAnalysisSettings->Fill(2.5,fR0);
2574 fAnalysisSettings->GetXaxis()->SetBinLabel(4,"Print RF results");
2575 fAnalysisSettings->Fill(3.5,fPrintFinalResults[0]);
2576 fAnalysisSettings->GetXaxis()->SetBinLabel(5,"Print RP results");
2577 fAnalysisSettings->Fill(4.5,fPrintFinalResults[1]);
2578 fAnalysisSettings->GetXaxis()->SetBinLabel(6,"Print POI results");
2579 fAnalysisSettings->Fill(5.5,fPrintFinalResults[2]);
2580 fHistList->Add(fAnalysisSettings);
2581
2582} // end of void AliFlowAnalysisWithCumulants::BookProfileHoldingSettings()
2583
2584//================================================================================================================
2585
2586void AliFlowAnalysisWithCumulants::BookCommonHistograms()
2587{
2588 // Book common control histograms and common histograms for final results.
2589
2590 // Common control histogram:
2591 TString commonHistsName = "AliFlowCommonHistGFC";
2592 fCommonHists = new AliFlowCommonHist(commonHistsName.Data());
2593 fHistList->Add(fCommonHists);
2594 // Common histograms for final results from 2nd order GFC:
2595 TString commonHistResults2ndOrderName = "AliFlowCommonHistResults2ndOrderGFC";
2596 fCommonHistsResults2nd = new AliFlowCommonHistResults(commonHistResults2ndOrderName.Data());
2597 fHistList->Add(fCommonHistsResults2nd);
2598 // Common histograms for final results from 4th order GFC:
2599 TString commonHistResults4thOrderName = "AliFlowCommonHistResults4thOrderGFC";
2600 fCommonHistsResults4th = new AliFlowCommonHistResults(commonHistResults4thOrderName.Data());
2601 fHistList->Add(fCommonHistsResults4th);
2602 // Common histograms for final results from 6th order GFC:
2603 TString commonHistResults6thOrderName = "AliFlowCommonHistResults6thOrderGFC";
2604 fCommonHistsResults6th = new AliFlowCommonHistResults(commonHistResults6thOrderName.Data());
2605 fHistList->Add(fCommonHistsResults6th);
2606 // Common histograms for final results from 8th order GFC:
2607 TString commonHistResults8thOrderName = "AliFlowCommonHistResults8thOrderGFC";
2608 fCommonHistsResults8th = new AliFlowCommonHistResults(commonHistResults8thOrderName.Data());
2609 fHistList->Add(fCommonHistsResults8th);
2610
2611} // end of void AliFlowAnalysisWithCumulants::BookCommonHistograms()
2612
2613//================================================================================================================
2614
2615void AliFlowAnalysisWithCumulants::CheckPointersUsedInMake()
2616{
2617 // Check pointers used in method Make().
2618
2619 if(!fCommonHists)
2620 {
2621 cout<<endl;
2622 cout<<" WARNING (GFC): !fCommonHists is NULL in CPUIM() !!!!"<<endl;
2623 cout<<endl;
2624 exit(0);
2625 }
2626 if(fUsePhiWeights && !fPhiWeights)
2627 {
2628 cout<<endl;
2629 cout<<" WARNING (GFC): !fPhiWeights is NULL in CPUIM() !!!!"<<endl;
2630 cout<<endl;
2631 exit(0);
2632 }
2633 if(fUsePtWeights && !fPtWeights)
2634 {
2635 cout<<endl;
2636 cout<<" WARNING (GFC): !fPtWeights is NULL in CPUIM() !!!!"<<endl;
2637 cout<<endl;
2638 exit(0);
2639 }
2640 if(fUseEtaWeights && !fEtaWeights)
2641 {
2642 cout<<endl;
2643 cout<<" WARNING (GFC): !fEtaWeights is NULL in CPUIM() !!!!"<<endl;
2644 cout<<endl;
2645 exit(0);
2646 }
2647 if(!fAverageOfSquaredWeight)
2648 {
2649 cout<<endl;
2650 cout<<" WARNING (GFC): !fAverageOfSquaredWeight is NULL in CPUIM() !!!!"<<endl;
2651 cout<<endl;
2652 exit(0);
2653 }
2654 if(!fReferenceFlowGenFun)
2655 {
2656 cout<<endl;
2657 cout<<" WARNING (GFC): !fReferenceFlowGenFun is NULL in CPUIM() !!!!"<<endl;
2658 cout<<endl;
2659 exit(0);
2660 }
2661 if(!fQvectorComponents)
2662 {
2663 cout<<endl;
2664 cout<<" WARNING (GFC): !fQvectorComponents is NULL in CPUIM() !!!!"<<endl;
2665 cout<<endl;
2666 exit(0);
2667 }
2668 if(!fGEBE)
2669 {
2670 cout<<endl;
2671 cout<<"WARNING (GFC): !fGEBE is NULL in CPUIM() !!!!"<<endl;
2672 cout<<endl;
2673 exit(0);
2674 }
2675
2676} // end of void AliFlowAnalysisWithCumulants::CheckPointersUsedInMake()
2677
2678//================================================================================================================
2679
2680void AliFlowAnalysisWithCumulants::CheckPointersUsedInFinish()
2681{
2682 // Check pointers used in method Finish().
2683
2684 if(!fAnalysisSettings)
2685 {
2686 cout<<endl;
2687 cout<<" WARNING (GFC): fAnalysisSettings is NULL in CPUIF() !!!!"<<endl;
2688 cout<<endl;
2689 exit(0);
2690 }
2691 if(!(fCommonHists && fCommonHists->GetHistMultRP()))
2692 {
2693 cout<<endl;
2694 cout<<" WARNING (GFC): (fCommonHists && fCommonHists->GetHistMultRP) is NULL in CPUIF() !!!!"<<endl;
2695 cout<<endl;
2696 exit(0);
2697 }
2698 if(!fReferenceFlowGenFun)
2699 {
2700 cout<<endl;
2701 cout<<" WARNING (GFC): fReferenceFlowGenFun is NULL in CPUIF() !!!!"<<endl;
2702 cout<<endl;
2703 exit(0);
2704 }
2705 if(!fReferenceFlowCumulants)
2706 {
2707 cout<<endl;
2708 cout<<" WARNING (GFC): fReferenceFlowCumulants is NULL in CPUIF() !!!!"<<endl;
2709 cout<<endl;
2710 exit(0);
2711 }
2712 if(!fQvectorComponents)
2713 {
2714 cout<<endl;
2715 cout<<" WARNING (GFC): fQvectorComponents is NULL in CPUIF() !!!!"<<endl;
2716 cout<<endl;
2717 exit(0);
fd46c3dd 2718 }
d6130938 2719 if(!fAverageOfSquaredWeight)
2720 {
2721 cout<<endl;
2722 cout<<" WARNING (GFC): fAverageOfSquaredWeight is NULL in CPUIF() !!!!"<<endl;
2723 cout<<endl;
2724 exit(0);
2725 }
2726 if(!(fCommonHistsResults2nd && fCommonHistsResults4th && fCommonHistsResults6th && fCommonHistsResults8th))
2727 {
2728 cout<<endl;
2729 cout<<" WARNING (GFC): fCommonHistsResults2nd && fCommonHistsResults4th && fCommonHistsResults6th && "<<endl;
2730 cout<<" fCommonHistsResults8th is NULL in CPUIF() !!!!"<<endl;
2731 cout<<endl;
2732 exit(0);
2733 }
2734 if(!fReferenceFlow)
2735 {
2736 cout<<endl;
2737 cout<<" WARNING (GFC): fReferenceFlow is NULL in CPUIF() !!!!"<<endl;
2738 cout<<endl;
2739 exit(0);
2740 }
2741 if(!fChi)
2742 {
2743 cout<<endl;
2744 cout<<" WARNING (GFC): fChi is NULL in CPUIF() !!!!"<<endl;
2745 cout<<endl;
2746 exit(0);
2747 }
2748 for(Int_t ri=0;ri<2;ri++)
2749 {
2750 for(Int_t rp=0;rp<2;rp++)
2751 {
2752 for(Int_t pe=0;pe<2;pe++)
2753 {
2754 if(!fDiffFlowGenFun[ri][rp][pe])
2755 {
2756 cout<<endl;
2757 cout<<" WARNING (GFC): "<<Form("fDiffFlowGenFun[%d][%d][%d]",ri,rp,pe)<<" is NULL in CPUIF() !!!!"<<endl;
2758 cout<<endl;
2759 exit(0);
2760 }
2761 }
2762 }
2763 }
2764 for(Int_t rp=0;rp<2;rp++)
2765 {
2766 for(Int_t pe=0;pe<2;pe++)
2767 {
2768 for(Int_t co=0;co<4;co++)
2769 {
2770 if(!fDiffFlowCumulants[rp][pe][co])
2771 {
2772 cout<<endl;
2773 cout<<" WARNING (GFC): "<<Form("fDiffFlowCumulants[%d][%d][%d]",rp,pe,co)<<" is NULL in CPUIF() !!!!"<<endl;
2774 cout<<endl;
2775 exit(0);
2776 }
2777 if(!fDiffFlow[rp][pe][co])
2778 {
2779 cout<<endl;
2780 cout<<" WARNING (GFC): "<<Form("fDiffFlow[%d][%d][%d]",rp,pe,co)<<" is NULL in CPUIF() !!!!"<<endl;
2781 cout<<endl;
2782 exit(0);
2783 }
2784 }
2785 }
2786 }
2787 for(Int_t rp=0;rp<2;rp++)
2788 {
2789 for(Int_t pe=0;pe<2;pe++)
2790 {
2791 if(!fNoOfParticlesInBin[rp][pe])
2792 {
2793 cout<<endl;
2794 cout<<" WARNING (GFC): "<<Form("fNoOfParticlesInBin[%d][%d]",rp,pe)<<" is NULL in CPUIF() !!!!"<<endl;
2795 cout<<endl;
2796 exit(0);
2797 }
2798 }
2799 }
2800
2801} // end of void AliFlowAnalysisWithCumulants::CheckPointersUsedInFinish()
fd46c3dd 2802
2803//================================================================================================================
2804
d6130938 2805void AliFlowAnalysisWithCumulants::AccessSettings()
2188af53 2806{
d6130938 2807 // Access the settings for analysis with Generating Function Cumulants.
2808
2809 fHarmonic = (Int_t)fAnalysisSettings->GetBinContent(1);
2810 fMultiple = (Int_t)fAnalysisSettings->GetBinContent(2);
2811 fR0 = (Double_t)fAnalysisSettings->GetBinContent(3);
2812 fPrintFinalResults[0] = (Bool_t)fAnalysisSettings->GetBinContent(4);
2813 fPrintFinalResults[1] = (Bool_t)fAnalysisSettings->GetBinContent(5);
2814 fPrintFinalResults[2] = (Bool_t)fAnalysisSettings->GetBinContent(6);
52021ae2 2815
d6130938 2816} // end of AliFlowAnalysisWithCumulants::AccessSettings()
f1d945a1 2817
1315fe58 2818//================================================================================================================
2819
2820void AliFlowAnalysisWithCumulants::WriteHistograms(TString* outputFileName)
2821{
d6130938 2822 // Store the final results in output .root file.
2823
1315fe58 2824 TFile *output = new TFile(outputFileName->Data(),"RECREATE");
0fe80f88 2825 //output->WriteObject(fHistList, "cobjGFC","SingleKey");
2826 fHistList->SetName("cobjGFC");
9455e15e 2827 fHistList->SetOwner(kTRUE);
0fe80f88 2828 fHistList->Write(fHistList->GetName(), TObject::kSingleKey);
1315fe58 2829 delete output;
1315fe58 2830
d6130938 2831} // end of void AliFlowAnalysisWithCumulants::WriteHistograms(TString* outputFileName)
1315fe58 2832
b0fda271 2833//================================================================================================================
2834
2835void AliFlowAnalysisWithCumulants::WriteHistograms(TString outputFileName)
2836{
d6130938 2837 // Store the final results in output .root file.
2838
b0fda271 2839 TFile *output = new TFile(outputFileName.Data(),"RECREATE");
0fe80f88 2840 //output->WriteObject(fHistList, "cobjGFC","SingleKey");
2841 fHistList->SetName("cobjGFC");
9455e15e 2842 fHistList->SetOwner(kTRUE);
0fe80f88 2843 fHistList->Write(fHistList->GetName(), TObject::kSingleKey);
b0fda271 2844 delete output;
b0fda271 2845
d6130938 2846} // end of void AliFlowAnalysisWithCumulants::WriteHistograms(TString outputFileName)
b0fda271 2847
ad87ae62 2848//================================================================================================================
2849
2850void AliFlowAnalysisWithCumulants::WriteHistograms(TDirectoryFile *outputFileName)
2851{
d6130938 2852 // Store the final results in output .root file.
2853
ad87ae62 2854 fHistList->SetName("cobjGFC");
2855 fHistList->SetOwner(kTRUE);
2856 outputFileName->Add(fHistList);
2857 outputFileName->Write(outputFileName->GetName(), TObject::kSingleKey);
d6130938 2858
2859} // end of void AliFlowAnalysisWithCumulants::WriteHistograms(TDirectoryFile *outputFileName)
f1d945a1 2860
ad87ae62 2861//================================================================================================================
f1d945a1 2862
2863
d6130938 2864
2865