PWG2/FLOW/AliFlowCommon/AliFlowAnalysisWithCumulants.cxx

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