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