PWG2/FLOW/AliFlowCommon/AliFlowAnalysisWithQCumulants.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 /**********************************
  17  * flow analysis with Q-cumulants *
  18  *                                *
  19  * author:  Ante Bilandzic        *
  20  *           (anteb@nikhef.nl)    *
  21  *********************************/
  22
  23 #define AliFlowAnalysisWithQCumulants_cxx
  24
  25 #include "Riostream.h"
  26 #include "AliFlowCommonConstants.h"
  27 #include "AliFlowCommonHist.h"
  28 #include "AliFlowCommonHistResults.h"
  29 #include "TChain.h"
  30 #include "TFile.h"
  31 #include "TList.h"
  32 #include "TGraph.h"
  33 #include "TParticle.h"
  34 #include "TRandom3.h"
  35 #include "TStyle.h"
  36 #include "TProfile.h"
  37 #include "TProfile2D.h"
  38 #include "TProfile3D.h"
  39 #include "TMath.h"
  40 #include "TArrow.h"
  41 #include "TPaveLabel.h"
  42 #include "TCanvas.h"
  43 #include "AliFlowEventSimple.h"
  44 #include "AliFlowTrackSimple.h"
  45 #include "AliFlowAnalysisWithQCumulants.h"
  46 #include "TArrayD.h"
  47 #include "TRandom.h"
  48
  49 class TH1;
  50 class TGraph;
  51 class TPave;
  52 class TLatex;
  53 class TMarker;
  54 class TRandom3;
  55 class TObjArray;
  56 class TList;
  57 class TCanvas;
  58 class TSystem;
  59 class TROOT;
  60 class AliFlowVector;
  61 class TVector;
  62
  63 //================================================================================================================
  64
  65 ClassImp(AliFlowAnalysisWithQCumulants)
  66
  67 AliFlowAnalysisWithQCumulants::AliFlowAnalysisWithQCumulants():
  68  fTrack(NULL),
  69  fHistList(NULL),
  70  fDiffFlowList(NULL),
  71  fWeightsList(NULL),
  72  fAvMultIntFlowQC(NULL),
  73  fQvectorComponents(NULL),
  74  fIntFlowResultsQC(NULL),
  75  fDiffFlowResults2ndOrderQC(NULL),
  76  fDiffFlowResults4thOrderQC(NULL),
  77  fCovariances(NULL),
  78  fQvectorForEachEventX(NULL),//to be removed
  79  fQvectorForEachEventY(NULL),//to be removed
  80  fQCorrelations(NULL),
  81  fWeightedQCorrelations(NULL),
  82  fQProduct(NULL),
  83  fDirectCorrelations(NULL),
  84  f2PerPtBin1n1nPOI(NULL),
  85  f4PerPtBin1n1n1n1nPOI(NULL),
  86  f2PerEtaBin1n1nPOI(NULL),
  87  f4PerEtaBin1n1n1n1nPOI(NULL),
  88  f2WPerPtBin1n1nPOI(NULL),
  89  f4WPerPtBin1n1n1n1nPOI(NULL),
  90  f2WPerEtaBin1n1nPOI(NULL),
  91  f4WPerEtaBin1n1n1n1nPOI(NULL),
  92  f2PerPtBin1n1nRP(NULL),
  93  f4PerPtBin1n1n1n1nRP(NULL),
  94  f2PerEtaBin1n1nRP(NULL),
  95  f4PerEtaBin1n1n1n1nRP(NULL),
  96  f2WPerPtBin1n1nRP(NULL),
  97  f4WPerPtBin1n1n1n1nRP(NULL),
  98  f2WPerEtaBin1n1nRP(NULL),
  99  f4WPerEtaBin1n1n1n1nRP(NULL),
 100  fCommonHists2nd(NULL),
 101  fCommonHists4th(NULL),
 102  fCommonHists6th(NULL),
 103  fCommonHists8th(NULL),
 104  fCommonHistsResults2nd(NULL),
 105  fCommonHistsResults4th(NULL),
 106  fCommonHistsResults6th(NULL),
 107  fCommonHistsResults8th(NULL),
 108  f2pDistribution(NULL),
 109  f4pDistribution(NULL),
 110  f6pDistribution(NULL),
 111  f8pDistribution(NULL),
 112  fnBinsPt(0),
 113  fPtMin(0),
 114  fPtMax(0),
 115  fnBinsEta(0),
 116  fEtaMin(0),
 117  fEtaMax(0),
 118  fEventCounter(0),
 119  fUsePhiWeights(kFALSE),
 120  fUsePtWeights(kFALSE),
 121  fUseEtaWeights(kFALSE)
 122 {
 123  //constructor
 124  fHistList = new TList();
 125  fDiffFlowList = new TList();
 126  fDiffFlowList->SetName("DifferentialFlow");
 127  fWeightsList = new TList();
 128  fWeightsList->SetName("Weights");
 129
 130  fnBinsPt = AliFlowCommonConstants::GetNbinsPt();
 131  fPtMin   = AliFlowCommonConstants::GetPtMin();
 132  fPtMax   = AliFlowCommonConstants::GetPtMax();
 133
 134  fnBinsEta = AliFlowCommonConstants::GetNbinsEta();
 135  fEtaMin   = AliFlowCommonConstants::GetEtaMin();
 136  fEtaMax   = AliFlowCommonConstants::GetEtaMax();
 137 }
 138
 139 AliFlowAnalysisWithQCumulants::~AliFlowAnalysisWithQCumulants()
 140 {
 141  //destructor
 142  delete fHistList;
 143  delete fDiffFlowList;
 144  delete fWeightsList;
 145 }
 146
 147 //================================================================================================================
 148
 149 void AliFlowAnalysisWithQCumulants::Init()
 150 {
 151  //various output histograms
 152  //avarage multiplicity
 153  fAvMultIntFlowQC = new TProfile("fAvMultIntFlowQC","Average Multiplicity",1,0,1,"s");
 154  fAvMultIntFlowQC->SetXTitle("");
 155  fAvMultIntFlowQC->SetYTitle("");
 156  fAvMultIntFlowQC->SetLabelSize(0.06);
 157  fAvMultIntFlowQC->SetMarkerStyle(25);
 158  fAvMultIntFlowQC->SetLabelOffset(0.01);
 159  (fAvMultIntFlowQC->GetXaxis())->SetBinLabel(1,"Average Multiplicity");
 160  fHistList->Add(fAvMultIntFlowQC);
 161
 162  //Q-vector stuff
 163  fQvectorComponents = new TProfile("fQvectorComponents","Avarage of Q-vector components",44,0.,44.,"s");
 164  fQvectorComponents->SetXTitle("");
 165  fQvectorComponents->SetYTitle("");
 166  //fHistList->Add(fQvectorComponents);
 167
 168  //final results for integrated flow from Q-cumulants
 169  fIntFlowResultsQC = new TH1D("fIntFlowResultsQC","Integrated Flow from Q-cumulants",4,0,4);
 170  //fIntFlowResults->SetXTitle("");
 171  //fIntFlowResultsQC->SetYTitle("IntegFALSrated Flow");
 172  fIntFlowResultsQC->SetLabelSize(0.06);
 173  //fIntFlowResultsQC->SetTickLength(1);
 174  fIntFlowResultsQC->SetMarkerStyle(25);
 175  (fIntFlowResultsQC->GetXaxis())->SetBinLabel(1,"v_{n}{2}");
 176  (fIntFlowResultsQC->GetXaxis())->SetBinLabel(2,"v_{n}{4}");
 177  (fIntFlowResultsQC->GetXaxis())->SetBinLabel(3,"v_{n}{6}");
 178  (fIntFlowResultsQC->GetXaxis())->SetBinLabel(4,"v_{n}{8}");
 179  fHistList->Add(fIntFlowResultsQC);
 180
 181  //final results for differential flow from 2nd order Q-cumulant
 182  fDiffFlowResults2ndOrderQC = new TH1D("fDiffFlowResults2ndOrderQC","Differential Flow from 2nd Order Q-cumulant",fnBinsPt,fPtMin,fPtMax);
 183  fDiffFlowResults2ndOrderQC->SetXTitle("p_{t} [GeV]");
 184  //fDiffFlowResults2ndOrderQC->SetYTitle("Differential Flow");
 185  fHistList->Add(fDiffFlowResults2ndOrderQC);
 186
 187  //final results for differential flow from 4th order Q-cumulant
 188  fDiffFlowResults4thOrderQC = new TH1D("fDiffFlowResults4thOrderQC","Differential Flow from 4th Order Q-cumulant",fnBinsPt,fPtMin,fPtMax);
 189  fDiffFlowResults4thOrderQC->SetXTitle("p_{t} [GeV]");
 190  //fDiffFlowResults4thOrderQC->SetYTitle("Differential Flow");
 191  fHistList->Add(fDiffFlowResults4thOrderQC);
 192
 193  //final results for covariances (1st bin: <2*4>-<2>*<4>, 2nd bin: <2*6>-<2>*<6>, ...)
 194  fCovariances = new TH1D("fCovariances","Covariances",6,0,6);
 195  //fCovariances->SetXTitle("");
 196  //fCovariances->SetYTitle("<covariance>");
 197  fCovariances->SetLabelSize(0.04);
 198  fCovariances->SetTickLength(1);
 199  fCovariances->SetMarkerStyle(25);
 200  (fCovariances->GetXaxis())->SetBinLabel(1,"Cov(2,4)");
 201  (fCovariances->GetXaxis())->SetBinLabel(2,"Cov(2,6)");
 202  (fCovariances->GetXaxis())->SetBinLabel(3,"Cov(2,8)");
 203  (fCovariances->GetXaxis())->SetBinLabel(4,"Cov(4,6)");
 204  (fCovariances->GetXaxis())->SetBinLabel(5,"Cov(4,8)");
 205  (fCovariances->GetXaxis())->SetBinLabel(6,"Cov(6,8)");
 206  fHistList->Add(fCovariances);
 207
 208  //xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
 209  //                   !!!! to be removed !!!!
 210  //profile containing the x-components of Q-vectors from all events
 211  fQvectorForEachEventX = new TProfile("fQvectorForEachEventX","x-components of Q-vectors",44000,1,44000,"s");
 212  fHistList->Add(fQvectorForEachEventX);
 213
 214  //profile containing the y-components of Q-vectors from all events
 215  fQvectorForEachEventY = new TProfile("fQvectorForEachEventY","y-components of Q-vectors",44000,1,44000,"s");
 216  fHistList->Add(fQvectorForEachEventY);
 217  //xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
 218
 219  //multi-particle correlations calculated from Q-vectors
 220  fQCorrelations = new TProfile("fQCorrelations","multi-particle correlations from Q-vectors",32,0,32,"s");
 221  //fQCorrelations->SetXTitle("correlations");
 222  //fQCorrelations->SetYTitle("");
 223  fQCorrelations->SetTickLength(-0.01,"Y");
 224  fQCorrelations->SetMarkerStyle(25);
 225  fQCorrelations->SetLabelSize(0.03);
 226  fQCorrelations->SetLabelOffset(0.01,"Y");
 227
 228  (fQCorrelations->GetXaxis())->SetBinLabel(1,"<<2>>_{n|n}");
 229  (fQCorrelations->GetXaxis())->SetBinLabel(2,"<<2>>_{2n|2n}");
 230  (fQCorrelations->GetXaxis())->SetBinLabel(3,"<<2>>_{3n|3n}");
 231  (fQCorrelations->GetXaxis())->SetBinLabel(4,"<<2>>_{4n|4n}");
 232
 233  (fQCorrelations->GetXaxis())->SetBinLabel(6,"<<3>>_{2n|n,n}");
 234  (fQCorrelations->GetXaxis())->SetBinLabel(7,"<<3>>_{3n|2n,n}");
 235  (fQCorrelations->GetXaxis())->SetBinLabel(8,"<<3>>_{4n|2n,2n}");
 236  (fQCorrelations->GetXaxis())->SetBinLabel(9,"<<3>>_{4n|3n,n}");
 237
 238  (fQCorrelations->GetXaxis())->SetBinLabel(11,"<<4>>_{n,n|n,n}");
 239  (fQCorrelations->GetXaxis())->SetBinLabel(12,"<<4>>_{2n,n|2n,n}");
 240  (fQCorrelations->GetXaxis())->SetBinLabel(13,"<<4>>_{2n,2n|2n,2n}");
 241  (fQCorrelations->GetXaxis())->SetBinLabel(14,"<<4>>_{3n|n,n,n}");
 242  (fQCorrelations->GetXaxis())->SetBinLabel(15,"<<4>>_{3n,n|3n,n}");
 243  (fQCorrelations->GetXaxis())->SetBinLabel(16,"<<4>>_{3n,n|2n,2n}");
 244  (fQCorrelations->GetXaxis())->SetBinLabel(17,"<<4>>_{4n|2n,n,n}");
 245
 246  (fQCorrelations->GetXaxis())->SetBinLabel(19,"<<5>>_{2n|n,n,n,n}");
 247  (fQCorrelations->GetXaxis())->SetBinLabel(20,"<<5>>_{2n,2n|2n,n,n}");
 248  (fQCorrelations->GetXaxis())->SetBinLabel(21,"<<5>>_{3n,n|2n,n,n}");
 249  (fQCorrelations->GetXaxis())->SetBinLabel(22,"<<5>>_{4n|n,n,n,n}");
 250
 251  (fQCorrelations->GetXaxis())->SetBinLabel(24,"<<6>>_{n,n,n|n,n,n}");
 252  (fQCorrelations->GetXaxis())->SetBinLabel(25,"<<6>>_{2n,n,n|2n,n,n}");
 253  (fQCorrelations->GetXaxis())->SetBinLabel(26,"<<6>>_{2n,2n|n,n,n,n}");
 254  (fQCorrelations->GetXaxis())->SetBinLabel(27,"<<6>>_{3n,n|n,n,n,n}");
 255
 256  (fQCorrelations->GetXaxis())->SetBinLabel(29,"<<7>>_{2n,n,n|n,n,n,n}");
 257
 258  (fQCorrelations->GetXaxis())->SetBinLabel(31,"<<8>>_{n,n,n,n|n,n,n,n}");
 259
 260  fHistList->Add(fQCorrelations);
 261
 262
 263
 264
 265  //weighted multi-particle correlations calculated from Q-vectors
 266  fWeightedQCorrelations = new TProfile("fWeightedQCorrelations","weighted multi-particle correlations from Q-vectors",100,0,100,"s");
 267  //fWeightedQCorrelations->SetXTitle("correlations");
 268  //fWeightedQCorrelations->SetYTitle("");
 269  fWeightedQCorrelations->SetTickLength(-0.01,"Y");
 270  fWeightedQCorrelations->SetMarkerStyle(25);
 271  fWeightedQCorrelations->SetLabelSize(0.03);
 272  fWeightedQCorrelations->SetLabelOffset(0.01,"Y");
 273
 274  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(1,"<w_{1}w_{2}cos(n(#phi_{1}-#phi_{2}))>");
 275  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(2,"<w_{1}^{2}w_{2}^{2}cos(2n(#phi_{1}-#phi_{2}))>");
 276  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(3,"<w_{1}^{3}w_{2}^{3}cos(3n(#phi_{1}-#phi_{2}))>");
 277  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(4,"<w_{1}^{4}w_{2}^{4}cos(4n(#phi_{1}-#phi_{2}))>");
 278  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(5,"<w_{1}^{3}w_{2}cos(n(#phi_{1}-#phi_{2}))>");
 279  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(6,"<w_{1}^{2}w_{2}w_{3}cos(n(#phi_{1}-#phi_{2}))>");
 280
 281  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(11,"<w_{1}w_{2}w_{3}^{2}cos(n(2#phi_{1}-#phi_{2}-#phi_{3}))>");
 282
 283  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(21,"<w_{1}w_{2}w_{3}w_{4}cos(n(#phi_{1}+#phi_{2}-#phi_{3}-#phi_{4}))>");
 284
 285  /*
 286  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(7,"<<3>>_{3n|2n,n}");
 287  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(8,"<<3>>_{4n|2n,2n}");
 288  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(9,"<<3>>_{4n|3n,n}");
 289
 290  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(11,"<<4>>_{n,n|n,n}");
 291  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(12,"<<4>>_{2n,n|2n,n}");
 292  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(13,"<<4>>_{2n,2n|2n,2n}");
 293  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(14,"<<4>>_{3n|n,n,n}");
 294  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(15,"<<4>>_{3n,n|3n,n}");
 295  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(16,"<<4>>_{3n,n|2n,2n}");
 296  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(17,"<<4>>_{4n|2n,n,n}");
 297
 298  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(19,"<<5>>_{2n|n,n,n,n}");
 299  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(20,"<<5>>_{2n,2n|2n,n,n}");
 300  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(21,"<<5>>_{3n,n|2n,n,n}");
 301  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(22,"<<5>>_{4n|n,n,n,n}");
 302
 303  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(24,"<<6>>_{n,n,n|n,n,n}");
 304  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(25,"<<6>>_{2n,n,n|2n,n,n}");
 305  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(26,"<<6>>_{2n,2n|n,n,n,n}");
 306  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(27,"<<6>>_{3n,n|n,n,n,n}");
 307
 308  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(29,"<<7>>_{2n,n,n|n,n,n,n}");
 309
 310  (fWeightedQCorrelations->GetXaxis())->SetBinLabel(31,"<<8>>_{n,n,n,n|n,n,n,n}");
 311  */
 312
 313  fHistList->Add(fWeightedQCorrelations);
 314
 315
 316
 317
 318  //average products
 319  fQProduct = new TProfile("fQProduct","average of products",6,0,6,"s");
 320  fQProduct->SetTickLength(-0.01,"Y");
 321  fQProduct->SetMarkerStyle(25);
 322  fQProduct->SetLabelSize(0.03);
 323  fQProduct->SetLabelOffset(0.01,"Y");
 324  (fQProduct->GetXaxis())->SetBinLabel(1,"<<2*4>>");
 325  (fQProduct->GetXaxis())->SetBinLabel(2,"<<2*6>>");
 326  (fQProduct->GetXaxis())->SetBinLabel(3,"<<2*8>>");
 327  (fQProduct->GetXaxis())->SetBinLabel(4,"<<4*6>>");
 328  (fQProduct->GetXaxis())->SetBinLabel(5,"<<4*8>>");
 329  (fQProduct->GetXaxis())->SetBinLabel(6,"<<6*8>>");
 330  fQProduct->SetXTitle("");
 331  fQProduct->SetYTitle("");
 332  fHistList->Add(fQProduct);
 333
 334  //weighted multi-particle correlations calculated with nested loops (0..100 integrated flow; 100..200 differential flow)
 335  fDirectCorrelations = new TProfile("fDirectCorrelations","multi-particle correlations with nested loops",200,0,200,"s");
 336  fDirectCorrelations->SetXTitle("");
 337  fDirectCorrelations->SetYTitle("correlations");
 338  fHistList->Add(fDirectCorrelations);
 339
 340  //f2PerPtBin1n1nRP
 341  f2PerPtBin1n1nRP = new TProfile("f2PerPtBin1n1nRP","<2'>_{n|n}",fnBinsPt,fPtMin,fPtMax,"s");
 342  f2PerPtBin1n1nRP->SetXTitle("p_{t} [GeV]");
 343  fDiffFlowList->Add(f2PerPtBin1n1nRP);
 344
 345  //f4PerPtBin1n1n1n1nRP
 346  f4PerPtBin1n1n1n1nRP = new TProfile("f4PerPtBin1n1n1n1nRP","<4'>_{n,n|n,n}",fnBinsPt,fPtMin,fPtMax,"s");
 347  f4PerPtBin1n1n1n1nRP->SetXTitle("p_{t} [GeV]");
 348  fDiffFlowList->Add(f4PerPtBin1n1n1n1nRP);
 349
 350  //f2PerEtaBin1n1nRP
 351  f2PerEtaBin1n1nRP = new TProfile("f2PerEtaBin1n1nRP","<2'>_{n|n}",fnBinsEta,fEtaMin,fEtaMax,"s");
 352  f2PerEtaBin1n1nRP->SetXTitle("#eta");
 353  fDiffFlowList->Add(f2PerEtaBin1n1nRP);
 354
 355  //f4PerEtaBin1n1n1n1nRP
 356  f4PerEtaBin1n1n1n1nRP = new TProfile("f4PerEtaBin1n1n1n1nRP","<4'>_{n,n|n,n}",fnBinsEta,fEtaMin,fEtaMax,"s");
 357  f4PerEtaBin1n1n1n1nRP->SetXTitle("#eta");
 358  fDiffFlowList->Add(f4PerEtaBin1n1n1n1nRP);
 359
 360  //f2PerPtBin1n1nPOI
 361  f2PerPtBin1n1nPOI = new TProfile("f2PerPtBin1n1nPOI","<2'>_{n|n}",fnBinsPt,fPtMin,fPtMax,"s");
 362  f2PerPtBin1n1nPOI->SetXTitle("#eta");
 363  fDiffFlowList->Add(f2PerPtBin1n1nPOI);
 364
 365  //f4PerPtBin1n1n1n1nPOI
 366  f4PerPtBin1n1n1n1nPOI = new TProfile("f4PerPtBin1n1n1n1nPOI","<4'>_{n,n|n,n}",fnBinsPt,fPtMin,fPtMax,"s");
 367  f4PerPtBin1n1n1n1nPOI->SetXTitle("p_{t} [GeV]");
 368  fDiffFlowList->Add(f4PerPtBin1n1n1n1nPOI);
 369
 370  //f2PerEtaBin1n1nPOI
 371  f2PerEtaBin1n1nPOI = new TProfile("f2PerEtaBin1n1nPOI","<2'>_{n|n}",fnBinsEta,fEtaMin,fEtaMax,"s");
 372  f2PerEtaBin1n1nPOI->SetXTitle("#eta");
 373  fDiffFlowList->Add(f2PerEtaBin1n1nPOI);
 374
 375  //f4PerEtaBin1n1n1n1nPOI
 376  f4PerEtaBin1n1n1n1nPOI = new TProfile("f4PerEtaBin1n1n1n1nPOI","<4'>_{n,n|n,n}",fnBinsEta,fEtaMin,fEtaMax,"s");
 377  f4PerEtaBin1n1n1n1nPOI->SetXTitle("#eta");
 378  fDiffFlowList->Add(f4PerEtaBin1n1n1n1nPOI);
 379
 380  //f2WPerPtBin1n1nPOI
 381  f2WPerPtBin1n1nPOI = new TProfile("f2WPerPtBin1n1nPOI","<2'>_{n|n}",fnBinsPt,fPtMin,fPtMax,"s");
 382  f2WPerPtBin1n1nPOI->SetXTitle("#pt");
 383  fDiffFlowList->Add(f2WPerPtBin1n1nPOI);
 384
 385  //f4WPerPtBin1n1n1n1nPOI
 386  f4WPerPtBin1n1n1n1nPOI = new TProfile("f4WPerPtBin1n1n1n1nPOI","<4'>_{n,n|n,n}",fnBinsPt,fPtMin,fPtMax,"s");
 387  f4WPerPtBin1n1n1n1nPOI->SetXTitle("#Pt");
 388  fDiffFlowList->Add(f4WPerPtBin1n1n1n1nPOI);
 389
 390  //f2WPerEtaBin1n1nPOI
 391  f2WPerEtaBin1n1nPOI = new TProfile("f2WPerEtaBin1n1nPOI","<2'>_{n|n}",fnBinsEta,fEtaMin,fEtaMax,"s");
 392  f2WPerEtaBin1n1nPOI->SetXTitle("#eta");
 393  fDiffFlowList->Add(f2WPerEtaBin1n1nPOI);
 394
 395  //f4WPerEtaBin1n1n1n1nPOI
 396  f4WPerEtaBin1n1n1n1nPOI = new TProfile("f4WPerEtaBin1n1n1n1nPOI","<4'>_{n,n|n,n}",fnBinsEta,fEtaMin,fEtaMax,"s");
 397  f4WPerEtaBin1n1n1n1nPOI->SetXTitle("#eta");
 398  fDiffFlowList->Add(f4WPerEtaBin1n1n1n1nPOI);
 399
 400  //f2WPerPtBin1n1nRP
 401  f2WPerPtBin1n1nRP = new TProfile("f2WPerPtBin1n1nRP","<2'>_{n|n}",fnBinsPt,fPtMin,fPtMax,"s");
 402  f2WPerPtBin1n1nRP->SetXTitle("#pt");
 403  fDiffFlowList->Add(f2WPerPtBin1n1nRP);
 404
 405  //f4WPerPtBin1n1n1n1nRP
 406  f4WPerPtBin1n1n1n1nRP = new TProfile("f4WPerPtBin1n1n1n1nRP","<4'>_{n,n|n,n}",fnBinsPt,fPtMin,fPtMax,"s");
 407  f4WPerPtBin1n1n1n1nRP->SetXTitle("#Pt");
 408  fDiffFlowList->Add(f4WPerPtBin1n1n1n1nRP);
 409
 410  //f2WPerEtaBin1n1nRP
 411  f2WPerEtaBin1n1nRP = new TProfile("f2WPerEtaBin1n1nRP","<2'>_{n|n}",fnBinsEta,fEtaMin,fEtaMax,"s");
 412  f2WPerEtaBin1n1nRP->SetXTitle("#eta");
 413  fDiffFlowList->Add(f2WPerEtaBin1n1nRP);
 414
 415  //f4WPerEtaBin1n1n1n1nRP
 416  f4WPerEtaBin1n1n1n1nRP = new TProfile("f4WPerEtaBin1n1n1n1nRP","<4'>_{n,n|n,n}",fnBinsEta,fEtaMin,fEtaMax,"s");
 417  f4WPerEtaBin1n1n1n1nRP->SetXTitle("#eta");
 418  fDiffFlowList->Add(f4WPerEtaBin1n1n1n1nRP);
 419
 420  //common control histogram (2nd order)
 421  fCommonHists2nd = new AliFlowCommonHist("AliFlowCommonHist2ndOrderQC");
 422  fHistList->Add(fCommonHists2nd);
 423
 424  //common control histogram (4th order)
 425  fCommonHists4th = new AliFlowCommonHist("AliFlowCommonHist4thOrderQC");
 426  fHistList->Add(fCommonHists4th);
 427
 428  //common control histogram (6th order)
 429  fCommonHists6th = new AliFlowCommonHist("AliFlowCommonHist6thOrderQC");
 430  fHistList->Add(fCommonHists6th);
 431
 432  //common control histogram (8th order)
 433  fCommonHists8th = new AliFlowCommonHist("AliFlowCommonHist8thOrderQC");
 434  fHistList->Add(fCommonHists8th);
 435
 436  //common histograms for final results (2nd order)
 437  fCommonHistsResults2nd = new AliFlowCommonHistResults("AliFlowCommonHistResults2ndOrderQC");
 438  fHistList->Add(fCommonHistsResults2nd);
 439
 440  //common histograms for final results (4th order)
 441  fCommonHistsResults4th = new AliFlowCommonHistResults("AliFlowCommonHistResults4thOrderQC");
 442  fHistList->Add(fCommonHistsResults4th);
 443
 444  //common histograms for final results (6th order)
 445  fCommonHistsResults6th = new AliFlowCommonHistResults("AliFlowCommonHistResults6thOrderQC");
 446  fHistList->Add(fCommonHistsResults6th);
 447
 448  //common histograms for final results (8th order)
 449  fCommonHistsResults8th = new AliFlowCommonHistResults("AliFlowCommonHistResults8thOrderQC");
 450  fHistList->Add(fCommonHistsResults8th);
 451
 452  //weighted <2>_{n|n} distribution
 453  f2pDistribution = new TH1D("f2pDistribution","<2>_{n|n} distribution",100000,-0.02,0.1);
 454  f2pDistribution->SetXTitle("<2>_{n|n}");
 455  f2pDistribution->SetYTitle("Counts");
 456  fHistList->Add(f2pDistribution);
 457
 458  //weighted <4>_{n,n|n,n} distribution
 459  f4pDistribution = new TH1D("f4pDistribution","<4>_{n,n|n,n} distribution",100000,-0.00025,0.002);
 460  f4pDistribution->SetXTitle("<4>_{n,n|n,n}");
 461  f4pDistribution->SetYTitle("Counts");
 462  fHistList->Add(f4pDistribution);
 463
 464  //weighted <6>_{n,n,n|n,n,n} distribution
 465  f6pDistribution = new TH1D("f6pDistribution","<6>_{n,n,n|n,n,n} distribution",100000,-0.000005,0.000025);
 466  f6pDistribution->SetXTitle("<6>_{n,n,n|n,n,n}");
 467  f6pDistribution->SetYTitle("Counts");
 468  fHistList->Add(f6pDistribution);
 469
 470  //weighted <8>_{n,n,n,n|n,n,n,n} distribution
 471  f8pDistribution = new TH1D("f8pDistribution","<8>_{n,n,n,n|n,n,n,n} distribution",100000,-0.000000001,0.00000001);
 472  f8pDistribution->SetXTitle("<8>_{n,n,n,n|n,n,n,n}");
 473  f8pDistribution->SetYTitle("Counts");
 474  fHistList->Add(f8pDistribution);
 475
 476  // add list fWeightsList with weights to the main list
 477  fHistList->Add(fWeightsList);
 478
 479  // add list fDiffFlowList with histograms and profiles needed for differential flow to the main list
 480  fHistList->Add(fDiffFlowList);
 481 }//end of Init()
 482
 483 //================================================================================================================
 484
 485 void AliFlowAnalysisWithQCumulants::Make(AliFlowEventSimple* anEvent)
 486 {
 487  // running over data
 488
 489  Int_t nPrim = anEvent->NumberOfTracks(); // nPrim = nRP + nPOI + rest
 490  Int_t nRP = anEvent->GetEventNSelTracksRP(); // nRP = number of particles used to determine the reaction plane
 491
 492  Int_t n = 2; // int flow harmonic (to be improved)
 493
 494  //needed for debugging: (to be improved - add explanation here)
 495  //Bool_t bNestedLoops=kFALSE;
 496  //if(!(bNestedLoops)||(nPrim>0&&nPrim<12))
 497  //{
 498  //if(nPrim>0&&nPrim<10)
 499  //{
 500
 501
 502
 503  //---------------------------------------------------------------------------------------------------------
 504  // non-weighted and weighted Q-vectors of an event built-up from RP particles evaluated in harmonics n, 2n, 3n and 4n:
 505  AliFlowVector afvQvector1n, afvQvector2n, afvQvector3n, afvQvector4n;
 506
 507  // non-weighted Q-vector in harmonic n:
 508  afvQvector1n.Set(0.,0.);
 509  afvQvector1n.SetMult(0);
 510  afvQvector1n = anEvent->GetQ(1*n);
 511
 512  // non-weighted Q-vector in harmonic 2n:
 513  afvQvector2n.Set(0.,0.);
 514  afvQvector2n.SetMult(0);
 515  afvQvector2n = anEvent->GetQ(2*n); // to be improved: weights
 516
 517  // non-weighted Q-vector in harmonic 3n:
 518  afvQvector3n.Set(0.,0.);
 519  afvQvector3n.SetMult(0);
 520  afvQvector3n = anEvent->GetQ(3*n); // to be improved: weights
 521
 522  // non-weighted Q-vector in harmonic 4n:
 523  afvQvector4n.Set(0.,0.);
 524  afvQvector4n.SetMult(0);
 525  afvQvector4n = anEvent->GetQ(4*n); // to be improved: weights
 526
 527
 528  //xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
 529  //                        !!!! to be removed !!!!
 530  fQvectorForEachEventX->Fill(1.*(++fEventCounter),afvQvector1n.X());
 531  fQvectorForEachEventY->Fill(1.*(fEventCounter),afvQvector1n.Y());
 532  //xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
 533
 534
 535
 536  //---------------------------------------------------------------------------------------------------------
 537
 538  //multiplicity of RP particles:
 539  Double_t dMult = afvQvector1n.GetMult(); // to be improved (name, this is actually weighted multiplicity)
 540
 541  fAvMultIntFlowQC->Fill(0.,dMult,1.); // to be removed (this info is also stored in one of control histograms)
 542
 543  //---------------------------------------------------------------------------------------------------------
 544  //
 545  //                                          *******************
 546  //                                          **** Q-vectors ****
 547  //                                          *******************
 548  //
 549  Double_t reQ2nQ1nstarQ1nstar = pow(afvQvector1n.X(),2.)*afvQvector2n.X()+2.*afvQvector1n.X()*afvQvector1n.Y()*afvQvector2n.Y()-pow(afvQvector1n.Y(),2.)*afvQvector2n.X();//Re[Q_{2n} Q_{n}^* Q_{n}^*]
 550  //Double_t imQ2nQ1nstarQ1nstar = pow(Qvector1n.X(),2.)*Qvector2n.Y()-2.*Qvector1n.X()*Qvector1n.Y()*Qvector2n.X()-pow(Qvector1n.Y(),2.)*Qvector2n.Y();//Im[Q_{2n} Q_{n}^* Q_{n}^*]
 551  Double_t reQ1nQ1nQ2nstar = reQ2nQ1nstarQ1nstar;//Re[Q_{n} Q_{n} Q_{2n}^*] = Re[Q_{2n} Q_{n}^* Q_{n}^*]
 552  Double_t reQ3nQ1nQ2nstarQ2nstar = (pow(afvQvector2n.X(),2.)-pow(afvQvector2n.Y(),2.))*(afvQvector3n.X()*afvQvector1n.X()-afvQvector3n.Y()*afvQvector1n.Y())+2.*afvQvector2n.X()*afvQvector2n.Y()*(afvQvector3n.X()*afvQvector1n.Y()+afvQvector3n.Y()*afvQvector1n.X());
 553  //Double_t imQ3nQ1nQ2nstarQ2nstar = calculate and implement this (deleteMe)
 554  Double_t reQ2nQ2nQ3nstarQ1nstar = reQ3nQ1nQ2nstarQ2nstar;
 555  Double_t reQ4nQ2nstarQ2nstar = pow(afvQvector2n.X(),2.)*afvQvector4n.X()+2.*afvQvector2n.X()*afvQvector2n.Y()*afvQvector4n.Y()-pow(afvQvector2n.Y(),2.)*afvQvector4n.X();//Re[Q_{4n} Q_{2n}^* Q_{2n}^*]
 556  //Double_t imQ4nQ2nstarQ2nstar = calculate and implement this (deleteMe)
 557  Double_t reQ2nQ2nQ4nstar = reQ4nQ2nstarQ2nstar;
 558  Double_t reQ4nQ3nstarQ1nstar = afvQvector4n.X()*(afvQvector3n.X()*afvQvector1n.X()-afvQvector3n.Y()*afvQvector1n.Y())+afvQvector4n.Y()*(afvQvector3n.X()*afvQvector1n.Y()+afvQvector3n.Y()*afvQvector1n.X());//Re[Q_{4n} Q_{3n}^* Q_{n}^*]
 559  Double_t reQ3nQ1nQ4nstar = reQ4nQ3nstarQ1nstar;//Re[Q_{3n} Q_{n} Q_{4n}^*] = Re[Q_{4n} Q_{3n}^* Q_{n}^*]
 560  //Double_t imQ4nQ3nstarQ1nstar = calculate and implement this (deleteMe)
 561  Double_t reQ3nQ2nstarQ1nstar = afvQvector3n.X()*afvQvector2n.X()*afvQvector1n.X()-afvQvector3n.X()*afvQvector2n.Y()*afvQvector1n.Y()+afvQvector3n.Y()*afvQvector2n.X()*afvQvector1n.Y()+afvQvector3n.Y()*afvQvector2n.Y()*afvQvector1n.X();//Re[Q_{3n} Q_{2n}^* Q_{n}^*]
 562  Double_t reQ2nQ1nQ3nstar = reQ3nQ2nstarQ1nstar;//Re[Q_{2n} Q_{n} Q_{3n}^*] = Re[Q_{3n} Q_{2n}^* Q_{n}^*]
 563  //Double_t imQ3nQ2nstarQ1nstar; //calculate and implement this (deleteMe)
 564  Double_t reQ3nQ1nstarQ1nstarQ1nstar = afvQvector3n.X()*pow(afvQvector1n.X(),3)-3.*afvQvector1n.X()*afvQvector3n.X()*pow(afvQvector1n.Y(),2)+3.*afvQvector1n.Y()*afvQvector3n.Y()*pow(afvQvector1n.X(),2)-afvQvector3n.Y()*pow(afvQvector1n.Y(),3);//Re[Q_{3n} Q_{n}^* Q_{n}^* Q_{n}^*]
 565  //Double_t imQ3nQ1nstarQ1nstarQ1nstar; //calculate and implement this (deleteMe)
 566  Double_t xQ2nQ1nQ2nstarQ1nstar = pow(afvQvector2n.Mod()*afvQvector1n.Mod(),2);//|Q_{2n}|^2 |Q_{n}|^2
 567  Double_t reQ4nQ2nstarQ1nstarQ1nstar = (afvQvector4n.X()*afvQvector2n.X()+afvQvector4n.Y()*afvQvector2n.Y())*(pow(afvQvector1n.X(),2)-pow(afvQvector1n.Y(),2))+2.*afvQvector1n.X()*afvQvector1n.Y()*(afvQvector4n.Y()*afvQvector2n.X()-afvQvector4n.X()*afvQvector2n.Y());//Re[Q_{4n} Q_{2n}^* Q_{n}^* Q_{n}^*]
 568  //Double_t imQ4nQ2nstarQ1nstarQ1nstar; //calculate and implement this (deleteMe)
 569  Double_t reQ2nQ1nQ1nstarQ1nstarQ1nstar = (afvQvector2n.X()*afvQvector1n.X()-afvQvector2n.Y()*afvQvector1n.Y())*(pow(afvQvector1n.X(),3)-3.*afvQvector1n.X()*pow(afvQvector1n.Y(),2))+(afvQvector2n.X()*afvQvector1n.Y()+afvQvector1n.X()*afvQvector2n.Y())*(3.*afvQvector1n.Y()*pow(afvQvector1n.X(),2)-pow(afvQvector1n.Y(),3));//Re[Q_{2n} Q_{n} Q_{n}^* Q_{n}^* Q_{n}^*]
 570  //Double_t imQ2nQ1nQ1nstarQ1nstarQ1nstar; //calculate and implement this (deleteMe)
 571  Double_t reQ2nQ2nQ2nstarQ1nstarQ1nstar = pow(afvQvector2n.Mod(),2.)*(afvQvector2n.X()*(pow(afvQvector1n.X(),2.)-pow(afvQvector1n.Y(),2.))+2.*afvQvector2n.Y()*afvQvector1n.X()*afvQvector1n.Y());//Re[Q_{2n} Q_{2n} Q_{2n}^* Q_{n}^* Q_{n}^*]
 572  //Double_t imQ2nQ2nQ2nstarQ1nstarQ1nstar = pow(Qvector2n.Mod(),2.)*(Qvector2n.Y()*(pow(Qvector1n.X(),2.)-pow(Qvector1n.Y(),2.))-2.*Qvector2n.X()*Qvector1n.X()*Qvector1n.Y());//Im[Q_{2n} Q_{2n} Q_{2n}^* Q_{n}^* Q_{n}^*]
 573  Double_t reQ4nQ1nstarQ1nstarQ1nstarQ1nstar = pow(afvQvector1n.X(),4.)*afvQvector4n.X()-6.*pow(afvQvector1n.X(),2.)*afvQvector4n.X()*pow(afvQvector1n.Y(),2.)+pow(afvQvector1n.Y(),4.)*afvQvector4n.X()+4.*pow(afvQvector1n.X(),3.)*afvQvector1n.Y()*afvQvector4n.Y()-4.*pow(afvQvector1n.Y(),3.)*afvQvector1n.X()*afvQvector4n.Y();//Re[Q_{4n} Q_{n}^* Q_{n}^* Q_{n}^* Q_{n}^*]
 574  //Double_t imQ4nQ1nstarQ1nstarQ1nstarQ1nstar = pow(Qvector1n.X(),4.)*Qvector4n.Y()-6.*pow(Qvector1n.X(),2.)*Qvector4n.Y()*pow(Qvector1n.Y(),2.)+pow(Qvector1n.Y(),4.)*Qvector4n.Y()+4.*pow(Qvector1n.Y(),3.)*Qvector1n.X()*Qvector4n.X()-4.*pow(Qvector1n.X(),3.)*Qvector1n.Y()*Qvector4n.X();//Im[Q_{4n} Q_{n}^* Q_{n}^* Q_{n}^* Q_{n}^*]
 575  Double_t reQ3nQ1nQ2nstarQ1nstarQ1nstar = pow(afvQvector1n.Mod(),2.)*(afvQvector1n.X()*afvQvector2n.X()*afvQvector3n.X()-afvQvector3n.X()*afvQvector1n.Y()*afvQvector2n.Y()+afvQvector2n.X()*afvQvector1n.Y()*afvQvector3n.Y()+afvQvector1n.X()*afvQvector2n.Y()*afvQvector3n.Y());//Re[Q_{3n} Q_{n} Q_{2n}^* Q_{n}^* Q_{n}^*]
 576  //Double_t imQ3nQ1nQ2nstarQ1nstarQ1nstar = pow(afvQvector1n.Mod(),2.)*(-afvQvector2n.X()*afvQvector3n.X()*afvQvector1n.Y()-afvQvector1n.X()*afvQvector3n.X()*afvQvector2n.Y()+afvQvector1n.X()*afvQvector2n.X()*afvQvector3n.Y()-afvQvector1n.Y()*afvQvector2n.Y()*afvQvector3n.Y());//Im[Q_{3n} Q_{n} Q_{2n}^* Q_{n}^* Q_{n}^*]
 577  Double_t reQ2nQ2nQ1nstarQ1nstarQ1nstarQ1nstar = (pow(afvQvector1n.X(),2.)*afvQvector2n.X()-2.*afvQvector1n.X()*afvQvector2n.X()*afvQvector1n.Y()-afvQvector2n.X()*pow(afvQvector1n.Y(),2.)+afvQvector2n.Y()*pow(afvQvector1n.X(),2.)+2.*afvQvector1n.X()*afvQvector1n.Y()*afvQvector2n.Y()-pow(afvQvector1n.Y(),2.)*afvQvector2n.Y())*(pow(afvQvector1n.X(),2.)*afvQvector2n.X()+2.*afvQvector1n.X()*afvQvector2n.X()*afvQvector1n.Y()-afvQvector2n.X()*pow(afvQvector1n.Y(),2.)-afvQvector2n.Y()*pow(afvQvector1n.X(),2.)+2.*afvQvector1n.X()*afvQvector1n.Y()*afvQvector2n.Y()+pow(afvQvector1n.Y(),2.)*afvQvector2n.Y());//Re[Q_{2n} Q_{2n} Q_{n}^* Q_{n}^* Q_{n}^* Q_{n}^*]
 578  //Double_t imQ2nQ2nQ1nstarQ1nstarQ1nstarQ1nstar = 2.*(pow(afvQvector1n.X(),2.)*afvQvector2n.X()-afvQvector2n.X()*pow(afvQvector1n.Y(),2.)+2.*afvQvector1n.X()*afvQvector1n.Y()*afvQvector2n.Y())*(pow(afvQvector1n.X(),2.)*afvQvector2n.Y()-2.*afvQvector1n.X()*afvQvector1n.Y()*afvQvector2n.X()-pow(afvQvector1n.Y(),2.)*afvQvector2n.Y());//Im[Q_{2n} Q_{2n} Q_{n}^* Q_{n}^* Q_{n}^* Q_{n}^*]
 579  Double_t reQ3nQ1nQ1nstarQ1nstarQ1nstarQ1nstar = pow(afvQvector1n.Mod(),2.)*(pow(afvQvector1n.X(),3.)*afvQvector3n.X()-3.*afvQvector1n.X()*afvQvector3n.X()*pow(afvQvector1n.Y(),2.)+3.*pow(afvQvector1n.X(),2.)*afvQvector1n.Y()*afvQvector3n.Y()-pow(afvQvector1n.Y(),3.)*afvQvector3n.Y());//Re[Q_{3n} Q_{n} Q_{n}^* Q_{n}^* Q_{n}^* Q_{n}^*]
 580  //Double_t imQ3nQ1nQ1nstarQ1nstarQ1nstarQ1nstar = pow(afvQvector1n.Mod(),2.)*(pow(afvQvector1n.Y(),3.)*afvQvector3n.X()-3.*afvQvector1n.Y()*afvQvector3n.X()*pow(afvQvector1n.X(),2.)-3.*pow(afvQvector1n.Y(),2.)*afvQvector1n.X()*afvQvector3n.Y()+pow(afvQvector1n.X(),3.)*afvQvector3n.Y());//Im[Q_{3n} Q_{n} Q_{n}^* Q_{n}^* Q_{n}^* Q_{n}^*]
 581  Double_t xQ2nQ1nQ1nQ2nstarQ1nstarQ1nstar = pow(afvQvector2n.Mod(),2.)*pow(afvQvector1n.Mod(),4.);//|Q_{2n}|^2 |Q_{n}|^4
 582  Double_t reQ2nQ1nQ1nQ1nstarQ1nstarQ1nstarQ1nstar = pow(afvQvector1n.Mod(),4.)*(pow(afvQvector1n.X(),2.)*afvQvector2n.X()-afvQvector2n.X()*pow(afvQvector1n.Y(),2.)+2.*afvQvector1n.X()*afvQvector1n.Y()*afvQvector2n.Y());//Re[Q_{2n} Q_{n} Q_{n} Q_{n}^* Q_{n}^* Q_{n}^* Q_{n}^*]
 583  //Double_t imQ2nQ1nQ1nQ1nstarQ1nstarQ1nstarQ1nstar = pow(afvQvector1n.Mod(),4.)*(pow(afvQvector1n.X(),2.)*afvQvector2n.Y()-afvQvector2n.Y()*pow(afvQvector1n.Y(),2.)-2.*afvQvector1n.X()*afvQvector2n.X()*afvQvector1n.Y());//Re[Q_{2n} Q_{n} Q_{n} Q_{n}^* Q_{n}^* Q_{n}^* Q_{n}^*]
 584  //---------------------------------------------------------------------------------------------------------
 585
 586  //---------------------------------------------------------------------------------------------------------
 587  //
 588  //                                        **************************************
 589  //                                        **** multi-particle correlations: ****
 590  //                                        **************************************
 591  //
 592  // Remark 1: multi-particle correlations calculated with Q-vectors are stored in fQCorrelations.
 593  // Remark 2: binning of fQCorrelations is organized as follows:
 594  //
 595  // 1st bin: <2>_{n|n} = two1n1n
 596  // 2nd bin: <2>_{2n|2n} = two2n2n
 597  // 3rd bin: <2>_{3n|3n} = two3n3n
 598  // 4th bin: <2>_{4n|4n} = two4n4n
 599  // 5th bin: --  EMPTY --
 600  // 6th bin: <3>_{2n|n,n} = three2n1n1n
 601  // 7th bin: <3>_{3n|2n,n} = three3n2n1n
 602  // 8th bin: <3>_{4n|2n,2n} = three4n2n2n
 603  // 9th bin: <3>_{4n|3n,n} = three4n3n1n
 604  //10th bin: --  EMPTY --
 605  //11th bin: <4>_{n,n|n,n} = four1n1n1n1n
 606  //12th bin: <4>_{2n,n|2n,n} = four2n1n2n1n
 607  //13th bin: <4>_{2n,2n|2n,2n} = four2n2n2n2n
 608  //14th bin: <4>_{3n|n,n,n} = four3n1n1n1n
 609  //15th bin: <4>_{3n,n|3n,n} = four3n1n3n1n
 610  //16th bin: <4>_{3n,n|2n,2n} = four3n1n2n2n
 611  //17th bin: <4>_{4n|2n,n,n} = four4n2n1n1n
 612  //18th bin: --  EMPTY --
 613  //19th bin: <5>_{2n|n,n,n,n} = five2n1n1n1n1n
 614  //20th bin: <5>_{2n,2n|2n,n,n} = five2n2n2n1n1n
 615  //21st bin: <5>_{3n,n|2n,n,n} = five3n1n2n1n1n
 616  //22nd bin: <5>_{4n|n,n,n,n} = five4n1n1n1n1n
 617  //23rd bin: --  EMPTY --
 618  //24th bin: <6>_{n,n,n|n,n,n} = six1n1n1n1n1n1n
 619  //25th bin: <6>_{2n,n,n|2n,n,n} = six2n1n1n2n1n1n
 620  //26th bin: <6>_{2n,2n|n,n,n,n} = six2n2n1n1n1n1n
 621  //27th bin: <6>_{3n,n|n,n,n,n} = six3n1n1n1n1n1n
 622  //28th bin: --  EMPTY --
 623  //29th bin: <7>_{2n,n,n|n,n,n,n} = seven2n1n1n1n1n1n1n
 624  //30th bin: --  EMPTY --
 625  //31st bin: <8>_{n,n,n,n|n,n,n,n} = eight1n1n1n1n1n1n1n1n
 626
 627
 628  // binning of fQProduct (all correlations are evaluated in harmonic n):
 629  // 1st bin: <2>*<4>
 630  // 2nd bin: <2>*<6>
 631  // 3rd bin: <2>*<8>
 632  // 4th bin: <4>*<6>
 633  // 5th bin: <4>*<8>
 634  // 6th bin: <6>*<8>
 635
 636  // 2-particle
 637  Double_t two1n1n = 0., two2n2n = 0., two3n3n = 0., two4n4n = 0.;
 638
 639  if(dMult>1)
 640  {
 641   //fill the common control histogram (2nd order):
 642   fCommonHists2nd->FillControlHistograms(anEvent);
 643
 644   two1n1n = (pow(afvQvector1n.Mod(),2.)-dMult)/(dMult*(dMult-1.)); // <2>_{n|n} = <cos(n*(phi1-phi2))>
 645   two2n2n = (pow(afvQvector2n.Mod(),2.)-dMult)/(dMult*(dMult-1.)); // <2>_{2n|2n} = <cos(2n*(phi1-phi2))>
 646   two3n3n = (pow(afvQvector3n.Mod(),2.)-dMult)/(dMult*(dMult-1.)); // <2>_{3n|3n} = <cos(3n*(phi1-phi2))>
 647   two4n4n = (pow(afvQvector4n.Mod(),2.)-dMult)/(dMult*(dMult-1.)); // <2>_{4n|4n} = <cos(4n*(phi1-phi2))>
 648
 649   fQCorrelations->Fill(0.,two1n1n,dMult*(dMult-1.));
 650   fQCorrelations->Fill(1.,two2n2n,dMult*(dMult-1.));
 651   fQCorrelations->Fill(2.,two3n3n,dMult*(dMult-1.));
 652   fQCorrelations->Fill(3.,two4n4n,dMult*(dMult-1.));
 653
 654   f2pDistribution->Fill(two1n1n,dMult*(dMult-1.));
 655  }
 656
 657  // 3-particle
 658  Double_t three2n1n1n=0., three3n2n1n=0., three4n2n2n=0., three4n3n1n=0.;
 659  if(dMult>2)
 660  {
 661   three2n1n1n = (reQ2nQ1nstarQ1nstar-2.*pow(afvQvector1n.Mod(),2.)-pow(afvQvector2n.Mod(),2.)+2.*dMult)/(dMult*(dMult-1.)*(dMult-2.)); //Re[<3>_{2n|n,n}] = Re[<3>_{n,n|2n}] = <cos(n*(2.*phi1-phi2-phi3))>
 662   three3n2n1n = (reQ3nQ2nstarQ1nstar-pow(afvQvector3n.Mod(),2.)-pow(afvQvector2n.Mod(),2.)-pow(afvQvector1n.Mod(),2.)+2.*dMult)/(dMult*(dMult-1.)*(dMult-2.)); //Re[<3>_{3n|2n,n}] = Re[<3>_{2n,n|3n}] = <cos(n*(3.*phi1-2.*phi2-phi3))>
 663   three4n2n2n = (reQ4nQ2nstarQ2nstar-2.*pow(afvQvector2n.Mod(),2.)-pow(afvQvector4n.Mod(),2.)+2.*dMult)/(dMult*(dMult-1.)*(dMult-2.)); //Re[<3>_{4n|2n,2n}] = Re[<3>_{2n,2n|4n}] = <cos(n*(4.*phi1-2.*phi2-2.*phi3))>
 664   three4n3n1n = (reQ4nQ3nstarQ1nstar-pow(afvQvector4n.Mod(),2.)-pow(afvQvector3n.Mod(),2.)-pow(afvQvector1n.Mod(),2.)+2.*dMult)/(dMult*(dMult-1.)*(dMult-2.)); //Re[<3>_{4n|3n,n}] = Re[<3>_{3n,n|4n}] = <cos(n*(4.*phi1-3.*phi2-phi3))>
 665
 666   fQCorrelations->Fill(5.,three2n1n1n,dMult*(dMult-1.)*(dMult-2.));
 667   fQCorrelations->Fill(6.,three3n2n1n,dMult*(dMult-1.)*(dMult-2.));
 668   fQCorrelations->Fill(7.,three4n2n2n,dMult*(dMult-1.)*(dMult-2.));
 669   fQCorrelations->Fill(8.,three4n3n1n,dMult*(dMult-1.)*(dMult-2.));
 670  }
 671
 672  //4-particle
 673  Double_t four1n1n1n1n=0., four2n2n2n2n=0., four2n1n2n1n=0., four3n1n1n1n=0., four4n2n1n1n=0., four3n1n2n2n=0., four3n1n3n1n=0.;
 674  if(dMult>3)
 675  {
 676   //fill the common control histogram (4th order):
 677   fCommonHists4th->FillControlHistograms(anEvent);
 678
 679   four1n1n1n1n = (2.*dMult*(dMult-3.)+pow(afvQvector1n.Mod(),4.)-4.*(dMult-2.)*pow(afvQvector1n.Mod(),2.)-2.*reQ2nQ1nstarQ1nstar+pow(afvQvector2n.Mod(),2.))/(dMult*(dMult-1)*(dMult-2.)*(dMult-3.));//<4>_{n,n|n,n}
 680   four2n2n2n2n = (2.*dMult*(dMult-3.)+pow(afvQvector2n.Mod(),4.)-4.*(dMult-2.)*pow(afvQvector2n.Mod(),2.)-2.*reQ4nQ2nstarQ2nstar+pow(afvQvector4n.Mod(),2.))/(dMult*(dMult-1)*(dMult-2.)*(dMult-3.));//<4>_{2n,2n|2n,2n}
 681   four2n1n2n1n = (xQ2nQ1nQ2nstarQ1nstar-2.*reQ3nQ2nstarQ1nstar-2.*reQ2nQ1nstarQ1nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.))-((dMult-5.)*pow(afvQvector1n.Mod(),2.)+(dMult-4.)*pow(afvQvector2n.Mod(),2.)-pow(afvQvector3n.Mod(),2.))/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.))+(dMult-6.)/((dMult-1.)*(dMult-2.)*(dMult-3.));//Re[<4>_{2n,n|2n,n}]
 682   four3n1n1n1n = (reQ3nQ1nstarQ1nstarQ1nstar-3.*reQ3nQ2nstarQ1nstar-3.*reQ2nQ1nstarQ1nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.))+(2.*pow(afvQvector3n.Mod(),2.)+3.*pow(afvQvector2n.Mod(),2.)+6.*pow(afvQvector1n.Mod(),2.)-6.*dMult)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.));//Re[<4>_{3n|n,n,n}]
 683   four4n2n1n1n = (reQ4nQ2nstarQ1nstarQ1nstar-2.*reQ4nQ3nstarQ1nstar-reQ4nQ2nstarQ2nstar-2.*reQ3nQ2nstarQ1nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.))-(reQ2nQ1nstarQ1nstar-2.*pow(afvQvector4n.Mod(),2.)-2.*pow(afvQvector3n.Mod(),2.)-3.*pow(afvQvector2n.Mod(),2.)-4.*pow(afvQvector1n.Mod(),2.))/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.))-(6.)/((dMult-1.)*(dMult-2.)*(dMult-3.));//Re[<4>_{4n|2n,n,n}]
 684   four3n1n2n2n = (reQ3nQ1nQ2nstarQ2nstar-reQ4nQ2nstarQ2nstar-reQ3nQ1nQ4nstar-2.*reQ3nQ2nstarQ1nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.))-(2.*reQ1nQ1nQ2nstar-pow(afvQvector4n.Mod(),2.)-2.*pow(afvQvector3n.Mod(),2.)-4.*pow(afvQvector2n.Mod(),2.)-4.*pow(afvQvector1n.Mod(),2.))/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.))-(6.)/((dMult-1.)*(dMult-2.)*(dMult-3.));//Re[<4>_{3n,n|2n,2n}]
 685   four3n1n3n1n = (pow(afvQvector3n.Mod(),2.)*pow(afvQvector1n.Mod(),2.)-2.*reQ4nQ3nstarQ1nstar-2.*reQ3nQ2nstarQ1nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.))+(pow(afvQvector4n.Mod(),2.)-(dMult-4.)*pow(afvQvector3n.Mod(),2.)+pow(afvQvector2n.Mod(),2.)-(dMult-4.)*pow(afvQvector1n.Mod(),2.))/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.))+(dMult-6.)/((dMult-1.)*(dMult-2.)*(dMult-3.));//<4>_{3n,n|3n,n}
 686   //four_3n1n3n1n = Q3nQ1nQ3nstarQ1nstar/(M*(M-1.)*(M-2.)*(M-3.))-(2.*three_3n2n1n+2.*three_4n3n1n)/(M-3.)-(two_4n4n+M*two_3n3n+two_2n2n+M*two_1n1n)/((M-2.)*(M-3.))-M/((M-1.)*(M-2.)*(M-3.));//<4>_{3n,n|3n,n}
 687
 688   fQCorrelations->Fill(10.,four1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.));
 689   fQCorrelations->Fill(11.,four2n1n2n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.));
 690   fQCorrelations->Fill(12.,four2n2n2n2n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.));
 691   fQCorrelations->Fill(13.,four3n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.));
 692   fQCorrelations->Fill(14.,four3n1n3n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.));
 693   fQCorrelations->Fill(15.,four3n1n2n2n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.));
 694   fQCorrelations->Fill(16.,four4n2n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.));
 695
 696   f4pDistribution->Fill(four1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.));
 697
 698   fQProduct->Fill(0.,two1n1n*four1n1n1n1n,dMult*(dMult-1.)*dMult*(dMult-1.)*(dMult-2.)*(dMult-3.));
 699  }
 700
 701  //5-particle
 702  Double_t five2n1n1n1n1n=0., five2n2n2n1n1n=0., five3n1n2n1n1n=0., five4n1n1n1n1n=0.;
 703  if(dMult>4)
 704  {
 705   five2n1n1n1n1n = (reQ2nQ1nQ1nstarQ1nstarQ1nstar-reQ3nQ1nstarQ1nstarQ1nstar+6.*reQ3nQ2nstarQ1nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))-(reQ2nQ1nQ3nstar+3.*(dMult-6.)*reQ2nQ1nstarQ1nstar+3.*reQ1nQ1nQ2nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))-(2.*pow(afvQvector3n.Mod(),2.)+3.*pow(afvQvector2n.Mod()*afvQvector1n.Mod(),2.)-3.*(dMult-4.)*pow(afvQvector2n.Mod(),2.))/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))-3.*(pow(afvQvector1n.Mod(),4.)-2.*(2*dMult-5.)*pow(afvQvector1n.Mod(),2.)+2.*dMult*(dMult-4.))/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.));//Re[<5>_{2n,n|n,n,n}]
 706
 707   five2n2n2n1n1n = (reQ2nQ2nQ2nstarQ1nstarQ1nstar-reQ4nQ2nstarQ1nstarQ1nstar-2.*reQ2nQ2nQ3nstarQ1nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))+2.*(reQ4nQ2nstarQ2nstar+4.*reQ3nQ2nstarQ1nstar+reQ3nQ1nQ4nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))+(reQ2nQ2nQ4nstar-2.*(dMult-5.)*reQ2nQ1nstarQ1nstar+2.*reQ1nQ1nQ2nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))-(2.*pow(afvQvector4n.Mod(),2.)+4.*pow(afvQvector3n.Mod(),2.)+1.*pow(afvQvector2n.Mod(),4.)-2.*(3.*dMult-10.)*pow(afvQvector2n.Mod(),2.))/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))-(4.*pow(afvQvector1n.Mod(),2.)*pow(afvQvector2n.Mod(),2.)-4.*(dMult-5.)*pow(afvQvector1n.Mod(),2.)+4.*dMult*(dMult-6.))/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.));//Re[<5>_{2n,2n|2n,n,n}]
 708
 709   //five_2n2n2n1n1n = reQ2nQ2nQ2nstarQ1nstarQ1nstar/(M*(M-1.)*(M-2.)*(M-3.)*(M-4.))-(4.*four_2n1n2n1n+2.*four_3n1n2n2n+1.*four_2n2n2n2n+four_4n2n1n1n)/(M-4.)-(2.*three_4n3n1n+three_4n2n2n+three_4n2n2n+2.*three_3n2n1n)/((M-3.)*(M-4.))-(4.*three_3n2n1n+(2.*M-1.)*three_2n1n1n+2.*three_2n1n1n)/((M-3.)*(M-4.))-(two_4n4n+2.*two_3n3n+4.*(M-1.)*two_2n2n+2.*(2.*M-1.)*two_1n1n)/((M-2.)*(M-3.)*(M-4.))-(2.*M-1.)/((M-1.)*(M-2.)*(M-3.)*(M-4.)); //OK!
 710
 711   five4n1n1n1n1n = (reQ4nQ1nstarQ1nstarQ1nstarQ1nstar-6.*reQ4nQ2nstarQ1nstarQ1nstar-4.*reQ3nQ1nstarQ1nstarQ1nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))+(8.*reQ4nQ3nstarQ1nstar+3.*reQ4nQ2nstarQ2nstar+12.*reQ3nQ2nstarQ1nstar+12.*reQ2nQ1nstarQ1nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))-(6.*pow(afvQvector4n.Mod(),2.)+8.*pow(afvQvector3n.Mod(),2.)+12.*pow(afvQvector2n.Mod(),2.)+24.*pow(afvQvector1n.Mod(),2.)-24.*dMult)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.));//Re[<5>_{4n|n,n,n,n}]
 712
 713   //five_4n1n1n1n1n = reQ4nQ1nstarQ1nstarQ1nstarQ1nstar/(M*(M-1.)*(M-2.)*(M-3.)*(M-4.)) -  (4.*four_3n1n1n1n+6.*four_4n2n1n1n)/(M-4.)  -  (6.*three_2n1n1n  + 12.*three_3n2n1n + 4.*three_4n3n1n + 3.*three_4n2n2n)/((M-3.)*(M-4.))  -  (4.*two_1n1n + 6.*two_2n2n + 4.*two_3n3n + 1.*two_4n4n)/((M-2.)*(M-3.)*(M-4.)) - 1./((M-1.)*(M-2.)*(M-3.)*(M-4.)); //OK!
 714
 715   five3n1n2n1n1n = (reQ3nQ1nQ2nstarQ1nstarQ1nstar-reQ4nQ2nstarQ1nstarQ1nstar-reQ3nQ1nstarQ1nstarQ1nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))-(reQ3nQ1nQ2nstarQ2nstar-3.*reQ4nQ3nstarQ1nstar-reQ4nQ2nstarQ2nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))-((2.*dMult-13.)*reQ3nQ2nstarQ1nstar-reQ3nQ1nQ4nstar-9.*reQ2nQ1nstarQ1nstar)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))-(2.*reQ1nQ1nQ2nstar+2.*pow(afvQvector4n.Mod(),2.)-2.*(dMult-5.)*pow(afvQvector3n.Mod(),2.)+2.*pow(afvQvector3n.Mod(),2.)*pow(afvQvector1n.Mod(),2.))/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))+(2.*(dMult-6.)*pow(afvQvector2n.Mod(),2.)-2.*pow(afvQvector2n.Mod(),2.)*pow(afvQvector1n.Mod(),2.)-pow(afvQvector1n.Mod(),4.)+2.*(3.*dMult-11.)*pow(afvQvector1n.Mod(),2.))/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.))-4.*(dMult-6.)/((dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.));//Re[<5>_{3n,n|2n,n,n}]
 716
 717   //five3n1n2n1n1n = reQ3nQ1nQ2nstarQ1nstarQ1nstar/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)) -  (four4n2n1n1n+four1n1n1n1n+four3n1n1n1n+2.*four2n1n2n1n+2.*three3n2n1n+2.*four3n1n3n1n+four3n1n2n2n)/(dMult-4.)  -   (2.*three4n3n1n+three4n2n2n+6.*three3n2n1n+three4n3n1n+2.*three3n2n1n+3.*three2n1n1n+2.*three2n1n1n+4.*two1n1n+2.*two2n2n+2.*two3n3n)/((dMult-3.)*(dMult-4.))  -  (5.*two1n1n + 4.*two2n2n + 3.*two3n3n + 1.*two4n4n + 2.)/((dMult-2.)*(dMult-3.)*(dMult-4.))  - 1./((dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)) ;//Re[<5>_{3n,n|2n,n,n}] //OK!
 718
 719   //five3n1n2n1n1n = reQ3nQ1nQ2nstarQ1nstarQ1nstar/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)) -  (four4n2n1n1n+four1n1n1n1n+four3n1n1n1n+2.*four2n1n2n1n+2.*four3n1n3n1n+four3n1n2n2n)/(dMult-4.)  -      (2.*three4n3n1n+three4n2n2n+2.*dMult*three3n2n1n+three4n3n1n+2.*three3n2n1n+3.*three2n1n1n+2.*three2n1n1n)/((dMult-3.)*(dMult-4.))  -  ((4.*dMult-3.)*two1n1n + 2.*dMult*two2n2n + (2.*dMult-1.)*two3n3n + two4n4n)/((dMult-2.)*(dMult-3.)*(dMult-4.))  - (2.*dMult-1.)/((dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)) ;//Re[<5>_{3n,n|2n,n,n}] //OK!
 720
 721   fQCorrelations->Fill(18.,five2n1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.));
 722   fQCorrelations->Fill(19.,five2n2n2n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.));
 723   fQCorrelations->Fill(20.,five3n1n2n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.));
 724   fQCorrelations->Fill(21.,five4n1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.));
 725  }
 726
 727  //xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
 728  //      !!!! to be removed: temporary fix !!!!
 729  if(dMult>1)
 730  {
 731   two1n1n = (pow(afvQvector1n.Mod(),2.)-dMult)/(dMult*(dMult-1.));
 732  }
 733  //xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
 734
 735  //6-particle
 736  Double_t six1n1n1n1n1n1n=0., six2n2n1n1n1n1n=0., six3n1n1n1n1n1n=0., six2n1n1n2n1n1n=0.;
 737  if(dMult>5)
 738  {
 739   //fill the common control histogram (6th order):
 740   fCommonHists6th->FillControlHistograms(anEvent);
 741
 742   six1n1n1n1n1n1n = (pow(afvQvector1n.Mod(),6.)+9.*xQ2nQ1nQ2nstarQ1nstar-6.*reQ2nQ1nQ1nstarQ1nstarQ1nstar)/(dMult*(dMult-1)*(dMult-2)*(dMult-3)*(dMult-4)*(dMult-5))+4.*(reQ3nQ1nstarQ1nstarQ1nstar-3.*reQ3nQ2nstarQ1nstar)/(dMult*(dMult-1)*(dMult-2)*(dMult-3)*(dMult-4)*(dMult-5))+2.*(9.*(dMult-4.)*reQ2nQ1nstarQ1nstar+2.*pow(afvQvector3n.Mod(),2.))/(dMult*(dMult-1)*(dMult-2)*(dMult-3)*(dMult-4)*(dMult-5))-9.*(pow(afvQvector1n.Mod(),4.)+pow(afvQvector2n.Mod(),2.))/(dMult*(dMult-1)*(dMult-2)*(dMult-3)*(dMult-5))+(18.*pow(afvQvector1n.Mod(),2.))/(dMult*(dMult-1)*(dMult-3)*(dMult-4))-(6.)/((dMult-1)*(dMult-2)*(dMult-3));//<6>_{n,n,n|n,n,n}
 743
 744   six2n1n1n2n1n1n = (xQ2nQ1nQ1nQ2nstarQ1nstarQ1nstar-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(2.*five2n2n2n1n1n+4.*five2n1n1n1n1n+4.*five3n1n2n1n1n+4.*four2n1n2n1n+1.*four1n1n1n1n)-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(4.*four1n1n1n1n+4.*two1n1n+2.*three2n1n1n+2.*three2n1n1n+4.*four3n1n1n1n+8.*three2n1n1n+2.*four4n2n1n1n+4.*four2n1n2n1n+2.*two2n2n+8.*four2n1n2n1n+4.*four3n1n3n1n+8.*three3n2n1n+4.*four3n1n2n2n+4.*four1n1n1n1n+4.*four2n1n2n1n+1.*four2n2n2n2n)-dMult*(dMult-1.)*(dMult-2.)*(2.*three2n1n1n+8.*two1n1n+4.*two1n1n+2.+4.*two1n1n+4.*three2n1n1n+2.*two2n2n+4.*three2n1n1n+8.*three3n2n1n+8.*two2n2n+4.*three4n3n1n+4.*two3n3n+4.*three3n2n1n+4.*two1n1n+8.*three2n1n1n+4.*two1n1n+4.*three3n2n1n+4.*three2n1n1n+2.*two2n2n+4.*three3n2n1n+2.*three4n2n2n)-dMult*(dMult-1.)*(4.*two1n1n+4.+4.*two1n1n+2.*two2n2n+1.+4.*two1n1n+4.*two2n2n+4.*two3n3n+   1.+2.*two2n2n+1.*two4n4n)-dMult)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.));//<6>_{2n,n,n|2n,n,n}
 745
 746   six2n2n1n1n1n1n = (reQ2nQ2nQ1nstarQ1nstarQ1nstarQ1nstar-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(five4n1n1n1n1n+8.*five2n1n1n1n1n+6.*five2n2n2n1n1n)-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(4.*four3n1n1n1n+6.*four4n2n1n1n+12.*three2n1n1n+12.*four1n1n1n1n+24.*four2n1n2n1n+4.*four3n1n2n2n+3.*four2n2n2n2n)-dMult*(dMult-1.)*(dMult-2.)*(6.*three2n1n1n+12.*three3n2n1n+4.*three4n3n1n+3.*three4n2n2n+8.*three2n1n1n+24.*two1n1n+12.*two2n2n+12.*three2n1n1n+8.*three3n2n1n+1.*three4n2n2n)-dMult*(dMult-1.)*(4.*two1n1n+6.*two2n2n+4.*two3n3n+1.*two4n4n+2.*two2n2n+8.*two1n1n+6.)-dMult)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.));//<6>_{2n,2n,n|n,n,n}
 747
 748   six3n1n1n1n1n1n = (reQ3nQ1nQ1nstarQ1nstarQ1nstarQ1nstar-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(five4n1n1n1n1n+4.*five2n1n1n1n1n+6.*five3n1n2n1n1n+4.*four3n1n1n1n)-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(4.*four3n1n1n1n+6.*four4n2n1n1n+6.*four1n1n1n1n+12.*three2n1n1n+12.*four2n1n2n1n+6.*four3n1n1n1n+12.*three3n2n1n+4.*four3n1n3n1n+3.*four3n1n2n2n)-dMult*(dMult-1.)*(dMult-2.)*(6.*three2n1n1n+12.*three3n2n1n+4.*three4n3n1n+3.*three4n2n2n+4.*two1n1n+12.*two1n1n+6.*three2n1n1n+12.*three2n1n1n+4.*three3n2n1n+12.*two2n2n+4.*three3n2n1n+4.*two3n3n+1.*three4n3n1n+6.*three3n2n1n)-dMult*(dMult-1.)*(4.*two1n1n+6.*two2n2n+4.*two3n3n+1.*two4n4n+1.*two1n1n+4.+6.*two1n1n+4.*two2n2n+1.*two3n3n)-dMult)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.));//<6>_{3n,n|n,n,n,n}
 749
 750   fQCorrelations->Fill(23.,six1n1n1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.));
 751   fQCorrelations->Fill(24.,six2n1n1n2n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.));
 752   fQCorrelations->Fill(25.,six2n2n1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.));
 753   fQCorrelations->Fill(26.,six3n1n1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.));
 754
 755   f6pDistribution->Fill(six1n1n1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.));
 756
 757   fQProduct->Fill(1.,two1n1n*six1n1n1n1n1n1n,dMult*(dMult-1.)*dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.));
 758   fQProduct->Fill(3.,four1n1n1n1n*six1n1n1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.));
 759  }
 760
 761  //7-particle
 762  Double_t seven2n1n1n1n1n1n1n=0.;
 763  if(dMult>6)
 764  {
 765   seven2n1n1n1n1n1n1n = (reQ2nQ1nQ1nQ1nstarQ1nstarQ1nstarQ1nstar-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.)*(2.*six3n1n1n1n1n1n+4.*six1n1n1n1n1n1n+1.*six2n2n1n1n1n1n+6.*six2n1n1n2n1n1n+8.*five2n1n1n1n1n)-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(1.*five4n1n1n1n1n +8.*five2n1n1n1n1n+8.*four3n1n1n1n+12.*five3n1n2n1n1n+4.*five2n1n1n1n1n+3.*five2n2n2n1n1n+6.*five2n2n2n1n1n+6.*four1n1n1n1n+24.*four1n1n1n1n+12.*five2n1n1n1n1n+12.*five2n1n1n1n1n+12.*three2n1n1n+24.*four2n1n2n1n+4.*five3n1n2n1n1n+4.*five2n1n1n1n1n)-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(4.*four3n1n1n1n+6.*four4n2n1n1n+12.*four1n1n1n1n+24.*three2n1n1n+24.*four2n1n2n1n+12.*four3n1n1n1n+24.*three3n2n1n+8.*four3n1n3n1n+6.*four3n1n2n2n+6.*three2n1n1n+12.*four1n1n1n1n+12.*four2n1n2n1n+6.*three2n1n1n+12.*four2n1n2n1n+4.*four3n1n2n2n+3.*four2n2n2n2n+4.*four1n1n1n1n+6.*three2n1n1n+24.*two1n1n+24.*four1n1n1n1n+4.*four3n1n1n1n+24.*two1n1n+24.*three2n1n1n+12.*two2n2n+24.*three2n1n1n+12.*four2n1n2n1n+8.*three3n2n1n+8.*four2n1n2n1n+1.*four4n2n1n1n)-dMult*(dMult-1.)*(dMult-2.)*(6.*three2n1n1n+1.*three2n1n1n+8.*two1n1n+12.*three3n2n1n+24.*two1n1n+12.*three2n1n1n+4.*three2n1n1n+8.*two1n1n+4.*three4n3n1n+24.*three2n1n1n+8.*three3n2n1n+12.*two1n1n+12.*two1n1n+3.*three4n2n2n+24.*two2n2n+6.*two2n2n+12.+12.*three3n2n1n+8.*two3n3n+12.*three2n1n1n+24.*two1n1n+4.*three3n2n1n+8.*three3n2n1n+2.*three4n3n1n+12.*two1n1n+8.*three2n1n1n+4.*three2n1n1n+2.*three3n2n1n+6.*two2n2n+8.*two2n2n+1.*three4n2n2n+4.*three3n2n1n+6.*three2n1n1n)-dMult*(dMult-1.)*(4.*two1n1n+2.*two1n1n+6.*two2n2n+8.+1.*two2n2n+4.*two3n3n+12.*two1n1n+4.*two1n1n+1.*two4n4n+8.*two2n2n+6.+2.*two3n3n+4.*two1n1n+1.*two2n2n)-dMult)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.)*(dMult-6.));
 766
 767   fQCorrelations->Fill(28.,seven2n1n1n1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.)*(dMult-6.));
 768  }
 769
 770  //8-particle
 771  Double_t eight1n1n1n1n1n1n1n1n=0.;
 772  if(dMult>7)
 773  {
 774   //fill the common control histogram (8th order):
 775   fCommonHists8th->FillControlHistograms(anEvent);
 776
 777   eight1n1n1n1n1n1n1n1n = (pow(afvQvector1n.Mod(),8)-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.)*(dMult-6.)*(12.*seven2n1n1n1n1n1n1n+16.*six1n1n1n1n1n1n)-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.)*(8.*six3n1n1n1n1n1n+48.*six1n1n1n1n1n1n+6.*six2n2n1n1n1n1n+96.*five2n1n1n1n1n+72.*four1n1n1n1n+36.*six2n1n1n2n1n1n)-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(2.*five4n1n1n1n1n+32.*five2n1n1n1n1n+36.*four1n1n1n1n+32.*four3n1n1n1n+48.*five2n1n1n1n1n+48.*five3n1n2n1n1n+144.*five2n1n1n1n1n+288.*four1n1n1n1n+36.*five2n2n2n1n1n+144.*three2n1n1n+96.*two1n1n+144.*four2n1n2n1n)-dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(8.*four3n1n1n1n+48.*four1n1n1n1n+12.*four4n2n1n1n+96.*four2n1n2n1n+96.*three2n1n1n+72.*three2n1n1n+144.*two1n1n+16.*four3n1n3n1n+48.*four3n1n1n1n+144.*four1n1n1n1n+72.*four1n1n1n1n+96.*three3n2n1n+24.*four3n1n2n2n+144.*four2n1n2n1n+288.*two1n1n+288.*three2n1n1n+9.*four2n2n2n2n+72.*two2n2n+24.)-dMult*(dMult-1.)*(dMult-2.)*(12.*three2n1n1n+16.*two1n1n+24.*three3n2n1n+48.*three2n1n1n+96.*two1n1n+8.*three4n3n1n+32.*three3n2n1n+96.*three2n1n1n+144.*two1n1n+6.*three4n2n2n+96.*two2n2n+36.*two2n2n+72.+48.*three3n2n1n+16.*two3n3n+72.*three2n1n1n+144.*two1n1n)-dMult*(dMult-1.)*(8.*two1n1n+12.*two2n2n+16.+8.*two3n3n+48.*two1n1n+1.*two4n4n+16.*two2n2n+18.)-dMult)/(dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.)*(dMult-6.)*(dMult-7.));
 778
 779   fQCorrelations->Fill(30.,eight1n1n1n1n1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.)*(dMult-6.)*(dMult-7.));
 780
 781   f8pDistribution->Fill(eight1n1n1n1n1n1n1n1n,dMult*(dMult-1.)*(dMult-2.)*(dMult-3.)*(dMult-4.)*(dMult-5.)*(dMult-6.)*(dMult-7.));
 782  }
 783  //---------------------------------------------------------------------------------------------------------
 784
 785
 786
 787
 788  //---------------------------------------------------------------------------------------------------------
 789  // weights:
 790  Bool_t useWeights = fUsePhiWeights||fUsePtWeights||fUseEtaWeights;
 791
 792  TH1F *phiWeights = NULL; // histogram with phi weights
 793  TH1D *ptWeights  = NULL; // histogram with pt weights
 794  TH1D *etaWeights = NULL; // histogram with eta weights
 795
 796  if(useWeights)
 797  {
 798   if(!fWeightsList)
 799   {
 800    cout<<" WARNING: fWeightsList is NULL pointer. "<<endl;
 801    exit(0);
 802   }
 803   if(fUsePhiWeights)
 804   {
 805    phiWeights = dynamic_cast<TH1F *>(fWeightsList->FindObject("phi_weights"));
 806    if(!phiWeights)
 807    {
 808     cout<<" WARNING: couldn't access the histogram with phi weights. "<<endl;
 809     exit(0);
 810    }
 811   }
 812   if(fUsePtWeights)
 813   {
 814    ptWeights = dynamic_cast<TH1D *>(fWeightsList->FindObject("pt_weights"));
 815    if(!ptWeights)
 816    {
 817     cout<<" WARNING: couldn't access the histogram with pt weights. "<<endl;
 818     exit(0);
 819    }
 820   }
 821   if(fUseEtaWeights)
 822   {
 823    etaWeights = dynamic_cast<TH1D *>(fWeightsList->FindObject("eta_weights"));
 824    if(!etaWeights)
 825    {
 826     cout<<" WARNING: couldn't access the histogram with eta weights. "<<endl;
 827     exit(0);
 828    }
 829   }
 830  }
 831
 832  Int_t nBinsPhi = 0;
 833  Double_t dBinWidthPt=0.;
 834  Double_t dBinWidthEta=0.;
 835
 836  if(fnBinsPt)
 837  {
 838   dBinWidthPt=(fPtMax-fPtMin)/fnBinsPt;
 839  }
 840  if(fnBinsEta)
 841  {
 842   dBinWidthEta=(fEtaMax-fEtaMin)/fnBinsEta;
 843  }
 844
 845  if(fWeightsList)
 846   {
 847    if(fUsePhiWeights)
 848    {
 849     if(phiWeights) nBinsPhi = phiWeights->GetNbinsX();
 850    }
 851    if(fUsePtWeights)
 852    {
 853     if(ptWeights)
 854     {
 855      Double_t dBinWidthPtW = ptWeights->GetBinWidth(1); // assuming that all bins have the same width
 856      if(dBinWidthPtW != dBinWidthPt)
 857      {
 858       cout<<" WARNING: dBinWidthPtW != dBinWidthPt in AFAWQC::Make()"<<endl;
 859       exit(0);
 860      }
 861      Double_t dPtMinW = (ptWeights->GetXaxis())->GetXmin();
 862      if(dPtMinW != fPtMin)
 863      {
 864       cout<<" WARNING: dPtMinW != fPtMin in AFAWQC::Make()"<<endl;
 865       exit(0);
 866      }
 867     }
 868    }
 869    if(fUseEtaWeights)
 870    {
 871     if(etaWeights)
 872     {
 873      Double_t dBinWidthEtaW = etaWeights->GetBinWidth(1); // assuming that all bins have the same width
 874      if(dBinWidthEtaW != dBinWidthEta)
 875      {
 876       cout<<" WARNING: dBinWidthEtaW != dBinWidthEta in AFAWQC::Make()"<<endl;
 877       exit(0);
 878      }
 879      Double_t dEtaMinW = (etaWeights->GetXaxis())->GetXmin();
 880      if(dEtaMinW != fEtaMin)
 881      {
 882       cout<<" WARNING: dEtaMinW != fEtaMin in AFAWQC::Make()"<<endl;
 883       exit(0);
 884      }
 885     }
 886    }
 887   } // end of if(weightsList)
 888
 889  Double_t wPhi = 1.; // phi weight
 890  Double_t wPt  = 1.; // pt weight
 891  Double_t wEta = 1.; // eta weight
 892
 893  Double_t dPhi = 0.;
 894  Double_t dPt  = 0.;
 895  Double_t dEta = 0.;
 896
 897  Double_t dQnkX[4][8] = {{0.}}; // sum_{i=1}^{M} w_i^k cos(n phi_i)
 898  Double_t dQnkY[4][8] = {{0.}}; // sum_{i=1}^{M} w_i^k sin(n phi_i)
 899  Double_t dSnk[4][8] = {{0.}};  //(sum_{i=1}^{M} w_i^k)^n
 900
 901   for(Int_t i=0;i<nPrim;i++) // loop over all particles
 902   {
 903    fTrack=anEvent->GetTrack(i);
 904    if(fTrack)
 905    {
 906     if(fTrack->InRPSelection())
 907     {
 908      dPhi = fTrack->Phi();
 909      dPt  = fTrack->Pt();
 910      dEta = fTrack->Eta();
 911
 912      // determine Phi weight:
 913      if(phiWeights && nBinsPhi)
 914      {
 915       wPhi = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*nBinsPhi/TMath::TwoPi())));
 916      }
 917      // determine pt weight:
 918      if(ptWeights && dBinWidthPt)
 919      {
 920       wPt = ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-fPtMin)/dBinWidthPt)));
 921      }
 922      // determine eta weight:
 923      if(etaWeights && dBinWidthEta)
 924      {
 925       wEta = etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-fEtaMin)/dBinWidthEta)));
 926      }
 927
 928      // (Q_{n,k})_x, (Q_{n,k})_y and S_{n,k}
 929      for(Int_t nn=0;nn<4;nn++)
 930      {
 931       for(Int_t k=0;k<8;k++)
 932       {
 933        dQnkX[nn][k]+=pow(wPhi*wPt*wEta,k+1)*TMath::Cos(2*(nn+1)*dPhi);
 934        dQnkY[nn][k]+=pow(wPhi*wPt*wEta,k+1)*TMath::Sin(2*(nn+1)*dPhi);
 935        dSnk[nn][k]+=pow(wPhi*wPt*wEta,k+1);
 936       }
 937      }
 938
 939     } // end of if (pTrack->InRPSelection())
 940    } // end of if (pTrack)
 941    else {cerr << "no particle!!!"<<endl;}
 942   } // loop over particles
 943
 944   for(Int_t nn=0;nn<4;nn++)
 945   {
 946    for(Int_t k=0;k<8;k++)
 947    {
 948     dSnk[nn][k]=pow(dSnk[nn][k],nn+1);
 949    }
 950   }
 951
 952  //..............................................................................................
 953  // Ms (introduced in order to simplify some Eqs. bellow)
 954  Double_t dM11 = dSnk[1][0]-dSnk[0][1]; // dM11 = sum_{i,j=1,i!=j}^M w_i w_j
 955  Double_t dM22 = dSnk[1][1]-dSnk[0][3]; // dM22 = sum_{i,j=1,i!=j}^M w_i^2 w_j^2
 956  Double_t dM33 = dSnk[1][2]-dSnk[0][5]; // dM33 = sum_{i,j=1,i!=j}^M w_i^3 w_j^3
 957  Double_t dM44 = dSnk[1][3]-dSnk[0][7]; // dM44 = sum_{i,j=1,i!=j}^M w_i^4 w_j^4
 958  Double_t dM31 = dSnk[0][2]*dSnk[0][0]-dSnk[0][3]; // dM31 = sum_{i,j=1,i!=j}^M w_i^3 w_j
 959  Double_t dM211 = dSnk[0][1]*dSnk[1][0]-2.*dSnk[0][2]*dSnk[0][0]-dSnk[1][1]+2.*dSnk[0][3]; // dM211 = sum_{i,j,k=1,i!=j!=k}^M w_i^2 w_j w_k
 960  Double_t dM1111 = dSnk[3][0]-6.*dM211-4.*dM31-3.*dM22-dSnk[0][3]; // dM1111 = sum_{i,j,k,l=1,i!=j!=k!=l}^M w_i w_j w_k w_l
 961  //..............................................................................................
 962
 963  //---------------------------------------------------------------------------------------------------------
 964  //
 965  //                                        ***********************************************
 966  //                                        **** weighted multi-particle correlations: ****
 967  //                                        ***********************************************
 968  //
 969  // Remark 1: weighted multi-particle correlations calculated with Q-vectors are stored in fWeightedQCorrelations.
 970  // Remark 2: binning of fWeightedQCorrelations is organized as follows:
 971
 972  // binning
 973  //..............................................................................................
 974  // 1st bin: weighted <2>_{n|n}   = <w1   w2   cos( n*(phi1-phi2))>
 975  // 2nd bin: weighted <2>_{2n|2n} = <w1^2 w2^2 cos(2n*(phi1-phi2))>
 976  // 3rd bin: weighted <2>_{3n|3n} = <w1^3 w2^3 cos(3n*(phi1-phi2))>
 977  // 4th bin: weighted <2>_{4n|4n} = <w1^4 w2^4 cos(4n*(phi1-phi2))>
 978  // 5th bin: weighted <2>_{n|n} = <w1^3 w2 cos(n*(phi1-phi2))>
 979  // 6th bin: weighted <2>_{n|n} = <w1 w2 w3^2 cos(n*(phi1-phi2))>
 980
 981  // 11th bin: weighted <3>_{2n|n,n} = <w1^2 w2 w3 cos(n*(2phi1-phi2-phi3))>
 982
 983  // 21st bin: weighted <4>_{n,n|n,n} = <w1 w2 w3 w4 cos(n*(phi1+phi2-phi3-phi4))>
 984  //..............................................................................................
 985
 986  //..............................................................................................
 987  // weighted 2-particle correlations:
 988  Double_t two1n1nW1W1=0., two2n2nW2W2=0., two3n3nW3W3=0., two4n4nW4W4=0., two1n1nW3W1=0., two1n1nW2W1W1=0.;
 989
 990  if(nRP>1) // nRP = number of particles used to determine the reaction plane
 991  {
 992   if(dM11 != 0)
 993   {
 994    two1n1nW1W1 = (pow(dQnkX[0][0],2)+pow(dQnkY[0][0],2)-dSnk[0][1])/dM11; // <2>_{n|n}=<w1 w2 cos(n*(phi1-phi2))>
 995    fWeightedQCorrelations->Fill(0.,two1n1nW1W1,dM11);
 996   }
 997   if(dM22 != 0)
 998   {
 999    two2n2nW2W2 = (pow(dQnkX[1][1],2)+pow(dQnkY[1][1],2)-dSnk[0][3])/dM22; // <2>_{2n|2n}=<w1^2 w2^2 cos(2n*(phi1-phi2))>
1000    fWeightedQCorrelations->Fill(1.,two2n2nW2W2,dM22);
1001   }
1002   if(dM33 != 0)
1003   {
1004    two3n3nW3W3 = (pow(dQnkX[2][2],2)+pow(dQnkY[2][2],2)-dSnk[0][5])/dM33; // <2>_{2n|2n}=<w1^3 w2^3 cos(3n*(phi1-phi2))>
1005    fWeightedQCorrelations->Fill(2.,two3n3nW3W3,dM33);
1006   }
1007   if(dM44 != 0)
1008   {
1009    two4n4nW4W4 = (pow(dQnkX[3][3],2)+pow(dQnkY[3][3],2)-dSnk[0][7])/dM44; // <2>_{2n|2n}=<w1^4 w2^4 cos(4n*(phi1-phi2))>
1010    fWeightedQCorrelations->Fill(3.,two4n4nW4W4,dM44);
1011   }
1012   if(dM31 != 0)
1013   {
1014    two1n1nW3W1 = (dQnkX[0][2]*dQnkX[0][0]+dQnkY[0][2]*dQnkY[0][0]-dSnk[0][3])/dM31; // <2>_{n|n}=<w1^3 w2 cos(n*(phi1-phi2))>
1015    fWeightedQCorrelations->Fill(4.,two1n1nW3W1,dM31);
1016   }
1017   if(dM211 != 0)
1018   {
1019    two1n1nW2W1W1 = (dSnk[0][1]*dM11*two1n1nW1W1-2.*dM31*two1n1nW3W1)/dM211; // <2>_{n|n}=<w1^2 w2 w3 cos(n*(phi1-phi2))>
1020    fWeightedQCorrelations->Fill(5.,two1n1nW2W1W1,dM211);
1021   }
1022  } // end of if(nRP>1)
1023  //..............................................................................................
1024
1025  //..............................................................................................
1026  // weighted 3-particle correlations:
1027  Double_t three2n1n1nW2W1W1=0.;
1028
1029  if(nRP>2) // nRP = number of particles used to determine the reaction plane
1030  {
1031   if(dM211 != 0)
1032   {
1033    three2n1n1nW2W1W1 = (pow(dQnkX[0][0],2.)*dQnkX[1][1]+2.*dQnkX[0][0]*dQnkY[0][0]*dQnkY[1][1]-pow(dQnkY[0][0],2.)*dQnkX[1][1]-2.*dM31*two1n1nW3W1-dM22*two2n2nW2W2-dSnk[0][3])/dM211;
1034    fWeightedQCorrelations->Fill(10.,three2n1n1nW2W1W1,dM211);
1035   }
1036  } // end of if(nRP>2)
1037  //..............................................................................................
1038
1039  //..............................................................................................
1040  // weighted 4-particle correlations:
1041  Double_t four1n1n1n1nW1W1W1W1=0.;
1042
1043  if(nRP>3) // nRP = number of particles used to determine the reaction plane
1044  {
1045   if(dM1111 != 0)
1046   {
1047    four1n1n1n1nW1W1W1W1 = (pow(pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.),2)-2.*dM211*three2n1n1nW2W1W1-4.*dM211*two1n1nW2W1W1-4.*dM31*two1n1nW3W1-dM22*two2n2nW2W2-2.*dM22-dSnk[0][3])/dM1111;
1048    fWeightedQCorrelations->Fill(20.,four1n1n1n1nW1W1W1W1,dM1111);
1049   }
1050  } // end of if(nRP>3)
1051  //..............................................................................................
1052
1053
1054
1055
1056  //---------------------------------------------------------------------------------------------------------
1057  //
1058  //                                        *******************************************************
1059  //                                        **** weighted reduced multi-particle correlations: ****
1060  //                                        *******************************************************
1061  //
1062  // pt POI
1063  TProfile *ptReq1nPrime = new TProfile("ptReq1nPrime","Re[q_{n}^{''}]",fnBinsPt,fPtMin,fPtMax,"s");
1064  TProfile *ptImq1nPrime = new TProfile("ptImq1nPrime","Im[q_{n}^{''}]",fnBinsPt,fPtMin,fPtMax,"s");
1065  TProfile *ptReq2nPrime = new TProfile("ptReq2nPrime","Re[q_{2n}^{''}]",fnBinsPt,fPtMin,fPtMax,"s");
1066  TProfile *ptImq2nPrime = new TProfile("ptImq2nPrime","Im[q_{2n}^{''}]",fnBinsPt,fPtMin,fPtMax,"s");
1067
1068  TProfile *ptReq1nPrimePrime = new TProfile("ptReq1nPrimePrime","Re[q_{n}^{''}]",fnBinsPt,fPtMin,fPtMax,"s");
1069  TProfile *ptImq1nPrimePrime = new TProfile("ptImq1nPrimePrime","Im[q_{n}^{''}]",fnBinsPt,fPtMin,fPtMax,"s");
1070  TProfile *ptReq2nPrimePrime = new TProfile("ptReq2nPrimePrime","Re[q_{2n}^{''}]",fnBinsPt,fPtMin,fPtMax,"s");
1071  TProfile *ptImq2nPrimePrime = new TProfile("ptImq2nPrimePrime","Im[q_{2n}^{''}]",fnBinsPt,fPtMin,fPtMax,"s");
1072
1073  TProfile *req1nW2PrimePrimePt = new TProfile("req1nW2PrimePrimePt","#sum_{i=1}^{m''} w_{i}^{2} cos(n(#phi_{i}))''",fnBinsPt,fPtMin,fPtMax,"s");
1074  TProfile *imq1nW2PrimePrimePt = new TProfile("imq1nW2PrimePrimePt","#sum_{i=1}^{m''} w_{i}^{2} sin(n(#phi_{i}))''",fnBinsPt,fPtMin,fPtMax,"s");
1075  TProfile *req2nW1PrimePrimePt = new TProfile("req2nW1PrimePrimePt","#sum_{i=1}^{m''} w_{i} cos(2n(#phi_{i}))''",fnBinsPt,fPtMin,fPtMax,"s");
1076  TProfile *imq2nW1PrimePrimePt = new TProfile("imq2nW1PrimePrimePt","#sum_{i=1}^{m''} w_{i} sin(2n(#phi_{i}))''",fnBinsPt,fPtMin,fPtMax,"s");
1077
1078  TProfile *sumOfW1upTomPrimePrimePt = new TProfile("sumOfW1upTomPrimePrimePt","#sum_{i=1}^{m''} w_{i}''",fnBinsPt,fPtMin,fPtMax,"s");
1079  TProfile *sumOfW2upTomPrimePrimePt = new TProfile("sumOfW2upTomPrimePrimePt","#sum_{i=1}^{m''} w_{i}^{2}''",fnBinsPt,fPtMin,fPtMax,"s");
1080  TProfile *sumOfW3upTomPrimePrimePt = new TProfile("sumOfW3upTomPrimePrimePt","#sum_{i=1}^{m''} w_{i}^{3}''",fnBinsPt,fPtMin,fPtMax,"s");
1081
1082  // eta POI
1083  TProfile *etaReq1nPrime = new TProfile("etaReq1nPrime","Re[q_{n}^{''}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1084  TProfile *etaImq1nPrime = new TProfile("etaImq1nPrime","Im[q_{n}^{''}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1085  TProfile *etaReq2nPrime = new TProfile("etaReq2nPrime","Re[q_{2n}^{''}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1086  TProfile *etaImq2nPrime = new TProfile("etaImq2nPrime","Im[q_{2n}^{''}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1087
1088  TProfile *etaReq1nPrimePrime = new TProfile("etaReq1nPrimePrime","Re[q_{n}^{''}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1089  TProfile *etaImq1nPrimePrime = new TProfile("etaImq1nPrimePrime","Im[q_{n}^{''}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1090  TProfile *etaReq2nPrimePrime = new TProfile("etaReq2nPrimePrime","Re[q_{2n}^{''}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1091  TProfile *etaImq2nPrimePrime = new TProfile("etaImq2nPrimePrime","Im[q_{2n}^{''}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1092
1093  TProfile *req1nW2PrimePrimeEta = new TProfile("req1nW2PrimePrimeEta","#sum_{i=1}^{m''} w_{i}^{2} cos(n(#phi_{i}))''",fnBinsEta,fEtaMin,fEtaMax,"s");
1094  TProfile *imq1nW2PrimePrimeEta = new TProfile("imq1nW2PrimePrimeEta","#sum_{i=1}^{m''} w_{i}^{2} sin(n(#phi_{i}))''",fnBinsEta,fEtaMin,fEtaMax,"s");
1095  TProfile *req2nW1PrimePrimeEta = new TProfile("req2nW1PrimePrimeEta","#sum_{i=1}^{m''} w_{i} cos(2n(#phi_{i}))''",fnBinsEta,fEtaMin,fEtaMax,"s");
1096  TProfile *imq2nW1PrimePrimeEta = new TProfile("imq2nW1PrimePrimeEta","#sum_{i=1}^{m''} w_{i} sin(2n(#phi_{i}))''",fnBinsEta,fEtaMin,fEtaMax,"s");
1097
1098  TProfile *sumOfW1upTomPrimePrimeEta = new TProfile("sumOfW1upTomPrimePrimeEta","#sum_{i=1}^{m''} w_{i}''",fnBinsEta,fEtaMin,fEtaMax,"s");
1099  TProfile *sumOfW2upTomPrimePrimeEta = new TProfile("sumOfW2upTomPrimePrimeEta","#sum_{i=1}^{m''} w_{i}^{2}''",fnBinsEta,fEtaMin,fEtaMax,"s");
1100  TProfile *sumOfW3upTomPrimePrimeEta = new TProfile("sumOfW3upTomPrimePrimeEta","#sum_{i=1}^{m''} w_{i}^{3}''",fnBinsEta,fEtaMin,fEtaMax,"s");
1101
1102  // pt RP
1103  TProfile *ptReq1n = new TProfile("ptReq1n","Re[q_{n}]",fnBinsPt,fPtMin,fPtMax,"s");
1104  TProfile *ptImq1n = new TProfile("ptImq1n","Im[q_{n}]",fnBinsPt,fPtMin,fPtMax,"s");
1105  TProfile *ptReq2n = new TProfile("ptReq2n","Re[q_{2n}]",fnBinsPt,fPtMin,fPtMax,"s");
1106  TProfile *ptImq2n = new TProfile("ptImq2n","Im[q_{2n}]",fnBinsPt,fPtMin,fPtMax,"s");
1107
1108  TProfile *req1nW2Pt = new TProfile("req1nW2Pt","#sum_{i=1}^{m} w_{i}^{2} cos(n(#phi_{i}))''",fnBinsPt,fPtMin,fPtMax,"s");
1109  TProfile *imq1nW2Pt = new TProfile("imq1nW2Pt","#sum_{i=1}^{m} w_{i}^{2} sin(n(#phi_{i}))''",fnBinsPt,fPtMin,fPtMax,"s");
1110  TProfile *req2nW1Pt = new TProfile("req2nW1Pt","#sum_{i=1}^{m} w_{i} cos(2n(#phi_{i}))''",fnBinsPt,fPtMin,fPtMax,"s");
1111  TProfile *imq2nW1Pt = new TProfile("imq2nW1Pt","#sum_{i=1}^{m} w_{i} sin(2n(#phi_{i}))''",fnBinsPt,fPtMin,fPtMax,"s");
1112
1113  TProfile *sumOfW1upTomPt = new TProfile("sumOfW1upTomPt","#sum_{i=1}^{m} w_{i}''",fnBinsPt,fPtMin,fPtMax,"s");
1114  TProfile *sumOfW2upTomPt = new TProfile("sumOfW2upTomPt","#sum_{i=1}^{m} w_{i}^{2}''",fnBinsPt,fPtMin,fPtMax,"s");
1115  TProfile *sumOfW3upTomPt = new TProfile("sumOfW3upTomPt","#sum_{i=1}^{m} w_{i}^{3}''",fnBinsPt,fPtMin,fPtMax,"s");
1116
1117  // eta RP
1118  TProfile *etaReq1n = new TProfile("etaReq1n","Re[q_{n}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1119  TProfile *etaImq1n = new TProfile("etaImq1n","Im[q_{n}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1120  TProfile *etaReq2n = new TProfile("etaReq2n","Re[q_{2n}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1121  TProfile *etaImq2n = new TProfile("etaImq2n","Im[q_{2n}]",fnBinsEta,fEtaMin,fEtaMax,"s");
1122
1123  TProfile *req1nW2Eta = new TProfile("req1nW2Eta","#sum_{i=1}^{m} w_{i}^{2} cos(n(#phi_{i}))''",fnBinsEta,fEtaMin,fEtaMax,"s");
1124  TProfile *imq1nW2Eta = new TProfile("imq1nW2Eta","#sum_{i=1}^{m} w_{i}^{2} sin(n(#phi_{i}))''",fnBinsEta,fEtaMin,fEtaMax,"s");
1125  TProfile *req2nW1Eta = new TProfile("req2nW1Eta","#sum_{i=1}^{m} w_{i} cos(2n(#phi_{i}))''",fnBinsEta,fEtaMin,fEtaMax,"s");
1126  TProfile *imq2nW1Eta = new TProfile("imq2nW1Eta","#sum_{i=1}^{m} w_{i} sin(2n(#phi_{i}))''",fnBinsEta,fEtaMin,fEtaMax,"s");
1127
1128  TProfile *sumOfW1upTomEta = new TProfile("sumOfW1upTomEta","#sum_{i=1}^{m} w_{i}''",fnBinsEta,fEtaMin,fEtaMax,"s");
1129  TProfile *sumOfW2upTomEta = new TProfile("sumOfW2upTomEta","#sum_{i=1}^{m} w_{i}^{2}''",fnBinsEta,fEtaMin,fEtaMax,"s");
1130  TProfile *sumOfW3upTomEta = new TProfile("sumOfW3upTomEta","#sum_{i=1}^{m} w_{i}^{3}''",fnBinsEta,fEtaMin,fEtaMax,"s");
1131
1132  for(Int_t i=0;i<nPrim;i++) // loop over all particles
1133  {
1134   fTrack=anEvent->GetTrack(i);
1135   if(fTrack)
1136   {
1137    if(fTrack->InPOISelection()) // checking if particle is POI
1138    {
1139     if(fTrack->InRPSelection()) // checking if particle is both POI and RP
1140     {
1141      // get azimuthal angle, momentum and pseudorapidity of a particle:
1142      dPhi = fTrack->Phi();
1143      dPt  = fTrack->Pt();
1144      dEta = fTrack->Eta();
1145      // phi weights:
1146      if(fUsePhiWeights)
1147      {
1148       nBinsPhi = phiWeights->GetNbinsX();
1149       if(nBinsPhi)
1150       {
1151        wPhi = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*nBinsPhi/TMath::TwoPi())));
1152       }
1153      }
1154      // pt weights:
1155      if(fUsePtWeights)
1156      {
1157       if(dBinWidthPt)
1158       {
1159        wPt = ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-fPtMin)/dBinWidthPt)));
1160       }
1161      }
1162      // eta weights:
1163      if(fUseEtaWeights)
1164      {
1165       if(dBinWidthEta)
1166       {
1167        wEta = etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-fEtaMin)/dBinWidthEta)));
1168       }
1169      }
1170      // pt:
1171      ptReq1nPrimePrime->Fill(dPt,cos(n*dPhi),1.);
1172      ptImq1nPrimePrime->Fill(dPt,sin(n*dPhi),1.);
1173      ptReq2nPrimePrime->Fill(dPt,cos(2.*n*dPhi),1.);
1174      ptImq2nPrimePrime->Fill(dPt,sin(2.*n*dPhi),1.);
1175      // weighted pt:
1176      req1nW2PrimePrimePt->Fill(dPt,cos(n*dPhi),pow(wPhi*wPt*wEta,2.));
1177      imq1nW2PrimePrimePt->Fill(dPt,sin(n*dPhi),pow(wPhi*wPt*wEta,2.));
1178      req2nW1PrimePrimePt->Fill(dPt,cos(2.*n*dPhi),wPhi*wPt*wEta);
1179      imq2nW1PrimePrimePt->Fill(dPt,sin(2.*n*dPhi),wPhi*wPt*wEta);
1180      sumOfW1upTomPrimePrimePt->Fill(dPt,wPhi*wPt*wEta,1.);
1181      sumOfW2upTomPrimePrimePt->Fill(dPt,pow(wPhi*wPt*wEta,2),1.);
1182      sumOfW3upTomPrimePrimePt->Fill(dPt,pow(wPhi*wPt*wEta,3),1.);
1183
1184      // eta:
1185      etaReq1nPrimePrime->Fill(dEta,cos(n*dPhi),1.);
1186      etaImq1nPrimePrime->Fill(dEta,sin(n*dPhi),1.);
1187      etaReq2nPrimePrime->Fill(dEta,cos(2.*n*dPhi),1.);
1188      etaImq2nPrimePrime->Fill(dEta,sin(2.*n*dPhi),1.);
1189      // weighted eta:
1190      req1nW2PrimePrimeEta->Fill(dEta,cos(n*dPhi),pow(wPhi*wPt*wEta,2.));
1191      imq1nW2PrimePrimeEta->Fill(dEta,sin(n*dPhi),pow(wPhi*wPt*wEta,2.));
1192      req2nW1PrimePrimeEta->Fill(dEta,cos(2.*n*dPhi),wPhi*wPt*wEta);
1193      imq2nW1PrimePrimeEta->Fill(dEta,sin(2.*n*dPhi),wPhi*wPt*wEta);
1194      sumOfW1upTomPrimePrimeEta->Fill(dEta,wPhi*wPt*wEta,1.);
1195      sumOfW2upTomPrimePrimeEta->Fill(dEta,pow(wPhi*wPt*wEta,2),1.);
1196      sumOfW3upTomPrimePrimeEta->Fill(dEta,pow(wPhi*wPt*wEta,3),1.);
1197
1198     }else if(!(fTrack->InRPSelection())) // checking if particles is POI and not RP
1199      {
1200       // get azimuthal angle, momentum and pseudorapidity of a particle:
1201       dPhi = fTrack->Phi();
1202       dPt  = fTrack->Pt();
1203       dEta = fTrack->Eta();
1204       // pt:
1205       ptReq1nPrime->Fill(dPt,cos(n*dPhi),1.);
1206       ptImq1nPrime->Fill(dPt,sin(n*dPhi),1.);
1207       ptReq2nPrime->Fill(dPt,cos(2.*n*dPhi),1.);
1208       ptImq2nPrime->Fill(dPt,sin(2.*n*dPhi),1.);
1209
1210       // eta:
1211       etaReq1nPrime->Fill(dEta,cos(n*dPhi),1.);
1212       etaImq1nPrime->Fill(dEta,sin(n*dPhi),1.);
1213       etaReq2nPrime->Fill(dEta,cos(2.*n*dPhi),1.);
1214       etaImq2nPrime->Fill(dEta,sin(2.*n*dPhi),1.);
1215
1216      } // end of else if(!(fTrack->InRPSelection())) // checking if particles is POI and not RP
1217    } // end of if(fTrack->InPOISelection()) // checking if particle is POI
1218
1219    if(fTrack->InRPSelection()) // checking if particles is only RP:
1220    {
1221     dPhi = fTrack->Phi();
1222     dPt  = fTrack->Pt();
1223     dEta = fTrack->Eta();
1224
1225     // phi weights:
1226     if(fUsePhiWeights)
1227     {
1228      nBinsPhi = phiWeights->GetNbinsX();
1229      if(nBinsPhi)
1230      {
1231       wPhi = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*nBinsPhi/TMath::TwoPi())));
1232      }
1233     }
1234     // pt weights:
1235     if(fUsePtWeights)
1236     {
1237      if(dBinWidthPt)
1238      {
1239       wPt = ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-fPtMin)/dBinWidthPt)));
1240      }
1241     }
1242     // eta weights:
1243     if(fUseEtaWeights)
1244     {
1245      if(dBinWidthEta)
1246      {
1247       wEta = etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-fEtaMin)/dBinWidthEta)));
1248      }
1249     }
1250     // pt:
1251     ptReq1n->Fill(dPt,cos(n*dPhi),1.);
1252     ptImq1n->Fill(dPt,sin(n*dPhi),1.);
1253     ptReq2n->Fill(dPt,cos(2.*n*dPhi),1.);
1254     ptImq2n->Fill(dPt,sin(2.*n*dPhi),1.);
1255     // weighted pt:
1256     req1nW2Pt->Fill(dPt,cos(n*dPhi),pow(wPhi*wPt*wEta,2.));
1257     imq1nW2Pt->Fill(dPt,sin(n*dPhi),pow(wPhi*wPt*wEta,2.));
1258     req2nW1Pt->Fill(dPt,cos(2.*n*dPhi),wPhi*wPt*wEta);
1259     imq2nW1Pt->Fill(dPt,sin(2.*n*dPhi),wPhi*wPt*wEta);
1260     sumOfW1upTomPt->Fill(dPt,wPhi*wPt*wEta,1.);
1261     sumOfW2upTomPt->Fill(dPt,pow(wPhi*wPt*wEta,2),1.);
1262     sumOfW3upTomPt->Fill(dPt,pow(wPhi*wPt*wEta,3),1.);
1263
1264     // eta:
1265     etaReq1n->Fill(dEta,cos(n*dPhi),1.);
1266     etaImq1n->Fill(dEta,sin(n*dPhi),1.);
1267     etaReq2n->Fill(dEta,cos(2.*n*dPhi),1.);
1268     etaImq2n->Fill(dEta,sin(2.*n*dPhi),1.);
1269     // weighted eta:
1270     req1nW2Eta->Fill(dEta,cos(n*dPhi),pow(wPhi*wPt*wEta,2.));
1271     imq1nW2Eta->Fill(dEta,sin(n*dPhi),pow(wPhi*wPt*wEta,2.));
1272     req2nW1Eta->Fill(dEta,cos(2.*n*dPhi),wPhi*wPt*wEta);
1273     imq2nW1Eta->Fill(dEta,sin(2.*n*dPhi),wPhi*wPt*wEta);
1274     sumOfW1upTomEta->Fill(dEta,wPhi*wPt*wEta,1.);
1275     sumOfW2upTomEta->Fill(dEta,pow(wPhi*wPt*wEta,2),1.);
1276     sumOfW3upTomEta->Fill(dEta,pow(wPhi*wPt*wEta,3),1.);
1277    } // end of if(fTrack->InRPSelection()) // checking if particles is only RP:
1278   } // end of if(fTrack}
1279  } // end of for(Int_t i=0;i<nPrim;i++)
1280
1281  //...........................................................................................................
1282  // PrimePrime Pt POI
1283  Double_t qxPrimePrimePtPOI=0.,qyPrimePrimePtPOI=0.,q2xPrimePrimePtPOIHere=0.,q2yPrimePrimePtPOIHere=0.;//add comments for these variable
1284  Double_t qxW2PrimePrimePtPOI=0.,qyW2PrimePrimePtPOI=0.,q2xW1PrimePrimePtPOI=0.,q2yW1PrimePrimePtPOI=0.;//add comments for these variable
1285  Double_t dS11mPrimePrimePtPOI=0.; // to be improved (name)
1286  Double_t dS12mPrimePrimePtPOI=0.; // to be improved (name)
1287  Double_t dS13mPrimePrimePtPOI=0.; // to be improved (name)
1288  Double_t mPrimePrimePtPOI=0.; // to be improved (name)
1289
1290  Double_t dM1pp11PtPOI=0.; // to be improved (name)
1291  Double_t dM0pp111PtPOI=0.; // to be improved (name)
1292  Double_t dM0pp12PtPOI=0.;
1293  Double_t dM2pp1PtPOI=0.;
1294  Double_t dM1pp2PtPOI=0.;
1295  Double_t dM0pp3PtPOI=0.;
1296
1297  // Prime Pt POI
1298  Double_t qxPrimePtPOI=0.,qyPrimePtPOI=0.;
1299  Double_t mPrimePtPOI=0.;
1300  Double_t dM0p111PtPOI=0.; // to be improved (name)
1301
1302  for(Int_t bin=1;bin<(fnBinsPt+1);bin++) // loop over pt-bins
1303  {
1304   // q'':
1305   qxPrimePrimePtPOI = (ptReq1nPrimePrime->GetBinContent(bin))*(ptReq1nPrimePrime->GetBinEntries(bin));
1306   qyPrimePrimePtPOI = (ptImq1nPrimePrime->GetBinContent(bin))*(ptImq1nPrimePrime->GetBinEntries(bin));
1307   q2xPrimePrimePtPOIHere = (ptReq2nPrimePrime->GetBinContent(bin))*(ptReq2nPrimePrime->GetBinEntries(bin));
1308   q2yPrimePrimePtPOIHere = (ptImq2nPrimePrime->GetBinContent(bin))*(ptImq2nPrimePrime->GetBinEntries(bin));
1309
1310   qxW2PrimePrimePtPOI = (req1nW2PrimePrimePt->GetBinContent(bin))*(req1nW2PrimePrimePt->GetBinEntries(bin));
1311   qyW2PrimePrimePtPOI = (imq1nW2PrimePrimePt->GetBinContent(bin))*(imq1nW2PrimePrimePt->GetBinEntries(bin));
1312
1313   q2xW1PrimePrimePtPOI = (req2nW1PrimePrimePt->GetBinContent(bin))*(req2nW1PrimePrimePt->GetBinEntries(bin));
1314   q2yW1PrimePrimePtPOI = (imq2nW1PrimePrimePt->GetBinContent(bin))*(imq2nW1PrimePrimePt->GetBinEntries(bin));
1315
1316   dS11mPrimePrimePtPOI = (sumOfW1upTomPrimePrimePt->GetBinContent(bin))*(sumOfW1upTomPrimePrimePt->GetBinEntries(bin));
1317   dS12mPrimePrimePtPOI = (sumOfW2upTomPrimePrimePt->GetBinContent(bin))*(sumOfW2upTomPrimePrimePt->GetBinEntries(bin));
1318   dS13mPrimePrimePtPOI = (sumOfW3upTomPrimePrimePt->GetBinContent(bin))*(sumOfW3upTomPrimePrimePt->GetBinEntries(bin));
1319
1320   mPrimePrimePtPOI = sumOfW1upTomPrimePrimePt->GetBinEntries(bin); // to be improved
1321
1322   dM1pp11PtPOI=dS11mPrimePrimePtPOI*(dSnk[1][0]-dSnk[0][1])-2.*dS12mPrimePrimePtPOI*dSnk[0][0]+2.*dS13mPrimePrimePtPOI;
1323   dM1pp2PtPOI=dS11mPrimePrimePtPOI*dSnk[0][1]-dS13mPrimePrimePtPOI;
1324   dM2pp1PtPOI=dS12mPrimePrimePtPOI*dSnk[0][0]-dS13mPrimePrimePtPOI;
1325   dM0pp3PtPOI=mPrimePrimePtPOI*dSnk[0][2]-dS13mPrimePrimePtPOI;
1326   dM0pp12PtPOI=mPrimePrimePtPOI*dSnk[0][0]*dSnk[0][1]-dM2pp1PtPOI-dM1pp2PtPOI-dM0pp3PtPOI-dS13mPrimePrimePtPOI;
1327
1328   // recursive formula (OK)
1329   // dM0pp111PtPOI=mPrimePrimePtPOI*dSnk[2][0]-3.*dM1pp11PtPOI-3.*dM0pp12PtPOI-3.*dM2pp1PtPOI-3.*dM1pp2PtPOI-dM0pp3PtPOI-dS13mPrimePrimePtPOI;
1330
1331   // direct formula (OK)
1332   dM0pp111PtPOI = mPrimePrimePtPOI*(dSnk[2][0]-3.*dSnk[0][0]*dSnk[0][1]+2.*dSnk[0][2])-3.*(dS11mPrimePrimePtPOI*(dSnk[1][0]-dSnk[0][1])+2.*(dS13mPrimePrimePtPOI-dS12mPrimePrimePtPOI*dSnk[0][0]));
1333
1334   // q':
1335   qxPrimePtPOI = (ptReq1nPrime->GetBinContent(bin))*(ptReq1nPrime->GetBinEntries(bin));
1336   qyPrimePtPOI = (ptImq1nPrime->GetBinContent(bin))*(ptImq1nPrime->GetBinEntries(bin));
1337
1338   mPrimePtPOI = ptReq1nPrime->GetBinEntries(bin); // to be improved
1339   dM0p111PtPOI=mPrimePtPOI*(dSnk[2][0]-3.*dSnk[0][1]*dSnk[0][0]+2.*dSnk[0][2]);
1340
1341   // 2-p the needed one
1342   Double_t two1n1nWPerPtBinPOI=0.;
1343   if((mPrimePrimePtPOI+mPrimePtPOI)*dSnk[0][0]-dS11mPrimePrimePtPOI>0)
1344   {
1345    two1n1nWPerPtBinPOI = (qxPrimePrimePtPOI*dQnkX[0][0]+qyPrimePrimePtPOI*dQnkY[0][0]+qxPrimePtPOI*dQnkX[0][0]+qyPrimePtPOI*dQnkY[0][0]-dS11mPrimePrimePtPOI)/((mPrimePrimePtPOI+mPrimePtPOI)*dSnk[0][0]-dS11mPrimePrimePtPOI);
1346
1347    f2WPerPtBin1n1nPOI->Fill(fPtMin+(bin-1)*dBinWidthPt,two1n1nWPerPtBinPOI,(mPrimePrimePtPOI+mPrimePtPOI)*dSnk[0][0]-dS11mPrimePrimePtPOI);
1348   }
1349
1350   // 2-p temporary one
1351   Double_t two1n1nW1ppW1W1PtPOI=0.; // <w1 w2 w3 cos(n(phi2-phi3))> // OK!!!
1352   if(dM1pp11PtPOI)
1353   {
1354    two1n1nW1ppW1W1PtPOI = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*dS11mPrimePrimePtPOI-2.*(qxW2PrimePrimePtPOI*dQnkX[0][0]+qyW2PrimePrimePtPOI*dQnkY[0][0]) - dS11mPrimePrimePtPOI*dSnk[0][1]+2.*dS13mPrimePrimePtPOI)/dM1pp11PtPOI; // CORRECT !!! <w1 w2 w3 cos(n(phi2-phi3))>
1355   }
1356
1357   // 2-p temporary one
1358   Double_t two1npp1nW1W2PtPOI=0.; // <w2 w3^2 cos(n(psi1-phi2))> // OK !!!
1359   if(dM0pp12PtPOI)
1360   {
1361    two1npp1nW1W2PtPOI = (dSnk[0][1]*(qxPrimePrimePtPOI*dQnkX[0][0]+qyPrimePrimePtPOI*dQnkY[0][0])-(qxW2PrimePrimePtPOI*dQnkX[0][0]+qyW2PrimePrimePtPOI*dQnkY[0][0])-dM1pp2PtPOI-(qxPrimePrimePtPOI*dQnkX[0][2]+qyPrimePrimePtPOI*dQnkY[0][2])+dS13mPrimePrimePtPOI)/dM0pp12PtPOI; // CORRECT !!! <w2 w3^2 cos(n(psi1-phi2))>
1362   }
1363
1364   // 2-p temporary one
1365   Double_t two1npp1nW2ppW1PtPOI=0.; // <w1^2 w2 cos(n(psi1-phi2))> // OK !!!
1366   if(dM2pp1PtPOI)
1367   {
1368    two1npp1nW2ppW1PtPOI = ((qxW2PrimePrimePtPOI*dQnkX[0][0]+qyW2PrimePrimePtPOI*dQnkY[0][0])-dS13mPrimePrimePtPOI)/dM2pp1PtPOI; // CORRECT !!! <w1^2 w2 cos(n(psi1-phi2))>
1369   }
1370
1371   // 2-p temporary one
1372   Double_t two2npp2nW1ppW2PtPOI=0.; // <w1 w2^2 cos(2n(psi1-phi2))> OK !!!
1373   if(dM1pp2PtPOI)
1374   {
1375    two2npp2nW1ppW2PtPOI = ((q2xW1PrimePrimePtPOI*dQnkX[1][1]+q2yW1PrimePrimePtPOI*dQnkY[1][1])-dS13mPrimePrimePtPOI)/dM1pp2PtPOI; // CORRECT !!! <w1 w2^2 cos(2n(psi1-phi2))>
1376   }
1377
1378   // 2-p temporary one
1379   Double_t two1npp1nW3PtPOI=0.; // <w2^3 cos(n(psi1-phi2))> // OK !!!
1380   if(dM0pp3PtPOI)
1381   {
1382    two1npp1nW3PtPOI = (qxPrimePrimePtPOI*dQnkX[0][2]+qyPrimePrimePtPOI*dQnkY[0][2]-dS13mPrimePrimePtPOI)/dM0pp3PtPOI; // CORRECT !!! <w2^3 cos(n(psi1-phi2))>
1383   }
1384
1385   // 3-p temporary one
1386   Double_t three2npp1n1nW1ppW1W1PtPOI=0.; // <w1 w2 w3 cos(n(2psi1-phi2-phi3))> // OK!!!
1387   if(dM1pp11PtPOI)
1388   {
1389    three2npp1n1nW1ppW1W1PtPOI = (q2xW1PrimePrimePtPOI*(dQnkX[0][0]*dQnkX[0][0]-dQnkY[0][0]*dQnkY[0][0])+2.*q2yW1PrimePrimePtPOI*dQnkX[0][0]*dQnkY[0][0]-2.*(qxW2PrimePrimePtPOI*dQnkX[0][0]+qyW2PrimePrimePtPOI*dQnkY[0][0])-(q2xW1PrimePrimePtPOI*dQnkX[1][1]+q2yW1PrimePrimePtPOI*dQnkY[1][1])+2.*dS13mPrimePrimePtPOI)/dM1pp11PtPOI; // CORRECT !!! <w1 w2 w3 cos(n(2psi1-phi2-phi3))>
1390   }
1391
1392   // 3-p temporary one
1393   Double_t three1npp1n2nW0ppW1W2PtPOI=0.; // <w2 w3^2 cos(n(psi1+phi2-2*phi3))> // OK!!!
1394   if(dM0pp12PtPOI)
1395   {
1396    three1npp1n2nW0ppW1W2PtPOI = (qxPrimePrimePtPOI*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])-qyPrimePrimePtPOI*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1])-(qxW2PrimePrimePtPOI*dQnkX[0][0]+qyW2PrimePrimePtPOI*dQnkY[0][0])-(q2xW1PrimePrimePtPOI*dQnkX[1][1]+q2yW1PrimePrimePtPOI*dQnkY[1][1])-(qxPrimePrimePtPOI*dQnkX[0][2]+qyPrimePrimePtPOI*dQnkY[0][2])+2.*dS13mPrimePrimePtPOI)/dM0pp12PtPOI; // CORRECT !!! <w2 w3^2 cos(n(psi1+phi2-2.*phi3))>
1397   }
1398
1399   /*
1400   // 4-p RP part
1401   Double_t four1npp1n1n1nW1W1W1=0.; // <w1 w2 w3 cos(n(psi1+phi1-phi2-phi3))>
1402   if(dM0pp111PtPOI)
1403   {
1404    four1npp1n1n1nW1W1W1 = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimePrimePtPOI*dQnkX[0][0]+qyPrimePrimePtPOI*dQnkY[0][0])-2.*dM1pp11PtPOI*two1n1nW1ppW1W1PtPOI-dM1pp11PtPOI*three2npp1n1nW1ppW1W1PtPOI-dM0pp12PtPOI*three1npp1n2nW0ppW1W2PtPOI-2.*dM0pp12PtPOI*two1npp1nW1W2PtPOI-3.*dM2pp1PtPOI*two1npp1nW2ppW1PtPOI-2.*dM1pp2PtPOI-dM1pp2PtPOI*two2npp2nW1ppW2PtPOI-dM0pp3PtPOI*two1npp1nW3PtPOI-dS13mPrimePrimePtPOI)/(dM0pp111PtPOI);
1405   }
1406   */
1407
1408   /*
1409   // 4-p POI part
1410   Double_t four1npp1n1n1nW1W1W1POI=0.;
1411   if(dM0p111PtPOI>0&&mPrimePtPOI>0&&nRP>0)
1412   {
1413    four1npp1n1n1nW1W1W1POI = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimePtPOI*dQnkX[0][0]+qyPrimePtPOI*dQnkY[0][0])-2.*dSnk[0][1]* (qxPrimePtPOI*dQnkX[0][0]+qyPrimePtPOI*dQnkY[0][0])+2.*(qxPrimePtPOI*dQnkX[0][2]+qyPrimePtPOI*dQnkY[0][2])-qxPrimePtPOI*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])+qyPrimePtPOI*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1]))/dM0p111PtPOI;
1414   }
1415   */
1416
1417   // recursive formula for 4-p RP and POI in all combinations (full, partial and no overlap) (OK)
1418   Double_t four1npp1n1n1nW1W1W1PtPOI=0.;
1419   if(dM0pp111PtPOI+dM0p111PtPOI)
1420   {
1421    four1npp1n1n1nW1W1W1PtPOI = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimePrimePtPOI*dQnkX[0][0]+qyPrimePrimePtPOI*dQnkY[0][0])-2.*dM1pp11PtPOI*two1n1nW1ppW1W1PtPOI-dM1pp11PtPOI*three2npp1n1nW1ppW1W1PtPOI-dM0pp12PtPOI*three1npp1n2nW0ppW1W2PtPOI-2.*dM0pp12PtPOI*two1npp1nW1W2PtPOI-3.*dM2pp1PtPOI*two1npp1nW2ppW1PtPOI-2.*dM1pp2PtPOI-dM1pp2PtPOI*two2npp2nW1ppW2PtPOI-dM0pp3PtPOI*two1npp1nW3PtPOI-dS13mPrimePrimePtPOI+(pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimePtPOI*dQnkX[0][0]+qyPrimePtPOI*dQnkY[0][0])-2.*dSnk[0][1]* (qxPrimePtPOI*dQnkX[0][0]+qyPrimePtPOI*dQnkY[0][0])+2.*(qxPrimePtPOI*dQnkX[0][2]+qyPrimePtPOI*dQnkY[0][2])-qxPrimePtPOI*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])+qyPrimePtPOI*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1]))/(dM0pp111PtPOI+dM0p111PtPOI);
1422
1423
1424  /*
1425  Double_t four1npp1n1n1nW1W1W1PtPOIb =
1426           ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimePrimePtPOI*dQnkX[0][0]+qyPrimePrimePtPOI*dQnkY[0][0])
1427           -  q2xW1PrimePrimePtPOI*(dQnkX[0][0]*dQnkX[0][0]-dQnkY[0][0]*dQnkY[0][0])+2.*q2yW1PrimePrimePtPOI*dQnkX[0][0]*dQnkY[0][0]
1428           -  qxPrimePrimePtPOI*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])-qyPrimePrimePtPOI*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1])
1429           -  2.*dSnk[0][1]*(qxPrimePrimePtPOI*dQnkX[0][0]+qyPrimePrimePtPOI*dQnkY[0][0])
1430           -  2.*(pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*dS11mPrimePrimePtPOI
1431           +  7.*(qxW2PrimePrimePtPOI*dQnkX[0][0]+qyW2PrimePrimePtPOI*dQnkY[0][0])
1432           -  1.*(qxW2PrimePrimePtPOI*dQnkX[0][0]+qyW2PrimePrimePtPOI*dQnkY[0][0])
1433           +  1.*(q2xW1PrimePrimePtPOI*dQnkX[1][1]+q2yW1PrimePrimePtPOI*dQnkY[1][1])
1434           +  2.*(qxPrimePrimePtPOI*dQnkX[0][2]+qyPrimePrimePtPOI*dQnkY[0][2])
1435           +  2.*dS11mPrimePrimePtPOI*dSnk[0][1]
1436           -  6.*dS13mPrimePrimePtPOI
1437           +  (pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimePtPOI*dQnkX[0][0]+qyPrimePtPOI*dQnkY[0][0])
1438           -  2.*dSnk[0][1]*(qxPrimePtPOI*dQnkX[0][0]+qyPrimePtPOI*dQnkY[0][0])
1439           -  qxPrimePtPOI*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])-qyPrimePtPOI*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1])
1440           +  2.*(qxPrimePrimePtPOI*dQnkX[0][2]+qyPrimePrimePtPOI*dQnkY[0][2]))/(dM0pp111PtPOI+dM0p111PtPOI);
1441
1442
1443
1444
1445
1446
1447
1448    cout<<endl;
1449    cout<<"ab"<<endl;
1450    cout<<four1npp1n1n1nW1W1W1PtPOI<<endl;
1451    cout<<four1npp1n1n1nW1W1W1PtPOIb<<endl;
1452    cout<<endl;
1453
1454    */
1455
1456    f4WPerPtBin1n1n1n1nPOI->Fill(fPtMin+(bin-1)*dBinWidthPt,four1npp1n1n1nW1W1W1PtPOI,dM0pp111PtPOI+dM0p111PtPOI);
1457   } // end of if(dM0pp111PtPOI+dM0p111PtPOI)
1458  } // for(Int_t bin=1;bin<(fnBinsPt+1);bin++) // loop over pt-bins
1459  //...........................................................................................................
1460
1461  //...........................................................................................................
1462  // PrimePrime Eta POI
1463  Double_t qxPrimePrimeEtaPOI=0.,qyPrimePrimeEtaPOI=0.,q2xPrimePrimeEtaPOIHere=0.,q2yPrimePrimeEtaPOIHere=0.;//add comments for these variable
1464  Double_t qxW2PrimePrimeEtaPOI=0.,qyW2PrimePrimeEtaPOI=0.,q2xW1PrimePrimeEtaPOI=0.,q2yW1PrimePrimeEtaPOI=0.;//add comments for these variable
1465  Double_t dS11mPrimePrimeEtaPOI=0.; // to be improved (name)
1466  Double_t dS12mPrimePrimeEtaPOI=0.; // to be improved (name)
1467  Double_t dS13mPrimePrimeEtaPOI=0.; // to be improved (name)
1468  Double_t mPrimePrimeEtaPOIHere=0.; // to be improved (name)
1469
1470  Double_t dM1pp11EtaPOI=0.; // to be improved (name)
1471  Double_t dM0pp111EtaPOI=0.; // to be improved (name)
1472  Double_t dM0pp12EtaPOI=0.;
1473  Double_t dM2pp1EtaPOI=0.;
1474  Double_t dM1pp2EtaPOI=0.;
1475  Double_t dM0pp3EtaPOI=0.;
1476
1477  // Prime Eta POI
1478  Double_t qxPrimeEtaPOI=0.,qyPrimeEtaPOI=0.;
1479  Double_t mPrimeEtaPOI=0.;
1480  Double_t dM0p111EtaPOI=0.; // to be improved (name)
1481
1482  for(Int_t bin=1;bin<(fnBinsEta+1);bin++) // loop over eta-bins
1483  {
1484   // q'':
1485   qxPrimePrimeEtaPOI = (etaReq1nPrimePrime->GetBinContent(bin))*(etaReq1nPrimePrime->GetBinEntries(bin));
1486   qyPrimePrimeEtaPOI = (etaImq1nPrimePrime->GetBinContent(bin))*(etaImq1nPrimePrime->GetBinEntries(bin));
1487   q2xPrimePrimeEtaPOIHere = (etaReq2nPrimePrime->GetBinContent(bin))*(etaReq2nPrimePrime->GetBinEntries(bin));
1488   q2yPrimePrimeEtaPOIHere = (etaImq2nPrimePrime->GetBinContent(bin))*(etaImq2nPrimePrime->GetBinEntries(bin));
1489
1490   qxW2PrimePrimeEtaPOI = (req1nW2PrimePrimeEta->GetBinContent(bin))*(req1nW2PrimePrimeEta->GetBinEntries(bin));
1491   qyW2PrimePrimeEtaPOI = (imq1nW2PrimePrimeEta->GetBinContent(bin))*(imq1nW2PrimePrimeEta->GetBinEntries(bin));
1492
1493   q2xW1PrimePrimeEtaPOI = (req2nW1PrimePrimeEta->GetBinContent(bin))*(req2nW1PrimePrimeEta->GetBinEntries(bin));
1494   q2yW1PrimePrimeEtaPOI = (imq2nW1PrimePrimeEta->GetBinContent(bin))*(imq2nW1PrimePrimeEta->GetBinEntries(bin));
1495
1496   dS11mPrimePrimeEtaPOI = (sumOfW1upTomPrimePrimeEta->GetBinContent(bin))*(sumOfW1upTomPrimePrimeEta->GetBinEntries(bin));
1497   dS12mPrimePrimeEtaPOI = (sumOfW2upTomPrimePrimeEta->GetBinContent(bin))*(sumOfW2upTomPrimePrimeEta->GetBinEntries(bin));
1498   dS13mPrimePrimeEtaPOI = (sumOfW3upTomPrimePrimeEta->GetBinContent(bin))*(sumOfW3upTomPrimePrimeEta->GetBinEntries(bin));
1499
1500   mPrimePrimeEtaPOIHere = sumOfW1upTomPrimePrimeEta->GetBinEntries(bin); // to be improved
1501
1502   dM1pp11EtaPOI=dS11mPrimePrimeEtaPOI*(dSnk[1][0]-dSnk[0][1])-2.*dS12mPrimePrimeEtaPOI*dSnk[0][0]+2.*dS13mPrimePrimeEtaPOI;
1503   dM1pp2EtaPOI=dS11mPrimePrimeEtaPOI*dSnk[0][1]-dS13mPrimePrimeEtaPOI;
1504   dM2pp1EtaPOI=dS12mPrimePrimeEtaPOI*dSnk[0][0]-dS13mPrimePrimeEtaPOI;
1505   dM0pp3EtaPOI=mPrimePrimeEtaPOIHere*dSnk[0][2]-dS13mPrimePrimeEtaPOI;
1506   dM0pp12EtaPOI=mPrimePrimeEtaPOIHere*dSnk[0][0]*dSnk[0][1]-dM2pp1EtaPOI-dM1pp2EtaPOI-dM0pp3EtaPOI-dS13mPrimePrimeEtaPOI;
1507   dM0pp111EtaPOI=mPrimePrimeEtaPOIHere*dSnk[2][0]-3.*dM1pp11EtaPOI-3.*dM0pp12EtaPOI-3.*dM2pp1EtaPOI-3.*dM1pp2EtaPOI-dM0pp3EtaPOI-dS13mPrimePrimeEtaPOI;
1508
1509   // q':
1510   qxPrimeEtaPOI = (etaReq1nPrime->GetBinContent(bin))*(etaReq1nPrime->GetBinEntries(bin));
1511   qyPrimeEtaPOI = (etaImq1nPrime->GetBinContent(bin))*(etaImq1nPrime->GetBinEntries(bin));
1512
1513   mPrimeEtaPOI = etaReq1nPrime->GetBinEntries(bin); // to be improved
1514   dM0p111EtaPOI=mPrimeEtaPOI*(dSnk[2][0]-3.*dSnk[0][1]*dSnk[0][0]+2.*dSnk[0][2]);
1515
1516   // 2-p the needed one
1517   Double_t two1n1nWPerEtaBinPOI=0.;
1518   if((mPrimePrimeEtaPOIHere+mPrimeEtaPOI)*dSnk[0][0]-dS11mPrimePrimeEtaPOI>0)
1519   {
1520    two1n1nWPerEtaBinPOI = (qxPrimePrimeEtaPOI*dQnkX[0][0]+qyPrimePrimeEtaPOI*dQnkY[0][0]+qxPrimeEtaPOI*dQnkX[0][0]+qyPrimeEtaPOI*dQnkY[0][0]-dS11mPrimePrimeEtaPOI)/((mPrimePrimeEtaPOIHere+mPrimeEtaPOI)*dSnk[0][0]-dS11mPrimePrimeEtaPOI);
1521
1522    f2WPerEtaBin1n1nPOI->Fill(fEtaMin+(bin-1)*dBinWidthEta,two1n1nWPerEtaBinPOI,(mPrimePrimeEtaPOIHere+mPrimeEtaPOI)*dSnk[0][0]-dS11mPrimePrimeEtaPOI); // <2'>_{n|n}
1523   }
1524
1525   // 2-p the temporary one
1526   Double_t two1n1nW1ppW1W1EtaPOI=0.; // <w1 w2 w3 cos(n(phi2-phi3))> // OK!!!
1527   if(dM1pp11EtaPOI)
1528   {
1529    two1n1nW1ppW1W1EtaPOI = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*dS11mPrimePrimeEtaPOI-2.*(qxW2PrimePrimeEtaPOI*dQnkX[0][0]+qyW2PrimePrimeEtaPOI*dQnkY[0][0])-dS11mPrimePrimeEtaPOI*dSnk[0][1]+2.*dS13mPrimePrimeEtaPOI)/dM1pp11EtaPOI; // CORRECT !!! <w1 w2 w3 cos(n(phi2-phi3))>
1530   }
1531
1532   // 2-p the temporary one
1533   Double_t two1npp1nW1W2EtaPOI=0.; // <w2 w3^2 cos(n(psi1-phi2))> // OK !!!
1534   if(dM0pp12EtaPOI)
1535   {
1536    two1npp1nW1W2EtaPOI = (dSnk[0][1]*(qxPrimePrimeEtaPOI*dQnkX[0][0]+qyPrimePrimeEtaPOI*dQnkY[0][0])-(qxW2PrimePrimeEtaPOI*dQnkX[0][0]+qyW2PrimePrimeEtaPOI*dQnkY[0][0])-dM1pp2EtaPOI-(qxPrimePrimeEtaPOI*dQnkX[0][2]+qyPrimePrimeEtaPOI*dQnkY[0][2])+dS13mPrimePrimeEtaPOI)/dM0pp12EtaPOI; // CORRECT !!! <w2 w3^2 cos(n(psi1-phi2))>
1537   }
1538
1539   // 2-p the temporary one
1540   Double_t two1npp1nW2ppW1EtaPOI=0.; // <w1^2 w2 cos(n(psi1-phi2))> // OK !!!
1541   if(dM2pp1EtaPOI)
1542   {
1543    two1npp1nW2ppW1EtaPOI = ((qxW2PrimePrimeEtaPOI*dQnkX[0][0]+qyW2PrimePrimeEtaPOI*dQnkY[0][0])-dS13mPrimePrimeEtaPOI)/dM2pp1EtaPOI; // CORRECT !!! <w1^2 w2 cos(n(psi1-phi2))>
1544   }
1545
1546   // 2-p the temporary one
1547   Double_t two2npp2nW1ppW2EtaPOI=0.; // <w1 w2^2 cos(2n(psi1-phi2))> OK !!!
1548   if(dM1pp2EtaPOI)
1549   {
1550    two2npp2nW1ppW2EtaPOI = ((q2xW1PrimePrimeEtaPOI*dQnkX[1][1]+q2yW1PrimePrimeEtaPOI*dQnkY[1][1])-dS13mPrimePrimeEtaPOI)/dM1pp2EtaPOI; // CORRECT !!! <w1 w2^2 cos(2n(psi1-phi2))>
1551   }
1552
1553   // 2-p the temporary one
1554   Double_t two1npp1nW3EtaPOI=0.; // <w2^3 cos(n(psi1-phi2))> // OK !!!
1555   if(dM0pp3EtaPOI)
1556   {
1557    two1npp1nW3EtaPOI = (qxPrimePrimeEtaPOI*dQnkX[0][2]+qyPrimePrimeEtaPOI*dQnkY[0][2]-dS13mPrimePrimeEtaPOI)/dM0pp3EtaPOI; // CORRECT !!! <w2^3 cos(n(psi1-phi2))>
1558   }
1559
1560   // 3-p the temporary one
1561   Double_t three2npp1n1nW1ppW1W1EtaPOI=0.; // <w1 w2 w3 cos(n(2psi1-phi2-phi3))> // OK!!!
1562   if(dM1pp11EtaPOI)
1563   {
1564    three2npp1n1nW1ppW1W1EtaPOI = (q2xW1PrimePrimeEtaPOI*(dQnkX[0][0]*dQnkX[0][0]-dQnkY[0][0]*dQnkY[0][0])+2.*q2yW1PrimePrimeEtaPOI*dQnkX[0][0]*dQnkY[0][0]-2.*(qxW2PrimePrimeEtaPOI*dQnkX[0][0]+qyW2PrimePrimeEtaPOI*dQnkY[0][0])-(q2xW1PrimePrimeEtaPOI*dQnkX[1][1]+q2yW1PrimePrimeEtaPOI*dQnkY[1][1])+2.*dS13mPrimePrimeEtaPOI)/dM1pp11EtaPOI; // CORRECT !!! <w1 w2 w3 cos(n(2psi1-phi2-phi3))>
1565   }
1566
1567   // 3-p the temporary one
1568   Double_t three1npp1n2nW0ppW1W2EtaPOI=0.; // <w2 w3^2 cos(n(psi1+phi2-2*phi3))> // OK!!!
1569   if(dM0pp12EtaPOI)
1570   {
1571    three1npp1n2nW0ppW1W2EtaPOI = (qxPrimePrimeEtaPOI*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])-qyPrimePrimeEtaPOI*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1])-(qxW2PrimePrimeEtaPOI*dQnkX[0][0]+qyW2PrimePrimeEtaPOI*dQnkY[0][0])-(q2xW1PrimePrimeEtaPOI*dQnkX[1][1]+q2yW1PrimePrimeEtaPOI*dQnkY[1][1])-(qxPrimePrimeEtaPOI*dQnkX[0][2]+qyPrimePrimeEtaPOI*dQnkY[0][2])+2.*dS13mPrimePrimeEtaPOI)/dM0pp12EtaPOI; // CORRECT !!! <w2 w3^2 cos(n(psi1+phi2-2.*phi3))>
1572   }
1573
1574   /*
1575   // 4-p RP part
1576   Double_t four1npp1n1n1nW1W1W1=0.; // <w1 w2 w3 cos(n(psi1+phi1-phi2-phi3))>
1577   if(dM0pp111EtaPOI)
1578   {
1579    four1npp1n1n1nW1W1W1 = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimePrimeEtaPOI*dQnkX[0][0]+qyPrimePrimeEtaPOI*dQnkY[0][0])-2.*dM1pp11EtaPOI*two1n1nW1ppW1W1EtaPOI-dM1pp11EtaPOI*three2npp1n1nW1ppW1W1EtaPOI-dM0pp12EtaPOI*three1npp1n2nW0ppW1W2EtaPOI-2.*dM0pp12EtaPOI*two1npp1nW1W2EtaPOI-3.*dM2pp1EtaPOI*two1npp1nW2ppW1EtaPOI-2.*dM1pp2EtaPOI-dM1pp2EtaPOI*two2npp2nW1ppW2EtaPOI-dM0pp3EtaPOI*two1npp1nW3EtaPOI-dS13mPrimePrimeEtaPOI)/(dM0pp111EtaPOI);
1580   }
1581   */
1582   /*
1583   // 4-p POI part
1584   Double_t four1npp1n1n1nW1W1W1POI=0.;
1585   if(dM0p111EtaPOI>0&&mPrimeEtaPOI>0&&nRP>0)
1586   {
1587    four1npp1n1n1nW1W1W1POI = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimeEtaPOI*dQnkX[0][0]+qyPrimeEtaPOI*dQnkY[0][0])-2.*dSnk[0][1]* (qxPrimeEtaPOI*dQnkX[0][0]+qyPrimeEtaPOI*dQnkY[0][0])+2.*(qxPrimeEtaPOI*dQnkX[0][2]+qyPrimeEtaPOI*dQnkY[0][2])-qxPrimeEtaPOI*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])+qyPrimeEtaPOI*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1]))/dM0p111EtaPOI;
1588   }
1589   */
1590
1591   // 4-p RP and POI in all combinations (full, partial and no overlap)
1592   Double_t four1npp1n1n1nW1W1W1EtaPOI=0.;
1593
1594   if(dM0pp111EtaPOI+dM0p111EtaPOI)
1595   {
1596    four1npp1n1n1nW1W1W1EtaPOI = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimePrimeEtaPOI*dQnkX[0][0]+qyPrimePrimeEtaPOI*dQnkY[0][0])-2.*dM1pp11EtaPOI*two1n1nW1ppW1W1EtaPOI-dM1pp11EtaPOI*three2npp1n1nW1ppW1W1EtaPOI-dM0pp12EtaPOI*three1npp1n2nW0ppW1W2EtaPOI-2.*dM0pp12EtaPOI*two1npp1nW1W2EtaPOI-3.*dM2pp1EtaPOI*two1npp1nW2ppW1EtaPOI-2.*dM1pp2EtaPOI-dM1pp2EtaPOI*two2npp2nW1ppW2EtaPOI-dM0pp3EtaPOI*two1npp1nW3EtaPOI-dS13mPrimePrimeEtaPOI+(pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimeEtaPOI*dQnkX[0][0]+qyPrimeEtaPOI*dQnkY[0][0])-2.*dSnk[0][1]*(qxPrimeEtaPOI*dQnkX[0][0]+qyPrimeEtaPOI*dQnkY[0][0])+2.*(qxPrimeEtaPOI*dQnkX[0][2]+qyPrimeEtaPOI*dQnkY[0][2])-qxPrimeEtaPOI*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])+qyPrimeEtaPOI*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1]))/(dM0pp111EtaPOI+dM0p111EtaPOI);
1597
1598    f4WPerEtaBin1n1n1n1nPOI->Fill(fEtaMin+(bin-1)*dBinWidthEta,four1npp1n1n1nW1W1W1EtaPOI,dM0pp111EtaPOI+dM0p111EtaPOI);
1599   }
1600  }
1601  //...........................................................................................................
1602
1603
1604
1605
1606
1607
1608
1609
1610
1611
1612
1613
1614
1615
1616  // RPs
1617  ptReq1nPrime->Reset(); // to be improved
1618  ptImq1nPrime->Reset(); // to be improved
1619  ptReq2nPrime->Reset(); // to be improved
1620  ptImq2nPrime->Reset(); // to be improved
1621
1622  etaReq1nPrime->Reset(); // to be improved
1623  etaImq1nPrime->Reset(); // to be improved
1624  etaReq2nPrime->Reset(); // to be improved
1625  etaImq2nPrime->Reset(); // to be improved
1626
1627
1628  //...........................................................................................................
1629  // PrimePrime Pt RP
1630  Double_t qxPtRP=0.,qyPtRP=0.,q2xPtRP=0.,q2yPtRP=0.;//add comments for these variable
1631  Double_t qxW2PtRP=0.,qyW2PtRP=0.,q2xW1PtRP=0.,q2yW1PtRP=0.;//add comments for these variable
1632  Double_t dS11mPtRP=0.; // to be improved (name)
1633  Double_t dS12mPtRP=0.; // to be improved (name)
1634  Double_t dS13mPtRP=0.; // to be improved (name)
1635  Double_t mPtRP=0.; // to be improved (name)
1636
1637  Double_t dM1pp11PtRP=0.; // to be improved (name)
1638  Double_t dM0pp111PtRP=0.; // to be improved (name)
1639  Double_t dM0pp12PtRP=0.;
1640  Double_t dM2pp1PtRP=0.;
1641  Double_t dM1pp2PtRP=0.;
1642  Double_t dM0pp3PtRP=0.;
1643
1644  // Prime Pt RP
1645  Double_t qxPrimePtRP=0.,qyPrimePtRP=0.;
1646  Double_t mPrimePtRP=0.;
1647  Double_t dM0p111PtRP=0.; // to be improved (name)
1648
1649  for(Int_t bin=1;bin<(fnBinsPt+1);bin++) // loop over pt-bins
1650  {
1651   // q'':
1652   qxPtRP = (ptReq1n->GetBinContent(bin))*(ptReq1n->GetBinEntries(bin));
1653   qyPtRP = (ptImq1n->GetBinContent(bin))*(ptImq1n->GetBinEntries(bin));
1654   q2xPtRP = (ptReq2n->GetBinContent(bin))*(ptReq2n->GetBinEntries(bin));
1655   q2yPtRP = (ptImq2n->GetBinContent(bin))*(ptImq2n->GetBinEntries(bin));
1656
1657   qxW2PtRP = (req1nW2Pt->GetBinContent(bin))*(req1nW2Pt->GetBinEntries(bin));
1658   qyW2PtRP = (imq1nW2Pt->GetBinContent(bin))*(imq1nW2Pt->GetBinEntries(bin));
1659
1660   q2xW1PtRP = (req2nW1Pt->GetBinContent(bin))*(req2nW1Pt->GetBinEntries(bin));
1661   q2yW1PtRP = (imq2nW1Pt->GetBinContent(bin))*(imq2nW1Pt->GetBinEntries(bin));
1662
1663   dS11mPtRP = (sumOfW1upTomPt->GetBinContent(bin))*(sumOfW1upTomPt->GetBinEntries(bin));
1664   dS12mPtRP = (sumOfW2upTomPt->GetBinContent(bin))*(sumOfW2upTomPt->GetBinEntries(bin));
1665   dS13mPtRP = (sumOfW3upTomPt->GetBinContent(bin))*(sumOfW3upTomPt->GetBinEntries(bin));
1666
1667   mPtRP = sumOfW1upTomPt->GetBinEntries(bin); // to be improved
1668
1669   dM1pp11PtRP=dS11mPtRP*(dSnk[1][0]-dSnk[0][1])-2.*dS12mPtRP*dSnk[0][0]+2.*dS13mPtRP;
1670   dM1pp2PtRP=dS11mPtRP*dSnk[0][1]-dS13mPtRP;
1671   dM2pp1PtRP=dS12mPtRP*dSnk[0][0]-dS13mPtRP;
1672   dM0pp3PtRP=mPtRP*dSnk[0][2]-dS13mPtRP;
1673   dM0pp12PtRP=mPtRP*dSnk[0][0]*dSnk[0][1]-dM2pp1PtRP-dM1pp2PtRP-dM0pp3PtRP-dS13mPtRP;
1674   dM0pp111PtRP=mPtRP*dSnk[2][0]-3.*dM1pp11PtRP-3.*dM0pp12PtRP-3.*dM2pp1PtRP-3.*dM1pp2PtRP-dM0pp3PtRP-dS13mPtRP;
1675
1676   // q':
1677   qxPrimePtRP = (ptReq1nPrime->GetBinContent(bin))*(ptReq1nPrime->GetBinEntries(bin));
1678   qyPrimePtRP = (ptImq1nPrime->GetBinContent(bin))*(ptImq1nPrime->GetBinEntries(bin));
1679
1680   mPrimePtRP = ptReq1nPrime->GetBinEntries(bin); // to be improved
1681   dM0p111PtRP=mPrimePtRP*(dSnk[2][0]-3.*dSnk[0][1]*dSnk[0][0]+2.*dSnk[0][2]);
1682
1683   // 2-p the needed one
1684   Double_t two1n1nWPerPtBinRP=0.;
1685   if((mPtRP+mPrimePtRP)*dSnk[0][0]-dS11mPtRP>0)
1686   {
1687    two1n1nWPerPtBinRP = (qxPtRP*dQnkX[0][0]+qyPtRP*dQnkY[0][0]+qxPrimePtRP*dQnkX[0][0]+qyPrimePtRP*dQnkY[0][0]-dS11mPtRP)/((mPtRP+mPrimePtRP)*dSnk[0][0]-dS11mPtRP);
1688    f2WPerPtBin1n1nRP->Fill(fPtMin+(bin-1)*dBinWidthPt,two1n1nWPerPtBinRP,(mPtRP+mPrimePtRP)*dSnk[0][0]-dS11mPtRP);
1689   }
1690
1691   // 2-p temporary one
1692   Double_t two1n1nW1ppW1W1PtRP=0.; // <w1 w2 w3 cos(n(phi2-phi3))> // OK!!!
1693   if(dM1pp11PtRP)
1694   {
1695    two1n1nW1ppW1W1PtRP = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*dS11mPtRP-2.*(qxW2PtRP*dQnkX[0][0]+qyW2PtRP*dQnkY[0][0])-dS11mPtRP*dSnk[0][1]+2.*dS13mPtRP)/dM1pp11PtRP; // CORRECT !!! <w1 w2 w3 cos(n(phi2-phi3))>
1696   }
1697
1698   // 2-p temporary one
1699   Double_t two1npp1nW1W2PtRP=0.; // <w2 w3^2 cos(n(psi1-phi2))> // OK !!!
1700   if(dM0pp12PtRP)
1701   {
1702    two1npp1nW1W2PtRP = (dSnk[0][1]*(qxPtRP*dQnkX[0][0]+qyPtRP*dQnkY[0][0])-(qxW2PtRP*dQnkX[0][0]+qyW2PtRP*dQnkY[0][0])-dM1pp2PtRP-(qxPtRP*dQnkX[0][2]+qyPtRP*dQnkY[0][2])+dS13mPtRP)/dM0pp12PtRP; // CORRECT !!! <w2 w3^2 cos(n(psi1-phi2))>
1703   }
1704
1705   // 2-p temporary one
1706   Double_t two1npp1nW2ppW1PtRP=0.; // <w1^2 w2 cos(n(psi1-phi2))> // OK !!!
1707   if(dM2pp1PtRP)
1708   {
1709    two1npp1nW2ppW1PtRP = ((qxW2PtRP*dQnkX[0][0]+qyW2PtRP*dQnkY[0][0])-dS13mPtRP)/dM2pp1PtRP; // CORRECT !!! <w1^2 w2 cos(n(psi1-phi2))>
1710   }
1711
1712   // 2-p temporary one
1713   Double_t two2npp2nW1ppW2PtRP=0.; // <w1 w2^2 cos(2n(psi1-phi2))> OK !!!
1714   if(dM1pp2PtRP)
1715   {
1716    two2npp2nW1ppW2PtRP = ((q2xW1PtRP*dQnkX[1][1]+q2yW1PtRP*dQnkY[1][1])-dS13mPtRP)/dM1pp2PtRP; // CORRECT !!! <w1 w2^2 cos(2n(psi1-phi2))>
1717   }
1718
1719   // 2-p temporary one
1720   Double_t two1npp1nW3PtRP=0.; // <w2^3 cos(n(psi1-phi2))> // OK !!!
1721   if(dM0pp3PtRP)
1722   {
1723    two1npp1nW3PtRP = (qxPtRP*dQnkX[0][2]+qyPtRP*dQnkY[0][2]-dS13mPtRP)/dM0pp3PtRP; // CORRECT !!! <w2^3 cos(n(psi1-phi2))>
1724   }
1725
1726   // 3-p temporary one
1727   Double_t three2npp1n1nW1ppW1W1PtRP=0.; // <w1 w2 w3 cos(n(2psi1-phi2-phi3))> // OK!!!
1728   if(dM1pp11PtRP)
1729   {
1730    three2npp1n1nW1ppW1W1PtRP = (q2xW1PtRP*(dQnkX[0][0]*dQnkX[0][0]-dQnkY[0][0]*dQnkY[0][0])+2.*q2yW1PtRP*dQnkX[0][0]*dQnkY[0][0]-2.*(qxW2PtRP*dQnkX[0][0]+qyW2PtRP*dQnkY[0][0])-(q2xW1PtRP*dQnkX[1][1]+q2yW1PtRP*dQnkY[1][1])+2.*dS13mPtRP)/dM1pp11PtRP; // CORRECT !!! <w1 w2 w3 cos(n(2psi1-phi2-phi3))>
1731   }
1732
1733   // 3-p temporary one
1734   Double_t three1npp1n2nW0ppW1W2PtRP=0.; // <w2 w3^2 cos(n(psi1+phi2-2*phi3))> // OK!!!
1735   if(dM0pp12PtRP)
1736   {
1737    three1npp1n2nW0ppW1W2PtRP = (qxPtRP*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])-qyPtRP*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1])-(qxW2PtRP*dQnkX[0][0]+qyW2PtRP*dQnkY[0][0])-(q2xW1PtRP*dQnkX[1][1]+q2yW1PtRP*dQnkY[1][1])-(qxPtRP*dQnkX[0][2]+qyPtRP*dQnkY[0][2])+2.*dS13mPtRP)/dM0pp12PtRP; // CORRECT !!! <w2 w3^2 cos(n(psi1+phi2-2.*phi3))>
1738   }
1739
1740   /*
1741   // 4-p RP part
1742   Double_t four1npp1n1n1nW1W1W1=0.; // <w1 w2 w3 cos(n(psi1+phi1-phi2-phi3))>
1743   if(dM0pp111PtRP)
1744   {
1745    four1npp1n1n1nW1W1W1 = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPtRP*dQnkX[0][0]+qyPtRP*dQnkY[0][0])-2.*dM1pp11PtRP*two1n1nW1ppW1W1PtRP-dM1pp11PtRP*three2npp1n1nW1ppW1W1PtRP-dM0pp12PtRP*three1npp1n2nW0ppW1W2PtRP-2.*dM0pp12PtRP*two1npp1nW1W2PtRP-3.*dM2pp1PtRP*two1npp1nW2ppW1PtRP-2.*dM1pp2PtRP-dM1pp2PtRP*two2npp2nW1ppW2PtRP-dM0pp3PtRP*two1npp1nW3PtRP-dS13mPtRP)/(dM0pp111PtRP);
1746   }
1747   */
1748
1749   /*
1750   // 4-p POI part
1751   Double_t four1npp1n1n1nW1W1W1RP=0.;
1752   if(dM0p111PtRP>0&&mPrimePtRP>0&&nRP>0)
1753   {
1754    four1npp1n1n1nW1W1W1RP = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimePtRP*dQnkX[0][0]+qyPrimePtRP*dQnkY[0][0])-2.*dSnk[0][1]* (qxPrimePtRP*dQnkX[0][0]+qyPrimePtRP*dQnkY[0][0])+2.*(qxPrimePtRP*dQnkX[0][2]+qyPrimePtRP*dQnkY[0][2])-qxPrimePtRP*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])+qyPrimePtRP*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1]))/dM0p111PtRP;
1755   }
1756   */
1757
1758   // 4-p RP and POI in all combinations (full, partial and no overlap)
1759   Double_t four1npp1n1n1nW1W1W1PtRP=0.;
1760   if(dM0pp111PtRP+dM0p111PtRP)
1761   {
1762   four1npp1n1n1nW1W1W1PtRP = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPtRP*dQnkX[0][0]+qyPtRP*dQnkY[0][0])-2.*dM1pp11PtRP*two1n1nW1ppW1W1PtRP-dM1pp11PtRP*three2npp1n1nW1ppW1W1PtRP-dM0pp12PtRP*three1npp1n2nW0ppW1W2PtRP-2.*dM0pp12PtRP*two1npp1nW1W2PtRP-3.*dM2pp1PtRP*two1npp1nW2ppW1PtRP-2.*dM1pp2PtRP-dM1pp2PtRP*two2npp2nW1ppW2PtRP-dM0pp3PtRP*two1npp1nW3PtRP-dS13mPtRP+(pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimePtRP*dQnkX[0][0]+qyPrimePtRP*dQnkY[0][0])-2.*dSnk[0][1]* (qxPrimePtRP*dQnkX[0][0]+qyPrimePtRP*dQnkY[0][0])+2.*(qxPrimePtRP*dQnkX[0][2]+qyPrimePtRP*dQnkY[0][2])-qxPrimePtRP*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])+qyPrimePtRP*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1]))/(dM0pp111PtRP+dM0p111PtRP);
1763
1764    f4WPerPtBin1n1n1n1nRP->Fill(fPtMin+(bin-1)*dBinWidthPt,four1npp1n1n1nW1W1W1PtRP,dM0pp111PtRP+dM0p111PtRP);
1765   } // end of if(dM0pp111PtRP+dM0p111PtRP)
1766  } // for(Int_t bin=1;bin<(fnBinsPt+1);bin++) // loop over pt-bins
1767
1768  delete ptReq1nPrime;
1769  delete ptImq1nPrime;
1770  delete ptReq2nPrime;
1771  delete ptImq2nPrime;
1772  delete ptReq1n;
1773  delete ptImq1n;
1774  delete ptReq2n;
1775  delete ptImq2n;
1776
1777  delete req1nW2Pt;
1778  delete imq1nW2Pt;
1779  delete req2nW1Pt;
1780  delete imq2nW1Pt;
1781  delete sumOfW1upTomPt;
1782  delete sumOfW2upTomPt;
1783  delete sumOfW3upTomPt;
1784  //...........................................................................................................
1785
1786  //...........................................................................................................
1787  // PrimePrime Eta RP
1788  Double_t qxEtaRP=0.,qyEtaRP=0.,q2xEtaRP=0.,q2yEtaRP=0.;//add comments for these variable
1789  Double_t qxW2EtaRP=0.,qyW2EtaRP=0.,q2xW1EtaRP=0.,q2yW1EtaRP=0.;//add comments for these variable
1790  Double_t dS11mEtaRP=0.; // to be improved (name)
1791  Double_t dS12mEtaRP=0.; // to be improved (name)
1792  Double_t dS13mEtaRP=0.; // to be improved (name)
1793  Double_t mEtaRPHere=0.; // to be improved (name)
1794
1795  Double_t dM1pp11EtaRP=0.; // to be improved (name)
1796  Double_t dM0pp111EtaRP=0.; // to be improved (name)
1797  Double_t dM0pp12EtaRP=0.;
1798  Double_t dM2pp1EtaRP=0.;
1799  Double_t dM1pp2EtaRP=0.;
1800  Double_t dM0pp3EtaRP=0.;
1801
1802  // Prime Eta RP
1803  Double_t qxPrimeEtaRPHere=0.,qyPrimeEtaRPHere=0.;
1804  Double_t mPrimeEtaRPHere=0.;
1805  Double_t dM0p111EtaRP=0.; // to be improved (name)
1806
1807  for(Int_t bin=1;bin<(fnBinsEta+1);bin++) // loop over eta-bins
1808  {
1809   // q'':
1810   qxEtaRP = (etaReq1n->GetBinContent(bin))*(etaReq1n->GetBinEntries(bin));
1811   qyEtaRP = (etaImq1n->GetBinContent(bin))*(etaImq1n->GetBinEntries(bin));
1812   q2xEtaRP = (etaReq2n->GetBinContent(bin))*(etaReq2n->GetBinEntries(bin));
1813   q2yEtaRP = (etaImq2n->GetBinContent(bin))*(etaImq2n->GetBinEntries(bin));
1814
1815   qxW2EtaRP = (req1nW2Eta->GetBinContent(bin))*(req1nW2Eta->GetBinEntries(bin));
1816   qyW2EtaRP = (imq1nW2Eta->GetBinContent(bin))*(imq1nW2Eta->GetBinEntries(bin));
1817
1818   q2xW1EtaRP = (req2nW1Eta->GetBinContent(bin))*(req2nW1Eta->GetBinEntries(bin));
1819   q2yW1EtaRP = (imq2nW1Eta->GetBinContent(bin))*(imq2nW1Eta->GetBinEntries(bin));
1820
1821   dS11mEtaRP = (sumOfW1upTomEta->GetBinContent(bin))*(sumOfW1upTomEta->GetBinEntries(bin));
1822   dS12mEtaRP = (sumOfW2upTomEta->GetBinContent(bin))*(sumOfW2upTomEta->GetBinEntries(bin));
1823   dS13mEtaRP = (sumOfW3upTomEta->GetBinContent(bin))*(sumOfW3upTomEta->GetBinEntries(bin));
1824
1825   mEtaRPHere = sumOfW1upTomEta->GetBinEntries(bin); // to be improved
1826
1827   dM1pp11EtaRP=dS11mEtaRP*(dSnk[1][0]-dSnk[0][1])-2.*dS12mEtaRP*dSnk[0][0]+2.*dS13mEtaRP;
1828   dM1pp2EtaRP=dS11mEtaRP*dSnk[0][1]-dS13mEtaRP;
1829   dM2pp1EtaRP=dS12mEtaRP*dSnk[0][0]-dS13mEtaRP;
1830   dM0pp3EtaRP=mEtaRPHere*dSnk[0][2]-dS13mEtaRP;
1831   dM0pp12EtaRP=mEtaRPHere*dSnk[0][0]*dSnk[0][1]-dM2pp1EtaRP-dM1pp2EtaRP-dM0pp3EtaRP-dS13mEtaRP;
1832   dM0pp111EtaRP=mEtaRPHere*dSnk[2][0]-3.*dM1pp11EtaRP-3.*dM0pp12EtaRP-3.*dM2pp1EtaRP-3.*dM1pp2EtaRP-dM0pp3EtaRP-dS13mEtaRP;
1833
1834   // q':
1835   qxPrimeEtaRPHere = (etaReq1nPrime->GetBinContent(bin))*(etaReq1nPrime->GetBinEntries(bin));
1836   qyPrimeEtaRPHere = (etaImq1nPrime->GetBinContent(bin))*(etaImq1nPrime->GetBinEntries(bin));
1837
1838   mPrimeEtaRPHere = etaReq1nPrime->GetBinEntries(bin); // to be improved
1839   dM0p111EtaRP=mPrimeEtaRPHere*(dSnk[2][0]-3.*dSnk[0][1]*dSnk[0][0]+2.*dSnk[0][2]);
1840
1841   // 2-p the needed one
1842   Double_t two1n1nWPerEtaBinRP=0.;
1843   if((mEtaRPHere+mPrimeEtaRPHere)*dSnk[0][0]-dS11mEtaRP>0)
1844   {
1845    two1n1nWPerEtaBinRP = (qxEtaRP*dQnkX[0][0]+qyEtaRP*dQnkY[0][0]+qxPrimeEtaRPHere*dQnkX[0][0]+qyPrimeEtaRPHere*dQnkY[0][0]-dS11mEtaRP)/((mEtaRPHere+mPrimeEtaRPHere)*dSnk[0][0]-dS11mEtaRP);
1846
1847    f2WPerEtaBin1n1nRP->Fill(fEtaMin+(bin-1)*dBinWidthEta,two1n1nWPerEtaBinRP,(mEtaRPHere+mPrimeEtaRPHere)*dSnk[0][0]-dS11mEtaRP); // <2'>_{n|n}
1848   }
1849
1850   // 2-p the temporary one
1851   Double_t two1n1nW1ppW1W1EtaRP=0.; // <w1 w2 w3 cos(n(phi2-phi3))> // OK!!!
1852   if(dM1pp11EtaRP)
1853   {
1854    two1n1nW1ppW1W1EtaRP = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*dS11mEtaRP-2.*(qxW2EtaRP*dQnkX[0][0]+qyW2EtaRP*dQnkY[0][0])-dS11mEtaRP*dSnk[0][1]+2.*dS13mEtaRP)/dM1pp11EtaRP; // CORRECT !!! <w1 w2 w3 cos(n(phi2-phi3))>
1855   }
1856
1857   // 2-p the temporary one
1858   Double_t two1npp1nW1W2EtaRP=0.; // <w2 w3^2 cos(n(psi1-phi2))> // OK !!!
1859   if(dM0pp12EtaRP)
1860   {
1861    two1npp1nW1W2EtaRP = (dSnk[0][1]*(qxEtaRP*dQnkX[0][0]+qyEtaRP*dQnkY[0][0])-(qxW2EtaRP*dQnkX[0][0]+qyW2EtaRP*dQnkY[0][0])-dM1pp2EtaRP-(qxEtaRP*dQnkX[0][2]+qyEtaRP*dQnkY[0][2])+dS13mEtaRP)/dM0pp12EtaRP; // CORRECT !!! <w2 w3^2 cos(n(psi1-phi2))>
1862   }
1863
1864   // 2-p the temporary one
1865   Double_t two1npp1nW2ppW1EtaRP=0.; // <w1^2 w2 cos(n(psi1-phi2))> // OK !!!
1866   if(dM2pp1EtaRP)
1867   {
1868    two1npp1nW2ppW1EtaRP = ((qxW2EtaRP*dQnkX[0][0]+qyW2EtaRP*dQnkY[0][0])-dS13mEtaRP)/dM2pp1EtaRP; // CORRECT !!! <w1^2 w2 cos(n(psi1-phi2))>
1869   }
1870
1871   // 2-p the temporary one
1872   Double_t two2npp2nW1ppW2EtaRP=0.; // <w1 w2^2 cos(2n(psi1-phi2))> OK !!!
1873   if(dM1pp2EtaRP)
1874   {
1875    two2npp2nW1ppW2EtaRP = ((q2xW1EtaRP*dQnkX[1][1]+q2yW1EtaRP*dQnkY[1][1])-dS13mEtaRP)/dM1pp2EtaRP; // CORRECT !!! <w1 w2^2 cos(2n(psi1-phi2))>
1876   }
1877
1878   // 2-p the temporary one
1879   Double_t two1npp1nW3EtaRP=0.; // <w2^3 cos(n(psi1-phi2))> // OK !!!
1880   if(dM0pp3EtaRP)
1881   {
1882    two1npp1nW3EtaRP = (qxEtaRP*dQnkX[0][2]+qyEtaRP*dQnkY[0][2]-dS13mEtaRP)/dM0pp3EtaRP; // CORRECT !!! <w2^3 cos(n(psi1-phi2))>
1883   }
1884
1885   // 3-p the temporary one
1886   Double_t three2npp1n1nW1ppW1W1EtaRP=0.; // <w1 w2 w3 cos(n(2psi1-phi2-phi3))> // OK!!!
1887   if(dM1pp11EtaRP)
1888   {
1889    three2npp1n1nW1ppW1W1EtaRP = (q2xW1EtaRP*(dQnkX[0][0]*dQnkX[0][0]-dQnkY[0][0]*dQnkY[0][0])+2.*q2yW1EtaRP*dQnkX[0][0]*dQnkY[0][0]-2.*(qxW2EtaRP*dQnkX[0][0]+qyW2EtaRP*dQnkY[0][0])-(q2xW1EtaRP*dQnkX[1][1]+q2yW1EtaRP*dQnkY[1][1])+2.*dS13mEtaRP)/dM1pp11EtaRP; // CORRECT !!! <w1 w2 w3 cos(n(2psi1-phi2-phi3))>
1890   }
1891
1892   // 3-p the temporary one
1893   Double_t three1npp1n2nW0ppW1W2EtaRP=0.; // <w2 w3^2 cos(n(psi1+phi2-2*phi3))> // OK!!!
1894   if(dM0pp12EtaRP)
1895   {
1896    three1npp1n2nW0ppW1W2EtaRP = (qxEtaRP*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])-qyEtaRP*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1])-(qxW2EtaRP*dQnkX[0][0]+qyW2EtaRP*dQnkY[0][0])-(q2xW1EtaRP*dQnkX[1][1]+q2yW1EtaRP*dQnkY[1][1])-(qxEtaRP*dQnkX[0][2]+qyEtaRP*dQnkY[0][2])+2.*dS13mEtaRP)/dM0pp12EtaRP; // CORRECT !!! <w2 w3^2 cos(n(psi1+phi2-2.*phi3))>
1897   }
1898
1899   /*
1900   // 4-p RP part
1901   Double_t four1npp1n1n1nW1W1W1=0.; // <w1 w2 w3 cos(n(psi1+phi1-phi2-phi3))>
1902   if(dM0pp111EtaRP)
1903   {
1904    four1npp1n1n1nW1W1W1 = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxEtaRP*dQnkX[0][0]+qyEtaRP*dQnkY[0][0])-2.*dM1pp11EtaRP*two1n1nW1ppW1W1EtaRP-dM1pp11EtaRP*three2npp1n1nW1ppW1W1EtaRP-dM0pp12EtaRP*three1npp1n2nW0ppW1W2EtaRP-2.*dM0pp12EtaRP*two1npp1nW1W2EtaRP-3.*dM2pp1EtaRP*two1npp1nW2ppW1EtaRP-2.*dM1pp2EtaRP-dM1pp2EtaRP*two2npp2nW1ppW2EtaRP-dM0pp3EtaRP*two1npp1nW3EtaRP-dS13mEtaRP)/(dM0pp111EtaRP);
1905   }
1906   */
1907   /*
1908   // 4-p POI part
1909   Double_t four1npp1n1n1nW1W1W1RP=0.;
1910   if(dM0p111EtaRP>0&&mPrimeEtaRPHere>0&&nRP>0)
1911   {
1912    four1npp1n1n1nW1W1W1RP = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimeEtaRPHere*dQnkX[0][0]+qyPrimeEtaRPHere*dQnkY[0][0])-2.*dSnk[0][1]* (qxPrimeEtaRPHere*dQnkX[0][0]+qyPrimeEtaRPHere*dQnkY[0][0])+2.*(qxPrimeEtaRPHere*dQnkX[0][2]+qyPrimeEtaRPHere*dQnkY[0][2])-qxPrimeEtaRPHere*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])+qyPrimeEtaRPHere*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1]))/dM0p111EtaRP;
1913   }
1914   */
1915
1916   // 4-p RP and POI in all combinations (full, partial and no overlap)
1917   Double_t four1npp1n1n1nW1W1W1EtaRP=0.;
1918
1919   if(dM0pp111EtaRP+dM0p111EtaRP)
1920   {
1921    four1npp1n1n1nW1W1W1EtaRP = ((pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxEtaRP*dQnkX[0][0]+qyEtaRP*dQnkY[0][0])-2.*dM1pp11EtaRP*two1n1nW1ppW1W1EtaRP-dM1pp11EtaRP*three2npp1n1nW1ppW1W1EtaRP-dM0pp12EtaRP*three1npp1n2nW0ppW1W2EtaRP-2.*dM0pp12EtaRP*two1npp1nW1W2EtaRP-3.*dM2pp1EtaRP*two1npp1nW2ppW1EtaRP-2.*dM1pp2EtaRP-dM1pp2EtaRP*two2npp2nW1ppW2EtaRP-dM0pp3EtaRP*two1npp1nW3EtaRP-dS13mEtaRP+(pow(dQnkX[0][0],2.)+pow(dQnkY[0][0],2.))*(qxPrimeEtaRPHere*dQnkX[0][0]+qyPrimeEtaRPHere*dQnkY[0][0])-2.*dSnk[0][1]*(qxPrimeEtaRPHere*dQnkX[0][0]+qyPrimeEtaRPHere*dQnkY[0][0])+2.*(qxPrimeEtaRPHere*dQnkX[0][2]+qyPrimeEtaRPHere*dQnkY[0][2])-qxPrimeEtaRPHere*(dQnkX[0][0]*dQnkX[1][1]+dQnkY[0][0]*dQnkY[1][1])+qyPrimeEtaRPHere*(dQnkY[0][0]*dQnkX[1][1]-dQnkX[0][0]*dQnkY[1][1]))/(dM0pp111EtaRP+dM0p111EtaRP);
1922
1923    f4WPerEtaBin1n1n1n1nRP->Fill(fEtaMin+(bin-1)*dBinWidthEta,four1npp1n1n1nW1W1W1EtaRP,dM0pp111EtaRP+dM0p111EtaRP);
1924   }
1925  }
1926
1927  delete etaReq1nPrime;
1928  delete etaImq1nPrime;
1929  delete etaReq2nPrime;
1930  delete etaImq2nPrime;
1931  delete etaReq1n;
1932  delete etaImq1n;
1933  delete etaReq2n;
1934  delete etaImq2n;
1935
1936  delete req1nW2Eta;
1937  delete imq1nW2Eta;
1938  delete req2nW1Eta;
1939  delete imq2nW1Eta;
1940  delete sumOfW1upTomEta;
1941  delete sumOfW2upTomEta;
1942  delete sumOfW3upTomEta;
1943  //...........................................................................................................
1944
1945
1946
1947
1948
1949
1950
1951
1952
1953
1954
1955
1956
1957
1958
1959
1960
1961
1962
1963
1964
1965
1966
1967
1968  delete ptReq1nPrimePrime;
1969  delete ptImq1nPrimePrime;
1970  delete ptReq2nPrimePrime;
1971  delete ptImq2nPrimePrime;
1972
1973  delete req1nW2PrimePrimePt;
1974  delete imq1nW2PrimePrimePt;
1975  delete req2nW1PrimePrimePt;
1976  delete imq2nW1PrimePrimePt;
1977
1978  delete sumOfW1upTomPrimePrimePt;
1979  delete sumOfW2upTomPrimePrimePt;
1980  delete sumOfW3upTomPrimePrimePt;
1981
1982  delete etaReq1nPrimePrime;
1983  delete etaImq1nPrimePrime;
1984  delete etaReq2nPrimePrime;
1985  delete etaImq2nPrimePrime;
1986
1987  delete req1nW2PrimePrimeEta;
1988  delete imq1nW2PrimePrimeEta;
1989  delete req2nW1PrimePrimeEta;
1990  delete imq2nW1PrimePrimeEta;
1991
1992  delete sumOfW1upTomPrimePrimeEta;
1993  delete sumOfW2upTomPrimePrimeEta;
1994  delete sumOfW3upTomPrimePrimeEta;
1995
1996
1997
1998
1999
2000
2001
2002
2003
2004
2005
2006
2007
2008
2009
2010
2011
2012
2013
2014
2015
2016
2017
2018
2019
2020
2021
2022
2023
2024
2025
2026
2027
2028
2029
2030
2031
2032
2033
2034
2035
2036
2037
2038
2039
2040
2041
2042
2043
2044
2045
2046
2047
2048
2049
2050
2051
2052
2053
2054
2055
2056
2057
2058
2059
2060
2061
2062
2063
2064
2065
2066
2067
2068
2069
2070
2071
2072
2073
2074
2075
2076
2077
2078
2079
2080
2081
2082
2083
2084
2085
2086
2087
2088
2089
2090
2091
2092
2093
2094
2095
2096
2097
2098
2099
2100
2101
2102
2103
2104
2105
2106
2107
2108
2109
2110
2111
2112
2113
2114
2115
2116
2117
2118
2119
2120
2121
2122
2123
2124
2125
2126
2127
2128
2129
2130
2131
2132
2133
2134
2135
2136
2137
2138
2139
2140  Bool_t nestedLoops = kFALSE;
2141
2142
2143
2144
2145  if(nestedLoops) // to be improved
2146  {
2147  //--------------------------------------------------------------------------------------------------------------------------------
2148  //
2149  //                                          **********************
2150  //                                          **** NESTED LOOPS ****
2151  //                                          **********************
2152  //
2153  // Remark 1: multi-particle correlations calculated with nested loops are stored in fDirectCorrelations.
2154  // Remark 2: binning of fDirectCorrelations: bins 0..100 - correlations needed for integrated flow; bins 100..200 - correlations needed for differential flow (taking as an example bin 0.5 < pt < 0.6)
2155  //
2156  // binning details of fDirectCorrelations (integrated flow):
2157  //..........................................................................
2158  // 1st bin: weighted <2>_{n|n}   = <w1   w2   cos( n*(phi1-phi2))>
2159  // 2nd bin: weighted <2>_{2n|2n} = <w1^2 w2^2 cos(2n*(phi1-phi2))>
2160  // 3rd bin: weighted <2>_{3n|3n} = <w1^3 w2^3 cos(3n*(phi1-phi2))>
2161  // 4th bin: weighted <2>_{4n|4n} = <w1^4 w2^4 cos(4n*(phi1-phi2))>
2162  // 5th bin: weighted <2>_{n|n} = <w1^3 w2 cos(n*(phi1-phi2))>
2163
2164  // 11th bin: weighted <3>_{2n|n,n} = <w1^2 w2 w3 cos(n*(2phi1-phi2-phi3))>
2165  //..........................................................................
2166
2167
2168  Double_t phi1=0., phi2=0., phi3=0., phi4=0.;
2169  // phi5=0., phi6=0., phi7=0., phi8=0.;
2170  Double_t wPhi1=1., wPhi2=1., wPhi3=1., wPhi4=1.;
2171  // wPhi5=1., wPhi6=1., wPhi7=1., wPhi8=1.;
2172
2173
2174
2175
2176
2177
2178  Double_t tempLoop = 0.;
2179  Int_t tempCounter = 0;
2180
2181
2182
2183
2184
2185  for(Int_t i1=0;i1<dMult;i1++)
2186  {
2187   fTrack=anEvent->GetTrack(i1);
2188   phi1=fTrack->Phi();
2189   if(phiWeights) wPhi1 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi1*nBinsPhi/TMath::TwoPi())));
2190   for(Int_t i2=0;i2<dMult;i2++)
2191   {
2192    if(i2==i1)continue;
2193    fTrack=anEvent->GetTrack(i2);
2194    phi2=fTrack->Phi();
2195    if(phiWeights) wPhi2 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi2*nBinsPhi/TMath::TwoPi())));
2196    // 2-p
2197    fDirectCorrelations->Fill(0.,cos(n*(phi1-phi2)),wPhi1*wPhi2);                  // <w1   w2   cos( n*(phi1-phi2))>
2198    fDirectCorrelations->Fill(1.,cos(2.*n*(phi1-phi2)),pow(wPhi1,2)*pow(wPhi2,2)); // <w1^2 w2^2 cos(2n*(phi1-phi2))>
2199    fDirectCorrelations->Fill(2.,cos(3.*n*(phi1-phi2)),pow(wPhi1,3)*pow(wPhi2,3)); // <w1^3 w2^3 cos(3n*(phi1-phi2))>
2200    fDirectCorrelations->Fill(3.,cos(4.*n*(phi1-phi2)),pow(wPhi1,4)*pow(wPhi2,4)); // <w1^4 w2^4 cos(4n*(phi1-phi2))>
2201    fDirectCorrelations->Fill(4.,cos(n*(phi1-phi2)),pow(wPhi1,3)*wPhi2); // <w1^3 w2 cos(n*(phi1-phi2))>
2202   }
2203  }
2204
2205
2206
2207
2208  /*
2209
2210  for(Int_t i1=0;i1<dMult;i1++)
2211  {
2212   fTrack=anEvent->GetTrack(i1);
2213   phi1=fTrack->Phi();
2214   if(phiWeights) wPhi1 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi1*nBinsPhi/TMath::TwoPi())));
2215   for(Int_t i2=0;i2<dMult;i2++)
2216   {
2217    if(i2==i1)continue;
2218    fTrack=anEvent->GetTrack(i2);
2219    phi2=fTrack->Phi();
2220    if(phiWeights) wPhi2 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi2*nBinsPhi/TMath::TwoPi())));
2221    for(Int_t i3=0;i3<dMult;i3++)
2222    {
2223     if(i3==i1||i3==i2)continue;
2224     fTrack=anEvent->GetTrack(i3);
2225     phi3=fTrack->Phi();
2226     if(phiWeights) wPhi3 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi3*nBinsPhi/TMath::TwoPi())));
2227     // 2-p
2228     fDirectCorrelations->Fill(5.,cos(n*(phi1-phi2)),wPhi1*wPhi2*pow(wPhi3,2)); // <w1 w2 w3^2 cos(n*(phi1-phi2))>
2229
2230     // 3-p
2231     fDirectCorrelations->Fill(10.,cos(2.*n*phi1-n*(phi2+phi3)),pow(wPhi1,2)*wPhi2*wPhi3); // <w1^2 w2 w3 cos(n*(2phi1-phi2-phi3))>
2232
2233
2234     fDirectCorrelations->Fill(6.,cos(3.*n*phi1-2.*n*phi2-n*phi3),pow(wPhi1,3)*pow(wPhi2,2)*wPhi3);           //<3>_{3n|2n,n}
2235     fDirectCorrelations->Fill(7.,cos(4.*n*phi1-2.*n*phi2-2.*n*phi3),pow(wPhi1,4)*pow(wPhi2,2)*pow(wPhi3,2)); //<3>_{4n|2n,2n}
2236     fDirectCorrelations->Fill(8.,cos(4.*n*phi1-3.*n*phi2-n*phi3),pow(wPhi1,4)*pow(wPhi2,3)*wPhi3);           //<3>_{4n|3n,n}
2237
2238
2239    }
2240   }
2241  }
2242  */
2243
2244
2245
2246
2247
2248
2249  //<4>_{n,n|n,n}, <4>_{2n,n|2n,n}, <4>_{2n,2n|2n,2n}, <4>_{3n|n,n,n}, <4>_{3n,n|3n,n}, <4>_{3n,n|2n,2n} and <4>_{4n|2n,n,n}
2250  for(Int_t i1=0;i1<dMult;i1++)
2251  {
2252   fTrack=anEvent->GetTrack(i1);
2253   phi1=fTrack->Phi();
2254   if(phiWeights) wPhi1 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi1*nBinsPhi/TMath::TwoPi())));
2255   for(Int_t i2=0;i2<dMult;i2++)
2256   {
2257    if(i2==i1)continue;
2258    fTrack=anEvent->GetTrack(i2);
2259    phi2=fTrack->Phi();
2260    if(phiWeights) wPhi2 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi2*nBinsPhi/TMath::TwoPi())));
2261    for(Int_t i3=0;i3<dMult;i3++)
2262    {
2263     if(i3==i1||i3==i2)continue;
2264     fTrack=anEvent->GetTrack(i3);
2265     phi3=fTrack->Phi();
2266     if(phiWeights) wPhi3 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi3*nBinsPhi/TMath::TwoPi())));
2267     for(Int_t i4=0;i4<dMult;i4++)
2268     {
2269      if(i4==i1||i4==i2||i4==i3)continue;
2270      fTrack=anEvent->GetTrack(i4);
2271      phi4=fTrack->Phi();
2272      if(phiWeights) wPhi4 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi4*nBinsPhi/TMath::TwoPi())));
2273      fDirectCorrelations->Fill(20.,cos(n*phi1+n*phi2-n*phi3-n*phi4),wPhi1*wPhi2*wPhi3*wPhi4); // <4>_{n,n|n,n} = <w1 w2 w3 w4 cos(n*(phi1+phi2-phi3-phi4))>
2274
2275
2276      fDirectCorrelations->Fill(11.,cos(2.*n*phi1+n*phi2-2.*n*phi3-n*phi4),wPhi1*wPhi2*wPhi3*wPhi4);      //<4>_{2n,n|2n,n}
2277      fDirectCorrelations->Fill(12.,cos(2.*n*phi1+2*n*phi2-2.*n*phi3-2.*n*phi4),wPhi1*wPhi2*wPhi3*wPhi4); //<4>_{2n,2n|2n,2n}
2278      fDirectCorrelations->Fill(13.,cos(3.*n*phi1-n*phi2-n*phi3-n*phi4),wPhi1*wPhi2*wPhi3*wPhi4);         //<4>_{3n|n,n,n}
2279      fDirectCorrelations->Fill(14.,cos(3.*n*phi1+n*phi2-3.*n*phi3-n*phi4),wPhi1*wPhi2*wPhi3*wPhi4);      //<4>_{3n,n|3n,n}
2280      fDirectCorrelations->Fill(15.,cos(3.*n*phi1+n*phi2-2.*n*phi3-2.*n*phi4),wPhi1*wPhi2*wPhi3*wPhi4);   //<4>_{3n,n|2n,2n}
2281      fDirectCorrelations->Fill(16.,cos(4.*n*phi1-2.*n*phi2-n*phi3-n*phi4),wPhi1*wPhi2*wPhi3*wPhi4);      //<4>_{4n|2n,n,n}
2282
2283     }
2284    }
2285   }
2286  }
2287
2288
2289
2290  /*
2291
2292  //<5>_{2n,n,n,n,n}, //<5>_{2n,2n|2n,n,n}, <5>_{3n,n|2n,n,n} and <5>_{4n|n,n,n,n}
2293  for(Int_t i1=0;i1<dMult;i1++)
2294  {
2295   //cout<<"i1 = "<<i1<<endl;
2296   fTrack=anEvent->GetTrack(i1);
2297   phi1=fTrack->Phi();
2298   if(phiWeights) wPhi1 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi1*nBinsPhi/TMath::TwoPi())));
2299   for(Int_t i2=0;i2<dMult;i2++)
2300   {
2301    if(i2==i1)continue;
2302    fTrack=anEvent->GetTrack(i2);
2303    phi2=fTrack->Phi();
2304    if(phiWeights) wPhi2 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi2*nBinsPhi/TMath::TwoPi())));
2305    for(Int_t i3=0;i3<dMult;i3++)
2306    {
2307     if(i3==i1||i3==i2)continue;
2308     fTrack=anEvent->GetTrack(i3);
2309     phi3=fTrack->Phi();
2310     if(phiWeights) wPhi3 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi3*nBinsPhi/TMath::TwoPi())));
2311     for(Int_t i4=0;i4<dMult;i4++)
2312     {
2313      if(i4==i1||i4==i2||i4==i3)continue;
2314      fTrack=anEvent->GetTrack(i4);
2315      phi4=fTrack->Phi();
2316      if(phiWeights) wPhi4 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi4*nBinsPhi/TMath::TwoPi())));
2317      for(Int_t i5=0;i5<dMult;i5++)
2318      {
2319       if(i5==i1||i5==i2||i5==i3||i5==i4)continue;
2320       fTrack=anEvent->GetTrack(i5);
2321       phi5=fTrack->Phi();
2322       if(phiWeights) wPhi5 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi5*nBinsPhi/TMath::TwoPi())));
2323       fDirectCorrelations->Fill(18.,cos(2.*n*phi1+n*phi2-n*phi3-n*phi4-n*phi5),wPhi1*wPhi2*wPhi3*wPhi4*wPhi5);       //<5>_{2n,n|n,n,n}
2324       fDirectCorrelations->Fill(19.,cos(2.*n*phi1+2.*n*phi2-2.*n*phi3-n*phi4-n*phi5),wPhi1*wPhi2*wPhi3*wPhi4*wPhi5); //<5>_{2n,2n|2n,n,n}
2325       fDirectCorrelations->Fill(20.,cos(3.*n*phi1+n*phi2-2.*n*phi3-n*phi4-n*phi5),wPhi1*wPhi2*wPhi3*wPhi4*wPhi5);    //<5>_{3n,n|2n,n,n}
2326       fDirectCorrelations->Fill(21.,cos(4.*n*phi1-n*phi2-n*phi3-n*phi4-n*phi5),wPhi1*wPhi2*wPhi3*wPhi4*wPhi5);       //<5>_{4n|n,n,n,n}
2327      }
2328     }
2329    }
2330   }
2331  }
2332
2333
2334  //<6>_{n,n,n,n,n,n}, <6>_{2n,n,n|2n,n,n}, <6>_{2n,2n|n,n,n,n} and <6>_{3n,n|n,n,n,n}
2335  for(Int_t i1=0;i1<dMult;i1++)
2336  {
2337   //cout<<"i1 = "<<i1<<endl;
2338   fTrack=anEvent->GetTrack(i1);
2339   phi1=fTrack->Phi();
2340   if(phiWeights) wPhi1 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi1*nBinsPhi/TMath::TwoPi())));
2341   for(Int_t i2=0;i2<dMult;i2++)
2342   {
2343    if(i2==i1)continue;
2344    fTrack=anEvent->GetTrack(i2);
2345    phi2=fTrack->Phi();
2346    if(phiWeights) wPhi2 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi2*nBinsPhi/TMath::TwoPi())));
2347    for(Int_t i3=0;i3<dMult;i3++)
2348    {
2349     if(i3==i1||i3==i2)continue;
2350     fTrack=anEvent->GetTrack(i3);
2351     phi3=fTrack->Phi();
2352     if(phiWeights) wPhi3 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi3*nBinsPhi/TMath::TwoPi())));
2353     for(Int_t i4=0;i4<dMult;i4++)
2354     {
2355      if(i4==i1||i4==i2||i4==i3)continue;
2356      fTrack=anEvent->GetTrack(i4);
2357      phi4=fTrack->Phi();
2358      if(phiWeights) wPhi4 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi4*nBinsPhi/TMath::TwoPi())));
2359      for(Int_t i5=0;i5<dMult;i5++)
2360      {
2361       if(i5==i1||i5==i2||i5==i3||i5==i4)continue;
2362       fTrack=anEvent->GetTrack(i5);
2363       phi5=fTrack->Phi();
2364       if(phiWeights) wPhi5 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi5*nBinsPhi/TMath::TwoPi())));
2365       for(Int_t i6=0;i6<dMult;i6++)
2366       {
2367        if(i6==i1||i6==i2||i6==i3||i6==i4||i6==i5)continue;
2368        fTrack=anEvent->GetTrack(i6);
2369        phi6=fTrack->Phi();
2370        if(phiWeights) wPhi6 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi6*nBinsPhi/TMath::TwoPi())));
2371        fDirectCorrelations->Fill(23.,cos(n*phi1+n*phi2+n*phi3-n*phi4-n*phi5-n*phi6),wPhi1*wPhi2*wPhi3*wPhi4*wPhi5*wPhi6);       //<6>_{n,n,n|n,n,n}
2372        fDirectCorrelations->Fill(24.,cos(2.*n*phi1+n*phi2+n*phi3-2.*n*phi4-n*phi5-n*phi6),wPhi1*wPhi2*wPhi3*wPhi4*wPhi5*wPhi6); //<6>_{2n,n,n|2n,n,n}
2373        fDirectCorrelations->Fill(25.,cos(2.*n*phi1+2.*n*phi2-n*phi3-n*phi4-n*phi5-n*phi6),wPhi1*wPhi2*wPhi3*wPhi4*wPhi5*wPhi6); //<6>_{2n,2n|n,n,n,n}
2374        fDirectCorrelations->Fill(26.,cos(3.*n*phi1+n*phi2-n*phi3-n*phi4-n*phi5-n*phi6),wPhi1*wPhi2*wPhi3*wPhi4*wPhi5*wPhi6);    //<6>_{3n,n|n,n,n,n}
2375       }
2376      }
2377     }
2378    }
2379   }
2380  }
2381
2382
2383
2384  //<7>_{2n,n,n|n,n,n,n}
2385  for(Int_t i1=0;i1<dMult;i1++)
2386  {
2387   //cout<<"i1 = "<<i1<<endl;
2388   fTrack=anEvent->GetTrack(i1);
2389   phi1=fTrack->Phi();
2390   if(phiWeights) wPhi1 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi1*nBinsPhi/TMath::TwoPi())));
2391   for(Int_t i2=0;i2<dMult;i2++)
2392   {
2393    if(i2==i1)continue;
2394    fTrack=anEvent->GetTrack(i2);
2395    phi2=fTrack->Phi();
2396    if(phiWeights) wPhi2 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi2*nBinsPhi/TMath::TwoPi())));
2397    for(Int_t i3=0;i3<dMult;i3++)
2398    {
2399     if(i3==i1||i3==i2)continue;
2400     fTrack=anEvent->GetTrack(i3);
2401     phi3=fTrack->Phi();
2402     if(phiWeights) wPhi3 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi3*nBinsPhi/TMath::TwoPi())));
2403     for(Int_t i4=0;i4<dMult;i4++)
2404     {
2405      if(i4==i1||i4==i2||i4==i3)continue;
2406      fTrack=anEvent->GetTrack(i4);
2407      phi4=fTrack->Phi();
2408      if(phiWeights) wPhi4 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi4*nBinsPhi/TMath::TwoPi())));
2409      for(Int_t i5=0;i5<dMult;i5++)
2410      {
2411       if(i5==i1||i5==i2||i5==i3||i5==i4)continue;
2412       fTrack=anEvent->GetTrack(i5);
2413       phi5=fTrack->Phi();
2414       if(phiWeights) wPhi5 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi5*nBinsPhi/TMath::TwoPi())));
2415       for(Int_t i6=0;i6<dMult;i6++)
2416       {
2417        if(i6==i1||i6==i2||i6==i3||i6==i4||i6==i5)continue;
2418        fTrack=anEvent->GetTrack(i6);
2419        phi6=fTrack->Phi();
2420        if(phiWeights) wPhi6 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi6*nBinsPhi/TMath::TwoPi())));
2421        for(Int_t i7=0;i7<dMult;i7++)
2422        {
2423         if(i7==i1||i7==i2||i7==i3||i7==i4||i7==i5||i7==i6)continue;
2424         fTrack=anEvent->GetTrack(i7);
2425         phi7=fTrack->Phi();
2426         if(phiWeights) wPhi7 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi7*nBinsPhi/TMath::TwoPi())));
2427         fDirectCorrelations->Fill(28.,cos(2.*n*phi1+n*phi2+n*phi3-n*phi4-n*phi5-n*phi6-n*phi7),wPhi1*wPhi2*wPhi3*wPhi4*wPhi5*wPhi6*wPhi7);//<7>_{2n,n,n|n,n,n,n}
2428        }
2429       }
2430      }
2431     }
2432    }
2433   }
2434  }
2435
2436
2437
2438  //<8>_{n,n,n,n|n,n,n,n}
2439  for(Int_t i1=0;i1<dMult;i1++)
2440  {
2441   cout<<"i1 = "<<i1<<endl;
2442   fTrack=anEvent->GetTrack(i1);
2443   phi1=fTrack->Phi();
2444   if(phiWeights) wPhi1 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi1*nBinsPhi/TMath::TwoPi())));
2445   for(Int_t i2=0;i2<dMult;i2++)
2446   {
2447    if(i2==i1)continue;
2448    fTrack=anEvent->GetTrack(i2);
2449    phi2=fTrack->Phi();
2450    if(phiWeights) wPhi2 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi2*nBinsPhi/TMath::TwoPi())));
2451    for(Int_t i3=0;i3<dMult;i3++)
2452    {
2453     if(i3==i1||i3==i2)continue;
2454     fTrack=anEvent->GetTrack(i3);
2455     phi3=fTrack->Phi();
2456     if(phiWeights) wPhi3 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi3*nBinsPhi/TMath::TwoPi())));
2457     for(Int_t i4=0;i4<dMult;i4++)
2458     {
2459      if(i4==i1||i4==i2||i4==i3)continue;
2460      fTrack=anEvent->GetTrack(i4);
2461      phi4=fTrack->Phi();
2462      if(phiWeights) wPhi4 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi4*nBinsPhi/TMath::TwoPi())));
2463      for(Int_t i5=0;i5<dMult;i5++)
2464      {
2465       if(i5==i1||i5==i2||i5==i3||i5==i4)continue;
2466       fTrack=anEvent->GetTrack(i5);
2467       phi5=fTrack->Phi();
2468       if(phiWeights) wPhi5 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi5*nBinsPhi/TMath::TwoPi())));
2469       for(Int_t i6=0;i6<dMult;i6++)
2470       {
2471        if(i6==i1||i6==i2||i6==i3||i6==i4||i6==i5)continue;
2472        fTrack=anEvent->GetTrack(i6);
2473        phi6=fTrack->Phi();
2474        if(phiWeights) wPhi6 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi6*nBinsPhi/TMath::TwoPi())));
2475        for(Int_t i7=0;i7<dMult;i7++)
2476        {
2477         if(i7==i1||i7==i2||i7==i3||i7==i4||i7==i5||i7==i6)continue;
2478         fTrack=anEvent->GetTrack(i7);
2479         phi7=fTrack->Phi();
2480         if(phiWeights) wPhi7 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi7*nBinsPhi/TMath::TwoPi())));
2481         for(Int_t i8=0;i8<dMult;i8++)
2482         {
2483          if(i8==i1||i8==i2||i8==i3||i8==i4||i8==i5||i8==i6||i8==i7)continue;
2484          fTrack=anEvent->GetTrack(i8);
2485          phi8=fTrack->Phi();
2486          if(phiWeights) wPhi8 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi8*nBinsPhi/TMath::TwoPi())));
2487          fDirectCorrelations->Fill(30.,cos(n*phi1+n*phi2+n*phi3+n*phi4-n*phi5-n*phi6-n*phi7-n*phi8),wPhi1*wPhi2*wPhi3*wPhi4*wPhi5*wPhi6*wPhi7*wPhi8);//<8>_{n,n,n,n|n,n,n,n}
2488         }
2489        }
2490       }
2491      }
2492     }
2493    }
2494   }
2495  }
2496
2497  */
2498
2499
2500
2501  // binning details of fDirectCorrelations (differential flow):
2502  //
2503  //101st bin: <2'>_{n|n}
2504
2505
2506
2507
2508  //<2'>_{n|n}
2509  for(Int_t i1=0;i1<nPrim;i1++)
2510  {
2511   fTrack=anEvent->GetTrack(i1);
2512   if(!((fTrack->Pt()>=0.5&&fTrack->Pt()<0.6)&&(fTrack->InPOISelection())))continue;//POI condition
2513   phi1=fTrack->Phi();
2514   if(phiWeights) wPhi1 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi1*nBinsPhi/TMath::TwoPi())));
2515   for(Int_t i2=0;i2<nPrim;i2++)
2516   {
2517    if(i2==i1)continue;
2518    fTrack=anEvent->GetTrack(i2);
2519    if(!(fTrack->InRPSelection()))continue;//RP condition
2520    phi2=fTrack->Phi();
2521    if(phiWeights) wPhi2 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi2*nBinsPhi/TMath::TwoPi())));
2522     //cout<<"1st = "<<i1<<"     "<< (anEvent->GetTrack(i1))->Eta() << " " << (anEvent->GetTrack(i1))->Pt()<<endl;
2523     //cout<<"2nd = "<<i2<<"     "<< (anEvent->GetTrack(i2))->Eta() << " " << (anEvent->GetTrack(i2))->Pt()<<endl;
2524    //fill the fDirectCorrelations:
2525    fDirectCorrelations->Fill(100.,cos(1.*n*(phi1-phi2)),wPhi2);//<2'>_{n,n}
2526    fDirectCorrelations->Fill(103.,cos(1.*n*(phi1-phi2)),pow(wPhi1,2)*wPhi2);//<2'>_{n,n}
2527    fDirectCorrelations->Fill(104.,cos(2.*n*(phi1-phi2)),wPhi1*pow(wPhi2,2));//<2'>_{n,n}
2528    fDirectCorrelations->Fill(105.,cos(1.*n*(phi1-phi2)),pow(wPhi2,3));//<2'>_{n,n}
2529
2530
2531
2532    fDirectCorrelations->Fill(41.,cos(2.*n*(phi1-phi2)),1);//<2'>_{2n,2n}
2533    fDirectCorrelations->Fill(42.,cos(3.*n*(phi1-phi2)),1);//<2'>_{3n,3n}
2534    fDirectCorrelations->Fill(43.,cos(4.*n*(phi1-phi2)),1);//<2'>_{4n,4n}
2535
2536   }//end of for(Int_t i2=0;i2<nPrim;i2++)
2537  }//end of for(Int_t i1=0;i1<nPrim;i1++)
2538
2539
2540
2541
2542
2543  /*
2544
2545  //<3'>_{2n|n,n}
2546  for(Int_t i1=0;i1<nPrim;i1++)
2547  {
2548   fTrack=anEvent->GetTrack(i1);
2549   if(!((fTrack->Pt()>=0.5&&fTrack->Pt()<0.6)&&(fTrack->InPOISelection())))continue;//POI condition
2550   phi1=fTrack->Phi();
2551   if(phiWeights) wPhi1 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi1*nBinsPhi/TMath::TwoPi())));
2552   for(Int_t i2=0;i2<nPrim;i2++)
2553   {
2554    if(i2==i1)continue;
2555    fTrack=anEvent->GetTrack(i2);
2556    if(!(fTrack->InRPSelection()))continue;//RP condition
2557    phi2=fTrack->Phi();
2558    if(phiWeights) wPhi2 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi2*nBinsPhi/TMath::TwoPi())));
2559    for(Int_t i3=0;i3<nPrim;i3++)
2560    {
2561     if(i3==i1||i3==i2)continue;
2562     fTrack=anEvent->GetTrack(i3);
2563     if(!(fTrack->InRPSelection()))continue;//RP condition
2564     phi3=fTrack->Phi();
2565     if(phiWeights) wPhi3 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi3*nBinsPhi/TMath::TwoPi())));
2566     //fill the fDirectCorrelations:
2567
2568     // 2-p
2569     fDirectCorrelations->Fill(101.,cos(n*(phi2-phi3)),wPhi1*wPhi2*wPhi3); // <w1 w2 w3 cos(n(phi2-phi3))>
2570     fDirectCorrelations->Fill(102.,cos(n*(phi1-phi3)),pow(wPhi2,2.)*wPhi3); // <w2^2 w3 cos(n(psi1-phi2))>
2571
2572     // 3-p
2573     fDirectCorrelations->Fill(110.,cos(n*(2.*phi1-phi2-phi3)),wPhi1*wPhi2*wPhi3); // <w1 w2 w3 cos(n(2psi1-phi2-phi3))>
2574     fDirectCorrelations->Fill(111.,cos(n*(phi1+phi2-2.*phi3)),wPhi2*pow(wPhi3,2.)); // <w2 w3^2 cos(n(psi1+phi2-2.*phi3))>
2575
2576
2577     //fDirectCorrelations->Fill(46.,cos(n*(phi1+phi2-2.*phi3)),1);//<3'>_{n,n|2n}
2578    }//end of for(Int_t i3=0;i3<nPrim;i3++)
2579   }//end of for(Int_t i2=0;i2<nPrim;i2++)
2580  }//end of for(Int_t i1=0;i1<nPrim;i1++)
2581
2582
2583  */
2584
2585
2586
2587
2588  //<4'>_{n,n|n,n}
2589  for(Int_t i1=0;i1<nPrim;i1++)
2590  {
2591   fTrack=anEvent->GetTrack(i1);
2592   if(!((fTrack->Pt()>=0.5&&fTrack->Pt()<0.6)&&(fTrack->InPOISelection())))continue;//POI condition
2593   tempCounter++;
2594   phi1=fTrack->Phi();
2595   if(phiWeights) wPhi1 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi1*nBinsPhi/TMath::TwoPi())));
2596   for(Int_t i2=0;i2<nPrim;i2++)
2597   {
2598    if(i2==i1)continue;
2599    fTrack=anEvent->GetTrack(i2);
2600    if(!(fTrack->InRPSelection()))continue;//RP condition
2601    phi2=fTrack->Phi();
2602    if(phiWeights) wPhi2 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi2*nBinsPhi/TMath::TwoPi())));
2603    for(Int_t i3=0;i3<nPrim;i3++)
2604    {
2605     if(i3==i1||i3==i2)continue;
2606     fTrack=anEvent->GetTrack(i3);
2607     if(!(fTrack->InRPSelection()))continue;//RP condition
2608     phi3=fTrack->Phi();
2609     if(phiWeights) wPhi3 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi3*nBinsPhi/TMath::TwoPi())));
2610     for(Int_t i4=0;i4<nPrim;i4++)
2611     {
2612      if(i4==i1||i4==i2||i4==i3)continue;
2613      fTrack=anEvent->GetTrack(i4);
2614      if(!(fTrack->InRPSelection()))continue;//RP condition
2615      phi4=fTrack->Phi();
2616      if(phiWeights) wPhi4 = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(phi4*nBinsPhi/TMath::TwoPi())));
2617      //fill the fDirectCorrelations:
2618      // 4-p
2619      fDirectCorrelations->Fill(120.,cos(n*(phi1+phi2-phi3-phi4)),wPhi2*wPhi3*wPhi4); // <w1 w2 w3 cos(n(psi1+phi1-phi2-phi3))>
2620
2621        tempLoop+=wPhi2*wPhi3*wPhi4;
2622
2623     }//end of for(Int_t i4=0;i4<nPrim;i4++)
2624    }//end of for(Int_t i3=0;i3<nPrim;i3++)
2625   }//end of for(Int_t i2=0;i2<nPrim;i2++)
2626  }//end of for(Int_t i1=0;i1<nPrim;i1++)
2627
2628
2629
2630
2631
2632
2633
2634
2635
2636
2637
2638
2639
2640
2641
2642
2643
2644
2645  /*
2646
2647
2648
2649
2650
2651
2652
2653  //<5'>_{2n,n|n,n,n}
2654  for(Int_t i1=0;i1<nPrim;i1++)
2655  {
2656   fTrack=anEvent->GetTrack(i1);
2657   if(!((fTrack->Pt()>=0.5&&fTrack->Pt()<0.6)&&(fTrack->InPOISelection())))continue;//POI condition
2658   phi1=fTrack->Phi();
2659   for(Int_t i2=0;i2<nPrim;i2++)
2660   {
2661    if(i2==i1)continue;
2662    fTrack=anEvent->GetTrack(i2);
2663    if(!(fTrack->InRPSelection()))continue;//RP condition
2664    phi2=fTrack->Phi();
2665    for(Int_t i3=0;i3<nPrim;i3++)
2666    {
2667     if(i3==i1||i3==i2)continue;
2668     fTrack=anEvent->GetTrack(i3);
2669     if(!(fTrack->InRPSelection()))continue;//RP condition
2670     phi3=fTrack->Phi();
2671     for(Int_t i4=0;i4<nPrim;i4++)
2672     {
2673      if(i4==i1||i4==i2||i4==i3)continue;
2674      fTrack=anEvent->GetTrack(i4);
2675      if(!(fTrack->InRPSelection()))continue;//RP condition
2676      phi4=fTrack->Phi();//
2677      for(Int_t i5=0;i5<nPrim;i5++)
2678      {
2679       if(i5==i1||i5==i2||i5==i3||i5==i4)continue;
2680       fTrack=anEvent->GetTrack(i5);
2681       if(!(fTrack->InRPSelection()))continue;//RP condition
2682       phi5=fTrack->Phi();
2683       //fill the fDirectCorrelations:if(bNestedLoops)
2684       fDirectCorrelations->Fill(55.,cos(2.*n*phi1+n*phi2-n*phi3-n*phi4-n*phi5),1);//<5'>_{2n,n|n,n,n}
2685      }//end of for(Int_t i5=0;i5<nPrim;i5++)
2686     }//end of for(Int_t i4=0;i4<nPrim;i4++)
2687    }//end of for(Int_t i3=0;i3<nPrim;i3++)
2688   }//end of for(Int_t i2=0;i2<nPrim;i2++)
2689  }//end of for(Int_t i1=0;i1<nPrim;i1++)
2690
2691  //<6'>_{n,n,n|n,n,n}
2692  for(Int_t i1=0;i1<nPrim;i1++)
2693  {
2694   fTrack=anEvent->GetTrack(i1);
2695   if(!((fTrack->Pt()>=0.5&&fTrack->Pt()<0.6)&&(fTrack->InPOISelection())))continue;//POI condition
2696   phi1=fTrack->Phi();
2697   for(Int_t i2=0;i2<nPrim;i2++)
2698   {
2699    if(i2==i1)continue;
2700    fTrack=anEvent->GetTrack(i2);
2701    if(!(fTrack->InRPSelection()))continue;//RP condition
2702    phi2=fTrack->Phi();
2703    for(Int_t i3=0;i3<nPrim;i3++)
2704    {
2705     if(i3==i1||i3==i2)continue;
2706     fTrack=anEvent->GetTrack(i3);
2707     if(!(fTrack->InRPSelection()))continue;//RP condition
2708     phi3=fTrack->Phi();
2709     for(Int_t i4=0;i4<nPrim;i4++)
2710     {
2711      if(i4==i1||i4==i2||i4==i3)continue;
2712      fTrack=anEvent->GetTrack(i4);
2713      if(!(fTrack->InRPSelection()))continue;//RP condition
2714      phi4=fTrack->Phi();
2715      for(Int_t i5=0;i5<nPrim;i5++)
2716      {
2717       if(i5==i1||i5==i2||i5==i3||i5==i4)continue;
2718       fTrack=anEvent->GetTrack(i5);
2719       if(!(fTrack->InRPSelection()))continue;//RP condition
2720       phi5=fTrack->Phi();
2721       for(Int_t i6=0;i6<nPrim;i6++)
2722       {
2723        if(i6==i1||i6==i2||i6==i3||i6==i4||i6==i5)continue;
2724        fTrack=anEvent->GetTrack(i6);
2725        if(!(fTrack->InRPSelection()))continue;//RP condition
2726        phi6=fTrack->Phi();
2727        //fill the fDirectCorrelations:
2728        fDirectCorrelations->Fill(60.,cos(n*(phi1+phi2+phi3-phi4-phi5-phi6)),1);//<6'>_{n,n,n|n,n,n}
2729       }//end of for(Int_t i6=0;i6<nPrim;i6++)
2730      }//end of for(Int_t i5=0;i5<nPrim;i5++)
2731     }//end of for(Int_t i4=0;i4<nPrim;i4++)
2732    }//end of for(Int_t i3=0;i3<nPrim;i3++)
2733   }//end of for(Int_t i2=0;i2<nPrim;i2++)
2734  }//end of for(Int_t i1=0;i1<nPrim;i1++)
2735
2736  //<7'>_{2n,n,n|n,n,n,n}
2737  for(Int_t i1=0;i1<nPrim;i1++)
2738  {
2739   fTrack=anEvent->GetTrack(i1);
2740   if(!((fTrack->Pt()>=0.5&&fTrack->Pt()<0.6)&&(fTrack->InPOISelection())))continue;//POI condition
2741   phi1=fTrack->Phi();
2742   for(Int_t i2=0;i2<nPrim;i2++)
2743   {
2744    if(i2==i1)continue;
2745    fTrack=anEvent->GetTrack(i2);
2746    if(!(fTrack->InRPSelection()))continue;//RP condition
2747    phi2=fTrack->Phi();
2748    for(Int_t i3=0;i3<nPrim;i3++)
2749    {
2750     if(i3==i1||i3==i2)continue;
2751     fTrack=anEvent->GetTrack(i3);
2752     if(!(fTrack->InRPSelection()))continue;//RP condition
2753     phi3=fTrack->Phi();
2754     for(Int_t i4=0;i4<nPrim;i4++)
2755     {
2756      if(i4==i1||i4==i2||i4==i3)continue;
2757      fTrack=anEvent->GetTrack(i4);
2758      if(!(fTrack->InRPSelection()))continue;//RP condition
2759      phi4=fTrack->Phi();
2760      for(Int_t i5=0;i5<nPrim;i5++)
2761      {
2762       if(i5==i1||i5==i2||i5==i3||i5==i4)continue;
2763       fTrack=anEvent->GetTrack(i5);
2764       if(!(fTrack->InRPSelection()))continue;//RP condition
2765       phi5=fTrack->Phi();
2766       for(Int_t i6=0;i6<nPrim;i6++)
2767       {
2768        if(i6==i1||i6==i2||i6==i3||i6==i4||i6==i5)continue;
2769        fTrack=anEvent->GetTrack(i6);
2770        if(!(fTrack->InRPSelection()))continue;//RP condition
2771        phi6=fTrack->Phi();
2772        for(Int_t i7=0;i7<nPrim;i7++)
2773        {
2774         if(i7==i1||i7==i2||i7==i3||i7==i4||i7==i5||i7==i6)continue;
2775         fTrack=anEvent->GetTrack(i7);
2776         if(!(fTrack->InRPSelection()))continue;//RP condition
2777         phi7=fTrack->Phi();
2778         //fill the fDirectCorrelations:
2779         fDirectCorrelations->Fill(65.,cos(2.*n*phi1+n*phi2+n*phi3-n*phi4-n*phi5-n*phi6-n*phi7),1);//<7'>_{2n,n,n|n,n,n,n}
2780        }//end of for(Int_t i7=0;i7<nPrim;i7++)
2781       }//end of for(Int_t i6=0;i6<nPrim;i6++)
2782      }//end of for(Int_t i5=0;i5<nPrim;i5++)
2783     }//end of for(Int_t i4=0;i4<nPrim;i4++)
2784    }//end of for(Int_t i3=0;i3<nPrim;i3++)
2785   }//end of for(Int_t i2=0;i2<nPrim;i2++)
2786  }//end of for(Int_t i1=0;i1<nPrim;i1++)
2787
2788  //<8'>_{n,n,n,n|n,n,n,n}
2789  for(Int_t i1=0;i1<nPrim;i1++)
2790  {
2791   fTrack=anEvent->GetTrack(i1);
2792   if(!((fTrack->Pt()>=0.5&&fTrack->Pt()<0.6)&&(fTrack->InPOISelection())))continue;//POI condition
2793   phi1=fTrack->Phi();
2794   for(Int_t i2=0;i2<nPrim;i2++)
2795   {
2796    if(i2==i1)continue;
2797    fTrack=anEvent->GetTrack(i2);
2798    if(!(fTrack->InRPSelection()))continue;//RP condition
2799    phi2=fTrack->Phi();
2800    for(Int_t i3=0;i3<nPrim;i3++)
2801    {
2802     if(i3==i1||i3==i2)continue;
2803     fTrack=anEvent->GetTrack(i3);
2804     if(!(fTrack->InRPSelection()))continue;//RP condition
2805     phi3=fTrack->Phi();
2806     for(Int_t i4=0;i4<nPrim;i4++)
2807     {
2808      if(i4==i1||i4==i2||i4==i3)continue;
2809      fTrack=anEvent->GetTrack(i4);
2810      if(!(fTrack->InRPSelection()))continue;//RP condition
2811      phi4=fTrack->Phi();
2812      for(Int_t i5=0;i5<nPrim;i5++)
2813      {
2814       if(i5==i1||i5==i2||i5==i3||i5==i4)continue;
2815       fTrack=anEvent->GetTrack(i5);
2816       if(!(fTrack->InRPSelection()))continue;//RP condition
2817       phi5=fTrack->Phi();
2818       for(Int_t i6=0;i6<nPrim;i6++)
2819       {
2820        if(i6==i1||i6==i2||i6==i3||i6==i4||i6==i5)continue;
2821        fTrack=anEvent->GetTrack(i6);
2822        if(!(fTrack->InRPSelection()))continue;//RP condition
2823        phi6=fTrack->Phi();
2824        for(Int_t i7=0;i7<nPrim;i7++)
2825        {
2826         if(i7==i1||i7==i2||i7==i3||i7==i4||i7==i5||i7==i6)continue;
2827         fTrack=anEvent->GetTrack(i7);
2828         if(!(fTrack->InRPSelection()))continue;//RP condition
2829         phi7=fTrack->Phi();
2830         for(Int_t i8=0;i8<nPrim;i8++)
2831         {
2832          if(i8==i1||i8==i2||i8==i3||i8==i4||i8==i5||i8==i6||i8==i7)continue;
2833          fTrack=anEvent->GetTrack(i8);
2834          if(!(fTrack->InRPSelection()))continue;//RP condition
2835          phi8=fTrack->Phi();
2836          //fill the fDirectCorrelations:
2837          fDirectCorrelations->Fill(70.,cos(n*(phi1+phi2+phi3+phi4-phi5-phi6-phi7-phi8)),1);//<8'>_{n,n,n,n|n,n,n,n}
2838         }//end of for(Int_t i8=0;i8<nPrim;i8++)
2839        }//end of for(Int_t i7=0;i7<nPrim;i7++)
2840       }//end of for(Int_t i6=0;i6<nPrim;i6++)
2841      }//end of for(Int_t i5=0;i5<nPrim;i5++)
2842     }//end of for(Int_t i4=0;i4<nPrim;i4++)
2843    }//end of for(Int_t i3=0;i3<nPrim;i3++)
2844   }//end of for(Int_t i2=0;i2<nPrim;i2++)
2845  }//end of for(Int_t i1=0;i1<nPrim;i1++)
2846
2847
2848
2849
2850
2851
2852
2853
2854
2855
2856  */
2857
2858
2859
2860  }// end of if (nestedLoops)
2861
2862
2863
2864
2865
2866
2867
2868
2869
2870
2871
2872
2873
2874
2875
2876
2877
2878
2879
2880
2881
2882
2883
2884
2885
2886
2887
2888
2889 //--------------------------------------------------------------------------------------------------------------------------------
2890
2891  //}//end of if(nPrim>0&&nPrim<12)
2892 }//end of Make()
2893
2894 //================================================================================================================
2895
2896 void AliFlowAnalysisWithQCumulants::Finish()
2897 {
2898  //calculate the final results
2899
2900  // harmonics
2901  Int_t n = 2; // to be improved
2902
2903  //---------------------------------------------------------------------------------------------------------
2904  // avarage multiplicity
2905  Double_t AvMPOI = (fCommonHists2nd->GetHistMultPOI())->GetMean(); // to be improved
2906  Double_t AvMRP  = (fCommonHists2nd->GetHistMultRP())->GetMean(); // to be improved
2907
2908  // number of events
2909  Double_t nEvtsPOI = (fCommonHists2nd->GetHistMultPOI())->GetEntries(); // to be improved
2910  Double_t nEvtsRP  = (fCommonHists2nd->GetHistMultRP())->GetEntries(); // to be improved
2911  //---------------------------------------------------------------------------------------------------------
2912
2913  //---------------------------------------------------------------------------------------------------------
2914  // 2-, 4-, 6- and 8-particle azimuthal correlation:
2915  Double_t two   = fQCorrelations->GetBinContent(1);  //<<2>>_{n|n}
2916  Double_t four  = fQCorrelations->GetBinContent(11); //<<4>>_{n,n|n,n}
2917  Double_t six   = fQCorrelations->GetBinContent(24); //<<6>>_{n,n,n|n,n,n}
2918  Double_t eight = fQCorrelations->GetBinContent(31); //<<8>>_{n,n,n,n|n,n,n,n}
2919
2920  // 2nd, 4th, 6th and 8th order Q-cumulant:
2921  Double_t secondOrderQCumulant = two; //c_n{2}
2922  Double_t fourthOrderQCumulant = four-2.*pow(two,2.); //c_n{4}
2923  Double_t sixthOrderQCumulant  = six-9.*two*four+12.*pow(two,3.); //c_n{6}
2924  Double_t eightOrderQCumulant  = eight-16.*two*six-18.*pow(four,2.)+144.*pow(two,2.)*four-144.*pow(two,4.); //c_n{8}
2925
2926  // "no-name" integrated flow estimates from Q-cumulants:
2927  Double_t vn2=0.,vn4=0.,vn6=0.,vn8=0.;
2928  // Double_t sd2=0.,sd4=0.,sd6=0.,sd8=0.;
2929  Double_t sd6=0.,sd8=0.; // to be improved/removed
2930  if(secondOrderQCumulant>0.)
2931  {
2932   vn2 = pow(secondOrderQCumulant,0.5); //v_n{2}
2933  }
2934  if(fourthOrderQCumulant<0.)
2935  {
2936   vn4 = pow(-fourthOrderQCumulant,1./4.); //v_n{4}
2937  }
2938  if(sixthOrderQCumulant>0.)
2939  {
2940   vn6 = pow((1./4.)*sixthOrderQCumulant,1./6.); //v_n{6}
2941  }
2942  if(eightOrderQCumulant<0.)
2943  {
2944   vn8 = pow((-1./33.)*eightOrderQCumulant,1./8.); //v_n{8}
2945  }
2946  //---------------------------------------------------------------------------------------------------------
2947
2948  //---------------------------------------------------------------------------------------------------------
2949  // weighted 2-, 4-, 6- and 8-particle azimuthal correlation:
2950  Double_t twoW   = fWeightedQCorrelations->GetBinContent(1);  //<<2>>_{n|n}
2951  Double_t fourW  = fWeightedQCorrelations->GetBinContent(21); //<<4>>_{n,n|n,n}
2952  // Double_t sixW   = fWeightedQCorrelations->GetBinContent(24); //<<6>>_{n,n,n|n,n,n}
2953  // Double_t eightW = fWeightedQCorrelations->GetBinContent(31); //<<8>>_{n,n,n,n|n,n,n,n}
2954
2955  // 2nd, 4th, 6th and 8th order weighted Q-cumulant:
2956  Double_t secondOrderQCumulantW = twoW; //c_n{2}
2957  Double_t fourthOrderQCumulantW = fourW-2.*pow(twoW,2.); //c_n{4}
2958  // Double_t sixthOrderQCumulantW  = sixW-9.*twoW*fourW+12.*pow(twoW,3.); //c_n{6}
2959  // Double_t eightOrderQCumulantW  = eightW-16.*twoW*sixW-18.*pow(fourW,2.)+144.*pow(twoW,2.)*fourW-144.*pow(twoW,4.); //c_n{8}
2960
2961  // "no-name" integrated flow estimates from weighted Q-cumulants:
2962  cout<<endl;
2963  cout<<"**************************************"<<endl;
2964  cout<<"**************************************"<<endl;
2965  cout<<"flow estimates from Q-cumulants :"<<endl;
2966  cout<<endl;
2967
2968  Double_t vn2W=0.,vn4W=0.;
2969  Double_t sd2W=0.,sd4W=0.;
2970  if(secondOrderQCumulantW>0.){
2971   vn2W = pow(secondOrderQCumulantW,0.5); // weighted v_n{2}
2972   //sd2W = 0.5*pow(secondOrderQCumulantW,-0.5)*secondOrderQCumulantErrorW; // to be improved (correct treatment of errors needed)
2973   cout<<" v_"<<n<<"{2} = "<<vn2W<<" +/- "<<sd2W<<endl;
2974   fIntFlowResultsQC->SetBinContent(1,vn2W);
2975   //fIntFlowResultsQC->SetBinError(1,sd2W);
2976   //common histograms:
2977   fCommonHistsResults2nd->FillIntegratedFlow(vn2W,sd2W);
2978   //fCommonHistsResults2nd->FillChi(vn2W*pow(AvM,0.5)); // to be removed
2979  }else{
2980   cout<<" v_"<<n<<"{2} = Im"<<endl;
2981  }
2982  if(fourW!=0. && fourthOrderQCumulantW<0.){
2983   vn4W = pow(-fourthOrderQCumulantW,1./4.); //v_n{4}
2984   //sd4W = 0.25*pow(-fourthOrderQCumulantW,-3./4.)*fourthOrderQCumulantErrorW; // to be improved (correct treatment of errors needed)
2985   cout<<" v_"<<n<<"{4} = "<<vn4W<<" +/- "<<sd4W<<endl;
2986   fIntFlowResultsQC->SetBinContent(2,vn4W);
2987   //fIntFlowResultsQC->SetBinError(2,sd4W);
2988   //common histograms:
2989   fCommonHistsResults4th->FillIntegratedFlow(vn4W,sd4W);
2990   //fCommonHistsResults4th->FillChi(vn4W*pow(AvM,0.5)); // to be removed
2991  }else{
2992   cout<<" v_"<<n<<"{4} = Im"<<endl;
2993  }
2994  // !!! to be improved (6th and 8th order are without weights) !!!
2995  if(six!=0. && sixthOrderQCumulant>0.){
2996   vn6 = pow((1./4.)*sixthOrderQCumulant,1./6.); //v_n{6}
2997   //sd6 = (1./6.)*pow(2.,-1./3.)*pow(sixthOrderQCumulant,-5./6.)*sixthOrderQCumulantError;
2998   cout<<" v_"<<n<<"{6} = "<<vn6<<" +/- "<<sd6<<endl;
2999   fIntFlowResultsQC->SetBinContent(3,vn6);
3000   //fIntFlowResultsQC->SetBinError(3,sd6);
3001   //common histograms:
3002   fCommonHistsResults6th->FillIntegratedFlow(vn6,sd6);
3003   //fCommonHistsResults6th->FillChi(vn6*pow(AvM,0.5));//to be removed
3004  }else{
3005   cout<<" v_"<<n<<"{6} = Im"<<endl;
3006  }
3007  if(eight!=0. && eightOrderQCumulant<0.){
3008   vn8 = pow((-1./33.)*eightOrderQCumulant,1./8.); //v_n{8}
3009   cout<<" v_"<<n<<"{8} = "<<vn8<<" +/- "<<sd8<<endl;
3010   fIntFlowResultsQC->SetBinContent(4,vn8);
3011   //fIntFlowResultsQC->SetBinError(4,sd8);
3012   //common histograms:
3013   fCommonHistsResults8th->FillIntegratedFlow(vn8,sd8);
3014   //fCommonHistsResults8th->FillChi(vn8*pow(AvM,0.5));//to be removed
3015  }else{
3016   cout<<" v_"<<n<<"{8} = Im"<<endl;
3017  }
3018  cout<<endl;
3019  cout<<"   nEvts = "<<nEvtsRP<<", AvM = "<<AvMRP<<endl; // to be improved
3020  cout<<"**************************************"<<endl;
3021  cout<<"**************************************"<<endl;
3022  cout<<endl;
3023 //---------------------------------------------------------------------------------------------------------
3024
3025 //---------------------------------------------------------------------------------------------------------
3026 // differential flow (POI)
3027 Int_t nBinsPtPOI = f2WPerPtBin1n1nPOI->GetNbinsX();
3028 Int_t nBinsEtaPOI = f4WPerPtBin1n1n1n1nPOI->GetNbinsX();
3029
3030 // Pt:
3031 Double_t secondOrderQCumulantDiffFlowPtPOIW = 0.;
3032 Double_t fourthOrderQCumulantDiffFlowPtPOIW = 0.;
3033
3034 Double_t dVn2ndPOIW=0.,dSd2ndPOIW=0.,dDiffvn2ndPOIW=0.,dYield2ndPOIW=0.,dSum2ndPOIW=0.;
3035 Double_t dVn4thPOIW=0.,dSd4thPOIW=0.,dDiffvn4thPOIW=0.,dYield4thPOIW=0.,dSum4thPOIW=0.;
3036
3037 for(Int_t bb=1;bb<nBinsPtPOI+1;bb++)
3038 {
3039  // QC{2}
3040  if(f2WPerPtBin1n1nPOI->GetBinEntries(bb)>0.&&vn2W!=0)
3041  {
3042   secondOrderQCumulantDiffFlowPtPOIW = f2WPerPtBin1n1nPOI->GetBinContent(bb); // with weights
3043   fDiffFlowResults2ndOrderQC->SetBinContent(bb,secondOrderQCumulantDiffFlowPtPOIW/vn2W);
3044   // common histogram:
3045   fCommonHistsResults2nd->FillDifferentialFlowPtPOI(bb,secondOrderQCumulantDiffFlowPtPOIW/vn2W, 0.); //to be improved (errors && bb or bb+1 ?)
3046   // -------------------------------------------------------------------
3047   // integrated flow (weighted, POI, Pt, 2nd order):
3048   dDiffvn2ndPOIW=(fCommonHistsResults2nd->GetHistDiffFlowPtPOI())->GetBinContent(bb);
3049   dYield2ndPOIW=(fCommonHists2nd->GetHistPtPOI())->GetBinContent(bb);
3050   dVn2ndPOIW+=dDiffvn2ndPOIW*dYield2ndPOIW;
3051   dSum2ndPOIW+=dYield2ndPOIW;
3052   // -------------------------------------------------------------------
3053  }
3054  // QC{4]
3055  if(f4WPerPtBin1n1n1n1nPOI->GetBinEntries(bb)>0.&&vn4W!=0.)
3056  {
3057   fourthOrderQCumulantDiffFlowPtPOIW = f4WPerPtBin1n1n1n1nPOI->GetBinContent(bb)-2.*f2WPerPtBin1n1nPOI->GetBinContent(bb)*pow(vn2W,2.); // with weights
3058   fDiffFlowResults4thOrderQC->SetBinContent(bb,-1.*fourthOrderQCumulantDiffFlowPtPOIW/pow(vn4W,3.));
3059   //common histogram:
3060   fCommonHistsResults4th->FillDifferentialFlowPtPOI(bb,-1.*fourthOrderQCumulantDiffFlowPtPOIW/pow(vn4W,3.), 0.); //to be improved (errors)
3061   // -------------------------------------------------------------------
3062   //integrated flow (POI, Pt, 4th order):
3063   dDiffvn4thPOIW=(fCommonHistsResults4th->GetHistDiffFlowPtPOI())->GetBinContent(bb);
3064   dYield4thPOIW=(fCommonHists4th->GetHistPtPOI())->GetBinContent(bb);
3065   dVn4thPOIW+=dDiffvn4thPOIW*dYield4thPOIW;
3066   dSum4thPOIW+=dYield4thPOIW;
3067   // -------------------------------------------------------------------
3068  }
3069 }
3070
3071 cout<<endl;
3072 cout<<"**************************************"<<endl;
3073 cout<<"**************************************"<<endl;
3074 cout<<"flow estimates from Q-cumulants (POI):"<<endl;
3075 cout<<endl;
3076 //storing the final results for integrated flow (POI):
3077 // QC{2}
3078 if(dSum2ndPOIW && fCommonHistsResults2nd)
3079 {
3080  dVn2ndPOIW/=dSum2ndPOIW;
3081  fCommonHistsResults2nd->FillIntegratedFlowPOI(dVn2ndPOIW,0.); // to be improved (errors)
3082  cout<<" v_"<<n<<"{2} = "<<dVn2ndPOIW<<" +/- "<<dSd2ndPOIW<<endl;
3083 }else
3084  {
3085   cout<<" v_"<<n<<"{2} = Im"<<endl;
3086  }
3087
3088 // QC{4}
3089 if(dSum4thPOIW && fCommonHistsResults4th)
3090 {
3091  dVn4thPOIW/=dSum4thPOIW;
3092  fCommonHistsResults4th->FillIntegratedFlowPOI(dVn4thPOIW,0.); // to be improved (errors)
3093  cout<<" v_"<<n<<"{4} = "<<dVn4thPOIW<<" +/- "<<dSd4thPOIW<<endl;
3094 }else
3095  {
3096   cout<<" v_"<<n<<"{4} = Im"<<endl;
3097  }
3098
3099 cout<<endl;
3100 cout<<"   nEvts = "<<nEvtsPOI<<", AvM = "<<AvMPOI<<endl;
3101 cout<<"**************************************"<<endl;
3102 cout<<"**************************************"<<endl;
3103 cout<<endl;
3104
3105 //Eta:
3106 Double_t secondOrderQCumulantDiffFlowEtaPOIW = 0.;
3107 Double_t fourthOrderQCumulantDiffFlowEtaPOIW = 0.;
3108
3109 for(Int_t bb=1;bb<nBinsEtaPOI+1;bb++)
3110 {
3111  if(f2WPerEtaBin1n1nPOI->GetBinEntries(bb)>0.&&vn2W!=0)
3112  {
3113   secondOrderQCumulantDiffFlowEtaPOIW = f2WPerEtaBin1n1nPOI->GetBinContent(bb); // with weights
3114   fDiffFlowResults2ndOrderQC->SetBinContent(bb,secondOrderQCumulantDiffFlowEtaPOIW/vn2W);
3115   //common histogram:
3116   fCommonHistsResults2nd->FillDifferentialFlowEtaPOI(bb,secondOrderQCumulantDiffFlowEtaPOIW/vn2W, 0.);//to be improved (errors)
3117  }
3118  if(f4WPerEtaBin1n1n1n1nPOI->GetBinEntries(bb)>0.&&vn4W!=0.)
3119  {
3120   fourthOrderQCumulantDiffFlowEtaPOIW = f4WPerEtaBin1n1n1n1nPOI->GetBinContent(bb)-2.*f2WPerEtaBin1n1nPOI->GetBinContent(bb)*pow(vn2W,2.); // with weights
3121   fDiffFlowResults4thOrderQC->SetBinContent(bb,-1.*fourthOrderQCumulantDiffFlowEtaPOIW/pow(vn4W,3.));
3122   //common histogram:
3123   fCommonHistsResults4th->FillDifferentialFlowEtaPOI(bb,-1.*fourthOrderQCumulantDiffFlowEtaPOIW/pow(vn4W,3.), 0.);//to be improved (errors)
3124  }
3125 }
3126 //------------------------------------------------------------
3127
3128
3129
3130
3131
3132  //---------------------------------------------------------------------------------------------------------
3133 // differential flow (RP)
3134 Int_t nBinsPtRP = f2WPerPtBin1n1nRP->GetNbinsX();
3135 Int_t nBinsEtaRP = f4WPerPtBin1n1n1n1nRP->GetNbinsX();
3136
3137 // Pt:
3138 Double_t secondOrderQCumulantDiffFlowPtRPW = 0.;
3139 Double_t fourthOrderQCumulantDiffFlowPtRPW = 0.;
3140
3141 Double_t dVn2ndRPW=0.,dSd2ndRPW=0.,dDiffvn2ndRPW=0.,dYield2ndRPW=0.,dSum2ndRPW=0.;
3142 Double_t dVn4thRPW=0.,dSd4thRPW=0.,dDiffvn4thRPW=0.,dYield4thRPW=0.,dSum4thRPW=0.;
3143
3144 for(Int_t bb=1;bb<nBinsPtRP+1;bb++)
3145 {
3146  // QC{2}
3147  if(f2WPerPtBin1n1nRP->GetBinEntries(bb)>0.&&vn2W!=0)
3148  {
3149   secondOrderQCumulantDiffFlowPtRPW = f2WPerPtBin1n1nRP->GetBinContent(bb); // with weights
3150   fDiffFlowResults2ndOrderQC->SetBinContent(bb,secondOrderQCumulantDiffFlowPtRPW/vn2W);
3151   // common histogram:
3152   fCommonHistsResults2nd->FillDifferentialFlowPtRP(bb,secondOrderQCumulantDiffFlowPtRPW/vn2W, 0.); //to be improved (errors && bb or bb+1 ?)
3153   // -------------------------------------------------------------------
3154   // integrated flow (weighted, RP, Pt, 2nd order):
3155   dDiffvn2ndRPW=(fCommonHistsResults2nd->GetHistDiffFlowPtRP())->GetBinContent(bb);
3156   dYield2ndRPW=(fCommonHists2nd->GetHistPtRP())->GetBinContent(bb);
3157   dVn2ndRPW+=dDiffvn2ndRPW*dYield2ndRPW;
3158   dSum2ndRPW+=dYield2ndRPW;
3159   // -------------------------------------------------------------------
3160  }
3161  // QC{4]
3162  if(f4WPerPtBin1n1n1n1nRP->GetBinEntries(bb)>0.&&vn4W!=0.)
3163  {
3164   fourthOrderQCumulantDiffFlowPtRPW = f4WPerPtBin1n1n1n1nRP->GetBinContent(bb)-2.*f2WPerPtBin1n1nRP->GetBinContent(bb)*pow(vn2W,2.); // with weights
3165   fDiffFlowResults4thOrderQC->SetBinContent(bb,-1.*fourthOrderQCumulantDiffFlowPtRPW/pow(vn4W,3.));
3166   //common histogram:
3167   fCommonHistsResults4th->FillDifferentialFlowPtRP(bb,-1.*fourthOrderQCumulantDiffFlowPtRPW/pow(vn4W,3.), 0.); //to be improved (errors)
3168   // -------------------------------------------------------------------
3169   //integrated flow (RP, Pt, 4th order):
3170   dDiffvn4thRPW=(fCommonHistsResults4th->GetHistDiffFlowPtRP())->GetBinContent(bb);
3171   dYield4thRPW=(fCommonHists4th->GetHistPtRP())->GetBinContent(bb);
3172   dVn4thRPW+=dDiffvn4thRPW*dYield4thRPW;
3173   dSum4thRPW+=dYield4thRPW;
3174   // -------------------------------------------------------------------
3175  }
3176 }
3177
3178 cout<<endl;
3179 cout<<"**************************************"<<endl;
3180 cout<<"**************************************"<<endl;
3181 cout<<"flow estimates from Q-cumulants (RP):"<<endl;
3182 cout<<endl;
3183 //storing the final results for integrated flow (RP):
3184 // QC{2}
3185 if(dSum2ndRPW && fCommonHistsResults2nd)
3186 {
3187  dVn2ndRPW/=dSum2ndRPW;
3188  fCommonHistsResults2nd->FillIntegratedFlowRP(dVn2ndRPW,0.); // to be improved (errors)
3189  cout<<" v_"<<n<<"{2} = "<<dVn2ndRPW<<" +/- "<<dSd2ndRPW<<endl;
3190 }else
3191  {
3192   cout<<" v_"<<n<<"{2} = Im"<<endl;
3193  }
3194
3195 // QC{4}
3196 if(dSum4thRPW && fCommonHistsResults4th)
3197 {
3198  dVn4thRPW/=dSum4thRPW;
3199  fCommonHistsResults4th->FillIntegratedFlowRP(dVn4thRPW,0.); // to be improved (errors)
3200  cout<<" v_"<<n<<"{4} = "<<dVn4thRPW<<" +/- "<<dSd4thRPW<<endl;
3201 }else
3202  {
3203   cout<<" v_"<<n<<"{4} = Im"<<endl;
3204  }
3205
3206 cout<<endl;
3207 cout<<"   nEvts = "<<nEvtsRP<<", AvM = "<<AvMRP<<endl;
3208 cout<<"**************************************"<<endl;
3209 cout<<"**************************************"<<endl;
3210 cout<<endl;
3211
3212 //Eta:
3213 Double_t secondOrderQCumulantDiffFlowEtaRPW = 0.;
3214 Double_t fourthOrderQCumulantDiffFlowEtaRPW = 0.;
3215
3216 for(Int_t bb=1;bb<nBinsEtaRP+1;bb++)
3217 {
3218  if(f2WPerEtaBin1n1nRP->GetBinEntries(bb)>0.&&vn2W!=0)
3219  {
3220   secondOrderQCumulantDiffFlowEtaRPW = f2WPerEtaBin1n1nRP->GetBinContent(bb); // with weights
3221   fDiffFlowResults2ndOrderQC->SetBinContent(bb,secondOrderQCumulantDiffFlowEtaRPW/vn2W);
3222   //common histogram:
3223   fCommonHistsResults2nd->FillDifferentialFlowEtaRP(bb,secondOrderQCumulantDiffFlowEtaRPW/vn2W, 0.);//to be improved (errors)
3224  }
3225  if(f4WPerEtaBin1n1n1n1nRP->GetBinEntries(bb)>0.&&vn4W!=0.)
3226  {
3227   fourthOrderQCumulantDiffFlowEtaRPW = f4WPerEtaBin1n1n1n1nRP->GetBinContent(bb)-2.*f2WPerEtaBin1n1nRP->GetBinContent(bb)*pow(vn2W,2.); // with weights
3228   fDiffFlowResults4thOrderQC->SetBinContent(bb,-1.*fourthOrderQCumulantDiffFlowEtaRPW/pow(vn4W,3.));
3229   //common histogram:
3230   fCommonHistsResults4th->FillDifferentialFlowEtaRP(bb,-1.*fourthOrderQCumulantDiffFlowEtaRPW/pow(vn4W,3.), 0.);//to be improved (errors)
3231  }
3232 }
3233 //------------------------------------------------------------
3234
3235
3236
3237
3238
3239
3240
3241
3242
3243
3244  Bool_t nestedLoops = kFALSE;
3245  if(nestedLoops)
3246  {
3247  //needed for direct correlations (obtained from nested loops)
3248  cout<<endl;
3249  cout<<endl;
3250  cout<<"   **** cross-checking the formulas ****"<<endl;
3251  cout<<"   ****     for integrated flow     ****"<<endl;
3252  cout<<"(selected only events for which 0 < M < 12 "<<endl;
3253  cout<<"  from dataset in /data/alice2/ante/AOD)   "<<endl;
3254
3255  cout<<endl;
3256  //cout<<"   nEvts = "<<nEvts<<", AvM = "<<AvM<<endl;
3257  cout<<endl;
3258  cout<<"<w1 w2 cos(n*(phi1-phi2))> from Q-vectors               = "<<fWeightedQCorrelations->GetBinContent(1)<<endl;
3259  cout<<"<w1 w2 cos(n*(phi1-phi2))> from nested loops            = "<<fDirectCorrelations->GetBinContent(1)<<endl;
3260  cout<<endl;
3261  cout<<"<w1^2 w2^2 cos(2n*(phi1-phi2))> from Q-vectors          = "<<fWeightedQCorrelations->GetBinContent(2)<<endl;
3262  cout<<"<w1^2 w2^2 cos(2n*(phi1-phi2))> from nested loops       = "<<fDirectCorrelations->GetBinContent(2)<<endl;
3263  cout<<endl;
3264  cout<<"<w1^3 w2^3 cos(3n*(phi1-phi2))> from Q-vectors          = "<<fWeightedQCorrelations->GetBinContent(3)<<endl;
3265  cout<<"<w1^3 w2^3 cos(3n*(phi1-phi2))> from nested loops       = "<<fDirectCorrelations->GetBinContent(3)<<endl;
3266  cout<<endl;
3267  cout<<"<w1^4 w2^4 cos(4n*(phi1-phi2))> from Q-vectors          = "<<fWeightedQCorrelations->GetBinContent(4)<<endl;
3268  cout<<"<w1^4 w2^4 cos(4n*(phi1-phi2))> from nested loops       = "<<fDirectCorrelations->GetBinContent(4)<<endl;
3269  cout<<endl;
3270  cout<<"<w1^3 w2 cos(n*(phi1-phi2))> from Q-vectors             = "<<fWeightedQCorrelations->GetBinContent(5)<<endl;
3271  cout<<"<w1^3 w2 cos(n*(phi1-phi2))> from nested loops          = "<<fDirectCorrelations->GetBinContent(5)<<endl;
3272  cout<<endl;
3273  cout<<"<w1 w2 w3^2 cos(n*(phi1-phi2))> from Q-vectors          = "<<fWeightedQCorrelations->GetBinContent(6)<<endl;
3274  cout<<"<w1 w2 w3^2 cos(n*(phi1-phi2))> from nested loops       = "<<fDirectCorrelations->GetBinContent(6)<<endl;
3275  cout<<endl;
3276  cout<<"<w1^2 w2 w3 cos(n*(2phi1-phi2-phi3))> from Q-vectors    = "<<fWeightedQCorrelations->GetBinContent(11)<<endl;
3277  cout<<"<w1^2 w2 w3 cos(n*(2phi1-phi2-phi3))> from nested loops = "<<fDirectCorrelations->GetBinContent(11)<<endl;
3278  cout<<endl;
3279
3280  cout<<"<w1 w2 w3 w4 cos(n*(phi1+phi2-phi3-phi4))> from Q-vectors    = "<<fWeightedQCorrelations->GetBinContent(21)<<endl;
3281  cout<<"<w1 w2 w3 w4 cos(n*(phi1+phi2-phi3-phi4))> from nested loops = "<<fDirectCorrelations->GetBinContent(21)<<endl;
3282  cout<<endl;;
3283
3284
3285  /*
3286  cout<<"<3>_{4n,2n,2n} from Q-vectors           = "<<fQCorrelations->GetBinContent(8)<<endl;
3287  cout<<"<3>_{4n,2n,2n} from nested loops        = "<<fDirectCorrelations->GetBinContent(8)<<endl;
3288  cout<<endl;
3289  cout<<"<3>_{4n,3n,n} from Q-vectors            = "<<fQCorrelations->GetBinContent(9)<<endl;
3290  cout<<"<3>_{4n,3n,n} from nested loops         = "<<fDirectCorrelations->GetBinContent(9)<<endl;
3291  cout<<endl;
3292
3293  cout<<"<4>_{2n,n|2n,n} from Q-vectors          = "<<fQCorrelations->GetBinContent(12)<<endl;
3294  cout<<"<4>_{2n,n|2n,n} from nested loops       = "<<fDirectCorrelations->GetBinContent(12)<<endl;
3295  cout<<endl;
3296  cout<<"<4>_{2n,2n|2n,2n} from Q-vectors        = "<<fQCorrelations->GetBinContent(13)<<endl;
3297  cout<<"<4>_{2n,2n|2n,2n} from nested loops     = "<<fDirectCorrelations->GetBinContent(13)<<endl;
3298  cout<<endl;
3299  cout<<"<4>_{3n|n,n,n} from Q-vectors           = "<<fQCorrelations->GetBinContent(14)<<endl;
3300  cout<<"<4>_{3n|n,n,n} from nested loops        = "<<fDirectCorrelations->GetBinContent(14)<<endl;
3301  cout<<endl;
3302  cout<<"<4>_{3n,n|3n,n} from Q-vectors          = "<<fQCorrelations->GetBinContent(15)<<endl;
3303  cout<<"<4>_{3n,n|3n,n} from nested loops       = "<<fDirectCorrelations->GetBinContent(15)<<endl;
3304  cout<<endl;
3305  cout<<"<4>_{3n,n|2n,2n} from Q-vectors         = "<<fQCorrelations->GetBinContent(16)<<endl;
3306  cout<<"<4>_{3n,n|2n,2n} from nested loops      = "<<fDirectCorrelations->GetBinContent(16)<<endl;
3307  cout<<endl;
3308  cout<<"<4>_{4n|2n,n,n} from Q-vectors          = "<<fQCorrelations->GetBinContent(17)<<endl;
3309  cout<<"<4>_{4n|2n,n,n} from nested loops       = "<<fDirectCorrelations->GetBinContent(17)<<endl;
3310  cout<<endl;
3311  cout<<"<5>_{2n,n|n,n,n} from Q-vectors         = "<<fQCorrelations->GetBinContent(19)<<endl;
3312  cout<<"<5>_{2n,n|n,n,n} from nested loops      = "<<fDirectCorrelations->GetBinContent(19)<<endl;
3313  cout<<endl;
3314  cout<<"<5>_{2n,2n|2n,n,n} from Q-vectors       = "<<fQCorrelations->GetBinContent(20)<<endl;
3315  cout<<"<5>_{2n,2n|2n,n,n} from nested loops    = "<<fDirectCorrelations->GetBinContent(20)<<endl;
3316  cout<<endl;
3317  cout<<"<5>_{3n,n|2n,n,n} from Q-vectors        = "<<fQCorrelations->GetBinContent(21)<<endl;
3318  cout<<"<5>_{3n,n|2n,n,n} from nested loops     = "<<fDirectCorrelations->GetBinContent(21)<<endl;
3319  cout<<endl;
3320  cout<<"<5>_{4n|n,n,n,n} from Q-vectors         = "<<fQCorrelations->GetBinContent(22)<<endl;
3321  cout<<"<5>_{4n|n,n,n,n} from nested loops      = "<<fDirectCorrelations->GetBinContent(22)<<endl;
3322  cout<<endl;
3323  cout<<"<6>_{n,n,n|n,n,n} from Q-vectors        = "<<fQCorrelations->GetBinContent(24)<<endl;
3324  cout<<"<6>_{n,n,n|n,n,n} from nested loops     = "<<fDirectCorrelations->GetBinContent(24)<<endl;
3325  cout<<endl;
3326  cout<<"<6>_{2n,n,n|2n,n,n} from Q-vectors      = "<<fQCorrelations->GetBinContent(25)<<endl;
3327  cout<<"<6>_{2n,n,n|2n,n,n} from nested loops   = "<<fDirectCorrelations->GetBinContent(25)<<endl;
3328  cout<<endl;
3329  cout<<"<6>_{2n,2n|n,n,n,n} from Q-vectors      = "<<fQCorrelations->GetBinContent(26)<<endl;
3330  cout<<"<6>_{2n,2n|n,n,n,n} from nested loops   = "<<fDirectCorrelations->GetBinContent(26)<<endl;
3331  cout<<endl;
3332  cout<<"<6>_{3n,n|n,n,n,n} from Q-vectors       = "<<fQCorrelations->GetBinContent(27)<<endl;
3333  cout<<"<6>_{3n,n|n,n,n,n} from nested loops    = "<<fDirectCorrelations->GetBinContent(27)<<endl;
3334  cout<<endl;
3335  cout<<"<7>_{2n,n,n|n,n,n,n} from Q-vectors     = "<<fQCorrelations->GetBinContent(29)<<endl;
3336  cout<<"<7>_{2n,n,n|n,n,n,n} from nested loops  = "<<fDirectCorrelations->GetBinContent(29)<<endl;
3337  cout<<endl;
3338  cout<<"<8>_{n,n,n,n|n,n,n,n} from Q-vectors    = "<<fQCorrelations->GetBinContent(31)<<endl;
3339  cout<<"<8>_{n,n,n,n|n,n,n,n} from nested loops = "<<fDirectCorrelations->GetBinContent(31)<<endl;
3340  cout<<endl;
3341
3342  */
3343
3344
3345
3346
3347
3348
3349
3350
3351
3352
3353
3354
3355
3356
3357
3358  cout<<"0.5 < Pt < 0.6 GeV"<<endl;
3359  cout<<endl;
3360  cout<<"<w2 cos(n(psi1-phi2))> from Q-vectors                    = "<<f2WPerPtBin1n1nPOI->GetBinContent(6)<<endl;
3361  //cout<<"<w2 cos(n(psi1-phi2))> from Q-vectors                    = "<<f2WPerPtBin1n1nRP->GetName()<<endl;
3362  cout<<"<w2 cos(n(psi1-phi2))> from Q-vectors                    = "<<fDirectCorrelations->GetBinContent(101)<<endl;
3363  cout<<endl;
3364  cout<<"<w2 w3 w4 cos(n(psi1+phi2-phi3-phi4))> from Q-vectors    = "<<f4WPerPtBin1n1n1n1nPOI->GetBinContent(6)<<endl;
3365  //cout<<"<w2 w3 w4 cos(n(psi1+phi2-phi3-phi4))> from Q-vectors    = "<<f4WPerPtBin1n1n1n1nRP->GetName()<<endl;
3366  cout<<"<w2 w3 w4 cos(n(psi1+phi2-phi3-phi4))> from nested loops = "<<fDirectCorrelations->GetBinContent(121)<<endl;
3367  cout<<endl;
3368
3369  }//end of if(nestedLoops)
3370
3371
3372
3373
3374
3375
3376
3377
3378   /*
3379  //xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
3380  //                   !!!! to be removed !!!!
3381
3382  TCanvas* qvectorPlot = new TCanvas("qvectorPlot","Q-vector Plot",1000,1000);
3383
3384  qvectorPlot->cd(1);
3385
3386  TH1D* style = new TH1D("style","Q-vectors",100,-244,244);
3387  (style->GetYaxis())->SetRangeUser(-244,244);
3388
3389  style->Draw();
3390
3391  Int_t nBins=fQvectorForEachEventX->GetNbinsX();
3392  Double_t qxxx=0.,qyyy=0.;
3393  //cout<<"nBins = "<<nBins<<endl;
3394  //cout<<fQvectorForEachEventX->GetBinEntries(4)<<endl;
3395  //cout<<fQvectorForEachEventY->GetBinEntries(4)<<endl;
3396
3397  for(Int_t b=1;b<nBins+1;b++)
3398  {
3399   if(fQvectorForEachEventX->GetBinEntries(b)==1 && fQvectorForEachEventY->GetBinEntries(b)==1)
3400   {
3401    qxxx=fQvectorForEachEventX->GetBinContent(b);
3402    qyyy=fQvectorForEachEventY->GetBinContent(b);
3403    //cout<<qxxx<<" "<<qyyy<<endl;
3404    TArrow *qvector = new TArrow(0.0,0.0,qxxx,qyyy,0.0144,"|>");
3405    qvector->SetAngle(40);
3406    qvector->SetLineWidth(2);
3407    qvector->Draw("");
3408   }
3409  }
3410
3411  //xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
3412  */
3413
3414 }
3415
3416 //================================================================================================================
3417
3418 void AliFlowAnalysisWithQCumulants::WriteHistograms(TString outputFileName)
3419 {
3420  //store the final results in output .root file
3421  TFile *output = new TFile(outputFileName.Data(),"RECREATE");
3422  output->WriteObject(fHistList, "cobjQC","SingleKey");
3423  delete output;
3424 }
3425
3426 //================================================================================================================
3427
3428
3429