PWG/FLOW/Base/AliFlowAnalysisWithNestedLoops.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  * Only in this class nested loops are used for flow analysis. *
  20  * Nested loops are used to evaluate:                          *
  21  *                                                             *
  22  *  a) Distribution of relative angle difference (phi1-phi2);  *
  23  *  b) Cross-check the results for mixed harmonics.            *
  24  *                                                             *
  25  *       Author: Ante Bilandzic (abilandzic@gmail.com)         *
  26  ***************************************************************/
  27
  28 #define AliFlowAnalysisWithNestedLoops_cxx
  29
  30 #include "Riostream.h"
  31 #include "AliFlowCommonConstants.h"
  32 #include "AliFlowCommonHist.h"
  33 #include "AliFlowCommonHistResults.h"
  34
  35 #include "TMath.h"
  36 #include "TFile.h"
  37 #include "TList.h"
  38 #include "TProfile.h"
  39
  40 #include "AliFlowEventSimple.h"
  41 #include "AliFlowTrackSimple.h"
  42 #include "AliFlowAnalysisWithNestedLoops.h"
  43
  44 class TH1;
  45 class TList;
  46
  47 using std::endl;
  48 using std::cout;
  49 ClassImp(AliFlowAnalysisWithNestedLoops)
  50
  51 //================================================================================================================
  52
  53 AliFlowAnalysisWithNestedLoops::AliFlowAnalysisWithNestedLoops():
  54 fHistList(NULL),
  55 fHistListName(NULL),
  56 fHarmonic(0),
  57 fAnalysisLabel(NULL),
  58 fAnalysisSettings(NULL),
  59 fOppositeChargesPOI(kFALSE),
  60 fEvaluateDifferential3pCorrelator(kFALSE),
  61 fPrintOnTheScreen(kTRUE),
  62 fCommonHists(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 fWeightsList(NULL),
  76 fUsePhiWeights(kFALSE),
  77 fUsePtWeights(kFALSE),
  78 fUseEtaWeights(kFALSE),
  79 fUseParticleWeights(NULL),
  80 fPhiWeights(NULL),
  81 fPtWeights(NULL),
  82 fEtaWeights(NULL),
  83 fListRAD(NULL),
  84 fEvaluateNestedLoopsForRAD(kTRUE),
  85 fRelativeAngleDistribution(NULL),
  86 fCharge(NULL),
  87 fListQC(NULL),
  88 fEvaluateNestedLoopsForQC(kFALSE),
  89 fListMH(NULL),
  90 fEvaluateNestedLoopsForMH(kFALSE),
  91 f3pCorrelatorPro(NULL),
  92 f5pCorrelatorPro(NULL)
  93 {
  94  // Constructor.
  95
  96  // Base list to hold all output objects:
  97  fHistList = new TList();
  98  fHistListName = new TString("cobjNL");
  99  fHistList->SetName(fHistListName->Data());
 100  fHistList->SetOwner(kTRUE);
 101
 102  // List to hold histograms with phi, pt and eta weights:
 103  fWeightsList = new TList();
 104
 105  // List to hold objects relevant for relative angle distributions:
 106  fListRAD = new TList();
 107
 108  // List holding objects relevant for debugging and cross-checking of Q-cumulants class:
 109  fListQC = new TList();
 110
 111  // List holding objects relevant for debugging and cross-checking of MH class:
 112  fListMH = new TList();
 113
 114  // Initialize all arrays:
 115  this->InitializeArraysForMH();
 116
 117 } // AliFlowAnalysisWithNestedLoops::AliFlowAnalysisWithNestedLoops()
 118
 119 //================================================================================================================
 120
 121 AliFlowAnalysisWithNestedLoops::~AliFlowAnalysisWithNestedLoops()
 122 {
 123  // Destructor.
 124
 125  delete fHistList;
 126
 127 } // end of AliFlowAnalysisWithNestedLoops::~AliFlowAnalysisWithNestedLoops()
 128
 129 //================================================================================================================
 130
 131 void AliFlowAnalysisWithNestedLoops::Init()
 132 {
 133  // Initialize and book all objects.
 134
 135  // a) Cross check if the user settings make sense before starting;
 136  // b) Access all common constants;
 137  // c) Book and nest all lists in the base list fHistList;
 138  // d) Book profile holding seetings for analysis with nested loops;
 139  // e) Book common control histograms;
 140  // f) Book all objects relevant for distributions;
 141  // g) Book and fill histograms to hold phi, pt and eta weights;
 142  // h) Store harmonic n.
 143
 144
 145  //save old value and prevent histograms from being added to directory
 146  //to avoid name clashes in case multiple analaysis objects are used
 147  //in an analysis
 148  Bool_t oldHistAddStatus = TH1::AddDirectoryStatus();
 149  TH1::AddDirectory(kFALSE);
 150
 151  TH1::SetDefaultSumw2();
 152
 153  this->CrossCheckSettings();
 154  this->AccessConstants();
 155  this->BookAndNestAllLists();
 156  this->BookAndFillProfileHoldingSettings();
 157  this->BookCommonHistograms();
 158  this->BookEverythingForRAD();
 159  this->BookEverythingForMH();
 160  this->BookAndFillWeightsHistograms();
 161  this->StoreHarmonic();
 162
 163  //restore old status
 164  TH1::AddDirectory(oldHistAddStatus);
 165 } // end of void AliFlowAnalysisWithNestedLoops::Init()
 166
 167 //================================================================================================================
 168
 169 void AliFlowAnalysisWithNestedLoops::Make(AliFlowEventSimple* anEvent)
 170 {
 171  // Running over data only in this method.
 172
 173  // a) Check all pointers used in this method;
 174  // b) Fill common control histograms;
 175  // c) Evaluate nested loops for relative angle distribution;
 176  // d) Evaluate nested loops for mixed harmonics.
 177
 178  this->CheckPointersUsedInMake();
 179  fCommonHists->FillControlHistograms(anEvent);
 180  if(fEvaluateNestedLoopsForRAD) this->EvaluateNestedLoopsForRAD(anEvent);
 181  if(fEvaluateNestedLoopsForMH) this->EvaluateNestedLoopsForMH(anEvent);
 182
 183 } // end of AliFlowAnalysisWithNestedLoops::Make(AliFlowEventSimple* anEvent)
 184
 185 //================================================================================================================
 186
 187 void AliFlowAnalysisWithNestedLoops::Finish()
 188 {
 189  // Calculate the final results.
 190
 191  // a) Access settings for analysis with mixed harmonics;
 192  // b) Check all pointers used in this method;
 193  // c) Print on the screen.
 194
 195  this->AccessSettings();
 196  this->CheckPointersUsedInFinish();
 197  if(fPrintOnTheScreen) this->PrintOnTheScreen();
 198
 199 } // end of AliFlowAnalysisWithNestedLoops::Finish()
 200
 201 //================================================================================================================
 202
 203 void AliFlowAnalysisWithNestedLoops::GetOutputHistograms(TList *outputListHistos)
 204 {
 205  // Get pointers to all objects saved in the output file.
 206
 207  // a) Get pointers for common control histograms.
 208  if(outputListHistos)
 209  {
 210   this->SetHistList(outputListHistos);
 211   if(!fHistList)
 212   {
 213    cout<<endl;
 214    cout<<" WARNING (NL): fHistList is NULL in AFAWNL::GOH() !!!!"<<endl;
 215    cout<<endl;
 216    exit(0);
 217   }
 218   this->GetPointersForBaseHistograms();
 219   this->GetPointersForCommonHistograms();
 220   this->AccessSettings();
 221   if(fEvaluateNestedLoopsForRAD) this->GetPointersForRAD();
 222   if(fEvaluateNestedLoopsForMH) this->GetPointersForMH();
 223  } else
 224    {
 225     cout<<endl;
 226     cout<<" WARNING (NL): outputListHistos is NULL in AFAWNL::GOH() !!!!"<<endl;
 227     cout<<endl;
 228     exit(0);
 229    }
 230
 231 } // end of void AliFlowAnalysisWithNestedLoops::GetOutputHistograms(TList *outputListHistos)
 232
 233 //================================================================================================================
 234
 235 void AliFlowAnalysisWithNestedLoops::GetPointersForBaseHistograms()
 236 {
 237  // Get pointers to base histograms.
 238
 239  TString analysisSettingsName = "fAnalysisSettings";
 240  TProfile *analysisSettings = dynamic_cast<TProfile*>(fHistList->FindObject(analysisSettingsName.Data()));
 241  if(analysisSettings)
 242  {
 243   this->SetAnalysisSettings(analysisSettings);
 244  } else
 245    {
 246     cout<<endl;
 247     cout<<" WARNING (NL): analysisSettings is NULL in AFAWNL::GPFBH() !!!!"<<endl;
 248     cout<<endl;
 249     exit(0);
 250    }
 251
 252 } // end of void AliFlowAnalysisWithNestedLoops::GetPointersForBaseHistograms()
 253
 254 //================================================================================================================
 255
 256 void AliFlowAnalysisWithNestedLoops::GetPointersForCommonHistograms()
 257 {
 258  // Get pointers to common control histograms.
 259
 260  TString commonHistsName = "AliFlowCommonHistNL";
 261  AliFlowCommonHist *commonHist = dynamic_cast<AliFlowCommonHist*>(fHistList->FindObject(commonHistsName.Data()));
 262  if(commonHist)
 263  {
 264   this->SetCommonHists(commonHist);
 265  } else
 266    {
 267     cout<<endl;
 268     cout<<" WARNING (NL): commonHist is NULL in AFAWNL::GPFCH() !!!!"<<endl;
 269     cout<<endl;
 270     exit(0);
 271    }
 272
 273 } // end of void AliFlowAnalysisWithNestedLoops::GetPointersForCommonHistograms()
 274
 275 //================================================================================================================
 276
 277 void AliFlowAnalysisWithNestedLoops::GetPointersForRAD()
 278 {
 279  // Get pointers to objects relevant for relative angle distributions.
 280
 281  TList *listRAD = NULL;
 282  listRAD = dynamic_cast<TList*>(fHistList->FindObject("Relative Angle Distribution"));
 283  if(!listRAD)
 284  {
 285   cout<<"WARNING: listRAD is NULL in AFAWNL::GPFRAD() !!!!"<<endl;
 286   exit(0);
 287  }
 288
 289  TString relativeAngleDistributionName = "fRelativeAngleDistribution";
 290  TH1D *relativeAngleDistribution = dynamic_cast<TH1D*>(listRAD->FindObject(relativeAngleDistributionName.Data()));
 291  if(relativeAngleDistribution)
 292  {
 293   this->SetRelativeAngleDistribution(relativeAngleDistribution);
 294  }
 295
 296  TString chargeName = "fCharge";
 297  TH1D *charge = dynamic_cast<TH1D*>(listRAD->FindObject(chargeName.Data()));
 298  if(charge)
 299  {
 300   this->SetCharge(charge);
 301  }
 302
 303 } // end of void AliFlowAnalysisWithNestedLoops::GetPointersForRAD()
 304
 305 //================================================================================================================
 306
 307 void AliFlowAnalysisWithNestedLoops::GetPointersForMH()
 308 {
 309  // Get pointers to objects evaluated with nested loops.
 310
 311  TList *listMH = NULL;
 312  listMH = dynamic_cast<TList*>(fHistList->FindObject("Mixed Harmonics"));
 313  if(!listMH)
 314  {
 315   cout<<"WARNING: listMH is NULL in AFAWNL::GPFMH() !!!!"<<endl;
 316   exit(0);
 317  }
 318
 319  TString s3pCorrelatorProName = "f3pCorrelatorPro";
 320  TProfile *p3pCorrelatorPro = dynamic_cast<TProfile*>(listMH->FindObject(s3pCorrelatorProName.Data()));
 321  if(p3pCorrelatorPro)
 322  {
 323   this->Set3pCorrelatorPro(p3pCorrelatorPro);
 324  }
 325  TString s5pCorrelatorProName = "f5pCorrelatorPro";
 326  TProfile *p5pCorrelatorPro = dynamic_cast<TProfile*>(listMH->FindObject(s5pCorrelatorProName.Data()));
 327  if(p5pCorrelatorPro)
 328  {
 329   this->Set5pCorrelatorPro(p5pCorrelatorPro);
 330  }
 331  if(!fEvaluateDifferential3pCorrelator){return;}
 332  TString psdFlag[2] = {"PtSum","PtDiff"};
 333  for(Int_t sd=0;sd<2;sd++)
 334  {
 335   TProfile *p3pCorrelatorVsPtSumDiffDirectPro = dynamic_cast<TProfile*>(listMH->FindObject(Form("f3pCorrelatorDirectVs%s",psdFlag[sd].Data())));
 336   if(p3pCorrelatorVsPtSumDiffDirectPro)
 337   {
 338    this->Set3pCorrelatorVsPtSumDiffDirectPro(p3pCorrelatorVsPtSumDiffDirectPro,sd);
 339   }
 340  } // end of for(Int_t sd=0;sd<2;sd++)
 341
 342 } // end of void AliFlowAnalysisWithNestedLoops::GetPointersForMH()
 343
 344 //================================================================================================================
 345
 346 void AliFlowAnalysisWithNestedLoops::WriteHistograms(TString outputFileName)
 347 {
 348  // Store the final results in output .root file.
 349  TFile *output = new TFile(outputFileName.Data(),"RECREATE");
 350  fHistList->Write(fHistList->GetName(),TObject::kSingleKey);
 351  delete output;
 352 }
 353
 354 //================================================================================================================
 355
 356 void AliFlowAnalysisWithNestedLoops::WriteHistograms(TDirectoryFile *outputFileName)
 357 {
 358  // Store the final results in output .root file.
 359  fHistList->SetName("cobjNL");
 360  fHistList->SetOwner(kTRUE);
 361  outputFileName->Add(fHistList);
 362  outputFileName->Write(outputFileName->GetName(),TObject::kSingleKey);
 363 }
 364
 365 //================================================================================================================
 366
 367 void AliFlowAnalysisWithNestedLoops::StoreHarmonic()
 368 {
 369  // Store harmonic n used in cos[n*(phi1+phi2-2phi3)] and cos[n*(psi1+psi2-2phi3)].
 370
 371  (fCommonHists->GetHarmonic())->Fill(0.5,fHarmonic);
 372
 373 } // end of void AliFlowAnalysisWithNestedLoops::StoreHarmonic()
 374
 375 //================================================================================================================
 376
 377 void AliFlowAnalysisWithNestedLoops::BookAndNestAllLists()
 378 {
 379  // Book and nest all list in base list fHistList.
 380
 381  // Weights:
 382  fWeightsList->SetName("Weights");
 383  fWeightsList->SetOwner(kTRUE);
 384  fHistList->Add(fWeightsList);
 385  // List for Relative Angle Distribution:
 386  fListRAD->SetName("Relative Angle Distribution");
 387  fListRAD->SetOwner(kTRUE);
 388  if(fEvaluateNestedLoopsForRAD) fHistList->Add(fListRAD);
 389  // List for Q-cumulants:
 390  fListQC->SetName("Q-cumulants");
 391  fListQC->SetOwner(kTRUE);
 392  if(fEvaluateNestedLoopsForQC) fHistList->Add(fListQC);
 393  // List for Mixed Harmonics:
 394  fListMH->SetName("Mixed Harmonics");
 395  fListMH->SetOwner(kTRUE);
 396  if(fEvaluateNestedLoopsForMH) fHistList->Add(fListMH);
 397
 398 } // end of void AliFlowAnalysisWithNestedLoops::BookAndNestAllLists()
 399
 400 //================================================================================================================
 401
 402 void AliFlowAnalysisWithNestedLoops::BookAndFillProfileHoldingSettings()
 403 {
 404  // Book profile to hold all analysis settings.
 405
 406  TString analysisSettingsName = "fAnalysisSettings";
 407  fAnalysisSettings = new TProfile(analysisSettingsName.Data(),"Settings for analysis with nested loops",7,0,7);
 408  fAnalysisSettings->GetXaxis()->SetLabelSize(0.035);
 409  fAnalysisSettings->GetXaxis()->SetBinLabel(1,"Nested loops for RAD?");
 410  fAnalysisSettings->Fill(0.5,(Int_t)fEvaluateNestedLoopsForRAD);
 411  fAnalysisSettings->GetXaxis()->SetBinLabel(2,"Nested loops for QC?");
 412  fAnalysisSettings->Fill(1.5,(Int_t)fEvaluateNestedLoopsForQC);
 413  fAnalysisSettings->GetXaxis()->SetBinLabel(3,"Nested loops for MH?");
 414  fAnalysisSettings->Fill(2.5,(Int_t)fEvaluateNestedLoopsForMH);
 415  fAnalysisSettings->GetXaxis()->SetBinLabel(4,"fHarmonic");
 416  fAnalysisSettings->Fill(3.5,(Int_t)fHarmonic);
 417  fAnalysisSettings->GetXaxis()->SetBinLabel(5,"Print on the screen?");
 418  fAnalysisSettings->Fill(4.5,(Int_t)fPrintOnTheScreen);
 419  fAnalysisSettings->GetXaxis()->SetBinLabel(6,"fOppositeChargesPOI");
 420  fAnalysisSettings->Fill(5.5,fOppositeChargesPOI);
 421  fAnalysisSettings->GetXaxis()->SetBinLabel(7,"fEvaluateDifferential3pCorrelator");
 422  fAnalysisSettings->Fill(6.5,fEvaluateDifferential3pCorrelator);
 423  fHistList->Add(fAnalysisSettings);
 424
 425 } // end of void AliFlowAnalysisWithNestedLoops::BookAndFillProfileHoldingSettings()
 426
 427 //================================================================================================================
 428
 429 void AliFlowAnalysisWithNestedLoops::BookCommonHistograms()
 430 {
 431  // Book common control histograms and common histograms for final results.
 432
 433  TString commonHistsName = "AliFlowCommonHistNL";
 434  fCommonHists = new AliFlowCommonHist(commonHistsName.Data());
 435  fHistList->Add(fCommonHists);
 436
 437 } // end of void AliFlowAnalysisWithNestedLoops::BookCommonHistograms()
 438
 439 //================================================================================================================
 440
 441 void AliFlowAnalysisWithNestedLoops::BookEverythingForRAD()
 442 {
 443  // Book all objects relevant calculation of relative angle distribution.
 444
 445  TString relativeAngleDistributionName = "fRelativeAngleDistribution";
 446  fRelativeAngleDistribution = new TH1D(relativeAngleDistributionName.Data(),"Relative angle distribution",720,-TMath::TwoPi(),TMath::TwoPi());
 447  fRelativeAngleDistribution->GetYaxis()->SetTitle("#frac{dN}{#Delta #phi}");
 448  fRelativeAngleDistribution->GetXaxis()->SetTitle("#Delta #phi");
 449  fListRAD->Add(fRelativeAngleDistribution);
 450
 451  TString chargeName = "fCharge";
 452  fCharge = new TH1D(chargeName.Data(),"Charges",3,0,3);
 453
 454  fCharge->GetXaxis()->SetBinLabel(1,"+");
 455  fCharge->GetXaxis()->SetBinLabel(2,"-");
 456  fCharge->GetXaxis()->SetBinLabel(3,"rest");
 457  fListRAD->Add(fCharge);
 458
 459 } // end fo void AliFlowAnalysisWithNestedLoops::BookEverythingForRAD()
 460
 461 //================================================================================================================
 462
 463 void AliFlowAnalysisWithNestedLoops::AccessConstants()
 464 {
 465  // Access needed common constants from AliFlowCommonConstants.
 466
 467  fnBinsPhi = AliFlowCommonConstants::GetMaster()->GetNbinsPhi();
 468  fPhiMin = AliFlowCommonConstants::GetMaster()->GetPhiMin();
 469  fPhiMax = AliFlowCommonConstants::GetMaster()->GetPhiMax();
 470  if(fnBinsPhi) fPhiBinWidth = (fPhiMax-fPhiMin)/fnBinsPhi;
 471  fnBinsPt = AliFlowCommonConstants::GetMaster()->GetNbinsPt();
 472  fPtMin = AliFlowCommonConstants::GetMaster()->GetPtMin();
 473  fPtMax = AliFlowCommonConstants::GetMaster()->GetPtMax();
 474  if(fnBinsPt) fPtBinWidth = (fPtMax-fPtMin)/fnBinsPt;
 475  fnBinsEta = AliFlowCommonConstants::GetMaster()->GetNbinsEta();
 476  fEtaMin = AliFlowCommonConstants::GetMaster()->GetEtaMin();
 477  fEtaMax = AliFlowCommonConstants::GetMaster()->GetEtaMax();
 478  if(fnBinsEta) fEtaBinWidth = (fEtaMax-fEtaMin)/fnBinsEta;
 479
 480 } // end of void AliFlowAnalysisWithNestedLoops::AccessConstants()
 481
 482 //================================================================================================================
 483
 484 void AliFlowAnalysisWithNestedLoops::CrossCheckSettings()
 485 {
 486  // Cross-check if the user settings make sense.
 487
 488  // ...
 489
 490 } // end of void AliFlowAnalysisWithNestedLoops::CrossCheckSettings()
 491
 492 //================================================================================================================
 493
 494 void AliFlowAnalysisWithNestedLoops::BookAndFillWeightsHistograms()
 495 {
 496  // Book and fill (by accessing file "weights.root") histograms which hold phi, pt and eta weights.
 497
 498  if(!fWeightsList)
 499  {
 500   cout<<"WARNING: fWeightsList is NULL in AFAWNL::BAFWH() !!!!"<<endl;
 501   exit(0);
 502  }
 503  // Profile to hold flags for weights:
 504  TString fUseParticleWeightsName = "fUseParticleWeightsNL";
 505  fUseParticleWeights = new TProfile(fUseParticleWeightsName.Data(),"0 = particle weight not used, 1 = particle weight used ",3,0,3);
 506  fUseParticleWeights->SetLabelSize(0.06);
 507  (fUseParticleWeights->GetXaxis())->SetBinLabel(1,"w_{#phi}");
 508  (fUseParticleWeights->GetXaxis())->SetBinLabel(2,"w_{p_{T}}");
 509  (fUseParticleWeights->GetXaxis())->SetBinLabel(3,"w_{#eta}");
 510  fUseParticleWeights->Fill(0.5,(Int_t)fUsePhiWeights);
 511  fUseParticleWeights->Fill(1.5,(Int_t)fUsePtWeights);
 512  fUseParticleWeights->Fill(2.5,(Int_t)fUseEtaWeights);
 513  fWeightsList->Add(fUseParticleWeights);
 514  // Phi-weights:
 515  if(fUsePhiWeights)
 516  {
 517   if(fWeightsList->FindObject("phi_weights"))
 518   {
 519    fPhiWeights = dynamic_cast<TH1F*>(fWeightsList->FindObject("phi_weights"));
 520    if (!fPhiWeights)
 521    {
 522      printf("WARNING: no phi weights object, bye!\n");
 523      exit(0);
 524    }
 525    if(TMath::Abs(fPhiWeights->GetBinWidth(1)-fPhiBinWidth)>pow(10.,-6.))
 526    {
 527     cout<<endl;
 528     cout<<"WARNING (NL): Inconsistent binning in histograms for phi-weights throughout the code."<<endl;
 529     cout<<endl;
 530     exit(0);
 531    }
 532   } else
 533     {
 534      cout<<"WARNING (NL): fWeightsList->FindObject(\"phi_weights\") is NULL in AFAWNL::BAFWH() !!!!"<<endl;
 535      exit(0);
 536     }
 537  } // end of if(fUsePhiWeights)
 538  // Pt-weights:
 539  if(fUsePtWeights)
 540  {
 541   if(fWeightsList->FindObject("pt_weights"))
 542   {
 543    fPtWeights = dynamic_cast<TH1D*>(fWeightsList->FindObject("pt_weights"));
 544    if (!fPtWeights)
 545    {
 546      printf("WARNING: no pt weights object, bye!\n");
 547      exit(0);
 548    }
 549    if(TMath::Abs(fPtWeights->GetBinWidth(1)-fPtBinWidth)>pow(10.,-6.))
 550    {
 551     cout<<endl;
 552     cout<<"WARNING (NL): Inconsistent binning in histograms for pt-weights throughout the code."<<endl;
 553     cout<<endl;
 554     exit(0);
 555    }
 556   } else
 557     {
 558      cout<<"WARNING (NL): fWeightsList->FindObject(\"pt_weights\") is NULL in AFAWNL::BAFWH() !!!!"<<endl;
 559      exit(0);
 560     }
 561  } // end of if(fUsePtWeights)
 562  // Eta-weights:
 563  if(fUseEtaWeights)
 564  {
 565   if(fWeightsList->FindObject("eta_weights"))
 566   {
 567    fEtaWeights = dynamic_cast<TH1D*>(fWeightsList->FindObject("eta_weights"));
 568    if (!fEtaWeights)
 569    {
 570      printf("WARNING: no eta weights object, bye!\n");
 571      exit(0);
 572    }
 573    if(TMath::Abs(fEtaWeights->GetBinWidth(1)-fEtaBinWidth)>pow(10.,-6.))
 574    {
 575     cout<<endl;
 576     cout<<"WARNING (NL): Inconsistent binning in histograms for eta-weights throughout the code."<<endl;
 577     cout<<endl;
 578     exit(0);
 579    }
 580   } else
 581     {
 582      cout<<"WARNING: fUseEtaWeights && fWeightsList->FindObject(\"eta_weights\") is NULL in AFAWNL::BAFWH() !!!!"<<endl;
 583      exit(0);
 584     }
 585  } // end of if(fUseEtaWeights)
 586
 587 } // end of AliFlowAnalysisWithNestedLoops::BookAndFillWeightsHistograms()
 588
 589 //================================================================================================================
 590
 591 void AliFlowAnalysisWithNestedLoops::CheckPointersUsedInMake()
 592 {
 593  // Check pointers used in method Make().
 594
 595  if(fEvaluateNestedLoopsForRAD) CheckPointersForRAD("Make");
 596  if(fEvaluateNestedLoopsForMH) CheckPointersForMH("Make");
 597
 598 } // end of AliFlowAnalysisWithNestedLoops::CheckPointersUsedInMake()
 599
 600 //================================================================================================================
 601
 602 void AliFlowAnalysisWithNestedLoops::CheckPointersUsedInFinish()
 603 {
 604  // Check pointers used in method Finish().
 605
 606  if(fEvaluateNestedLoopsForRAD) CheckPointersForRAD("Finish");
 607  if(fEvaluateNestedLoopsForMH) CheckPointersForMH("Finish");
 608
 609 } // end of AliFlowAnalysisWithNestedLoops::CheckPointersUsedInFinish()
 610
 611 //================================================================================================================
 612
 613 void AliFlowAnalysisWithNestedLoops::CheckPointersForRAD(TString where)
 614 {
 615  // Check pointers relevant for calculation of relative angle distribution.
 616
 617  if(!fRelativeAngleDistribution)
 618  {
 619   cout<<endl;
 620   cout<<" WARNING (NL): fRelativeAngleDistribution is NULL in "<<where.Data()<<"() !!!!"<<endl;
 621   cout<<endl;
 622   exit(0);
 623  }
 624
 625  if(strcmp(where.Data(),"Make") == 0)
 626  {
 627   // Check pointers used only in method Make():
 628   // ...
 629  }
 630  else if(strcmp(where.Data(),"Finish") == 0)
 631  {
 632   // Check pointers used only in method Finish():
 633   // ...
 634  }
 635
 636 } // end of void AliFlowAnalysisWithNestedLoops::CheckPointersForRAD(TString where)
 637
 638 //================================================================================================================
 639
 640 void AliFlowAnalysisWithNestedLoops::CheckPointersForMH(TString where)
 641 {
 642  // Check pointers relevant for calculation of mixed harmonics.
 643
 644  if(strcmp(where.Data(),"Make") == 0)
 645  {
 646   // Check pointers used only in method Make():
 647   if(!f3pCorrelatorPro)
 648   {
 649    cout<<endl;
 650    cout<<" WARNING (NL): f3pCorrelatorPro is NULL in Make() !!!!"<<endl;
 651    cout<<endl;
 652    exit(0);
 653   }
 654   if(!f5pCorrelatorPro)
 655   {
 656    cout<<endl;
 657    cout<<" WARNING (NL): f5pCorrelatorPro is NULL in Make() !!!!"<<endl;
 658    cout<<endl;
 659    exit(0);
 660   }
 661   if(!fEvaluateDifferential3pCorrelator){return;}
 662   for(Int_t sd=0;sd<2;sd++)
 663   {
 664    if(!(f3pCorrelatorVsPtSumDiffDirectPro[sd]))
 665    {
 666     cout<<endl;
 667     cout<<" WARNING (NL): "<<Form("f3pCorrelatorVsPtSumDiffDirectPro[%d]",sd)<<" is NULL in Make() !!!!"<<endl;
 668     cout<<endl;
 669     exit(0);
 670    }
 671   } // end of for(Int_t sd=0;sd<2;sd++)
 672  } // if(strcmp(where.Data(),"Make") == 0)
 673  else if(strcmp(where.Data(),"Finish") == 0)
 674  {
 675   // Check pointers used only in method Finish():
 676   if(!f3pCorrelatorPro)
 677   {
 678    cout<<endl;
 679    cout<<" WARNING (NL): f3pCorrelatorPro is NULL in Finish() !!!!"<<endl;
 680    cout<<endl;
 681    exit(0);
 682   }
 683   if(!f5pCorrelatorPro)
 684   {
 685    cout<<endl;
 686    cout<<" WARNING (NL): f5pCorrelatorPro is NULL in Finish() !!!!"<<endl;
 687    cout<<endl;
 688    exit(0);
 689   }
 690  } // else if(strcmp(where.Data(),"Finish") == 0)
 691
 692 } // end of void AliFlowAnalysisWithNestedLoops::CheckPointersForMH(TString where)
 693
 694 //================================================================================================================
 695
 696 void AliFlowAnalysisWithNestedLoops::AccessSettings()
 697 {
 698  // Access the settings for analysis.
 699
 700  fEvaluateNestedLoopsForRAD = (Bool_t)fAnalysisSettings->GetBinContent(1);
 701  fEvaluateNestedLoopsForQC = (Bool_t)fAnalysisSettings->GetBinContent(2);
 702  fEvaluateNestedLoopsForMH = (Bool_t)fAnalysisSettings->GetBinContent(3);
 703  fHarmonic = (Int_t)fAnalysisSettings->GetBinContent(4);
 704  fPrintOnTheScreen = (Bool_t)fAnalysisSettings->GetBinContent(5);
 705  fOppositeChargesPOI = (Bool_t)fAnalysisSettings->GetBinContent(6);
 706  fEvaluateDifferential3pCorrelator = (Bool_t)fAnalysisSettings->GetBinContent(7);
 707
 708 } // end of AliFlowAnalysisWithNestedLoops::AccessSettings()
 709
 710 //================================================================================================================
 711
 712 void AliFlowAnalysisWithNestedLoops::InitializeArraysForMH()
 713 {
 714  // Initialize arrays mixed harmonics calculations.
 715
 716  for(Int_t sd=0;sd<2;sd++) // sum or difference
 717  {
 718   f3pCorrelatorVsPtSumDiffDirectPro[sd] = NULL;
 719  }
 720
 721 } // end of AliFlowAnalysisWithNestedLoops::InitializeArraysForMH()
 722
 723 //================================================================================================================
 724
 725 void AliFlowAnalysisWithNestedLoops::BookEverythingForMH()
 726 {
 727  // Book all objects relevant for mixed harmonics.
 728
 729  if(fEvaluateNestedLoopsForMH)
 730  {
 731   // 3-p:
 732   TString s3pCorrelatorProName = "f3pCorrelatorPro";
 733   f3pCorrelatorPro = new TProfile(s3pCorrelatorProName.Data(),"",1,0,1);
 734   f3pCorrelatorPro->SetStats(kFALSE);
 735   f3pCorrelatorPro->GetXaxis()->SetLabelOffset(0.01);
 736   f3pCorrelatorPro->GetXaxis()->SetLabelSize(0.05);
 737   if(fHarmonic == 1)
 738   {
 739    f3pCorrelatorPro->GetXaxis()->SetBinLabel(1,"#LT#LTcos(#phi_{1}+#phi_{2}-2#phi_{3})#GT#GT");
 740   } else
 741     {
 742      f3pCorrelatorPro->GetXaxis()->SetBinLabel(1,Form("#LT#LTcos[%i(#phi_{1}+#phi_{2}-2#phi_{3})]#GT#GT",fHarmonic));
 743     }
 744   fListMH->Add(f3pCorrelatorPro);
 745   // 5-p:
 746   TString s5pCorrelatorProName = "f5pCorrelatorPro";
 747   f5pCorrelatorPro = new TProfile(s5pCorrelatorProName.Data(),"",1,0,1);
 748   f5pCorrelatorPro->SetStats(kFALSE);
 749   f5pCorrelatorPro->GetXaxis()->SetLabelOffset(0.01);
 750   f5pCorrelatorPro->GetXaxis()->SetLabelSize(0.05);
 751   if(fHarmonic == 1)
 752   {
 753    f5pCorrelatorPro->GetXaxis()->SetBinLabel(1,"#LT#LTcos(2#phi_{1}+2#phi_{2}+2#phi_{3}-3#phi_{4}-3#phi_{5})#GT#GT");
 754   } else
 755     {
 756      f5pCorrelatorPro->GetXaxis()->SetBinLabel(1,Form("#LT#LTcos[%i(2#phi_{1}+2#phi_{2}+2#phi_{3}-3#phi_{4}-3#phi_{5})]#GT#GT",fHarmonic));
 757     }
 758   fListMH->Add(f5pCorrelatorPro);
 759   if(fEvaluateDifferential3pCorrelator)
 760   {
 761    TString psdFlag[2] = {"PtSum","PtDiff"};
 762    TString psdTitleFlag[2] = {"(p_{t,1}+p_{t,2})/2","#left|p_{t,1}-p_{t,2}#right|"};
 763    //TString s3pCorrelatorVsPtSumDiffDirectProName = "f3pCorrelatorVsPtSumDiffDirectPro";
 764    for(Int_t sd=0;sd<2;sd++)
 765    {
 766     // to be improved: hardwired ,fnBinsPt,0.,fPtMax):
 767     f3pCorrelatorVsPtSumDiffDirectPro[sd] = new TProfile(Form("f3pCorrelatorDirectVs%s",psdFlag[sd].Data()),"",fnBinsPt,0.,fPtMax);
 768     //f3pCorrelatorVsPtSumDiffDirectPro[sd]->SetLabelSize(0.05);
 769     //f3pCorrelatorVsPtSumDiffDirectPro[sd]->SetMarkerStyle(25);
 770     f3pCorrelatorVsPtSumDiffDirectPro[sd]->GetXaxis()->SetTitle(psdTitleFlag[sd].Data());
 771     fListMH->Add(f3pCorrelatorVsPtSumDiffDirectPro[sd]);
 772    }
 773   } // end of if(fEvaluateDifferential3pCorrelator)
 774  } // end of if(fEvaluateNestedLoopsForMH)
 775
 776 } // end of AliFlowAnalysisWithNestedLoops::BookEverythingForMH()
 777
 778 //================================================================================================================
 779
 780 void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForRAD(AliFlowEventSimple *anEvent)
 781 {
 782  // Evaluate nested loops needed for calculation of relative angle distribution.
 783
 784  Double_t dPhi1=0., dPhi2=0.; // azimuthal angles in the laboratory frame
 785  AliFlowTrackSimple *aftsTrack = NULL; // simple track
 786
 787  // Loop over data and store for each distinct pair phi1-phi2 in fRelativeAngleDistribution:
 788  Int_t nPrim = anEvent->NumberOfTracks();  // nPrim = total number of primary tracks, i.e. nPrim = nRP + nPOI + rest, where:
 789                                            // nRP   = # of particles used to determine the reaction plane ("Reference Particles");
 790                                            // nPOI  = # of particles of interest for a detailed flow analysis ("Particles of Interest");
 791                                            // rest  = # of particles which are not niether RPs nor POIs.
 792  // Start nested loops over data:
 793  for(Int_t i=0;i<nPrim;i++)
 794  {
 795   aftsTrack=anEvent->GetTrack(i);
 796   if(aftsTrack)
 797   {
 798    if(!aftsTrack->InRPSelection()) continue; // consider only tracks which are RPs
 799    dPhi1 = aftsTrack->Phi();
 800    for(Int_t j=0;j<nPrim;j++)
 801    {
 802     if(j==i) continue; // eliminating trivial contribution from autocorrelation
 803     aftsTrack=anEvent->GetTrack(j);
 804     if(aftsTrack)
 805     {
 806      if(!aftsTrack->InRPSelection()) continue; // consider only tracks which are RPs
 807      dPhi2 = aftsTrack->Phi();
 808      // Fill the histogram:
 809      fRelativeAngleDistribution->Fill(dPhi1-dPhi2);
 810     }
 811    } // end of for(Int_t j=0;j<nPrim;j++)
 812   } else // to if(aftsTrack)
 813     {
 814      cout<<endl;
 815      cout<<" WARNING (NL): No particle! (i.e. aftsTrack is a NULL pointer in AFAWNL::Make().)"<<endl;
 816      cout<<endl;
 817     }
 818  } // end of for(Int_t i=0;i<nPrim;i++)
 819
 820 } // end of void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForRAD(AliFlowEventSimple *anEvent)
 821
 822 //================================================================================================================
 823
 824 void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForMH(AliFlowEventSimple *anEvent)
 825 {
 826  // Evaluate nested loops needed for mixed harmonics.
 827  // Remark: phi labels the azimuthal angle of RP particle and psi labels azimuthal angle of POI particle.
 828
 829  Int_t nPrim = anEvent->NumberOfTracks();
 830  Int_t nRP = anEvent->GetEventNSelTracksRP();
 831  AliFlowTrackSimple *aftsTrack = NULL;
 832  Double_t phi1=0.,phi2=0.,phi3=0.,phi4=0.,phi5=0.; // azimuthal angles of RPs
 833  Double_t psi1=0.,psi2=0.; // azimuthal angles of POIs
 834  Int_t charge1=0,charge2=0; // charge of POIs
 835  Double_t pt1=0.,pt2=0.; // transverse momenta of POIs
 836  Int_t n = fHarmonic;
 837
 838  // Evaluting correlator cos[n(phi1+phi2-2*phi3)] with three nested loops:
 839  if(nRP>=3)
 840  {
 841   for(Int_t i1=0;i1<nPrim;i1++)
 842   {
 843    aftsTrack=anEvent->GetTrack(i1);
 844    if(!(aftsTrack->InRPSelection())) continue;
 845    phi1 = aftsTrack->Phi();
 846    for(Int_t i2=0;i2<nPrim;i2++)
 847    {
 848     if(i2==i1) continue;
 849     aftsTrack = anEvent->GetTrack(i2);
 850     if(!(aftsTrack->InRPSelection())) continue;
 851     phi2 = aftsTrack->Phi();
 852     for(Int_t i3=0;i3<nPrim;i3++)
 853     {
 854      if(i3==i1||i3==i2) continue;
 855      aftsTrack=anEvent->GetTrack(i3);
 856      if(!(aftsTrack->InRPSelection())) continue;
 857      phi3 = aftsTrack->Phi();
 858      f3pCorrelatorPro->Fill(0.5,cos(n*(phi1+phi2-2.*phi3)),1.);
 859     } // end of for(Int_t i3=0;i3<nPrim;i3++)
 860    } // end of for(Int_t i2=0;i2<nPrim;i2++)
 861   } // end of for(Int_t i1=0;i1<nPrim;i1++)
 862  }
 863
 864  // Evaluting correlator cos[n(2*phi1+2*phi2+2*phi3-3*phi4-3*phi5)] with five nested loops:
 865  if(nRP>=5)
 866  {
 867   for(Int_t i1=0;i1<nPrim;i1++)
 868   {
 869    aftsTrack=anEvent->GetTrack(i1);
 870    if(!(aftsTrack->InRPSelection())) continue;
 871    phi1=aftsTrack->Phi();
 872    for(Int_t i2=0;i2<nPrim;i2++)
 873    {
 874     if(i2==i1)continue;
 875     aftsTrack=anEvent->GetTrack(i2);
 876     if(!(aftsTrack->InRPSelection())) continue;
 877     phi2=aftsTrack->Phi();
 878     for(Int_t i3=0;i3<nPrim;i3++)
 879     {
 880      if(i3==i1||i3==i2)continue;
 881      aftsTrack=anEvent->GetTrack(i3);
 882      if(!(aftsTrack->InRPSelection())) continue;
 883      phi3=aftsTrack->Phi();
 884      for(Int_t i4=0;i4<nPrim;i4++)
 885      {
 886       if(i4==i1||i4==i2||i4==i3)continue;
 887       aftsTrack=anEvent->GetTrack(i4);
 888       if(!(aftsTrack->InRPSelection())) continue;
 889       phi4=aftsTrack->Phi();
 890       for(Int_t i5=0;i5<nPrim;i5++)
 891       {
 892        if(i5==i1||i5==i2||i5==i3||i5==i4)continue;
 893        aftsTrack=anEvent->GetTrack(i5);
 894        if(!(aftsTrack->InRPSelection())) continue;
 895        phi5=aftsTrack->Phi();
 896        f5pCorrelatorPro->Fill(0.5,cos(2.*n*phi1+2.*n*phi2+2.*n*phi3-3.*n*phi4-3.*n*phi5),1.);
 897       } // end of for(Int_t i5=0;i5<nPrim;i5++)
 898      } // end of for(Int_t i4=0;i4<nPrim;i4++)
 899     } // end of for(Int_t i3=0;i3<nPrim;i3++)
 900    } // end of for(Int_t i2=0;i2<nPrim;i2++)
 901   } // end of for(Int_t i1=0;i1<nPrim;i1++)
 902  } // end of if(nPrim>=5)
 903
 904  // Evaluting correlator cos[n(psi1+psi2-2*phi3)] with three nested loops:
 905  if(!fEvaluateDifferential3pCorrelator){return;}
 906  for(Int_t i1=0;i1<nPrim;i1++)
 907  {
 908   aftsTrack=anEvent->GetTrack(i1);
 909   if(!(aftsTrack->InPOISelection())){continue;}
 910   psi1 = aftsTrack->Phi();
 911   pt1 = aftsTrack->Pt();
 912   charge1 = aftsTrack->Charge();
 913   for(Int_t i2=0;i2<nPrim;i2++)
 914   {
 915    if(i2==i1){continue;}
 916    aftsTrack = anEvent->GetTrack(i2);
 917    if(!(aftsTrack->InPOISelection())){continue;}
 918    psi2 = aftsTrack->Phi();
 919    pt2 = aftsTrack->Pt();
 920    charge2 = aftsTrack->Charge();
 921    if(fOppositeChargesPOI && charge1 == charge2){continue;}
 922    for(Int_t i3=0;i3<nPrim;i3++)
 923    {
 924     if(i3==i1||i3==i2){continue;}
 925     aftsTrack=anEvent->GetTrack(i3);
 926     if(!(aftsTrack->InRPSelection())){continue;}
 927     phi3 = aftsTrack->Phi();
 928     // Evaluate and store differential correlator cos[n(psi1+psi2-2*phi3)]:
 929     Double_t ptSum = (pt1+pt2)/2.;
 930     Double_t ptDiff = TMath::Abs(pt1-pt2);
 931     Double_t diff3pCorrelator = TMath::Cos(n*(psi1+psi2-2.*phi3));
 932     f3pCorrelatorVsPtSumDiffDirectPro[0]->Fill(ptSum,diff3pCorrelator,1.);
 933     f3pCorrelatorVsPtSumDiffDirectPro[1]->Fill(ptDiff,diff3pCorrelator,1.);
 934    } // end of for(Int_t i3=0;i3<nPrim;i3++)
 935   } // end of for(Int_t i2=0;i2<nPrim;i2++)
 936  } // end of for(Int_t i1=0;i1<nPrim;i1++)
 937
 938 } // end of void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForMH(AliFlowEventSimple *anEvent)
 939
 940 //================================================================================================================
 941
 942 void AliFlowAnalysisWithNestedLoops::PrintOnTheScreen()
 943 {
 944  // Print on the screen.
 945
 946  cout<<endl;
 947  cout<<"********************************************************************"<<endl;
 948  cout<<"********************************************************************"<<endl;
 949  cout<<"                  Nested Loops                 "<<endl;
 950  cout<<endl;
 951
 952  if(fEvaluateNestedLoopsForRAD)
 953  {
 954   cout<<"  Evaluated for relative angle distribution."<<endl;
 955  }
 956
 957  if(fEvaluateNestedLoopsForMH)
 958  {
 959   cout<<"  Evaluated for mixed harmonics:"<<endl;
 960   if(fHarmonic != 1)
 961   {
 962    cout<<"   cos["<<fHarmonic<<"(phi1+phi2-2phi3)] = "<<f3pCorrelatorPro->GetBinContent(1)<<" +/- " <<f3pCorrelatorPro->GetBinError(1)<<endl;
 963    cout<<"   cos["<<fHarmonic<<"(2phi1+2phi2+2phi3-3phi4-3phi5)] = "<<f5pCorrelatorPro->GetBinContent(1)<<" +/- " <<f5pCorrelatorPro->GetBinError(1)<<endl;
 964   } else
 965     {
 966      cout<<"   cos(phi1+phi2-2phi3) = "<<f3pCorrelatorPro->GetBinContent(1)<<" +/- " <<f3pCorrelatorPro->GetBinError(1)<<endl;
 967      cout<<"   cos[(2phi1+2phi2+2phi3-3phi4-3phi5)] = "<<f5pCorrelatorPro->GetBinContent(1)<<" +/- " <<f5pCorrelatorPro->GetBinError(1)<<endl;
 968     }
 969   if(fEvaluateDifferential3pCorrelator)
 970   {
 971    cout<<"  Evaluated also for differential 3p correlator "<<endl;
 972    if(fHarmonic != 1)
 973    {
 974     cout<<"   cos["<<fHarmonic<<"(psi1+psi2-2phi3)]. "<<endl;
 975    } else
 976      {
 977       cout<<"   cos(psi1+psi2-2phi3). "<<endl;
 978      }
 979   } // end of if(fEvaluateDifferential3pCorrelator)
 980  } // end of if(fEvaluateNestedLoopsForMH)
 981
 982  if(!fEvaluateNestedLoopsForRAD && !fEvaluateNestedLoopsForMH)
 983  {
 984   cout<<"  Not evaluated."<<endl;
 985  }
 986  cout<<endl;
 987  cout<<"********************************************************************"<<endl;
 988  cout<<"********************************************************************"<<endl;
 989  cout<<endl;
 990
 991 } // end of void AliFlowAnalysisWithNestedLoops::PrintOnTheScreen()
 992
 993