PWG2/FLOW/AliFlowCommon/AliFlowAnalysisWithCumulants.cxx

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