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