+ // transfer 2D profile into 2D histogram:
+ // to be improved (see in documentation if there is a method to transfer values from 2D profile into 2D histogram)
+ for(Int_t ci=0;ci<4;ci++)
+ {
+ for(Int_t p=1;p<=fnBinsPt;p++)
+ {
+ for(Int_t e=1;e<=fnBinsEta;e++)
+ {
+ Double_t correlation = fCorrelationsPro[t][pW][eW][ci]->GetBinContent(fCorrelationsPro[t][pW][eW][ci]->GetBin(p,e));
+ Double_t spread = fCorrelationsPro[t][pW][eW][ci]->GetBinError(fCorrelationsPro[t][pW][eW][ci]->GetBin(p,e));
+ Double_t nEvts = fNonEmptyBins2D[t]->GetBinContent(fNonEmptyBins2D[t]->GetBin(p,e));
+ Double_t error = 0.;
+ fFinalCorrelations2D[t][pW][eW][ci]->SetBinContent(fFinalCorrelations2D[t][pW][eW][ci]->GetBin(p,e),correlation);
+ if(nEvts>0)
+ {
+ error = spread/pow(nEvts,0.5);
+ fFinalCorrelations2D[t][pW][eW][ci]->SetBinError(fFinalCorrelations2D[t][pW][eW][ci]->GetBin(p,e),error);
+ }
+ } // end of for(Int_t e=1;e<=fnBinsEta;e++)
+ } // end of for(Int_t p=1;p<=fnBinsPt;p++)
+ } // end of for(Int_t ci=0;ci<4;ci++)
+
+ // transfer 1D profile into 1D histogram (pt):
+ // to be improved (see in documentation if there is a method to transfer values from 1D profile into 1D histogram)
+ for(Int_t ci=0;ci<4;ci++)
+ {
+ for(Int_t p=1;p<=fnBinsPt;p++)
+ {
+ if(profile[0][ci])
+ {
+ Double_t correlation = profile[0][ci]->GetBinContent(p);
+ Double_t spread = profile[0][ci]->GetBinError(p);
+ Double_t nEvts = fNonEmptyBins1D[t][0]->GetBinContent(p);
+ Double_t error = 0.;
+ fFinalCorrelations1D[t][pW][eW][0][ci]->SetBinContent(p,correlation);
+ if(nEvts>0)
+ {
+ error = spread/pow(nEvts,0.5);
+ fFinalCorrelations1D[t][pW][eW][0][ci]->SetBinError(p,error);
+ }
+ }
+ } // end of for(Int_t p=1;p<=fnBinsPt;p++)
+ } // end of for(Int_t ci=0;ci<4;ci++)
+
+ // transfer 1D profile into 1D histogram (eta):
+ // to be improved (see in documentation if there is a method to transfer values from 1D profile into 1D histogram)
+ for(Int_t ci=0;ci<4;ci++)
+ {
+ for(Int_t e=1;e<=fnBinsEta;e++)
+ {
+ if(profile[1][ci])
+ {
+ Double_t correlation = profile[1][ci]->GetBinContent(e);
+ fFinalCorrelations1D[t][pW][eW][1][ci]->SetBinContent(e,correlation);
+ }
+ } // end of for(Int_t e=1;e<=fnBinsEta;e++)
+ } // end of for(Int_t ci=0;ci<4;ci++)
+
+} // end of void AliFlowAnalysisWithQCumulants::FinalizeCorrelationsForDiffFlow(TString type, Bool_t useParticleWeights, TString eventWeights)
+*/
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::CalculateDiffFlowCumulants(TString type, TString ptOrEta)
+{
+ // calcualate cumulants for differential flow from measured correlations
+ // Remark: cumulants calculated here are NOT corrected for non-uniform acceptance. This correction is applied in the method ...
+ // to be improved (description)
+
+ Int_t typeFlag = -1;
+ Int_t ptEtaFlag = -1;
+
+ if(type == "RP")
+ {
+ typeFlag = 0;
+ } else if(type == "POI")
+ {
+ typeFlag = 1;
+ }
+
+ if(ptOrEta == "Pt")
+ {
+ ptEtaFlag = 0;
+ } else if(ptOrEta == "Eta")
+ {
+ ptEtaFlag = 1;
+ }
+
+ // shortcuts:
+ Int_t t = typeFlag;
+ Int_t pe = ptEtaFlag;
+
+ // common:
+ Int_t nBinsPtEta[2] = {fnBinsPt,fnBinsEta};
+
+ // correlation <<2>>:
+ Double_t two = fIntFlowCorrelationsHist->GetBinContent(1);
+
+ // 1D:
+ for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+ {
+ // reduced correlations:
+ Double_t twoPrime = fDiffFlowCorrelationsHist[t][pe][0]->GetBinContent(b); // <<2'>>(pt)
+ Double_t fourPrime = fDiffFlowCorrelationsHist[t][pe][1]->GetBinContent(b); // <<4'>>(pt)
+ // final statistical error of reduced correlations:
+ //Double_t twoPrimeError = fFinalCorrelations1D[t][pW][eW][0][0]->GetBinError(p);
+ // QC{2'}:
+ Double_t qc2Prime = twoPrime; // QC{2'}
+ //Double_t qc2PrimeError = twoPrimeError; // final stat. error of QC{2'}
+ fDiffFlowCumulants[t][pe][0]->SetBinContent(b,qc2Prime);
+ //fFinalCumulantsPt[t][pW][eW][nua][0]->SetBinError(p,qc2PrimeError);
+ // QC{4'}:
+ Double_t qc4Prime = fourPrime - 2.*twoPrime*two; // QC{4'} = <<4'>> - 2*<<2'>><<2>>
+ fDiffFlowCumulants[t][pe][1]->SetBinContent(b,qc4Prime);
+ } // end of for(Int_t p=1;p<=fnBinsPt;p++)
+
+
+ /*
+ // 2D (pt,eta):
+ // to be improved (see documentation if I can do all this without looping)
+ for(Int_t p=1;p<=fnBinsPt;p++)
+ {
+ for(Int_t e=1;e<=fnBinsEta;e++)
+ {
+ // reduced correlations:
+ Double_t twoPrime = fFinalCorrelations2D[t][pW][eW][0]->GetBinContent(fFinalCorrelations2D[t][pW][eW][0]->GetBin(p,e)); // <<2'>>(pt,eta)
+ Double_t fourPrime = fFinalCorrelations2D[t][pW][eW][1]->GetBinContent(fFinalCorrelations2D[t][pW][eW][1]->GetBin(p,e)); // <<4'>>(pt,eta)
+ for(Int_t nua=0;nua<2;nua++)
+ {
+ // QC{2'}:
+ Double_t qc2Prime = twoPrime; // QC{2'} = <<2'>>
+ fFinalCumulants2D[t][pW][eW][nua][0]->SetBinContent(fFinalCumulants2D[t][pW][eW][nua][0]->GetBin(p,e),qc2Prime);
+ // QC{4'}:
+ Double_t qc4Prime = fourPrime - 2.*twoPrime*two; // QC{4'} = <<4'>> - 2*<<2'>><<2>>
+ fFinalCumulants2D[t][pW][eW][nua][1]->SetBinContent(fFinalCumulants2D[t][pW][eW][nua][1]->GetBin(p,e),qc4Prime);
+ } // end of for(Int_t nua=0;nua<2;nua++)
+ } // end of for(Int_t e=1;e<=fnBinsEta;e++)
+ } // end of for(Int_t p=1;p<=fnBinsPt;p++)
+ */
+
+} // end of void AliFlowAnalysisWithQCumulants::CalculateDiffFlowCumulants(TString type, Bool_t useParticleWeights, TString eventWeights);
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::CalculateFinalResultsForRPandPOIIntegratedFlow(TString type)
+{
+ // calculate final results for integrated flow of RPs and POIs
+
+ Int_t typeFlag = -1;
+
+ if(type == "RP")
+ {
+ typeFlag = 0;
+ } else if(type == "POI")
+ {
+ typeFlag = 1;
+ } else
+ {
+ cout<<"WARNING: type must be either RP or POI in AFAWQC::CDF() !!!!"<<endl;
+ exit(0);
+ }
+
+ // shortcuts:
+ Int_t t = typeFlag;
+
+ // pt yield:
+ TH1F *yield2ndPt = NULL;
+ TH1F *yield4thPt = NULL;
+ TH1F *yield6thPt = NULL;
+ TH1F *yield8thPt = NULL;
+
+ if(type == "POI")
+ {
+ yield2ndPt = (TH1F*)(fCommonHists2nd->GetHistPtPOI())->Clone();
+ yield4thPt = (TH1F*)(fCommonHists4th->GetHistPtPOI())->Clone();
+ yield6thPt = (TH1F*)(fCommonHists6th->GetHistPtPOI())->Clone();
+ yield8thPt = (TH1F*)(fCommonHists8th->GetHistPtPOI())->Clone();
+ }
+ else if(type == "RP")
+ {
+ yield2ndPt = (TH1F*)(fCommonHists2nd->GetHistPtRP())->Clone();
+ yield4thPt = (TH1F*)(fCommonHists4th->GetHistPtRP())->Clone();
+ yield6thPt = (TH1F*)(fCommonHists6th->GetHistPtRP())->Clone();
+ yield8thPt = (TH1F*)(fCommonHists8th->GetHistPtRP())->Clone();
+ }
+
+ Int_t nBinsPt = yield2ndPt->GetNbinsX();
+
+ TH1D *flow2ndPt = NULL;
+ TH1D *flow4thPt = NULL;
+ TH1D *flow6thPt = NULL;
+ TH1D *flow8thPt = NULL;
+
+ // to be improved (hardwired pt index)
+ flow2ndPt = (TH1D*)fDiffFlow[t][0][0]->Clone();
+ flow4thPt = (TH1D*)fDiffFlow[t][0][1]->Clone();
+ flow6thPt = (TH1D*)fDiffFlow[t][0][2]->Clone();
+ flow8thPt = (TH1D*)fDiffFlow[t][0][3]->Clone();
+
+ Double_t dvn2nd = 0., dvn4th = 0., dvn6th = 0., dvn8th = 0.; // differential flow
+ Double_t dErrvn2nd = 0., dErrvn4th = 0., dErrvn6th = 0., dErrvn8th = 0.; // error on differential flow
+
+ Double_t dVn2nd = 0., dVn4th = 0., dVn6th = 0., dVn8th = 0.; // integrated flow
+ Double_t dErrVn2nd = 0., dErrVn4th = 0., dErrVn6th = 0., dErrVn8th = 0.; // error on integrated flow
+
+ Double_t dYield2nd = 0., dYield4th = 0., dYield6th = 0., dYield8th = 0.; // pt yield
+ Double_t dSum2nd = 0., dSum4th = 0., dSum6th = 0., dSum8th = 0.; // needed for normalizing integrated flow
+
+ // looping over pt bins:
+ for(Int_t p=1;p<nBinsPt+1;p++)
+ {
+ dvn2nd = flow2ndPt->GetBinContent(p);
+ dvn4th = flow4thPt->GetBinContent(p);
+ dvn6th = flow6thPt->GetBinContent(p);
+ dvn8th = flow8thPt->GetBinContent(p);
+
+ dErrvn2nd = flow2ndPt->GetBinError(p);
+ dErrvn4th = flow4thPt->GetBinError(p);
+ dErrvn6th = flow6thPt->GetBinError(p);
+ dErrvn8th = flow8thPt->GetBinError(p);
+
+ dYield2nd = yield2ndPt->GetBinContent(p);
+ dYield4th = yield4thPt->GetBinContent(p);
+ dYield6th = yield6thPt->GetBinContent(p);
+ dYield8th = yield8thPt->GetBinContent(p);
+
+ dVn2nd += dvn2nd*dYield2nd;
+ dVn4th += dvn4th*dYield4th;
+ dVn6th += dvn6th*dYield6th;
+ dVn8th += dvn8th*dYield8th;
+
+ dSum2nd += dYield2nd;
+ dSum4th += dYield4th;
+ dSum6th += dYield6th;
+ dSum8th += dYield8th;
+
+ dErrVn2nd += dYield2nd*dYield2nd*dErrvn2nd*dErrvn2nd; // ro be improved (check this relation)
+ dErrVn4th += dYield4th*dYield4th*dErrvn4th*dErrvn4th;
+ dErrVn6th += dYield6th*dYield6th*dErrvn6th*dErrvn6th;
+ dErrVn8th += dYield8th*dYield8th*dErrvn8th*dErrvn8th;
+
+ } // end of for(Int_t p=1;p<nBinsPt+1;p++)
+
+ // normalizing the results for integrated flow:
+ if(dSum2nd)
+ {
+ dVn2nd /= dSum2nd;
+ dErrVn2nd /= (dSum2nd*dSum2nd);
+ dErrVn2nd = TMath::Sqrt(dErrVn2nd);
+ }
+ if(dSum4th)
+ {
+ dVn4th /= dSum4th;
+ dErrVn4th /= (dSum4th*dSum4th);
+ dErrVn4th = TMath::Sqrt(dErrVn4th);
+ }
+ //if(dSum6th) dVn6th/=dSum6th;
+ //if(dSum8th) dVn8th/=dSum8th;
+
+ // storing the results for integrated flow in common histos: (to be improved: new method for this?)
+ if(type == "POI")
+ {
+ fCommonHistsResults2nd->FillIntegratedFlowPOI(dVn2nd,dErrVn2nd);
+ fCommonHistsResults4th->FillIntegratedFlowPOI(dVn4th,dErrVn4th);
+ fCommonHistsResults6th->FillIntegratedFlowPOI(dVn6th,0.); // to be improved (errors)
+ fCommonHistsResults8th->FillIntegratedFlowPOI(dVn8th,0.); // to be improved (errors)
+ }
+ else if (type == "RP")
+ {
+ fCommonHistsResults2nd->FillIntegratedFlowRP(dVn2nd,dErrVn2nd);
+ fCommonHistsResults4th->FillIntegratedFlowRP(dVn4th,dErrVn4th);
+ fCommonHistsResults6th->FillIntegratedFlowRP(dVn6th,0.); // to be improved (errors)
+ fCommonHistsResults8th->FillIntegratedFlowRP(dVn8th,0.); // to be improved (errors)
+ }
+
+ delete flow2ndPt;
+ delete flow4thPt;
+ //delete flow6thPt;
+ //delete flow8thPt;
+
+ delete yield2ndPt;
+ delete yield4thPt;
+ delete yield6thPt;
+ delete yield8thPt;
+
+} // end of AliFlowAnalysisWithQCumulants::CalculateFinalResultsForRPandPOIIntegratedFlow(TString type)
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::InitializeArraysForDistributions()
+{
+ // Initialize all arrays used for distributions.
+
+ // a) Initialize arrays of histograms used to hold distributions of correlations;
+ // b) Initialize array to hold min and max values of correlations.
+
+ // a) Initialize arrays of histograms used to hold distributions of correlations:
+ for(Int_t di=0;di<4;di++) // distribution index
+ {
+ fDistributions[di] = NULL;
+ }
+
+ // b) Initialize default min and max values of correlations:
+ // (Remark: The default values bellow were chosen for v2=5% and M=500)
+ fMinValueOfCorrelation[0] = -0.01; // <2>_min
+ fMaxValueOfCorrelation[0] = 0.04; // <2>_max
+ fMinValueOfCorrelation[1] = -0.00002; // <4>_min
+ fMaxValueOfCorrelation[1] = 0.00015; // <4>_max
+ fMinValueOfCorrelation[2] = -0.0000003; // <6>_min
+ fMaxValueOfCorrelation[2] = 0.0000006; // <6>_max
+ fMinValueOfCorrelation[3] = -0.000000006; // <8>_min
+ fMaxValueOfCorrelation[3] = 0.000000003; // <8>_max
+
+} // end of void AliFlowAnalysisWithQCumulants::InitializeArraysForDistributions()
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::BookEverythingForDistributions()
+{
+ // a) Book profile to hold all flags for distributions of correlations;
+ // b) Book all histograms to hold distributions of correlations.
+
+ TString correlationIndex[4] = {"<2>","<4>","<6>","<8>"}; // to be improved (should I promote this to data members?)
+
+ // a) Book profile to hold all flags for distributions of correlations:
+ TString distributionsFlagsName = "fDistributionsFlags";
+ distributionsFlagsName += fAnalysisLabel->Data();
+ fDistributionsFlags = new TProfile(distributionsFlagsName.Data(),"Flags for Distributions of Correlations",9,0,9);
+ fDistributionsFlags->SetTickLength(-0.01,"Y");
+ fDistributionsFlags->SetMarkerStyle(25);
+ fDistributionsFlags->SetLabelSize(0.05);
+ fDistributionsFlags->SetLabelOffset(0.02,"Y");
+ fDistributionsFlags->GetXaxis()->SetBinLabel(1,"Store or not?");
+ fDistributionsFlags->GetXaxis()->SetBinLabel(2,"<2>_{min}");
+ fDistributionsFlags->GetXaxis()->SetBinLabel(3,"<2>_{max}");
+ fDistributionsFlags->GetXaxis()->SetBinLabel(4,"<4>_{min}");
+ fDistributionsFlags->GetXaxis()->SetBinLabel(5,"<4>_{max}");
+ fDistributionsFlags->GetXaxis()->SetBinLabel(6,"<6>_{min}");
+ fDistributionsFlags->GetXaxis()->SetBinLabel(7,"<6>_{max}");
+ fDistributionsFlags->GetXaxis()->SetBinLabel(8,"<8>_{min}");
+ fDistributionsFlags->GetXaxis()->SetBinLabel(9,"<8>_{max}");
+ fDistributionsList->Add(fDistributionsFlags);
+
+ // b) Book all histograms to hold distributions of correlations.
+ if(fStoreDistributions)
+ {
+ TString distributionsName = "fDistributions";
+ distributionsName += fAnalysisLabel->Data();
+ for(Int_t di=0;di<4;di++) // distribution index
+ {
+ fDistributions[di] = new TH1D(Form("Distribution of %s",correlationIndex[di].Data()),Form("Distribution of %s",correlationIndex[di].Data()),10000,fMinValueOfCorrelation[di],fMaxValueOfCorrelation[di]);
+ fDistributions[di]->SetXTitle(correlationIndex[di].Data());
+ fDistributionsList->Add(fDistributions[di]);
+ } // end of for(Int_t di=0;di<4;di++) // distribution index
+ } // end of if(fStoreDistributions)
+
+} // end of void AliFlowAnalysisWithQCumulants::BookEverythingForDistributions()
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::StoreFlagsForDistributions()
+{
+ // Store all flags for distributiuons of correlations in profile fDistributionsFlags.
+
+ if(!fDistributionsFlags)
+ {
+ cout<<"WARNING: fDistributionsFlags is NULL in AFAWQC::SDF() !!!!"<<endl;
+ exit(0);
+ }
+
+ fDistributionsFlags->Fill(0.5,(Int_t)fStoreDistributions); // histos with distributions of correlations stored or not in the output file
+ // store min and max values of correlations:
+ for(Int_t di=0;di<4;di++) // distribution index
+ {
+ fDistributionsFlags->Fill(1.5+2.*(Double_t)di,fMinValueOfCorrelation[di]);
+ fDistributionsFlags->Fill(2.5+2.*(Double_t)di,fMaxValueOfCorrelation[di]);
+ }
+
+} // end of void AliFlowAnalysisWithQCumulants::StoreFlagsForDistributions()
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::StoreDistributionsOfCorrelations()
+{
+ // Store distributions of correlations.
+
+ if(!(fIntFlowCorrelationsEBE && fIntFlowEventWeightsForCorrelationsEBE))
+ {
+ cout<<"WARNING: fIntFlowCorrelationsEBE && fIntFlowEventWeightsForCorrelationsEBE"<<endl;
+ cout<<" is NULL in AFAWQC::SDOC() !!!!"<<endl;
+ exit(0);
+ }
+
+ for(Int_t di=0;di<4;di++) // distribution index
+ {
+ if(!fDistributions[di])
+ {
+ cout<<"WARNING: fDistributions[di] is NULL in AFAWQC::SDOC() !!!!"<<endl;
+ cout<<"di = "<<di<<endl;
+ exit(0);
+ } else
+ {
+ fDistributions[di]->Fill(fIntFlowCorrelationsEBE->GetBinContent(di+1),fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(di+1));
+ }
+ } // end of for(Int_t di=0;di<4;di++) // distribution index
+
+} // end of void AliFlowAnalysisWithQCumulants::StoreDistributionsOfCorrelations()
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::BookAndNestAllLists()
+{
+ // Book and nest all lists nested in the base list fHistList.
+ // a) Book and nest lists for integrated flow;
+ // b) Book and nest lists for differential flow;
+ // c) Book and nest list for particle weights;
+ // d) Book and nest list for distributions;
+ // e) Book and nest list for nested loops;
+
+ // a) Book and nest all lists for integrated flow:
+ // base list for integrated flow:
+ fIntFlowList = new TList();
+ fIntFlowList->SetName("Integrated Flow");
+ fIntFlowList->SetOwner(kTRUE);
+ fHistList->Add(fIntFlowList);
+ // list holding profiles:
+ fIntFlowProfiles = new TList();
+ fIntFlowProfiles->SetName("Profiles");
+ fIntFlowProfiles->SetOwner(kTRUE);
+ fIntFlowList->Add(fIntFlowProfiles);
+ // list holding histograms with results:
+ fIntFlowResults = new TList();
+ fIntFlowResults->SetName("Results");
+ fIntFlowResults->SetOwner(kTRUE);
+ fIntFlowList->Add(fIntFlowResults);
+
+ // b) Book and nest lists for differential flow;
+ fDiffFlowList = new TList();
+ fDiffFlowList->SetName("Differential Flow");
+ fDiffFlowList->SetOwner(kTRUE);
+ fHistList->Add(fDiffFlowList);
+ // list holding profiles:
+ fDiffFlowProfiles = new TList();
+ fDiffFlowProfiles->SetName("Profiles");
+ fDiffFlowProfiles->SetOwner(kTRUE);
+ fDiffFlowList->Add(fDiffFlowProfiles);
+ // list holding histograms with results:
+ fDiffFlowResults = new TList();
+ fDiffFlowResults->SetName("Results");
+ fDiffFlowResults->SetOwner(kTRUE);
+ fDiffFlowList->Add(fDiffFlowResults);
+ // flags used for naming nested lists in list fDiffFlowProfiles and fDiffFlowResults:
+ TList list;
+ list.SetOwner(kTRUE);
+ TString typeFlag[2] = {"RP","POI"};
+ TString ptEtaFlag[2] = {"p_{T}","#eta"};
+ TString powerFlag[2] = {"linear","quadratic"};
+ // nested lists in fDiffFlowProfiles (~/Differential Flow/Profiles):
+ for(Int_t t=0;t<2;t++) // type: RP or POI
+ {
+ for(Int_t pe=0;pe<2;pe++) // pt or eta
+ {
+ // list holding profiles with correlations:
+ fDiffFlowCorrelationsProList[t][pe] = (TList*)list.Clone();
+ fDiffFlowCorrelationsProList[t][pe]->SetName(Form("Profiles with correlations (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data()));
+ fDiffFlowProfiles->Add(fDiffFlowCorrelationsProList[t][pe]);
+ // list holding profiles with products of correlations:
+ fDiffFlowProductOfCorrelationsProList[t][pe] = (TList*)list.Clone();
+ fDiffFlowProductOfCorrelationsProList[t][pe]->SetName(Form("Profiles with products of correlations (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data()));
+ fDiffFlowProfiles->Add(fDiffFlowProductOfCorrelationsProList[t][pe]);
+ // list holding profiles with corrections:
+ fDiffFlowCorrectionsProList[t][pe] = (TList*)list.Clone();
+ fDiffFlowCorrectionsProList[t][pe]->SetName(Form("Profiles with correction terms for NUA (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data()));
+ fDiffFlowProfiles->Add(fDiffFlowCorrectionsProList[t][pe]);
+ } // end of for(Int_t pe=0;pe<2;pe++) // pt or eta
+ } // end of for(Int_t t=0;t<2;t++) // type: RP or POI
+ // nested lists in fDiffFlowResults (~/Differential Flow/Results):
+ for(Int_t t=0;t<2;t++) // type: RP or POI
+ {
+ for(Int_t pe=0;pe<2;pe++) // pt or eta
+ {
+ // list holding histograms with correlations:
+ fDiffFlowCorrelationsHistList[t][pe] = (TList*)list.Clone();
+ fDiffFlowCorrelationsHistList[t][pe]->SetName(Form("Correlations (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data()));
+ fDiffFlowResults->Add(fDiffFlowCorrelationsHistList[t][pe]);
+ // list holding histograms with corrections:
+ fDiffFlowCorrectionsHistList[t][pe] = (TList*)list.Clone();
+ fDiffFlowCorrectionsHistList[t][pe]->SetName(Form("Histograms with correction terms for NUA (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data()));
+ fDiffFlowResults->Add(fDiffFlowCorrectionsHistList[t][pe]);
+ for(Int_t power=0;power<2;power++)
+ {
+ // list holding histograms with sums of event weights:
+ fDiffFlowSumOfEventWeightsHistList[t][pe][power] = (TList*)list.Clone();
+ fDiffFlowSumOfEventWeightsHistList[t][pe][power]->SetName(Form("Sum of %s event weights (%s, %s)",powerFlag[power].Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data()));
+ fDiffFlowResults->Add(fDiffFlowSumOfEventWeightsHistList[t][pe][power]);
+ } // end of for(Int_t power=0;power<2;power++)
+ // list holding histograms with sums of products of event weights:
+ fDiffFlowSumOfProductOfEventWeightsHistList[t][pe] = (TList*)list.Clone();
+ fDiffFlowSumOfProductOfEventWeightsHistList[t][pe]->SetName(Form("Sum of products of event weights (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data()));
+ fDiffFlowResults->Add(fDiffFlowSumOfProductOfEventWeightsHistList[t][pe]);
+ // list holding histograms with covariances of correlations:
+ fDiffFlowCovariancesHistList[t][pe] = (TList*)list.Clone();
+ fDiffFlowCovariancesHistList[t][pe]->SetName(Form("Covariances of correlations (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data()));
+ fDiffFlowResults->Add(fDiffFlowCovariancesHistList[t][pe]);
+ // list holding histograms with differential Q-cumulants:
+ fDiffFlowCumulantsHistList[t][pe] = (TList*)list.Clone();
+ fDiffFlowCumulantsHistList[t][pe]->SetName(Form("Differential Q-cumulants (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data()));
+ fDiffFlowResults->Add(fDiffFlowCumulantsHistList[t][pe]);
+ // list holding histograms with differential flow estimates from Q-cumulants:
+ fDiffFlowHistList[t][pe] = (TList*)list.Clone();
+ fDiffFlowHistList[t][pe]->SetName(Form("Differential flow (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data()));
+ fDiffFlowResults->Add(fDiffFlowHistList[t][pe]);
+ } // end of for(Int_t pe=0;pe<2;pe++) // pt or eta
+ } // end of for(Int_t t=0;t<2;t++) // type: RP or POI
+
+ // c) Book and nest list for particle weights:
+ fWeightsList->SetName("Weights");
+ fWeightsList->SetOwner(kTRUE);
+ fHistList->Add(fWeightsList);
+
+ // d) Book and nest list for distributions:
+ fDistributionsList = new TList();
+ fDistributionsList->SetName("Distributions");
+ fDistributionsList->SetOwner(kTRUE);
+ fHistList->Add(fDistributionsList);
+
+ // e) Book and nest list for nested loops:
+ fNestedLoopsList = new TList();
+ fNestedLoopsList->SetName("Nested Loops");
+ fNestedLoopsList->SetOwner(kTRUE);
+ fHistList->Add(fNestedLoopsList);
+
+} // end of void AliFlowAnalysisWithQCumulants::BookAndNestAllLists()
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::FillCommonHistResultsDiffFlow(TString type)
+{
+ // fill common result histograms for differential flow
+
+ Int_t typeFlag = -1;
+ //Int_t ptEtaFlag = -1;
+
+ if(type == "RP")
+ {
+ typeFlag = 0;
+ } else if(type == "POI")
+ {
+ typeFlag = 1;
+ }
+
+ // shortcuts:
+ Int_t t = typeFlag;
+ //Int_t pe = ptEtaFlag;
+
+ // to be improved (implement protection here)
+
+ if(!(fCommonHistsResults2nd && fCommonHistsResults4th && fCommonHistsResults6th && fCommonHistsResults8th))
+ {
+ cout<<"WARNING: fCommonHistsResults2nd && fCommonHistsResults4th && fCommonHistsResults6th && fCommonHistsResults8th"<<endl;
+ cout<<" is NULL in AFAWQC::FCHRIF() !!!!"<<endl;
+ exit(0);
+ }
+
+ // pt:
+ for(Int_t p=1;p<=fnBinsPt;p++)
+ {
+ Double_t v2 = fDiffFlow[t][0][0]->GetBinContent(p);
+ Double_t v4 = fDiffFlow[t][0][1]->GetBinContent(p);
+ Double_t v6 = fDiffFlow[t][0][2]->GetBinContent(p);
+ Double_t v8 = fDiffFlow[t][0][3]->GetBinContent(p);
+
+ Double_t v2Error = fDiffFlow[t][0][0]->GetBinError(p);
+ Double_t v4Error = fDiffFlow[t][0][1]->GetBinError(p);
+ //Double_t v6Error = fFinalFlow1D[t][pW][nua][0][2]->GetBinError(p);
+ //Double_t v8Error = fFinalFlow1D[t][pW][nua][0][3]->GetBinError(p);
+
+ if(type == "RP")
+ {
+ fCommonHistsResults2nd->FillDifferentialFlowPtRP(p,v2,v2Error);
+ fCommonHistsResults4th->FillDifferentialFlowPtRP(p,v4,v4Error);
+ fCommonHistsResults6th->FillDifferentialFlowPtRP(p,v6,0.);
+ fCommonHistsResults8th->FillDifferentialFlowPtRP(p,v8,0.);
+ } else if(type == "POI")
+ {
+ fCommonHistsResults2nd->FillDifferentialFlowPtPOI(p,v2,v2Error);
+ fCommonHistsResults4th->FillDifferentialFlowPtPOI(p,v4,v4Error);
+ fCommonHistsResults6th->FillDifferentialFlowPtPOI(p,v6,0.);
+ fCommonHistsResults8th->FillDifferentialFlowPtPOI(p,v8,0.);
+ }
+ } // end of for(Int_t p=1;p<=fnBinsPt;p++)
+
+ // eta:
+ for(Int_t e=1;e<=fnBinsEta;e++)
+ {
+ Double_t v2 = fDiffFlow[t][1][0]->GetBinContent(e);
+ Double_t v4 = fDiffFlow[t][1][1]->GetBinContent(e);
+ Double_t v6 = fDiffFlow[t][1][2]->GetBinContent(e);
+ Double_t v8 = fDiffFlow[t][1][3]->GetBinContent(e);
+
+ Double_t v2Error = fDiffFlow[t][1][0]->GetBinError(e);
+ Double_t v4Error = fDiffFlow[t][1][1]->GetBinError(e);
+ //Double_t v6Error = fDiffFlow[t][1][2]->GetBinError(e);
+ //Double_t v8Error = fDiffFlow[t][1][3]->GetBinError(e);
+
+ if(type == "RP")
+ {
+ fCommonHistsResults2nd->FillDifferentialFlowEtaRP(e,v2,v2Error);
+ fCommonHistsResults4th->FillDifferentialFlowEtaRP(e,v4,v4Error);
+ fCommonHistsResults6th->FillDifferentialFlowEtaRP(e,v6,0.);
+ fCommonHistsResults8th->FillDifferentialFlowEtaRP(e,v8,0.);
+ } else if(type == "POI")
+ {
+ fCommonHistsResults2nd->FillDifferentialFlowEtaPOI(e,v2,v2Error);
+ fCommonHistsResults4th->FillDifferentialFlowEtaPOI(e,v4,v4Error);
+ fCommonHistsResults6th->FillDifferentialFlowEtaPOI(e,v6,0.);
+ fCommonHistsResults8th->FillDifferentialFlowEtaPOI(e,v8,0.);
+ }
+ } // end of for(Int_t e=1;e<=fnBinsEta;e++)
+
+} // end of void AliFlowAnalysisWithQCumulants::FillCommonHistResultsDiffFlow(TString type, Bool_t useParticleWeights, TString eventWeights, Bool_t correctedForNUA)
+
+//================================================================================================================================
+
+void AliFlowAnalysisWithQCumulants::AccessConstants()
+{
+ // Access needed common constants from AliFlowCommonConstants.
+
+ fnBinsPhi = AliFlowCommonConstants::GetMaster()->GetNbinsPhi();
+ fPhiMin = AliFlowCommonConstants::GetMaster()->GetPhiMin();
+ fPhiMax = AliFlowCommonConstants::GetMaster()->GetPhiMax();
+ if(fnBinsPhi) fPhiBinWidth = (fPhiMax-fPhiMin)/fnBinsPhi;
+ fnBinsPt = AliFlowCommonConstants::GetMaster()->GetNbinsPt();
+ fPtMin = AliFlowCommonConstants::GetMaster()->GetPtMin();
+ fPtMax = AliFlowCommonConstants::GetMaster()->GetPtMax();
+ if(fnBinsPt) fPtBinWidth = (fPtMax-fPtMin)/fnBinsPt;
+ fnBinsEta = AliFlowCommonConstants::GetMaster()->GetNbinsEta();
+ fEtaMin = AliFlowCommonConstants::GetMaster()->GetEtaMin();
+ fEtaMax = AliFlowCommonConstants::GetMaster()->GetEtaMax();
+ if(fnBinsEta) fEtaBinWidth = (fEtaMax-fEtaMin)/fnBinsEta;
+
+} // end of void AliFlowAnalysisWithQCumulants::AccessConstants()
+
+//================================================================================================================================
+
+void AliFlowAnalysisWithQCumulants::CrossCheckSettings()
+{
+ // a) Cross check if the choice for multiplicity weights make sense;
+
+ // a) Cross check if the choice for multiplicity weights make sense:
+ if(strcmp(fMultiplicityWeight->Data(),"combinations") &&
+ strcmp(fMultiplicityWeight->Data(),"unit") &&
+ strcmp(fMultiplicityWeight->Data(),"multiplicity"))
+ {
+ cout<<"WARNING (QC): Multiplicity weight can be either \"combinations\", \"unit\""<<endl;
+ cout<<" or \"multiplicity\". Certainly not \""<<fMultiplicityWeight->Data()<<"\"."<<endl;
+ exit(0);
+ }
+
+} // end of void AliFlowAnalysisWithQCumulants::CrossCheckSettings()
+
+//================================================================================================================================
+
+void AliFlowAnalysisWithQCumulants::CalculateIntFlowSumOfEventWeights()
+{
+ // Calculate sum of linear and quadratic event weights for correlations.
+
+ // multiplicity:
+ Double_t dMult = (*fSMpk)(0,0);
+
+ for(Int_t p=0;p<2;p++) // power-1
+ {
+ for(Int_t ci=0;ci<4;ci++) // correlation index
+ {
+ fIntFlowSumOfEventWeights[p]->Fill(ci+0.5,pow(fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(ci+1),p+1));
+ if(fCalculateCumulantsVsM)
+ {
+ fIntFlowSumOfEventWeightsVsM[ci][p]->Fill(dMult+0.5,pow(fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(ci+1),p+1)); // to be improved: dMult => sum of weights?
+ }
+ }
+ }
+
+} // end of void AliFlowAnalysisWithQCumulants::CalculateIntFlowSumOfEventWeights()
+
+//================================================================================================================================
+
+void AliFlowAnalysisWithQCumulants::CalculateIntFlowSumOfEventWeightsNUA()
+{
+ // Calculate sum of linear and quadratic event weights for NUA terms.
+
+ for(Int_t sc=0;sc<2;sc++) // sin or cos terms
+ {
+ for(Int_t p=0;p<2;p++) // power-1
+ {
+ for(Int_t ci=0;ci<4;ci++) // nua term index
+ {
+ fIntFlowSumOfEventWeightsNUA[sc][p]->Fill(ci+0.5,pow(fIntFlowEventWeightForCorrectionTermsForNUAEBE[sc]->GetBinContent(ci+1),p+1));
+ }
+ }
+ }
+
+} // end of void AliFlowAnalysisWithQCumulants::CalculateIntFlowSumOfEventWeightsNUA()
+
+//================================================================================================================================
+
+void AliFlowAnalysisWithQCumulants::CalculateIntFlowSumOfProductOfEventWeights()
+{
+ // Calculate sum of product of event weights for correlations.
+
+ // multiplicity:
+ Double_t dMult = (*fSMpk)(0,0);
+
+ Int_t counter = 0;
+
+ for(Int_t ci1=1;ci1<4;ci1++)
+ {
+ for(Int_t ci2=ci1+1;ci2<=4;ci2++)
+ {
+ fIntFlowSumOfProductOfEventWeights->Fill(0.5+counter,
+ fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(ci1)*
+ fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(ci2));
+ if(fCalculateCumulantsVsM)
+ {
+ fIntFlowSumOfProductOfEventWeightsVsM[counter]->Fill(dMult+0.5, // to be improved: dMult => sum of weights?
+ fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(ci1)*
+ fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(ci2));
+ } // end of if(fCalculateCumulantsVsM)
+ counter++;
+ }
+ }
+
+} // end of void AliFlowAnalysisWithQCumulants::CalculateIntFlowSumOfProductOfEventWeights()
+
+//================================================================================================================================
+
+void AliFlowAnalysisWithQCumulants::CalculateIntFlowSumOfProductOfEventWeightsNUA()
+{
+ // Calculate sum of product of event weights for NUA terms.
+
+ // w_{<2>} * w_{<cos(#phi)>}:
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(0.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(1));
+ // w_{<2>} * w_{<sin(#phi)>}:
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(1.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(1));
+ // w_{<cos(#phi)> * w_{<sin(#phi)>}:
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(2.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(1));
+ // w_{<2>} * w{<cos(phi1+phi2)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(3.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(2));
+ // w_{<2>} * w{<sin(phi1+phi2)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(4.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(2));
+ // w_{<2>} * w{<cos(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(5.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(3));
+ // w_{<2>} * w{<sin(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(6.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(3));
+ // w_{<4>} * w{<cos(phi1)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(7.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(2)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(1));
+ // w_{<4>} * w{<sin(phi1)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(8.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(2)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(1));
+ // w_{<4>} * w{<cos(phi1+phi2)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(9.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(2)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(2));
+ // w_{<4>} * w{<sin(phi1+phi2)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(10.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(2)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(2));
+ // w_{<4>} * w{<cos(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(11.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(2)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(3));
+ // w_{<4>} * w{<sin(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(12.5,fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(2)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(3));
+ // w_{<cos(phi1)>} * w{<cos(phi1+phi2)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(13.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(2));
+ // w_{<cos(phi1)>} * w{<sin(phi1+phi2)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(14.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(2));
+ // w_{<cos(phi1)>} * w{<cos(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(15.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(3));
+ // w_{<cos(phi1)>} * w{<sin(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(16.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(3));
+ // w_{<sin(phi1)>} * w{<cos(phi1+phi2)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(17.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(2));
+ // w_{<sin(phi1)>} * w{<sin(phi1+phi2)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(18.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(2));
+ // w_{<sin(phi1)>} * w{<cos(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(19.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(3));
+ // w_{<sin(phi1)>} * w{<sin(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(20.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(1)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(3));
+ // w_{<cos(phi1+phi2)>} * w{<sin(phi1+phi2))>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(21.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(2)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(2));
+ // w_{<cos(phi1+phi2)>} * w{<cos(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(22.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(2)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(3));
+ // w_{<cos(phi1+phi2)>} * w{<sin(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(23.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(2)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(3));
+ // w_{<sin(phi1+phi2)>} * w{<cos(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(24.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(2)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(3));
+ // w_{<sin(phi1+phi2)>} * w{<sin(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(25.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(2)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(3));
+ // w_{<cos(phi1-phi2-phi3)>} * w{<sin(phi1-phi2-phi3)>}
+ fIntFlowSumOfProductOfEventWeightsNUA->Fill(26.5,fIntFlowEventWeightForCorrectionTermsForNUAEBE[1]->GetBinContent(3)*
+ fIntFlowEventWeightForCorrectionTermsForNUAEBE[0]->GetBinContent(3));
+
+} // end of void AliFlowAnalysisWithQCumulants::CalculateIntFlowIntFlowSumOfProductOfEventWeightsNUA()
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::CalculateDiffFlowCorrelations(TString type, TString ptOrEta)
+{
+ // calculate reduced correlations for RPs or POIs in pt or eta bins
+
+ // multiplicity:
+ Double_t dMult = (*fSMpk)(0,0);
+
+ // real and imaginary parts of non-weighted Q-vectors evaluated in harmonics n, 2n, 3n and 4n:
+ Double_t dReQ1n = (*fReQ)(0,0);
+ Double_t dReQ2n = (*fReQ)(1,0);
+ //Double_t dReQ3n = (*fReQ)(2,0);
+ //Double_t dReQ4n = (*fReQ)(3,0);
+ Double_t dImQ1n = (*fImQ)(0,0);
+ Double_t dImQ2n = (*fImQ)(1,0);
+ //Double_t dImQ3n = (*fImQ)(2,0);
+ //Double_t dImQ4n = (*fImQ)(3,0);
+
+ // reduced correlations are stored in fDiffFlowCorrelationsPro[0=RP,1=POI][0=pt,1=eta][correlation index]. Correlation index runs as follows:
+ //
+ // 0: <<2'>>
+ // 1: <<4'>>
+ // 2: <<6'>>
+ // 3: <<8'>>
+
+ Int_t t = -1; // type flag
+ Int_t pe = -1; // ptEta flag
+
+ if(type == "RP")
+ {
+ t = 0;
+ } else if(type == "POI")
+ {
+ t = 1;
+ }
+
+ if(ptOrEta == "Pt")
+ {
+ pe = 0;
+ } else if(ptOrEta == "Eta")
+ {
+ pe = 1;
+ }
+
+ Int_t nBinsPtEta[2] = {fnBinsPt,fnBinsEta};
+ Double_t minPtEta[2] = {fPtMin,fEtaMin};
+ //Double_t maxPtEta[2] = {fPtMax,fEtaMax};
+ Double_t binWidthPtEta[2] = {fPtBinWidth,fEtaBinWidth};
+
+ // looping over all bins and calculating reduced correlations:
+ for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+ {
+ // real and imaginary parts of p_{m*n,0} (non-weighted Q-vector evaluated for POIs in particular pt or eta bin):
+ Double_t p1n0kRe = 0.;
+ Double_t p1n0kIm = 0.;
+
+ // number of POIs in particular pt or eta bin:
+ Double_t mp = 0.;
+
+ // real and imaginary parts of q_{m*n,0} (non-weighted Q-vector evaluated for particles which are both RPs and POIs in particular pt or eta bin):
+ Double_t q1n0kRe = 0.;
+ Double_t q1n0kIm = 0.;
+ Double_t q2n0kRe = 0.;
+ Double_t q2n0kIm = 0.;
+
+ // number of particles which are both RPs and POIs in particular pt or eta bin:
+ Double_t mq = 0.;
+
+ if(type == "POI")
+ {
+ // q_{m*n,0}:
+ q1n0kRe = fReRPQ1dEBE[2][pe][0][0]->GetBinContent(fReRPQ1dEBE[2][pe][0][0]->GetBin(b))
+ * fReRPQ1dEBE[2][pe][0][0]->GetBinEntries(fReRPQ1dEBE[2][pe][0][0]->GetBin(b));
+ q1n0kIm = fImRPQ1dEBE[2][pe][0][0]->GetBinContent(fImRPQ1dEBE[2][pe][0][0]->GetBin(b))
+ * fImRPQ1dEBE[2][pe][0][0]->GetBinEntries(fImRPQ1dEBE[2][pe][0][0]->GetBin(b));
+ q2n0kRe = fReRPQ1dEBE[2][pe][1][0]->GetBinContent(fReRPQ1dEBE[2][pe][1][0]->GetBin(b))
+ * fReRPQ1dEBE[2][pe][1][0]->GetBinEntries(fReRPQ1dEBE[2][pe][1][0]->GetBin(b));
+ q2n0kIm = fImRPQ1dEBE[2][pe][1][0]->GetBinContent(fImRPQ1dEBE[2][pe][1][0]->GetBin(b))
+ * fImRPQ1dEBE[2][pe][1][0]->GetBinEntries(fImRPQ1dEBE[2][pe][1][0]->GetBin(b));
+
+ mq = fReRPQ1dEBE[2][pe][0][0]->GetBinEntries(fReRPQ1dEBE[2][pe][0][0]->GetBin(b)); // to be improved (cross-checked by accessing other profiles here)
+ }
+ else if(type == "RP")
+ {
+ // q_{m*n,0}:
+ q1n0kRe = fReRPQ1dEBE[0][pe][0][0]->GetBinContent(fReRPQ1dEBE[0][pe][0][0]->GetBin(b))
+ * fReRPQ1dEBE[0][pe][0][0]->GetBinEntries(fReRPQ1dEBE[0][pe][0][0]->GetBin(b));
+ q1n0kIm = fImRPQ1dEBE[0][pe][0][0]->GetBinContent(fImRPQ1dEBE[0][pe][0][0]->GetBin(b))
+ * fImRPQ1dEBE[0][pe][0][0]->GetBinEntries(fImRPQ1dEBE[0][pe][0][0]->GetBin(b));
+ q2n0kRe = fReRPQ1dEBE[0][pe][1][0]->GetBinContent(fReRPQ1dEBE[0][pe][1][0]->GetBin(b))
+ * fReRPQ1dEBE[0][pe][1][0]->GetBinEntries(fReRPQ1dEBE[0][pe][1][0]->GetBin(b));
+ q2n0kIm = fImRPQ1dEBE[0][pe][1][0]->GetBinContent(fImRPQ1dEBE[0][pe][1][0]->GetBin(b))
+ * fImRPQ1dEBE[0][pe][1][0]->GetBinEntries(fImRPQ1dEBE[0][pe][1][0]->GetBin(b));
+
+ mq = fReRPQ1dEBE[0][pe][0][0]->GetBinEntries(fReRPQ1dEBE[0][pe][0][0]->GetBin(b)); // to be improved (cross-checked by accessing other profiles here)
+ }
+
+ if(type == "POI")
+ {
+ // p_{m*n,0}:
+ p1n0kRe = fReRPQ1dEBE[1][pe][0][0]->GetBinContent(fReRPQ1dEBE[1][pe][0][0]->GetBin(b))
+ * fReRPQ1dEBE[1][pe][0][0]->GetBinEntries(fReRPQ1dEBE[1][pe][0][0]->GetBin(b));
+ p1n0kIm = fImRPQ1dEBE[1][pe][0][0]->GetBinContent(fImRPQ1dEBE[1][pe][0][0]->GetBin(b))
+ * fImRPQ1dEBE[1][pe][0][0]->GetBinEntries(fImRPQ1dEBE[1][pe][0][0]->GetBin(b));
+
+ mp = fReRPQ1dEBE[1][pe][0][0]->GetBinEntries(fReRPQ1dEBE[1][pe][0][0]->GetBin(b)); // to be improved (cross-checked by accessing other profiles here)
+
+ t = 1; // typeFlag = RP or POI
+ }
+ else if(type == "RP")
+ {
+ // p_{m*n,0} = q_{m*n,0}:
+ p1n0kRe = q1n0kRe;
+ p1n0kIm = q1n0kIm;
+
+ mp = mq;
+
+ t = 0; // typeFlag = RP or POI
+ }
+
+ // 2'-particle correlation for particular (pt,eta) bin:
+ Double_t two1n1nPtEta = 0.;
+ if(mp*dMult-mq)
+ {
+ two1n1nPtEta = (p1n0kRe*dReQ1n+p1n0kIm*dImQ1n-mq)
+ / (mp*dMult-mq);
+
+ if(type == "POI") // to be improved (I do not this if)
+ {
+ // fill profile to get <<2'>> for POIs
+ fDiffFlowCorrelationsPro[1][pe][0]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],two1n1nPtEta,mp*dMult-mq);
+ // histogram to store <2'> for POIs e-b-e (needed in some other methods):
+ fDiffFlowCorrelationsEBE[1][pe][0]->SetBinContent(b,two1n1nPtEta);
+ fDiffFlowEventWeightsForCorrelationsEBE[1][pe][0]->SetBinContent(b,mp*dMult-mq);
+ }
+ else if(type == "RP") // to be improved (I do not this if)
+ {
+ // profile to get <<2'>> for RPs:
+ fDiffFlowCorrelationsPro[0][pe][0]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],two1n1nPtEta,mp*dMult-mq);
+ // histogram to store <2'> for RPs e-b-e (needed in some other methods):
+ fDiffFlowCorrelationsEBE[0][pe][0]->SetBinContent(b,two1n1nPtEta);
+ fDiffFlowEventWeightsForCorrelationsEBE[0][pe][0]->SetBinContent(b,mp*dMult-mq);
+ }
+ } // end of if(mp*dMult-mq)
+
+ // 4'-particle correlation:
+ Double_t four1n1n1n1nPtEta = 0.;
+ if((mp-mq)*dMult*(dMult-1.)*(dMult-2.)
+ + mq*(dMult-1.)*(dMult-2.)*(dMult-3.)) // to be improved (introduce a new variable for this expression)
+ {
+ four1n1n1n1nPtEta = ((pow(dReQ1n,2.)+pow(dImQ1n,2.))*(p1n0kRe*dReQ1n+p1n0kIm*dImQ1n)
+ - q2n0kRe*(pow(dReQ1n,2.)-pow(dImQ1n,2.))
+ - 2.*q2n0kIm*dReQ1n*dImQ1n
+ - p1n0kRe*(dReQ1n*dReQ2n+dImQ1n*dImQ2n)
+ + p1n0kIm*(dImQ1n*dReQ2n-dReQ1n*dImQ2n)
+ - 2.*dMult*(p1n0kRe*dReQ1n+p1n0kIm*dImQ1n)
+ - 2.*(pow(dReQ1n,2.)+pow(dImQ1n,2.))*mq
+ + 6.*(q1n0kRe*dReQ1n+q1n0kIm*dImQ1n)
+ + 1.*(q2n0kRe*dReQ2n+q2n0kIm*dImQ2n)
+ + 2.*(p1n0kRe*dReQ1n+p1n0kIm*dImQ1n)
+ + 2.*mq*dMult
+ - 6.*mq)
+ / ((mp-mq)*dMult*(dMult-1.)*(dMult-2.)
+ + mq*(dMult-1.)*(dMult-2.)*(dMult-3.));
+
+ if(type == "POI")
+ {
+ // profile to get <<4'>> for POIs:
+ fDiffFlowCorrelationsPro[1][pe][1]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],four1n1n1n1nPtEta,
+ (mp-mq)*dMult*(dMult-1.)*(dMult-2.)
+ + mq*(dMult-1.)*(dMult-2.)*(dMult-3.));
+ // histogram to store <4'> for POIs e-b-e (needed in some other methods):
+ fDiffFlowCorrelationsEBE[1][pe][1]->SetBinContent(b,four1n1n1n1nPtEta);
+ fDiffFlowEventWeightsForCorrelationsEBE[1][pe][1]->SetBinContent(b,(mp-mq)*dMult*(dMult-1.)*(dMult-2.)
+ + mq*(dMult-1.)*(dMult-2.)*(dMult-3.));
+ }
+ else if(type == "RP")
+ {
+ // profile to get <<4'>> for RPs:
+ fDiffFlowCorrelationsPro[0][pe][1]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],four1n1n1n1nPtEta,
+ (mp-mq)*dMult*(dMult-1.)*(dMult-2.)
+ + mq*(dMult-1.)*(dMult-2.)*(dMult-3.));
+ // histogram to store <4'> for RPs e-b-e (needed in some other methods):
+ fDiffFlowCorrelationsEBE[0][pe][1]->SetBinContent(b,four1n1n1n1nPtEta);
+ fDiffFlowEventWeightsForCorrelationsEBE[0][pe][1]->SetBinContent(b,(mp-mq)*dMult*(dMult-1.)*(dMult-2.)
+ + mq*(dMult-1.)*(dMult-2.)*(dMult-3.));
+ }
+ } // end of if((mp-mq)*dMult*(dMult-1.)*(dMult-2.)
+ // +mq*(dMult-1.)*(dMult-2.)*(dMult-3.))
+
+ } // end of for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+
+
+} // end of void AliFlowAnalysisWithQCumulants::CalculateDiffFlowCorrelations(TString type, TString ptOrEta);
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::CalculateDiffFlowSumOfEventWeights(TString type, TString ptOrEta)
+{
+ // Calculate sums of various event weights for reduced correlations.
+ // (These quantitites are needed in expressions for unbiased estimators relevant for the statistical errors.)
+
+ Int_t typeFlag = -1;
+ Int_t ptEtaFlag = -1;
+
+ if(type == "RP")
+ {
+ typeFlag = 0;
+ } else if(type == "POI")
+ {
+ typeFlag = 1;
+ }
+
+ if(ptOrEta == "Pt")
+ {
+ ptEtaFlag = 0;
+ } else if(ptOrEta == "Eta")
+ {
+ ptEtaFlag = 1;
+ }
+
+ // shortcuts:
+ Int_t t = typeFlag;
+ Int_t pe = ptEtaFlag;
+
+ // binning:
+ Int_t nBinsPtEta[2] = {fnBinsPt,fnBinsEta};
+ Double_t minPtEta[2] = {fPtMin,fEtaMin};
+ //Double_t maxPtEta[2] = {fPtMax,fEtaMax};
+ Double_t binWidthPtEta[2] = {fPtBinWidth,fEtaBinWidth};
+
+ for(Int_t rpq=0;rpq<3;rpq++)
+ {
+ for(Int_t m=0;m<4;m++)
+ {
+ for(Int_t k=0;k<9;k++)
+ {
+ if(!fReRPQ1dEBE[rpq][pe][m][k])
+ {
+ cout<<"WARNING: fReRPQ1dEBE[rpq][pe][m][k] is NULL in AFAWQC::CSAPOEWFDF() !!!!"<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"rpq = "<<rpq<<endl;
+ cout<<"m = "<<m<<endl;
+ cout<<"k = "<<k<<endl;
+ exit(0);
+ }
+ }
+ }
+ }
+
+ // multiplicities:
+ Double_t dMult = (*fSMpk)(0,0); // total event multiplicity
+ //Double_t mr = 0.; // number of RPs in particular pt or eta bin
+ Double_t mp = 0.; // number of POIs in particular pt or eta bin
+ Double_t mq = 0.; // number of particles which are both RPs and POIs in particular pt or eta bin
+
+ // event weights for reduced correlations:
+ Double_t dw2 = 0.; // event weight for <2'>
+ Double_t dw4 = 0.; // event weight for <4'>
+ //Double_t dw6 = 0.; // event weight for <6'>
+ //Double_t dw8 = 0.; // event weight for <8'>
+
+ // looping over bins:
+ for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+ {
+ if(type == "RP")
+ {
+ mq = fReRPQ1dEBE[0][pe][0][0]->GetBinEntries(b);
+ mp = mq; // trick to use the very same Eqs. bellow both for RP's and POI's diff. flow
+ } else if(type == "POI")
+ {
+ mp = fReRPQ1dEBE[1][pe][0][0]->GetBinEntries(b);
+ mq = fReRPQ1dEBE[2][pe][0][0]->GetBinEntries(b);
+ }
+
+ // event weight for <2'>:
+ dw2 = mp*dMult-mq;
+ fDiffFlowSumOfEventWeights[t][pe][0][0]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw2);
+ fDiffFlowSumOfEventWeights[t][pe][1][0]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],pow(dw2,2.));
+
+ // event weight for <4'>:
+ dw4 = (mp-mq)*dMult*(dMult-1.)*(dMult-2.)
+ + mq*(dMult-1.)*(dMult-2.)*(dMult-3.);
+ fDiffFlowSumOfEventWeights[t][pe][0][1]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw4);
+ fDiffFlowSumOfEventWeights[t][pe][1][1]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],pow(dw4,2.));
+
+ // event weight for <6'>:
+ //dw6 = ...;
+ //fDiffFlowSumOfEventWeights[t][pe][0][2]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw6);
+ //fDiffFlowSumOfEventWeights[t][pe][t][1][2]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],pow(dw6,2.));
+
+ // event weight for <8'>:
+ //dw8 = ...;
+ //fDiffFlowSumOfEventWeights[t][pe][0][3]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw8);
+ //fDiffFlowSumOfEventWeights[t][pe][1][3]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],pow(dw8,2.));
+ } // end of for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+
+} // end of void AliFlowAnalysisWithQCumulants::CalculateDiffFlowSumOfEventWeights()
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::CalculateDiffFlowSumOfProductOfEventWeights(TString type, TString ptOrEta)
+{
+ // Calculate sum of products of various event weights for both types of correlations (the ones for int. and diff. flow).
+ // (These quantitites are needed in expressions for unbiased estimators relevant for the statistical errors.)
+ //
+ // Important: To fill fDiffFlowSumOfProductOfEventWeights[][][][] use bellow table (i,j) with following constraints:
+ // 1.) i<j
+ // 2.) do not store terms which DO NOT include reduced correlations;
+ // Table:
+ // [0=<2>,1=<2'>,2=<4>,3=<4'>,4=<6>,5=<6'>,6=<8>,7=<8'>] x [0=<2>,1=<2'>,2=<4>,3=<4'>,4=<6>,5=<6'>,6=<8>,7=<8'>]
+
+ Int_t typeFlag = -1;
+ Int_t ptEtaFlag = -1;
+
+ if(type == "RP")
+ {
+ typeFlag = 0;
+ } else if(type == "POI")
+ {
+ typeFlag = 1;
+ }
+
+ if(ptOrEta == "Pt")
+ {
+ ptEtaFlag = 0;
+ } else if(ptOrEta == "Eta")
+ {
+ ptEtaFlag = 1;
+ }
+
+ // shortcuts:
+ Int_t t = typeFlag;
+ Int_t pe = ptEtaFlag;
+
+ // binning:
+ Int_t nBinsPtEta[2] = {fnBinsPt,fnBinsEta};
+ Double_t minPtEta[2] = {fPtMin,fEtaMin};
+ //Double_t maxPtEta[2] = {fPtMax,fEtaMax};
+ Double_t binWidthPtEta[2] = {fPtBinWidth,fEtaBinWidth};
+
+ // protection:
+ for(Int_t rpq=0;rpq<3;rpq++)
+ {
+ for(Int_t m=0;m<4;m++)
+ {
+ for(Int_t k=0;k<9;k++)
+ {
+ if(!fReRPQ1dEBE[rpq][pe][m][k])
+ {
+ cout<<"WARNING: fReRPQ1dEBE[rpq][pe][m][k] is NULL in AFAWQC::CSAPOEWFDF() !!!!"<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"rpq = "<<rpq<<endl;
+ cout<<"m = "<<m<<endl;
+ cout<<"k = "<<k<<endl;
+ exit(0);
+ }
+ }
+ }
+ }
+
+ // multiplicities:
+ Double_t dMult = (*fSMpk)(0,0); // total event multiplicity
+ //Double_t mr = 0.; // number of RPs in particular pt or eta bin
+ Double_t mp = 0.; // number of POIs in particular pt or eta bin
+ Double_t mq = 0.; // number of particles which are both RPs and POIs in particular pt or eta bin
+
+ // event weights for correlations:
+ Double_t dW2 = dMult*(dMult-1); // event weight for <2>
+ Double_t dW4 = dMult*(dMult-1)*(dMult-2)*(dMult-3); // event weight for <4>
+ Double_t dW6 = dMult*(dMult-1)*(dMult-2)*(dMult-3)*(dMult-4)*(dMult-5); // event weight for <6>
+ Double_t dW8 = dMult*(dMult-1)*(dMult-2)*(dMult-3)*(dMult-4)*(dMult-5)*(dMult-6)*(dMult-7); // event weight for <8>
+
+ // event weights for reduced correlations:
+ Double_t dw2 = 0.; // event weight for <2'>
+ Double_t dw4 = 0.; // event weight for <4'>
+ //Double_t dw6 = 0.; // event weight for <6'>
+ //Double_t dw8 = 0.; // event weight for <8'>
+
+ // looping over bins:
+ for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+ {
+ if(type == "RP")
+ {
+ mq = fReRPQ1dEBE[0][pe][0][0]->GetBinEntries(b);
+ mp = mq; // trick to use the very same Eqs. bellow both for RP's and POI's diff. flow
+ } else if(type == "POI")
+ {
+ mp = fReRPQ1dEBE[1][pe][0][0]->GetBinEntries(b);
+ mq = fReRPQ1dEBE[2][pe][0][0]->GetBinEntries(b);
+ }
+
+ // event weight for <2'>:
+ dw2 = mp*dMult-mq;
+ fDiffFlowSumOfProductOfEventWeights[t][pe][0][1]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dW2*dw2); // storing product of even weights for <2> and <2'>
+ fDiffFlowSumOfProductOfEventWeights[t][pe][1][2]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw2*dW4); // storing product of even weights for <4> and <2'>
+ fDiffFlowSumOfProductOfEventWeights[t][pe][1][4]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw2*dW6); // storing product of even weights for <6> and <2'>
+ fDiffFlowSumOfProductOfEventWeights[t][pe][1][6]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw2*dW8); // storing product of even weights for <8> and <2'>
+
+ // event weight for <4'>:
+ dw4 = (mp-mq)*dMult*(dMult-1.)*(dMult-2.)
+ + mq*(dMult-1.)*(dMult-2.)*(dMult-3.);
+ fDiffFlowSumOfProductOfEventWeights[t][pe][0][3]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dW2*dw4); // storing product of even weights for <2> and <4'>
+ fDiffFlowSumOfProductOfEventWeights[t][pe][1][3]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw2*dw4); // storing product of even weights for <2'> and <4'>
+ fDiffFlowSumOfProductOfEventWeights[t][pe][2][3]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dW4*dw4); // storing product of even weights for <4> and <4'>
+ fDiffFlowSumOfProductOfEventWeights[t][pe][3][4]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw4*dW6); // storing product of even weights for <6> and <4'>
+ fDiffFlowSumOfProductOfEventWeights[t][pe][3][6]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw4*dW8); // storing product of even weights for <8> and <4'>
+
+ // event weight for <6'>:
+ //dw6 = ...;
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][0][5]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dW2*dw6); // storing product of even weights for <2> and <6'>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][1][5]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw2*dw6); // storing product of even weights for <2'> and <6'>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][2][5]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dW4*dw6); // storing product of even weights for <4> and <6'>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][3][5]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw4*dw6); // storing product of even weights for <4'> and <6'>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][4][5]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dW6*dw6); // storing product of even weights for <6> and <6'>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][5][6]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw6*dW8); // storing product of even weights for <6'> and <8>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][5][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw6*dw8); // storing product of even weights for <6'> and <8'>
+
+ // event weight for <8'>:
+ //dw8 = ...;
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][0][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dW2*dw8); // storing product of even weights for <2> and <8'>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][1][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw2*dw8); // storing product of even weights for <2'> and <8'>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][2][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dW4*dw8); // storing product of even weights for <4> and <8'>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][3][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw4*dw8); // storing product of even weights for <4'> and <8'>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][4][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dW6*dw8); // storing product of even weights for <6> and <8'>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][5][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dw6*dw8); // storing product of even weights for <6'> and <8'>
+ //fDiffFlowSumOfProductOfEventWeights[t][pe][6][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],dW8*dw8); // storing product of even weights for <8> and <8'>
+
+ // Table:
+ // [0=<2>,1=<2'>,2=<4>,3=<4'>,4=<6>,5=<6'>,6=<8>,7=<8'>] x [0=<2>,1=<2'>,2=<4>,3=<4'>,4=<6>,5=<6'>,6=<8>,7=<8'>]
+
+ } // end of for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+
+
+
+} // end of void AliFlowAnalysisWithQCumulants::CalculateDiffFlowSumOfProductOfEventWeights(TString type, TString ptOrEta)
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::FinalizeReducedCorrelations(TString type, TString ptOrEta)
+{
+ // Transfer profiles into histograms and calculate statistical errors correctly.
+
+ Int_t typeFlag = -1;
+ Int_t ptEtaFlag = -1;
+
+ if(type == "RP")
+ {
+ typeFlag = 0;
+ } else if(type == "POI")
+ {
+ typeFlag = 1;
+ }
+
+ if(ptOrEta == "Pt")
+ {
+ ptEtaFlag = 0;
+ } else if(ptOrEta == "Eta")
+ {
+ ptEtaFlag = 1;
+ }
+
+ // shortcuts:
+ Int_t t = typeFlag;
+ Int_t pe = ptEtaFlag;
+
+ for(Int_t rci=0;rci<4;rci++)
+ {
+ if(!fDiffFlowCorrelationsPro[t][pe][rci])
+ {
+ cout<<"WARNING: fDiffFlowCorrelationsPro[t][pe][rci] is NULL in AFAWQC::FRC() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"rci = "<<rci<<endl;
+ exit(0);
+ }
+ for(Int_t power=0;power<2;power++)
+ {
+ if(!fDiffFlowSumOfEventWeights[t][pe][power][rci])
+ {
+ cout<<"WARNING: fDiffFlowSumOfEventWeights[t][pe][power][rci] is NULL in AFAWQC::FRC() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"power = "<<power<<endl;
+ cout<<"rci = "<<rci<<endl;
+ exit(0);
+ }
+ } // end of for(Int_t power=0;power<2;power++)
+ } // end of for(Int_t rci=0;rci<4;rci++)
+
+ // common:
+ Int_t nBinsPtEta[2] = {fnBinsPt,fnBinsEta};
+
+ // transfer 1D profile into 1D histogram:
+ Double_t correlation = 0.;
+ Double_t spread = 0.;
+ Double_t sumOfWeights = 0.; // sum of weights for particular reduced correlations for particular pt or eta bin
+ Double_t sumOfSquaredWeights = 0.; // sum of squared weights for particular reduced correlations for particular pt or eta bin
+ Double_t error = 0.; // error = termA * spread * termB
+ // termA = (sqrt(sumOfSquaredWeights)/sumOfWeights)
+ // termB = 1/pow(1-termA^2,0.5)
+ Double_t termA = 0.;
+ Double_t termB = 0.;
+ for(Int_t rci=0;rci<4;rci++) // index of reduced correlation
+ {
+ for(Int_t b=1;b<=nBinsPtEta[pe];b++) // number of pt or eta bins
+ {
+ correlation = fDiffFlowCorrelationsPro[t][pe][rci]->GetBinContent(b);
+ spread = fDiffFlowCorrelationsPro[t][pe][rci]->GetBinError(b);
+ sumOfWeights = fDiffFlowSumOfEventWeights[t][pe][0][rci]->GetBinContent(b);
+ sumOfSquaredWeights = fDiffFlowSumOfEventWeights[t][pe][1][rci]->GetBinContent(b);
+ if(sumOfWeights) termA = (pow(sumOfSquaredWeights,0.5)/sumOfWeights);
+ if(1.-pow(termA,2.)>0.) termB = 1./pow(1.-pow(termA,2.),0.5);
+ error = termA*spread*termB; // final error (unbiased estimator for standard deviation)
+ fDiffFlowCorrelationsHist[t][pe][rci]->SetBinContent(b,correlation);
+ fDiffFlowCorrelationsHist[t][pe][rci]->SetBinError(b,error);
+ } // end of for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+ } // end of for(Int_t rci=0;rci<4;rci++)
+
+} // end of void AliFlowAnalysisWithQCumulants::FinalizeReducedCorrelations(TString type, TString ptOrEta)
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::CalculateDiffFlowProductOfCorrelations(TString type, TString ptOrEta)
+{
+ // store products: <2><2'>, <2><4'>, <2><6'>, <2><8'>, <2'><4>,
+ // <2'><4'>, <2'><6>, <2'><6'>, <2'><8>, <2'><8'>,
+ // <4><4'>, <4><6'>, <4><8'>, <4'><6>, <4'><6'>,
+ // <4'><8>, <4'><8'>, <6><6'>, <6><8'>, <6'><8>,
+ // <6'><8'>, <8><8'>.
+
+ Int_t typeFlag = -1;
+ Int_t ptEtaFlag = -1;
+
+ if(type == "RP")
+ {
+ typeFlag = 0;
+ } else if(type == "POI")
+ {
+ typeFlag = 1;
+ }
+
+ if(ptOrEta == "Pt")
+ {
+ ptEtaFlag = 0;
+ } else if(ptOrEta == "Eta")
+ {
+ ptEtaFlag = 1;
+ }
+
+ // shortcuts:
+ Int_t t = typeFlag;
+ Int_t pe = ptEtaFlag;
+
+ // common:
+ Int_t nBinsPtEta[2] = {fnBinsPt,fnBinsEta};
+ Double_t minPtEta[2] = {fPtMin,fEtaMin};
+ Double_t binWidthPtEta[2] = {fPtBinWidth,fEtaBinWidth};
+
+ // protections // to be improved (add protection for all pointers in this method)
+ if(!fIntFlowCorrelationsEBE)
+ {
+ cout<<"WARNING: fIntFlowCorrelationsEBE is NULL in AFAWQC::CDFPOC() !!!!"<<endl;
+ exit(0);
+ }
+
+ /*
+ Double_t dMult = (*fSMpk)(0,0); // multiplicity (number of particles used to determine the reaction plane)
+ //Double_t mr = 0.; // number of RPs in particular pt or eta bin
+ Double_t mp = 0.; // number of POIs in particular pt or eta bin
+ Double_t mq = 0.; // number of particles which are both RPs and POIs in particular pt or eta bin
+ */
+
+ // e-b-e correlations:
+ Double_t twoEBE = fIntFlowCorrelationsEBE->GetBinContent(1); // <2>
+ Double_t fourEBE = fIntFlowCorrelationsEBE->GetBinContent(2); // <4>
+ Double_t sixEBE = fIntFlowCorrelationsEBE->GetBinContent(3); // <6>
+ Double_t eightEBE = fIntFlowCorrelationsEBE->GetBinContent(4); // <8>
+
+ // event weights for correlations:
+ Double_t dW2 = fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(1); // event weight for <2>
+ Double_t dW4 = fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(2); // event weight for <4>
+ Double_t dW6 = fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(3); // event weight for <6>
+ Double_t dW8 = fIntFlowEventWeightsForCorrelationsEBE->GetBinContent(4); // event weight for <8>
+
+ // e-b-e reduced correlations:
+ Double_t twoReducedEBE = 0.; // <2'>
+ Double_t fourReducedEBE = 0.; // <4'>
+ Double_t sixReducedEBE = 0.; // <6'>
+ Double_t eightReducedEBE = 0.; // <8'>
+
+ // event weights for reduced correlations:
+ Double_t dw2 = 0.; // event weight for <2'>
+ Double_t dw4 = 0.; // event weight for <4'>
+ //Double_t dw6 = 0.; // event weight for <6'>
+ //Double_t dw8 = 0.; // event weight for <8'>
+
+ // looping over bins:
+ for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+ {
+ // e-b-e reduced correlations:
+ twoReducedEBE = fDiffFlowCorrelationsEBE[t][pe][0]->GetBinContent(b);
+ fourReducedEBE = fDiffFlowCorrelationsEBE[t][pe][1]->GetBinContent(b);
+ sixReducedEBE = fDiffFlowCorrelationsEBE[t][pe][2]->GetBinContent(b);
+ eightReducedEBE = fDiffFlowCorrelationsEBE[t][pe][3]->GetBinContent(b);
+
+ /*
+ // to be improved (I should not do this here again)
+ if(type == "RP")
+ {
+ mq = fReRPQ1dEBE[0][pe][0][0]->GetBinEntries(b);
+ mp = mq; // trick to use the very same Eqs. bellow both for RP's and POI's diff. flow
+ } else if(type == "POI")
+ {
+ mp = fReRPQ1dEBE[1][pe][0][0]->GetBinEntries(b);
+ mq = fReRPQ1dEBE[2][pe][0][0]->GetBinEntries(b);
+ }
+
+ // event weights for reduced correlations:
+ dw2 = mp*dMult-mq; // weight for <2'>
+ dw4 = (mp-mq)*dMult*(dMult-1.)*(dMult-2.)
+ + mq*(dMult-1.)*(dMult-2.)*(dMult-3.); // weight for <4'>
+ //dw6 = ...
+ //dw8 = ...
+
+ */
+
+ dw2 = fDiffFlowEventWeightsForCorrelationsEBE[t][pe][0]->GetBinContent(b);
+ dw4 = fDiffFlowEventWeightsForCorrelationsEBE[t][pe][1]->GetBinContent(b);
+
+ // storing all products:
+ fDiffFlowProductOfCorrelationsPro[t][pe][0][1]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],twoEBE*twoReducedEBE,dW2*dw2); // storing <2><2'>
+ fDiffFlowProductOfCorrelationsPro[t][pe][1][2]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],fourEBE*twoReducedEBE,dW4*dw2); // storing <4><2'>
+ fDiffFlowProductOfCorrelationsPro[t][pe][1][4]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],sixEBE*twoReducedEBE,dW6*dw2); // storing <6><2'>
+ fDiffFlowProductOfCorrelationsPro[t][pe][1][6]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],eightEBE*twoReducedEBE,dW8*dw2); // storing <8><2'>
+
+ // event weight for <4'>:
+ fDiffFlowProductOfCorrelationsPro[t][pe][0][3]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],twoEBE*fourReducedEBE,dW2*dw4); // storing <2><4'>
+ fDiffFlowProductOfCorrelationsPro[t][pe][1][3]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],twoReducedEBE*fourReducedEBE,dw2*dw4); // storing <2'><4'>
+ fDiffFlowProductOfCorrelationsPro[t][pe][2][3]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],fourEBE*fourReducedEBE,dW4*dw4); // storing <4><4'>
+ fDiffFlowProductOfCorrelationsPro[t][pe][3][4]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],sixEBE*fourReducedEBE,dW6*dw4); // storing <6><4'>
+ fDiffFlowProductOfCorrelationsPro[t][pe][3][6]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],eightEBE*fourReducedEBE,dW8*dw4); // storing <8><4'>
+
+ // event weight for <6'>:
+ //dw6 = ...;
+ //fDiffFlowProductOfCorrelationsPro[t][pe][0][5]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],twoEBE*sixReducedEBE,dW2*dw6); // storing <2><6'>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][1][5]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],twoReducedEBE*sixReducedEBE,dw2*dw6); // storing <2'><6'>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][2][5]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],fourEBE*sixReducedEBE,dW4*dw6); // storing <4><6'>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][3][5]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],fourReducedEBE*sixReducedEBE,dw4*dw6); // storing <4'><6'>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][4][5]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],sixEBE*sixReducedEBE,dW6*dw6); // storing <6><6'>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][5][6]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],sixReducedEBE*eightEBE,dw6*dW8); // storing <6'><8>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][5][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],sixReducedEBE*eightReducedEBE,dw6*dw8); // storing <6'><8'>
+
+ // event weight for <8'>:
+ //dw8 = ...;
+ //fDiffFlowProductOfCorrelationsPro[t][pe][0][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],twoEBE*eightReducedEBE,dW2*dw8); // storing <2><8'>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][1][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],twoReducedEBE*eightReducedEBE,dw2*dw8); // storing <2'><8'>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][2][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],fourEBE*eightReducedEBE,dW4*dw8); // storing <4><8'>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][3][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],fourReducedEBE*eightReducedEBE,dw4*dw8); // storing <4'><8'>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][4][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],sixEBE*eightReducedEBE,dW6*dw8); // storing <6><8'>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][5][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],sixReducedEBE*eightReducedEBE,dw6*dw8); // storing <6'><8'>
+ //fDiffFlowProductOfCorrelationsPro[t][pe][6][7]->Fill(minPtEta[pe]+(b-1)*binWidthPtEta[pe],eightEBE*eightReducedEBE,dW8*dw8); // storing <8><8'>
+ } // end of for(Int_t b=1;b<=nBinsPtEta[pe];b++
+
+} // end of void AliFlowAnalysisWithQCumulants::CalculateDiffFlowProductOfCorrelations(TString type, TString ptOrEta)
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::CalculateDiffFlowCovariances(TString type, TString ptOrEta) // to be improved (reimplemented)
+{
+ // a) Calculate unbiased estimators Cov(<2>,<2'>), Cov(<2>,<4'>), Cov(<4>,<2'>), Cov(<4>,<4'>) and Cov(<2'>,<4'>)
+ // for covariances V(<2>,<2'>), V(<2>,<4'>), V(<4>,<2'>), V(<4>,<4'>) and V(<2'>,<4'>).
+ // b) Store in histogram fDiffFlowCovariances[t][pe][index] for instance the following:
+ //
+ // Cov(<2>,<2'>) * (sum_{i=1}^{N} w_{<2>}_i w_{<2'>}_i )/[(sum_{i=1}^{N} w_{<2>}_i) * (sum_{j=1}^{N} w_{<2'>}_j)]
+ //
+ // where N is the number of events, w_{<2>} is event weight for <2> and w_{<2'>} is event weight for <2'>.
+ // c) Binning of fDiffFlowCovariances[t][pe][index] is organized as follows:
+ //
+ // 1st bin: Cov(<2>,<2'>) * (sum_{i=1}^{N} w_{<2>}_i w_{<2'>}_i )/[(sum_{i=1}^{N} w_{<2>}_i) * (sum_{j=1}^{N} w_{<2'>}_j)]
+ // 2nd bin: Cov(<2>,<4'>) * (sum_{i=1}^{N} w_{<2>}_i w_{<4'>}_i )/[(sum_{i=1}^{N} w_{<2>}_i) * (sum_{j=1}^{N} w_{<4'>}_j)]
+ // 3rd bin: Cov(<4>,<2'>) * (sum_{i=1}^{N} w_{<4>}_i w_{<2'>}_i )/[(sum_{i=1}^{N} w_{<4>}_i) * (sum_{j=1}^{N} w_{<2'>}_j)]
+ // 4th bin: Cov(<4>,<4'>) * (sum_{i=1}^{N} w_{<4>}_i w_{<4'>}_i )/[(sum_{i=1}^{N} w_{<4>}_i) * (sum_{j=1}^{N} w_{<4'>}_j)]
+ // 5th bin: Cov(<2'>,<4'>) * (sum_{i=1}^{N} w_{<2'>}_i w_{<4'>}_i )/[(sum_{i=1}^{N} w_{<2'>}_i) * (sum_{j=1}^{N} w_{<4'>}_j)]
+ // ...
+
+ Int_t typeFlag = -1;
+ Int_t ptEtaFlag = -1;
+
+ if(type == "RP")
+ {
+ typeFlag = 0;
+ } else if(type == "POI")
+ {
+ typeFlag = 1;
+ }
+
+ if(ptOrEta == "Pt")
+ {
+ ptEtaFlag = 0;
+ } else if(ptOrEta == "Eta")
+ {
+ ptEtaFlag = 1;
+ }
+
+ // shortcuts:
+ Int_t t = typeFlag;
+ Int_t pe = ptEtaFlag;
+
+ // common:
+ Int_t nBinsPtEta[2] = {fnBinsPt,fnBinsEta};
+ //Double_t minPtEta[2] = {fPtMin,fEtaMin};
+ //Double_t maxPtEta[2] = {fPtMax,fEtaMax};
+ //Double_t binWidthPtEta[2] = {fPtBinWidth,fEtaBinWidth};
+
+ // average correlations:
+ Double_t two = fIntFlowCorrelationsHist->GetBinContent(1); // <<2>>
+ Double_t four = fIntFlowCorrelationsHist->GetBinContent(2); // <<4>>
+ //Double_t six = fIntFlowCorrelationsHist->GetBinContent(3); // <<6>>
+ //Double_t eight = fIntFlowCorrelationsHist->GetBinContent(4); // <<8>>
+
+ // sum of weights for correlation:
+ Double_t sumOfWeightsForTwo = fIntFlowSumOfEventWeights[0]->GetBinContent(1); // sum_{i=1}^{N} w_{<2>}
+ Double_t sumOfWeightsForFour = fIntFlowSumOfEventWeights[0]->GetBinContent(2); // sum_{i=1}^{N} w_{<4>}
+ //Double_t sumOfWeightsForSix = fIntFlowSumOfEventWeights[0]->GetBinContent(3); // sum_{i=1}^{N} w_{<6>}
+ //Double_t sumOfWeightsForEight = fIntFlowSumOfEventWeights[0]->GetBinContent(4); // sum_{i=1}^{N} w_{<8>}
+
+ // average reduced correlations:
+ Double_t twoReduced = 0.; // <<2'>>
+ Double_t fourReduced = 0.; // <<4'>>
+ //Double_t sixReduced = 0.; // <<6'>>
+ //Double_t eightReduced = 0.; // <<8'>>
+
+ // sum of weights for reduced correlation:
+ Double_t sumOfWeightsForTwoReduced = 0.; // sum_{i=1}^{N} w_{<2'>}
+ Double_t sumOfWeightsForFourReduced = 0.; // sum_{i=1}^{N} w_{<4'>}
+ //Double_t sumOfWeightsForSixReduced = 0.; // sum_{i=1}^{N} w_{<6'>}
+ //Double_t sumOfWeightsForEightReduced = 0.; // sum_{i=1}^{N} w_{<8'>}
+
+ // product of weights for reduced correlation:
+ Double_t productOfWeightsForTwoTwoReduced = 0.; // sum_{i=1}^{N} w_{<2>}w_{<2'>}
+ Double_t productOfWeightsForTwoFourReduced = 0.; // sum_{i=1}^{N} w_{<2>}w_{<4'>}
+ Double_t productOfWeightsForFourTwoReduced = 0.; // sum_{i=1}^{N} w_{<4>}w_{<2'>}
+ Double_t productOfWeightsForFourFourReduced = 0.; // sum_{i=1}^{N} w_{<4>}w_{<4'>}
+ Double_t productOfWeightsForTwoReducedFourReduced = 0.; // sum_{i=1}^{N} w_{<2'>}w_{<4'>}
+ // ...
+
+ // products for differential flow:
+ Double_t twoTwoReduced = 0; // <<2><2'>>
+ Double_t twoFourReduced = 0; // <<2><4'>>
+ Double_t fourTwoReduced = 0; // <<4><2'>>
+ Double_t fourFourReduced = 0; // <<4><4'>>
+ Double_t twoReducedFourReduced = 0; // <<2'><4'>>
+
+ // denominators in the expressions for the unbiased estimators for covariances:
+ // denominator = 1 - term1/(term2*term3)
+ // prefactor = term1/(term2*term3)
+ Double_t denominator = 0.;
+ Double_t prefactor = 0.;
+ Double_t term1 = 0.;
+ Double_t term2 = 0.;
+ Double_t term3 = 0.;
+
+ // unbiased estimators for covariances for differential flow:
+ Double_t covTwoTwoReduced = 0.; // Cov(<2>,<2'>)
+ Double_t wCovTwoTwoReduced = 0.; // Cov(<2>,<2'>) * prefactor(w_{<2>},w_{<2'>})
+ Double_t covTwoFourReduced = 0.; // Cov(<2>,<4'>)
+ Double_t wCovTwoFourReduced = 0.; // Cov(<2>,<4'>) * prefactor(w_{<2>},w_{<4'>})
+ Double_t covFourTwoReduced = 0.; // Cov(<4>,<2'>)
+ Double_t wCovFourTwoReduced = 0.; // Cov(<4>,<2'>) * prefactor(w_{<4>},w_{<2'>})
+ Double_t covFourFourReduced = 0.; // Cov(<4>,<4'>)
+ Double_t wCovFourFourReduced = 0.; // Cov(<4>,<4'>) * prefactor(w_{<4>},w_{<4'>})
+ Double_t covTwoReducedFourReduced = 0.; // Cov(<2'>,<4'>)
+ Double_t wCovTwoReducedFourReduced = 0.; // Cov(<2'>,<4'>) * prefactor(w_{<2'>},w_{<4'>})
+
+ for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+ {
+ // average reduced corelations:
+ twoReduced = fDiffFlowCorrelationsHist[t][pe][0]->GetBinContent(b);
+ fourReduced = fDiffFlowCorrelationsHist[t][pe][1]->GetBinContent(b);
+ // average products:
+ twoTwoReduced = fDiffFlowProductOfCorrelationsPro[t][pe][0][1]->GetBinContent(b);
+ twoFourReduced = fDiffFlowProductOfCorrelationsPro[t][pe][0][3]->GetBinContent(b);
+ fourTwoReduced = fDiffFlowProductOfCorrelationsPro[t][pe][1][2]->GetBinContent(b);
+ fourFourReduced = fDiffFlowProductOfCorrelationsPro[t][pe][2][3]->GetBinContent(b);
+ twoReducedFourReduced = fDiffFlowProductOfCorrelationsPro[t][pe][1][3]->GetBinContent(b);
+ // sum of weights for reduced correlations:
+ sumOfWeightsForTwoReduced = fDiffFlowSumOfEventWeights[t][pe][0][0]->GetBinContent(b);
+ sumOfWeightsForFourReduced = fDiffFlowSumOfEventWeights[t][pe][0][1]->GetBinContent(b);
+ // products of weights for correlations:
+ productOfWeightsForTwoTwoReduced = fDiffFlowSumOfProductOfEventWeights[t][pe][0][1]->GetBinContent(b);
+ productOfWeightsForTwoFourReduced = fDiffFlowSumOfProductOfEventWeights[t][pe][0][3]->GetBinContent(b);
+ productOfWeightsForFourTwoReduced = fDiffFlowSumOfProductOfEventWeights[t][pe][1][2]->GetBinContent(b);
+ productOfWeightsForFourFourReduced = fDiffFlowSumOfProductOfEventWeights[t][pe][2][3]->GetBinContent(b);
+ productOfWeightsForTwoReducedFourReduced = fDiffFlowSumOfProductOfEventWeights[t][pe][1][3]->GetBinContent(b);
+ // denominator for the unbiased estimator for covariances: 1 - term1/(term2*term3)
+ // prefactor (multiplies Cov's) = term1/(term2*term3)
+ // <2>,<2'>:
+ term1 = productOfWeightsForTwoTwoReduced;
+ term2 = sumOfWeightsForTwo;
+ term3 = sumOfWeightsForTwoReduced;
+ if(term2*term3>0.)
+ {
+ denominator = 1.-term1/(term2*term3);
+ prefactor = term1/(term2*term3);
+ if(TMath::Abs(denominator)>1e-6)
+ {
+ covTwoTwoReduced = (twoTwoReduced-two*twoReduced)/denominator;
+ wCovTwoTwoReduced = covTwoTwoReduced*prefactor;
+ fDiffFlowCovariances[t][pe][0]->SetBinContent(b,wCovTwoTwoReduced);
+ }
+ }
+ // <2>,<4'>:
+ term1 = productOfWeightsForTwoFourReduced;
+ term2 = sumOfWeightsForTwo;
+ term3 = sumOfWeightsForFourReduced;
+ if(term2*term3>0.)
+ {
+ denominator = 1.-term1/(term2*term3);
+ prefactor = term1/(term2*term3);
+ if(TMath::Abs(denominator)>1e-6)
+ {
+ covTwoFourReduced = (twoFourReduced-two*fourReduced)/denominator;
+ wCovTwoFourReduced = covTwoFourReduced*prefactor;
+ fDiffFlowCovariances[t][pe][1]->SetBinContent(b,wCovTwoFourReduced);
+ }
+ }
+ // <4>,<2'>:
+ term1 = productOfWeightsForFourTwoReduced;
+ term2 = sumOfWeightsForFour;
+ term3 = sumOfWeightsForTwoReduced;
+ if(term2*term3>0.)
+ {
+ denominator = 1.-term1/(term2*term3);
+ prefactor = term1/(term2*term3);
+ if(TMath::Abs(denominator)>1e-6)
+ {
+ covFourTwoReduced = (fourTwoReduced-four*twoReduced)/denominator;
+ wCovFourTwoReduced = covFourTwoReduced*prefactor;
+ fDiffFlowCovariances[t][pe][2]->SetBinContent(b,wCovFourTwoReduced);
+ }
+ }
+ // <4>,<4'>:
+ term1 = productOfWeightsForFourFourReduced;
+ term2 = sumOfWeightsForFour;
+ term3 = sumOfWeightsForFourReduced;
+ if(term2*term3>0.)
+ {
+ denominator = 1.-term1/(term2*term3);
+ prefactor = term1/(term2*term3);
+ if(TMath::Abs(denominator)>1e-6)
+ {
+ covFourFourReduced = (fourFourReduced-four*fourReduced)/denominator;
+ wCovFourFourReduced = covFourFourReduced*prefactor;
+ fDiffFlowCovariances[t][pe][3]->SetBinContent(b,wCovFourFourReduced);
+ }
+ }
+ // <2'>,<4'>:
+ term1 = productOfWeightsForTwoReducedFourReduced;
+ term2 = sumOfWeightsForTwoReduced;
+ term3 = sumOfWeightsForFourReduced;
+ if(term2*term3>0.)
+ {
+ denominator = 1.-term1/(term2*term3);
+ prefactor = term1/(term2*term3);
+ if(TMath::Abs(denominator)>1e-6)
+ {
+ covTwoReducedFourReduced = (twoReducedFourReduced-twoReduced*fourReduced)/denominator;
+ wCovTwoReducedFourReduced = covTwoReducedFourReduced*prefactor;
+ fDiffFlowCovariances[t][pe][4]->SetBinContent(b,wCovTwoReducedFourReduced);
+ }
+ }
+ } // end of for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+
+} // end of void AliFlowAnalysisWithQCumulants::CalculateDiffFlowCovariances(TString type, TString ptOrEta)
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::CalculateDiffFlow(TString type, TString ptOrEta)
+{
+ // calculate differential flow from differential cumulants and previously obtained integrated flow: (to be improved: description)
+
+ Int_t typeFlag = -1;
+ Int_t ptEtaFlag = -1;
+
+ if(type == "RP")
+ {
+ typeFlag = 0;
+ } else if(type == "POI")
+ {
+ typeFlag = 1;
+ }
+
+ if(ptOrEta == "Pt")
+ {
+ ptEtaFlag = 0;
+ } else if(ptOrEta == "Eta")
+ {
+ ptEtaFlag = 1;
+ }
+
+ // shortcuts:
+ Int_t t = typeFlag;
+ Int_t pe = ptEtaFlag;
+
+ // common:
+ Int_t nBinsPtEta[2] = {fnBinsPt,fnBinsEta};
+
+ // correlations:
+ Double_t two = fIntFlowCorrelationsHist->GetBinContent(1); // <<2>>
+ Double_t four = fIntFlowCorrelationsHist->GetBinContent(2); // <<4>>
+
+ // statistical errors of correlations:
+ Double_t twoError = fIntFlowCorrelationsHist->GetBinError(1);
+ Double_t fourError = fIntFlowCorrelationsHist->GetBinError(2);
+
+ // reduced correlations:
+ Double_t twoReduced = 0.; // <<2'>>
+ Double_t fourReduced = 0.; // <<4'>>
+
+ // statistical errors of reduced correlations:
+ Double_t twoReducedError = 0.;
+ Double_t fourReducedError = 0.;
+
+ // covariances:
+ Double_t wCovTwoFour = fIntFlowCovariances->GetBinContent(1);// // Cov(<2>,<4>) * prefactor(<2>,<4>)
+ Double_t wCovTwoTwoReduced = 0.; // Cov(<2>,<2'>) * prefactor(<2>,<2'>)
+ Double_t wCovTwoFourReduced = 0.; // Cov(<2>,<4'>) * prefactor(<2>,<4'>)
+ Double_t wCovFourTwoReduced = 0.; // Cov(<4>,<2'>) * prefactor(<4>,<2'>)
+ Double_t wCovFourFourReduced = 0.; // Cov(<4>,<4'>) * prefactor(<4>,<4'>)
+ Double_t wCovTwoReducedFourReduced = 0.; // Cov(<2'>,<4'>) * prefactor(<2'>,<4'>)
+
+ // differential flow:
+ Double_t v2Prime = 0.; // v'{2}
+ Double_t v4Prime = 0.; // v'{4}
+
+ // statistical error of differential flow:
+ Double_t v2PrimeError = 0.;
+ Double_t v4PrimeError = 0.;
+
+ // squared statistical error of differential flow:
+ Double_t v2PrimeErrorSquared = 0.;
+ Double_t v4PrimeErrorSquared = 0.;
+
+ // loop over pt or eta bins:
+ for(Int_t b=1;b<=nBinsPtEta[pe];b++)
+ {
+ // reduced correlations and statistical errors:
+ twoReduced = fDiffFlowCorrelationsHist[t][pe][0]->GetBinContent(b);
+ twoReducedError = fDiffFlowCorrelationsHist[t][pe][0]->GetBinError(b);
+ fourReduced = fDiffFlowCorrelationsHist[t][pe][1]->GetBinContent(b);
+ fourReducedError = fDiffFlowCorrelationsHist[t][pe][1]->GetBinError(b);
+ // covariances:
+ wCovTwoTwoReduced = fDiffFlowCovariances[t][pe][0]->GetBinContent(b);
+ wCovTwoFourReduced = fDiffFlowCovariances[t][pe][1]->GetBinContent(b);
+ wCovFourTwoReduced = fDiffFlowCovariances[t][pe][2]->GetBinContent(b);
+ wCovFourFourReduced = fDiffFlowCovariances[t][pe][3]->GetBinContent(b);
+ wCovTwoReducedFourReduced = fDiffFlowCovariances[t][pe][4]->GetBinContent(b);
+ // differential flow:
+ // v'{2}:
+ if(two>0.)
+ {
+ v2Prime = twoReduced/pow(two,0.5);
+ v2PrimeErrorSquared = (1./4.)*pow(two,-3.)*
+ (pow(twoReduced,2.)*pow(twoError,2.)
+ + 4.*pow(two,2.)*pow(twoReducedError,2.)
+ - 4.*two*twoReduced*wCovTwoTwoReduced);
+
+
+ if(v2PrimeErrorSquared>0.) v2PrimeError = pow(v2PrimeErrorSquared,0.5);
+ fDiffFlow[t][pe][0]->SetBinContent(b,v2Prime);
+ if(TMath::Abs(v2Prime)>1.e-44)fDiffFlow[t][pe][0]->SetBinError(b,v2PrimeError);
+ }
+ // differential flow:
+ // v'{4}
+ if(2.*pow(two,2.)-four > 0.)
+ {
+ v4Prime = (2.*two*twoReduced-fourReduced)/pow(2.*pow(two,2.)-four,3./4.);
+ v4PrimeErrorSquared = pow(2.*pow(two,2.)-four,-7./2.)*
+ (pow(2.*pow(two,2.)*twoReduced-3.*two*fourReduced+2.*four*twoReduced,2.)*pow(twoError,2.)
+ + (9./16.)*pow(2.*two*twoReduced-fourReduced,2.)*pow(fourError,2.)
+ + 4.*pow(two,2.)*pow(2.*pow(two,2.)-four,2.)*pow(twoReducedError,2.)
+ + pow(2.*pow(two,2.)-four,2.)*pow(fourReducedError,2.)
+ - (3./2.)*(2.*two*twoReduced-fourReduced)
+ * (2.*pow(two,2.)*twoReduced-3.*two*fourReduced+2.*four*twoReduced)*wCovTwoFour
+ - 4.*two*(2.*pow(two,2.)-four)
+ * (2.*pow(two,2.)*twoReduced-3.*two*fourReduced+2.*four*twoReduced)*wCovTwoTwoReduced
+ + 2.*(2.*pow(two,2.)-four)
+ * (2.*pow(two,2.)*twoReduced-3.*two*fourReduced+2.*four*twoReduced)*wCovTwoFourReduced
+ + 3.*two*(2.*pow(two,2.)-four)*(2.*two*twoReduced-fourReduced)*wCovFourTwoReduced
+ - (3./2.)*(2.*pow(two,2.)-four)*(2.*two*twoReduced-fourReduced)*wCovFourFourReduced
+ - 4.*two*pow(2.*pow(two,2.)-four,2.)*wCovTwoReducedFourReduced);
+ if(v4PrimeErrorSquared>0.) v4PrimeError = pow(v4PrimeErrorSquared,0.5);
+ fDiffFlow[t][pe][1]->SetBinContent(b,v4Prime);
+ if(TMath::Abs(v4Prime)>1.e-44)fDiffFlow[t][pe][1]->SetBinError(b,v4PrimeError);
+ }
+
+ } // end of for(Int_t b=1;b<=fnBinsPtEta[pe];b++)
+
+
+
+
+ /*
+ // 2D:
+ for(Int_t nua=0;nua<2;nua++)
+ {
+ for(Int_t p=1;p<=fnBinsPt;p++)
+ {
+ for(Int_t e=1;e<=fnBinsEta;e++)
+ {
+ // differential cumulants:
+ Double_t qc2Prime = fFinalCumulants2D[t][pW][eW][nua][0]->GetBinContent(fFinalCumulants2D[t][pW][eW][nua][0]->GetBin(p,e)); // QC{2'}
+ Double_t qc4Prime = fFinalCumulants2D[t][pW][eW][nua][1]->GetBinContent(fFinalCumulants2D[t][pW][eW][nua][1]->GetBin(p,e)); // QC{4'}
+ // differential flow:
+ Double_t v2Prime = 0.;
+ Double_t v4Prime = 0.;
+ if(v2)
+ {
+ v2Prime = qc2Prime/v2;
+ fFinalFlow2D[t][pW][eW][nua][0]->SetBinContent(fFinalFlow2D[t][pW][eW][nua][0]->GetBin(p,e),v2Prime);
+ }
+ if(v4)
+ {
+ v4Prime = -qc4Prime/pow(v4,3.);
+ fFinalFlow2D[t][pW][eW][nua][1]->SetBinContent(fFinalFlow2D[t][pW][eW][nua][1]->GetBin(p,e),v4Prime);
+ }
+ } // end of for(Int_t e=1;e<=fnBinsEta;e++)
+ } // end of for(Int_t p=1;p<=fnBinsPt;p++)
+ } // end of for(Int_t nua=0;nua<2;nua++)
+ */
+
+} // end of AliFlowAnalysisWithQCumulants::CalculateDiffFlow(TString type, Bool_t useParticleWeights)
+
+//================================================================================================================================
+
+void AliFlowAnalysisWithQCumulants::StoreIntFlowFlags()
+{
+ // a) Store all flags for integrated flow in profile fIntFlowFlags.
+
+ if(!fIntFlowFlags)
+ {
+ cout<<"WARNING: fIntFlowFlags is NULL in AFAWQC::SFFIF() !!!!"<<endl;
+ exit(0);
+ }
+
+ // particle weights used or not:
+ fIntFlowFlags->Fill(0.5,(Int_t)fUsePhiWeights||fUsePtWeights||fUseEtaWeights);
+ // which event weights were used:
+ if(strcmp(fMultiplicityWeight->Data(),"combinations"))
+ {
+ fIntFlowFlags->Fill(1.5,0); // 0 = "combinations" (default)
+ } else if(strcmp(fMultiplicityWeight->Data(),"unit"))
+ {
+ fIntFlowFlags->Fill(1.5,1); // 1 = "unit"
+ } else if(strcmp(fMultiplicityWeight->Data(),"multiplicity"))
+ {
+ fIntFlowFlags->Fill(1.5,2); // 2 = "multiplicity"
+ }
+ // corrected for non-uniform acceptance or not:
+ fIntFlowFlags->Fill(2.5,(Int_t)fApplyCorrectionForNUA);
+ fIntFlowFlags->Fill(3.5,(Int_t)fPrintFinalResults[0]);
+ fIntFlowFlags->Fill(4.5,(Int_t)fPrintFinalResults[1]);
+ fIntFlowFlags->Fill(5.5,(Int_t)fPrintFinalResults[2]);
+ fIntFlowFlags->Fill(6.5,(Int_t)fPrintFinalResults[3]);
+ fIntFlowFlags->Fill(7.5,(Int_t)fApplyCorrectionForNUAVsM);
+ fIntFlowFlags->Fill(8.5,(Int_t)fPropagateErrorAlsoFromNUATerms);
+ fIntFlowFlags->Fill(9.5,(Int_t)fCalculateCumulantsVsM);
+ fIntFlowFlags->Fill(10.5,(Int_t)fMinimumBiasReferenceFlow);
+
+} // end of void AliFlowAnalysisWithQCumulants::StoreIntFlowFlags()
+
+//================================================================================================================================
+
+void AliFlowAnalysisWithQCumulants::StoreDiffFlowFlags()
+{
+ // Store all flags for differential flow in the profile fDiffFlowFlags.
+
+ if(!fDiffFlowFlags)
+ {
+ cout<<"WARNING: fDiffFlowFlags is NULL in AFAWQC::SFFDF() !!!!"<<endl;
+ exit(0);
+ }
+
+ fDiffFlowFlags->Fill(0.5,fUsePhiWeights||fUsePtWeights||fUseEtaWeights); // particle weights used or not
+ //fDiffFlowFlags->Fill(1.5,""); // which event weight was used? // to be improved
+ fDiffFlowFlags->Fill(2.5,fApplyCorrectionForNUA); // corrected for non-uniform acceptance or not
+ fDiffFlowFlags->Fill(3.5,fCalculate2DFlow); // calculate also 2D differential flow in (pt,eta) or not
+
+} // end of void AliFlowAnalysisWithQCumulants::StoreDiffFlowFlags()
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::GetPointersForCommonHistograms()
+{
+ // Access all pointers to common control and common result histograms and profiles.
+
+ TString commonHistsName = "AliFlowCommonHistQC";
+ commonHistsName += fAnalysisLabel->Data();
+ AliFlowCommonHist *commonHist = dynamic_cast<AliFlowCommonHist*>(fHistList->FindObject(commonHistsName.Data()));
+ if(commonHist) this->SetCommonHists(commonHist);
+ TString commonHists2ndOrderName = "AliFlowCommonHist2ndOrderQC";
+ commonHists2ndOrderName += fAnalysisLabel->Data();
+ AliFlowCommonHist *commonHist2nd = dynamic_cast<AliFlowCommonHist*>(fHistList->FindObject(commonHists2ndOrderName.Data()));
+ if(commonHist2nd) this->SetCommonHists2nd(commonHist2nd);
+ TString commonHists4thOrderName = "AliFlowCommonHist4thOrderQC";
+ commonHists4thOrderName += fAnalysisLabel->Data();
+ AliFlowCommonHist *commonHist4th = dynamic_cast<AliFlowCommonHist*>(fHistList->FindObject(commonHists4thOrderName.Data()));
+ if(commonHist4th) this->SetCommonHists4th(commonHist4th);
+ TString commonHists6thOrderName = "AliFlowCommonHist6thOrderQC";
+ commonHists6thOrderName += fAnalysisLabel->Data();
+ AliFlowCommonHist *commonHist6th = dynamic_cast<AliFlowCommonHist*>(fHistList->FindObject(commonHists6thOrderName.Data()));
+ if(commonHist6th) this->SetCommonHists6th(commonHist6th);
+ TString commonHists8thOrderName = "AliFlowCommonHist8thOrderQC";
+ commonHists8thOrderName += fAnalysisLabel->Data();
+ AliFlowCommonHist *commonHist8th = dynamic_cast<AliFlowCommonHist*>(fHistList->FindObject(commonHists8thOrderName.Data()));
+ if(commonHist8th) this->SetCommonHists8th(commonHist8th);
+ TString commonHistResults2ndOrderName = "AliFlowCommonHistResults2ndOrderQC";
+ commonHistResults2ndOrderName += fAnalysisLabel->Data();
+ AliFlowCommonHistResults *commonHistRes2nd = dynamic_cast<AliFlowCommonHistResults*> (fHistList->FindObject(commonHistResults2ndOrderName.Data()));
+ if(commonHistRes2nd) this->SetCommonHistsResults2nd(commonHistRes2nd);
+ TString commonHistResults4thOrderName = "AliFlowCommonHistResults4thOrderQC";
+ commonHistResults4thOrderName += fAnalysisLabel->Data();
+ AliFlowCommonHistResults *commonHistRes4th = dynamic_cast<AliFlowCommonHistResults*>
+ (fHistList->FindObject(commonHistResults4thOrderName.Data()));
+ if(commonHistRes4th) this->SetCommonHistsResults4th(commonHistRes4th);
+ TString commonHistResults6thOrderName = "AliFlowCommonHistResults6thOrderQC";
+ commonHistResults6thOrderName += fAnalysisLabel->Data();
+ AliFlowCommonHistResults *commonHistRes6th = dynamic_cast<AliFlowCommonHistResults*>
+ (fHistList->FindObject(commonHistResults6thOrderName.Data()));
+ if(commonHistRes6th) this->SetCommonHistsResults6th(commonHistRes6th);
+ TString commonHistResults8thOrderName = "AliFlowCommonHistResults8thOrderQC";
+ commonHistResults8thOrderName += fAnalysisLabel->Data();
+ AliFlowCommonHistResults *commonHistRes8th = dynamic_cast<AliFlowCommonHistResults*>
+ (fHistList->FindObject(commonHistResults8thOrderName.Data()));
+ if(commonHistRes8th) this->SetCommonHistsResults8th(commonHistRes8th);
+
+} // end of void AliFlowAnalysisWithQCumulants::GetPointersForCommonHistograms()
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::GetPointersForParticleWeightsHistograms()
+{
+ // Get pointers for histograms with particle weights.
+
+ TList *weightsList = dynamic_cast<TList*>(fHistList->FindObject("Weights"));
+ if(weightsList) this->SetWeightsList(weightsList);
+ TString fUseParticleWeightsName = "fUseParticleWeightsQC"; // to be improved (hirdwired label QC)
+ fUseParticleWeightsName += fAnalysisLabel->Data();
+ TProfile *useParticleWeights = dynamic_cast<TProfile*>(weightsList->FindObject(fUseParticleWeightsName.Data()));
+ if(useParticleWeights)
+ {
+ this->SetUseParticleWeights(useParticleWeights);
+ fUsePhiWeights = (Int_t)fUseParticleWeights->GetBinContent(1);
+ fUsePtWeights = (Int_t)fUseParticleWeights->GetBinContent(2);
+ fUseEtaWeights = (Int_t)fUseParticleWeights->GetBinContent(3);
+ }
+} // end of void AliFlowAnalysisWithQCumulants::GetPointersForParticleWeightsHistograms();
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::GetPointersForIntFlowHistograms()
+{
+ // Get pointers for histograms and profiles relevant for integrated flow:
+ // a) Get pointer to base list for integrated flow holding profile fIntFlowFlags and lists fIntFlowProfiles and fIntFlowResults.
+ // b) Get pointer to profile fIntFlowFlags holding all flags for integrated flow.
+ // c) Get pointer to list fIntFlowProfiles and pointers to all objects that she holds.
+ // d) Get pointer to list fIntFlowResults and pointers to all objects that she holds.
+
+ TString sinCosFlag[2] = {"sin","cos"}; // to be improved (should I promote this to data member?)
+ TString powerFlag[2] = {"linear","quadratic"}; // to be improved (should I promote this to data member?)
+ TString correlationFlag[4] = {"<<2>>","<<4>>","<<6>>","<<8>>"}; // to be improved (should I promote this to data member?)
+
+ // a) Get pointer to base list for integrated flow holding profile fIntFlowFlags and lists fIntFlowProfiles and fIntFlowResults:
+ TList *intFlowList = NULL;
+ intFlowList = dynamic_cast<TList*>(fHistList->FindObject("Integrated Flow"));
+ if(!intFlowList)
+ {
+ cout<<"WARNING: intFlowList is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ exit(0);
+ }
+
+ // b) Get pointer to profile fIntFlowFlags holding all flags for integrated flow:
+ TString intFlowFlagsName = "fIntFlowFlags";
+ intFlowFlagsName += fAnalysisLabel->Data();
+ TProfile *intFlowFlags = dynamic_cast<TProfile*>(intFlowList->FindObject(intFlowFlagsName.Data()));
+ if(intFlowFlags)
+ {
+ this->SetIntFlowFlags(intFlowFlags);
+ fApplyCorrectionForNUA = (Bool_t)intFlowFlags->GetBinContent(3);
+ fApplyCorrectionForNUAVsM = (Bool_t)intFlowFlags->GetBinContent(8);
+ fCalculateCumulantsVsM = (Bool_t)intFlowFlags->GetBinContent(10);
+ } else
+ {
+ cout<<"WARNING: intFlowFlags is NULL in FAWQC::GPFIFH() !!!!"<<endl;
+ }
+
+ // c) Get pointer to list fIntFlowProfiles and pointers to all objects that she holds:
+ TList *intFlowProfiles = NULL;
+ intFlowProfiles = dynamic_cast<TList*>(intFlowList->FindObject("Profiles"));
+ if(intFlowProfiles)
+ {
+ // average multiplicities:
+ TString avMultiplicityName = "fAvMultiplicity";
+ avMultiplicityName += fAnalysisLabel->Data();
+ TProfile *avMultiplicity = dynamic_cast<TProfile*>(intFlowProfiles->FindObject(avMultiplicityName.Data()));
+ if(avMultiplicity)
+ {
+ this->SetAvMultiplicity(avMultiplicity);
+ } else
+ {
+ cout<<"WARNING: avMultiplicity is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // average correlations <<2>>, <<4>>, <<6>> and <<8>> (with wrong errors!):
+ TString intFlowCorrelationsProName = "fIntFlowCorrelationsPro";
+ intFlowCorrelationsProName += fAnalysisLabel->Data();
+ TProfile *intFlowCorrelationsPro = dynamic_cast<TProfile*>(intFlowProfiles->FindObject(intFlowCorrelationsProName.Data()));
+ if(intFlowCorrelationsPro)
+ {
+ this->SetIntFlowCorrelationsPro(intFlowCorrelationsPro);
+ } else
+ {
+ cout<<"WARNING: intFlowCorrelationsPro is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // average correlations <<2>>, <<4>>, <<6>> and <<8>> versus multiplicity for all events (error is wrong here):
+ if(fCalculateCumulantsVsM)
+ {
+ TString intFlowCorrelationsVsMProName = "fIntFlowCorrelationsVsMPro";
+ intFlowCorrelationsVsMProName += fAnalysisLabel->Data();
+ for(Int_t ci=0;ci<4;ci++) // correlation index
+ {
+ TProfile *intFlowCorrelationsVsMPro = dynamic_cast<TProfile*>
+ (intFlowProfiles->FindObject(Form("%s, %s",intFlowCorrelationsVsMProName.Data(),correlationFlag[ci].Data())));
+ if(intFlowCorrelationsVsMPro)
+ {
+ this->SetIntFlowCorrelationsVsMPro(intFlowCorrelationsVsMPro,ci);
+ } else
+ {
+ cout<<"WARNING: "<<Form("intFlowCorrelationsVsMPro[%d]",ci)<<" is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of for(Int_t ci=0;ci<4;ci++) // correlation index
+ } // end of if(fCalculateCumulantsVsM)
+ // average all correlations for integrated flow (with wrong errors!):
+ TString intFlowCorrelationsAllProName = "fIntFlowCorrelationsAllPro";
+ intFlowCorrelationsAllProName += fAnalysisLabel->Data();
+ TProfile *intFlowCorrelationsAllPro = dynamic_cast<TProfile*>(intFlowProfiles->FindObject(intFlowCorrelationsAllProName.Data()));
+ if(intFlowCorrelationsAllPro)
+ {
+ this->SetIntFlowCorrelationsAllPro(intFlowCorrelationsAllPro);
+ } else
+ {
+ cout<<"WARNING: intFlowCorrelationsAllPro is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // average extra correlations for integrated flow (which appear only when particle weights are used):
+ // (to be improved: Weak point in implementation, I am assuming here that method GetPointersForParticleWeightsHistograms() was called)
+ if(fUsePhiWeights||fUsePtWeights||fUseEtaWeights)
+ {
+ TString intFlowExtraCorrelationsProName = "fIntFlowExtraCorrelationsPro";
+ intFlowExtraCorrelationsProName += fAnalysisLabel->Data();
+ TProfile *intFlowExtraCorrelationsPro = dynamic_cast<TProfile*>(intFlowProfiles->FindObject(intFlowExtraCorrelationsProName.Data()));
+ if(intFlowExtraCorrelationsPro)
+ {
+ this->SetIntFlowExtraCorrelationsPro(intFlowExtraCorrelationsPro);
+ } else
+ {
+ cout<<"WARNING: intFlowExtraCorrelationsPro is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of if(fUsePhiWeights||fUsePtWeights||fUseEtaWeights)
+ // average products of correlations <2>, <4>, <6> and <8>:
+ TString intFlowProductOfCorrelationsProName = "fIntFlowProductOfCorrelationsPro";
+ intFlowProductOfCorrelationsProName += fAnalysisLabel->Data();
+ TProfile *intFlowProductOfCorrelationsPro = dynamic_cast<TProfile*>(intFlowProfiles->FindObject(intFlowProductOfCorrelationsProName.Data()));
+ if(intFlowProductOfCorrelationsPro)
+ {
+ this->SetIntFlowProductOfCorrelationsPro(intFlowProductOfCorrelationsPro);
+ } else
+ {
+ cout<<"WARNING: intFlowProductOfCorrelationsPro is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // average product of correlations <2>, <4>, <6> and <8> versus multiplicity
+ // [0=<<2><4>>,1=<<2><6>>,2=<<2><8>>,3=<<4><6>>,4=<<4><8>>,5=<<6><8>>]
+ if(fCalculateCumulantsVsM)
+ {
+ TString intFlowProductOfCorrelationsVsMProName = "fIntFlowProductOfCorrelationsVsMPro";
+ intFlowProductOfCorrelationsVsMProName += fAnalysisLabel->Data();
+ TString productFlag[6] = {"<<2><4>>","<<2><6>>","<<2><8>>","<<4><6>>","<<4><8>>","<<6><8>>"};
+ for(Int_t pi=0;pi<6;pi++)
+ {
+ TProfile *intFlowProductOfCorrelationsVsMPro = dynamic_cast<TProfile*>(intFlowProfiles->FindObject(Form("%s, %s",intFlowProductOfCorrelationsVsMProName.Data(),productFlag[pi].Data())));
+ if(intFlowProductOfCorrelationsVsMPro)
+ {
+ this->SetIntFlowProductOfCorrelationsVsMPro(intFlowProductOfCorrelationsVsMPro,pi);
+ } else
+ {
+ cout<<"WARNING: "<<Form("intFlowProductOfCorrelationsVsMPro[%d]",pi)<<" is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of for(Int_t pi=0;pi<6;pi++)
+ } // end of if(fCalculateCumulantsVsM)
+ // average correction terms for non-uniform acceptance (with wrong errors!):
+ for(Int_t sc=0;sc<2;sc++)
+ {
+ TString intFlowCorrectionTermsForNUAProName = "fIntFlowCorrectionTermsForNUAPro";
+ intFlowCorrectionTermsForNUAProName += fAnalysisLabel->Data();
+ TProfile *intFlowCorrectionTermsForNUAPro = dynamic_cast<TProfile*>(intFlowProfiles->FindObject((Form("%s: %s terms",intFlowCorrectionTermsForNUAProName.Data(),sinCosFlag[sc].Data()))));
+ if(intFlowCorrectionTermsForNUAPro)
+ {
+ this->SetIntFlowCorrectionTermsForNUAPro(intFlowCorrectionTermsForNUAPro,sc);
+ } else
+ {
+ cout<<"WARNING: intFlowCorrectionTermsForNUAPro is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ cout<<"sc = "<<sc<<endl;
+ }
+ // versus multiplicity:
+ if(fCalculateCumulantsVsM)
+ {
+ TString correctionTermFlag[4] = {"(n(phi1))","(n(phi1+phi2))","(n(phi1-phi2-phi3))","(n(2phi1-phi2))"}; // to be improved - hardwired 4
+ TString intFlowCorrectionTermsForNUAVsMProName = "fIntFlowCorrectionTermsForNUAVsMPro";
+ intFlowCorrectionTermsForNUAVsMProName += fAnalysisLabel->Data();
+ for(Int_t ci=0;ci<4;ci++) // correction term index (to be improved - hardwired 4)
+ {
+ TProfile *intFlowCorrectionTermsForNUAVsMPro = dynamic_cast<TProfile*>(intFlowProfiles->FindObject(Form("%s: #LT#LT%s%s#GT#GT",intFlowCorrectionTermsForNUAVsMProName.Data(),sinCosFlag[sc].Data(),correctionTermFlag[ci].Data())));
+ if(intFlowCorrectionTermsForNUAVsMPro)
+ {
+ this->SetIntFlowCorrectionTermsForNUAVsMPro(intFlowCorrectionTermsForNUAVsMPro,sc,ci);
+ } else
+ {
+ cout<<"WARNING: intFlowCorrectionTermsForNUAVsMPro is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ cout<<"sc = "<<sc<<endl;
+ cout<<"ci = "<<ci<<endl;
+ }
+ } // end of for(Int_t ci=0;ci<4;ci++) // correction term index (to be improved - hardwired 4)
+ } // end of if(fCalculateCumulantsVsM)
+ } // end of for(Int_t sc=0;sc<2;sc++)
+ // average products of correction terms for NUA:
+ TString intFlowProductOfCorrectionTermsForNUAProName = "fIntFlowProductOfCorrectionTermsForNUAPro";
+ intFlowProductOfCorrectionTermsForNUAProName += fAnalysisLabel->Data();
+ TProfile *intFlowProductOfCorrectionTermsForNUAPro = dynamic_cast<TProfile*>(intFlowProfiles->FindObject(intFlowProductOfCorrectionTermsForNUAProName.Data()));
+ if(intFlowProductOfCorrectionTermsForNUAPro)
+ {
+ this->SetIntFlowProductOfCorrectionTermsForNUAPro(intFlowProductOfCorrectionTermsForNUAPro);
+ } else
+ {
+ cout<<"WARNING: intFlowProductOfCorrectionTermsForNUAPro is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } else // to if(intFlowProfiles)
+ {
+ cout<<"WARNING: intFlowProfiles is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+
+ // d) Get pointer to list fIntFlowResults and pointers to all objects that she holds.
+ TList *intFlowResults = NULL;
+ intFlowResults = dynamic_cast<TList*>(intFlowList->FindObject("Results"));
+ if(intFlowResults)
+ {
+ // average correlations <<2>>, <<4>>, <<6>> and <<8>> (with correct errors!):
+ TString intFlowCorrelationsHistName = "fIntFlowCorrelationsHist";
+ intFlowCorrelationsHistName += fAnalysisLabel->Data();
+ TH1D *intFlowCorrelationsHist = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowCorrelationsHistName.Data()));
+ if(intFlowCorrelationsHist)
+ {
+ this->SetIntFlowCorrelationsHist(intFlowCorrelationsHist);
+ } else
+ {
+ cout<<"WARNING: intFlowCorrelationsHist is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // average correlations <<2>>, <<4>>, <<6>> and <<8>> (with correct errors!) vs M:
+ if(fCalculateCumulantsVsM)
+ {
+ TString intFlowCorrelationsVsMHistName = "fIntFlowCorrelationsVsMHist";
+ intFlowCorrelationsVsMHistName += fAnalysisLabel->Data();
+ for(Int_t ci=0;ci<4;ci++) // correlation index
+ {
+ TH1D *intFlowCorrelationsVsMHist = dynamic_cast<TH1D*>
+ (intFlowResults->FindObject(Form("%s, %s",intFlowCorrelationsVsMHistName.Data(),correlationFlag[ci].Data())));
+ if(intFlowCorrelationsVsMHist)
+ {
+ this->SetIntFlowCorrelationsVsMHist(intFlowCorrelationsVsMHist,ci);
+ } else
+ {
+ cout<<"WARNING: "<<Form("intFlowCorrelationsVsMHist[%d]",ci)<<" is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of for(Int_t ci=0;ci<4;ci++) // correlation index
+ } // end of if(fCalculateCumulantsVsM)
+ // average all correlations for integrated flow (with correct errors!):
+ TString intFlowCorrelationsAllHistName = "fIntFlowCorrelationsAllHist";
+ intFlowCorrelationsAllHistName += fAnalysisLabel->Data();
+ TH1D *intFlowCorrelationsAllHist = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowCorrelationsAllHistName.Data()));
+ if(intFlowCorrelationsAllHist)
+ {
+ this->SetIntFlowCorrelationsAllHist(intFlowCorrelationsAllHist);
+ } else
+ {
+ cout<<"WARNING: intFlowCorrelationsAllHist is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // average correction terms for non-uniform acceptance (with correct errors!):
+ TString intFlowCorrectionTermsForNUAHistName = "fIntFlowCorrectionTermsForNUAHist";
+ intFlowCorrectionTermsForNUAHistName += fAnalysisLabel->Data();
+ for(Int_t sc=0;sc<2;sc++)
+ {
+ TH1D *intFlowCorrectionTermsForNUAHist = dynamic_cast<TH1D*>(intFlowResults->FindObject((Form("%s: %s terms",intFlowCorrectionTermsForNUAHistName.Data(),sinCosFlag[sc].Data()))));
+ if(intFlowCorrectionTermsForNUAHist)
+ {
+ this->SetIntFlowCorrectionTermsForNUAHist(intFlowCorrectionTermsForNUAHist,sc);
+ } else
+ {
+ cout<<"WARNING: intFlowCorrectionTermsForNUAHist is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ cout<<"sc = "<<sc<<endl;
+ }
+ } // end of for(Int_t sc=0;sc<2;sc++)
+ // covariances (multiplied with weight dependent prefactor):
+ TString intFlowCovariancesName = "fIntFlowCovariances";
+ intFlowCovariancesName += fAnalysisLabel->Data();
+ TH1D *intFlowCovariances = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowCovariancesName.Data()));
+ if(intFlowCovariances)
+ {
+ this->SetIntFlowCovariances(intFlowCovariances);
+ } else
+ {
+ cout<<"WARNING: intFlowCovariances is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // sum of linear and quadratic event weights for <2>, <4>, <6> and <8>:
+ TString intFlowSumOfEventWeightsName = "fIntFlowSumOfEventWeights";
+ intFlowSumOfEventWeightsName += fAnalysisLabel->Data();
+ for(Int_t power=0;power<2;power++)
+ {
+ TH1D *intFlowSumOfEventWeights = dynamic_cast<TH1D*>(intFlowResults->FindObject(Form("%s: %s",intFlowSumOfEventWeightsName.Data(),powerFlag[power].Data())));
+ if(intFlowSumOfEventWeights)
+ {
+ this->SetIntFlowSumOfEventWeights(intFlowSumOfEventWeights,power);
+ } else
+ {
+ cout<<"WARNING: intFlowSumOfEventWeights is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ cout<<"power = "<<power<<endl;
+ }
+ } // end of for(Int_t power=0;power<2;power++)
+ // sum of products of event weights for correlations <2>, <4>, <6> and <8>:
+ TString intFlowSumOfProductOfEventWeightsName = "fIntFlowSumOfProductOfEventWeights";
+ intFlowSumOfProductOfEventWeightsName += fAnalysisLabel->Data();
+ TH1D *intFlowSumOfProductOfEventWeights = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowSumOfProductOfEventWeightsName.Data()));
+ if(intFlowSumOfProductOfEventWeights)
+ {
+ this->SetIntFlowSumOfProductOfEventWeights(intFlowSumOfProductOfEventWeights);
+ } else
+ {
+ cout<<"WARNING: intFlowSumOfProductOfEventWeights is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // final result for covariances of correlations (multiplied with weight dependent prefactor) versus M
+ // [0=Cov(2,4),1=Cov(2,6),2=Cov(2,8),3=Cov(4,6),4=Cov(4,8),5=Cov(6,8)]:
+ if(fCalculateCumulantsVsM)
+ {
+ TString intFlowCovariancesVsMName = "fIntFlowCovariancesVsM";
+ intFlowCovariancesVsMName += fAnalysisLabel->Data();
+ TString covarianceFlag[6] = {"Cov(<2>,<4>)","Cov(<2>,<6>)","Cov(<2>,<8>)","Cov(<4>,<6>)","Cov(<4>,<8>)","Cov(<6>,<8>)"};
+ for(Int_t ci=0;ci<6;ci++)
+ {
+ TH1D *intFlowCovariancesVsM = dynamic_cast<TH1D*>(intFlowResults->FindObject(Form("%s, %s",intFlowCovariancesVsMName.Data(),covarianceFlag[ci].Data())));
+ if(intFlowCovariancesVsM)
+ {
+ this->SetIntFlowCovariancesVsM(intFlowCovariancesVsM,ci);
+ } else
+ {
+ cout<<"WARNING: "<<Form("intFlowCovariancesVsM[%d]",ci)<<" is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of for(Int_t ci=0;ci<6;ci++)
+ } // end of if(fCalculateCumulantsVsM)
+ // sum of linear and quadratic event weights for <2>, <4>, <6> and <8> versus multiplicity
+ // [0=sum{w_{<2>}},1=sum{w_{<4>}},2=sum{w_{<6>}},3=sum{w_{<8>}}][0=linear 1,1=quadratic]:
+ if(fCalculateCumulantsVsM)
+ {
+ TString intFlowSumOfEventWeightsVsMName = "fIntFlowSumOfEventWeightsVsM";
+ intFlowSumOfEventWeightsVsMName += fAnalysisLabel->Data();
+ TString sumFlag[2][4] = {{"#sum_{i=1}^{N} w_{<2>}","#sum_{i=1}^{N} w_{<4>}","#sum_{i=1}^{N} w_{<6>}","#sum_{i=1}^{N} w_{<8>}"},
+ {"#sum_{i=1}^{N} w_{<2>}^{2}","#sum_{i=1}^{N} w_{<4>}^{2}","#sum_{i=1}^{N} w_{<6>}^{2}","#sum_{i=1}^{N} w_{<8>}^{2}"}};
+ for(Int_t si=0;si<4;si++)
+ {
+ for(Int_t power=0;power<2;power++)
+ {
+ TH1D *intFlowSumOfEventWeightsVsM = dynamic_cast<TH1D*>(intFlowResults->FindObject(Form("%s, %s",intFlowSumOfEventWeightsVsMName.Data(),sumFlag[power][si].Data())));
+ if(intFlowSumOfEventWeightsVsM)
+ {
+ this->SetIntFlowSumOfEventWeightsVsM(intFlowSumOfEventWeightsVsM,si,power);
+ } else
+ {
+ cout<<"WARNING: "<<Form("intFlowSumOfEventWeightsVsM[%d][%d]",si,power)<<" is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of for(Int_t power=0;power<2;power++)
+ } // end of for(Int_t si=0;si<4;si++)
+ } // end of if(fCalculateCumulantsVsM)
+ // sum of products of event weights for correlations <2>, <4>, <6> and <8> vs M
+ // [0=sum{w_{<2>}w_{<4>}},1=sum{w_{<2>}w_{<6>}},2=sum{w_{<2>}w_{<8>}},
+ // 3=sum{w_{<4>}w_{<6>}},4=sum{w_{<4>}w_{<8>}},5=sum{w_{<6>}w_{<8>}}]:
+ if(fCalculateCumulantsVsM)
+ {
+ TString intFlowSumOfProductOfEventWeightsVsMName = "fIntFlowSumOfProductOfEventWeightsVsM";
+ intFlowSumOfProductOfEventWeightsVsMName += fAnalysisLabel->Data();
+ TString sopowFlag[6] = {"#sum_{i=1}^{N} w_{<2>} w_{<4>}","#sum_{i=1}^{N} w_{<2>} w_{<6>}","#sum_{i=1}^{N} w_{<2>} w_{<8>}",
+ "#sum_{i=1}^{N} w_{<4>} w_{<6>}","#sum_{i=1}^{N} w_{<4>} w_{<8>}","#sum_{i=1}^{N} w_{<6>} w_{<8>}"};
+ for(Int_t pi=0;pi<6;pi++)
+ {
+ TH1D *intFlowSumOfProductOfEventWeightsVsM = dynamic_cast<TH1D*>(intFlowResults->FindObject(Form("%s, %s",intFlowSumOfProductOfEventWeightsVsMName.Data(),sopowFlag[pi].Data())));
+ if(intFlowSumOfProductOfEventWeightsVsM)
+ {
+ this->SetIntFlowSumOfProductOfEventWeightsVsM(intFlowSumOfProductOfEventWeightsVsM,pi);
+ } else
+ {
+ cout<<"WARNING: "<<Form("intFlowSumOfProductOfEventWeightsVsM[%d]",pi)<<" is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of for(Int_t pi=0;pi<6;pi++)
+ } // end of if(fCalculateCumulantsVsM)
+ // covariances for NUA (multiplied with weight dependent prefactor):
+ TString intFlowCovariancesNUAName = "fIntFlowCovariancesNUA";
+ intFlowCovariancesNUAName += fAnalysisLabel->Data();
+ TH1D *intFlowCovariancesNUA = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowCovariancesNUAName.Data()));
+ if(intFlowCovariancesNUA)
+ {
+ this->SetIntFlowCovariancesNUA(intFlowCovariancesNUA);
+ } else
+ {
+ cout<<"WARNING: intFlowCovariancesNUA is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // sum of linear and quadratic event weights NUA terms:
+ TString intFlowSumOfEventWeightsNUAName = "fIntFlowSumOfEventWeightsNUA";
+ intFlowSumOfEventWeightsNUAName += fAnalysisLabel->Data();
+ for(Int_t sc=0;sc<2;sc++)
+ {
+ for(Int_t power=0;power<2;power++)
+ {
+ TH1D *intFlowSumOfEventWeightsNUA = dynamic_cast<TH1D*>(intFlowResults->FindObject(Form("%s: %s, %s",intFlowSumOfEventWeightsNUAName.Data(),powerFlag[power].Data(),sinCosFlag[sc].Data())));
+ if(intFlowSumOfEventWeightsNUA)
+ {
+ this->SetIntFlowSumOfEventWeightsNUA(intFlowSumOfEventWeightsNUA,sc,power);
+ } else
+ {
+ cout<<"WARNING: intFlowSumOfEventWeightsNUA is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ cout<<"sc = "<<sc<<endl;
+ cout<<"power = "<<power<<endl;
+ }
+ } // end of for(Int_t power=0;power<2;power++)
+ } // end of for(Int_t sc=0;sc<2;sc++)
+ // sum of products of event weights for NUA terms:
+ TString intFlowSumOfProductOfEventWeightsNUAName = "fIntFlowSumOfProductOfEventWeightsNUA";
+ intFlowSumOfProductOfEventWeightsNUAName += fAnalysisLabel->Data();
+ TH1D *intFlowSumOfProductOfEventWeightsNUA = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowSumOfProductOfEventWeightsNUAName.Data()));
+ if(intFlowSumOfProductOfEventWeightsNUA)
+ {
+ this->SetIntFlowSumOfProductOfEventWeightsNUA(intFlowSumOfProductOfEventWeightsNUA);
+ } else
+ {
+ cout<<"WARNING: intFlowSumOfProductOfEventWeightsNUA is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // Final results for reference Q-cumulants:
+ TString intFlowQcumulantsName = "fIntFlowQcumulants";
+ intFlowQcumulantsName += fAnalysisLabel->Data();
+ TH1D *intFlowQcumulants = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowQcumulantsName.Data()));
+ if(intFlowQcumulants)
+ {
+ this->SetIntFlowQcumulants(intFlowQcumulants);
+ } else
+ {
+ cout<<"WARNING: intFlowQcumulants is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // Final results for reference Q-cumulants rebinned in M:
+ if(fCalculateCumulantsVsM)
+ {
+ TString intFlowQcumulantsRebinnedInMName = "fIntFlowQcumulantsRebinnedInM";
+ intFlowQcumulantsRebinnedInMName += fAnalysisLabel->Data();
+ TH1D *intFlowQcumulantsRebinnedInM = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowQcumulantsRebinnedInMName.Data()));
+ if(intFlowQcumulantsRebinnedInM)
+ {
+ this->SetIntFlowQcumulantsRebinnedInM(intFlowQcumulantsRebinnedInM);
+ } else
+ {
+ cout<<"WARNING: intFlowQcumulantsRebinnedInM is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of if(fCalculateCumulantsVsM)
+ // Ratio between error squared: with/without non-isotropic terms:
+ TString intFlowQcumulantsErrorSquaredRatioName = "fIntFlowQcumulantsErrorSquaredRatio";
+ intFlowQcumulantsErrorSquaredRatioName += fAnalysisLabel->Data();
+ TH1D *intFlowQcumulantsErrorSquaredRatio = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowQcumulantsErrorSquaredRatioName.Data()));
+ if(intFlowQcumulantsErrorSquaredRatio)
+ {
+ this->SetIntFlowQcumulantsErrorSquaredRatio(intFlowQcumulantsErrorSquaredRatio);
+ } else
+ {
+ cout<<" WARNING: intntFlowQcumulantsErrorSquaredRatio is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // final results for integrated Q-cumulants versus multiplicity:
+ TString cumulantFlag[4] = {"QC{2}","QC{4}","QC{6}","QC{8}"};
+ if(fCalculateCumulantsVsM)
+ {
+ TString intFlowQcumulantsVsMName = "fIntFlowQcumulantsVsM";
+ intFlowQcumulantsVsMName += fAnalysisLabel->Data();
+ for(Int_t co=0;co<4;co++) // cumulant order
+ {
+ TH1D *intFlowQcumulantsVsM = dynamic_cast<TH1D*>
+ (intFlowResults->FindObject(Form("%s, %s",intFlowQcumulantsVsMName.Data(),cumulantFlag[co].Data())));
+ if(intFlowQcumulantsVsM)
+ {
+ this->SetIntFlowQcumulantsVsM(intFlowQcumulantsVsM,co);
+ } else
+ {
+ cout<<"WARNING: "<<Form("intFlowQcumulantsVsM[%d]",co)<<" is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of for(Int_t co=0;co<4;co++) // cumulant order
+ } // end of if(fCalculateCumulantsVsM)
+ // Final reference flow estimates from Q-cumulants:
+ TString intFlowName = "fIntFlow";
+ intFlowName += fAnalysisLabel->Data();
+ TH1D *intFlow = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowName.Data()));
+ if(intFlow)
+ {
+ this->SetIntFlow(intFlow);
+ } else
+ {
+ cout<<"WARNING: intFlow is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // Final reference flow estimates from Q-cumulants vs M rebinned in M:
+ if(fCalculateCumulantsVsM)
+ {
+ TString intFlowRebinnedInMName = "fIntFlowRebinnedInM";
+ intFlowRebinnedInMName += fAnalysisLabel->Data();
+ TH1D *intFlowRebinnedInM = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowRebinnedInMName.Data()));
+ if(intFlowRebinnedInM)
+ {
+ this->SetIntFlowRebinnedInM(intFlowRebinnedInM);
+ } else
+ {
+ cout<<"WARNING: intFlowRebinnedInM is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of if(fCalculateCumulantsVsM)
+ // integrated flow from Q-cumulants versus multiplicity:
+ if(fCalculateCumulantsVsM)
+ {
+ TString intFlowVsMName = "fIntFlowVsM";
+ intFlowVsMName += fAnalysisLabel->Data();
+ TString flowFlag[4] = {"v_{2}{2,QC}","v_{2}{4,QC}","v_{2}{6,QC}","v_{2}{8,QC}"}; // to be improved (harwired harmonic)
+ for(Int_t co=0;co<4;co++) // cumulant order
+ {
+ TH1D *intFlowVsM = dynamic_cast<TH1D*>
+ (intFlowResults->FindObject(Form("%s, %s",intFlowVsMName.Data(),flowFlag[co].Data())));
+ if(intFlowVsM)
+ {
+ this->SetIntFlowVsM(intFlowVsM,co);
+ } else
+ {
+ cout<<"WARNING: "<<Form("intFlowVsM[%d]",co)<<" is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of for(Int_t co=0;co<4;co++) // cumulant order
+ } // end of if(fCalculateCumulantsVsM)
+ // quantifying detector effects effects to correlations:
+ TString intFlowDetectorBiasName = "fIntFlowDetectorBias";
+ intFlowDetectorBiasName += fAnalysisLabel->Data();
+ TH1D *intFlowDetectorBias = dynamic_cast<TH1D*>(intFlowResults->FindObject(intFlowDetectorBiasName.Data()));
+ if(intFlowDetectorBias)
+ {
+ this->SetIntFlowDetectorBias(intFlowDetectorBias);
+ } else
+ {
+ cout<<"WARNING: intFlowDetectorBias is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ // quantifying detector effects effects to correlations vs multiplicity:
+ if(fApplyCorrectionForNUAVsM && fCalculateCumulantsVsM)
+ {
+ TString intFlowDetectorBiasVsMName = "fIntFlowDetectorBiasVsM";
+ intFlowDetectorBiasVsMName += fAnalysisLabel->Data();
+ for(Int_t ci=0;ci<4;ci++) // correlation index
+ {
+ TH1D *intFlowDetectorBiasVsM = dynamic_cast<TH1D*>
+ (intFlowResults->FindObject(Form("%s for %s",intFlowDetectorBiasVsMName.Data(),cumulantFlag[ci].Data())));
+ if(intFlowDetectorBiasVsM)
+ {
+ this->SetIntFlowDetectorBiasVsM(intFlowDetectorBiasVsM,ci);
+ } else
+ {
+ cout<<"WARNING: "<<Form("intFlowDetectorBiasVsM[%d]",ci)<<" is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+ } // end of for(Int_t ci=0;ci<4;ci++) // correlation index
+ } // end of if(ApplyCorrectionForNUAVsM)
+ } else // to if(intFlowResults)
+ {
+ cout<<"WARNING: intFlowResults is NULL in AFAWQC::GPFIFH() !!!!"<<endl;
+ }
+
+} // end of void AliFlowAnalysisWithQCumulants::GetPointersForIntFlowHistograms()
+
+//================================================================================================================================
+
+void AliFlowAnalysisWithQCumulants::GetPointersForDiffFlowHistograms()
+{
+ // Get pointer to all objects relevant for differential flow.
+ // a) Define flags locally (to be improved: should I promote flags to data members?);
+ // b) Get pointer to base list for differential flow fDiffFlowList and nested lists fDiffFlowListProfiles and fDiffFlowListResults;
+ // c) Get pointer to profile fDiffFlowFlags holding all flags for differential flow;
+ // d) Get pointers to all nested lists in fDiffFlowListProfiles and to profiles which they hold;
+ // e) Get pointers to all nested lists in fDiffFlowListResults and to histograms which they hold.
+
+ // a) Define flags locally (to be improved: should I promote flags to data members?):
+ TString typeFlag[2] = {"RP","POI"};
+ TString ptEtaFlag[2] = {"p_{T}","#eta"};
+ TString powerFlag[2] = {"linear","quadratic"};
+ TString sinCosFlag[2] = {"sin","cos"};
+ TString differentialCumulantIndex[4] = {"QC{2'}","QC{4'}","QC{6'}","QC{8'}"};
+ TString differentialFlowIndex[4] = {"v'{2}","v'{4}","v'{6}","v'{8}"};
+ TString reducedCorrelationIndex[4] = {"<2'>","<4'>","<6'>","<8'>"};
+ TString mixedCorrelationIndex[8] = {"<2>","<2'>","<4>","<4'>","<6>","<6'>","<8>","<8'>"};
+ TString covarianceName[5] = {"Cov(<2>,<2'>)","Cov(<2>,<4'>)","Cov(<4>,<2'>)","Cov(<4>,<4'>)","Cov(<2'>,<4'>)"};
+
+ // b) Get pointer to base list for differential flow fDiffFlowList and nested lists fDiffFlowListProfiles and fDiffFlowListResults:
+ TList *diffFlowList = NULL;
+ diffFlowList = dynamic_cast<TList*>(fHistList->FindObject("Differential Flow"));
+ if(!diffFlowList)
+ {
+ cout<<"WARNING: diffFlowList is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ exit(0);
+ }
+ // list holding nested lists containing profiles:
+ TList *diffFlowListProfiles = NULL;
+ diffFlowListProfiles = dynamic_cast<TList*>(diffFlowList->FindObject("Profiles"));
+ if(!diffFlowListProfiles)
+ {
+ cout<<"WARNING: diffFlowListProfiles is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ exit(0);
+ }
+ // list holding nested lists containing 2D and 1D histograms with final results:
+ TList *diffFlowListResults = NULL;
+ diffFlowListResults = dynamic_cast<TList*>(diffFlowList->FindObject("Results"));
+ if(!diffFlowListResults)
+ {
+ cout<<"WARNING: diffFlowListResults is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ exit(0);
+ }
+
+ // c) Get pointer to profile holding all flags for differential flow;
+ TString diffFlowFlagsName = "fDiffFlowFlags";
+ diffFlowFlagsName += fAnalysisLabel->Data();
+ TProfile *diffFlowFlags = dynamic_cast<TProfile*>(diffFlowList->FindObject(diffFlowFlagsName.Data()));
+ Bool_t bCalculate2DFlow = kFALSE;
+ if(diffFlowFlags)
+ {
+ this->SetDiffFlowFlags(diffFlowFlags);
+ bCalculate2DFlow = (Int_t)diffFlowFlags->GetBinContent(4);
+ this->SetCalculate2DFlow(bCalculate2DFlow); // to be improved (shoul I call this setter somewhere else?)
+ }
+
+ // d) Get pointers to all nested lists in fDiffFlowListProfiles and to profiles which they hold;
+ // correlations:
+ TList *diffFlowCorrelationsProList[2][2] = {{NULL}};
+ TString diffFlowCorrelationsProName = "fDiffFlowCorrelationsPro";
+ diffFlowCorrelationsProName += fAnalysisLabel->Data();
+ TProfile *diffFlowCorrelationsPro[2][2][4] = {{{NULL}}};
+ // products of correlations:
+ TList *diffFlowProductOfCorrelationsProList[2][2] = {{NULL}};
+ TString diffFlowProductOfCorrelationsProName = "fDiffFlowProductOfCorrelationsPro";
+ diffFlowProductOfCorrelationsProName += fAnalysisLabel->Data();
+ TProfile *diffFlowProductOfCorrelationsPro[2][2][8][8] = {{{{NULL}}}};
+ // corrections:
+ TList *diffFlowCorrectionsProList[2][2] = {{NULL}};
+ TString diffFlowCorrectionTermsForNUAProName = "fDiffFlowCorrectionTermsForNUAPro";
+ diffFlowCorrectionTermsForNUAProName += fAnalysisLabel->Data();
+ TProfile *diffFlowCorrectionTermsForNUAPro[2][2][2][10] = {{{{NULL}}}};
+ for(Int_t t=0;t<2;t++)
+ {
+ for(Int_t pe=0;pe<2;pe++)
+ {
+ diffFlowCorrelationsProList[t][pe] = dynamic_cast<TList*>(diffFlowListProfiles->FindObject(Form("Profiles with correlations (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data())));
+ if(!diffFlowCorrelationsProList[t][pe])
+ {
+ cout<<"WARNING: diffFlowCorrelationsProList[t][pe] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ exit(0);
+ }
+ for(Int_t ci=0;ci<4;ci++) // correlation index
+ {
+ diffFlowCorrelationsPro[t][pe][ci] = dynamic_cast<TProfile*>(diffFlowCorrelationsProList[t][pe]->FindObject(Form("%s, %s, %s, %s",diffFlowCorrelationsProName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),reducedCorrelationIndex[ci].Data())));
+ if(diffFlowCorrelationsPro[t][pe][ci])
+ {
+ this->SetDiffFlowCorrelationsPro(diffFlowCorrelationsPro[t][pe][ci],t,pe,ci);
+ } else
+ {
+ cout<<"WARNING: diffFlowCorrelationsPro[t][pe][ci] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"ci = "<<ci<<endl;
+ }
+ } // end of for(Int_t ci=0;ci<4;ci++) // correlation index
+ // products of correlations:
+ diffFlowProductOfCorrelationsProList[t][pe] = dynamic_cast<TList*>(diffFlowListProfiles->FindObject(Form("Profiles with products of correlations (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data())));
+ if(!diffFlowProductOfCorrelationsProList[t][pe])
+ {
+ cout<<"WARNING: ddiffFlowProductOfCorrelationsProList[t][pe] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ exit(0);
+ }
+ for(Int_t mci1=0;mci1<8;mci1++) // mixed correlation index
+ {
+ for(Int_t mci2=mci1+1;mci2<8;mci2++) // mixed correlation index
+ {
+ diffFlowProductOfCorrelationsPro[t][pe][mci1][mci2] = dynamic_cast<TProfile*>(diffFlowProductOfCorrelationsProList[t][pe]->FindObject(Form("%s, %s, %s, %s, %s",diffFlowProductOfCorrelationsProName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),mixedCorrelationIndex[mci1].Data(),mixedCorrelationIndex[mci2].Data())));
+ if(diffFlowProductOfCorrelationsPro[t][pe][mci1][mci2])
+ {
+ this->SetDiffFlowProductOfCorrelationsPro(diffFlowProductOfCorrelationsPro[t][pe][mci1][mci2],t,pe,mci1,mci2);
+ } else
+ {
+ cout<<"WARNING: diffFlowCorrelationsPro[t][pe][ci] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"mci1 = "<<mci1<<endl;
+ cout<<"mci2 = "<<mci2<<endl;
+ }
+ if(mci1%2 == 0) mci2++; // products which DO NOT include reduced correlations are not stored here
+ } // end of for(Int_t mci2=mci1+1;mci2<8;mci2++) // mixed correlation index
+ } // end of for(Int_t mci1=0;mci1<8;mci1++) // mixed correlation index
+ // corrections:
+ diffFlowCorrectionsProList[t][pe] = dynamic_cast<TList*>(diffFlowListProfiles->FindObject(Form("Profiles with correction terms for NUA (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data())));
+ if(!diffFlowCorrectionsProList[t][pe])
+ {
+ cout<<"WARNING: diffFlowCorrectionsProList[t][pe] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ exit(0);
+ }
+ // correction terms for NUA:
+ for(Int_t sc=0;sc<2;sc++) // sin or cos
+ {
+ for(Int_t cti=0;cti<9;cti++) // correction term index
+ {
+ diffFlowCorrectionTermsForNUAPro[t][pe][sc][cti] = dynamic_cast<TProfile*>(diffFlowCorrectionsProList[t][pe]->FindObject(Form("%s, %s, %s, %s, cti = %d",diffFlowCorrectionTermsForNUAProName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),sinCosFlag[sc].Data(),cti+1)));
+ if(diffFlowCorrectionTermsForNUAPro[t][pe][sc][cti])
+ {
+ this->SetDiffFlowCorrectionTermsForNUAPro(diffFlowCorrectionTermsForNUAPro[t][pe][sc][cti],t,pe,sc,cti);
+ } else
+ {
+ cout<<"WARNING: diffFlowCorrectionTermsForNUAPro[t][pe][sc][cti] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"sc = "<<sc<<endl;
+ cout<<"cti = "<<cti<<endl;
+ }
+ } // end of for(Int_t cti=0;cti<9;cti++) // correction term index
+ } // end of for(Int_t sc=0;sc<2;sc++) // sin or cos
+ // ...
+ } // end of for(Int_t pe=0;pe<2;pe++)
+ } // end of for(Int_t t=0;t<2;t++)
+
+ // e) Get pointers to all nested lists in fDiffFlowListResults and to histograms which they hold.
+ // reduced correlations:
+ TList *diffFlowCorrelationsHistList[2][2] = {{NULL}};
+ TString diffFlowCorrelationsHistName = "fDiffFlowCorrelationsHist";
+ diffFlowCorrelationsHistName += fAnalysisLabel->Data();
+ TH1D *diffFlowCorrelationsHist[2][2][4] = {{{NULL}}};
+ // corrections for NUA:
+ TList *diffFlowCorrectionsHistList[2][2] = {{NULL}};
+ TString diffFlowCorrectionTermsForNUAHistName = "fDiffFlowCorrectionTermsForNUAHist";
+ diffFlowCorrectionTermsForNUAHistName += fAnalysisLabel->Data();
+ TH1D *diffFlowCorrectionTermsForNUAHist[2][2][2][10] = {{{{NULL}}}};
+ // differential Q-cumulants:
+ TList *diffFlowCumulantsHistList[2][2] = {{NULL}};
+ TString diffFlowCumulantsName = "fDiffFlowCumulants";
+ diffFlowCumulantsName += fAnalysisLabel->Data();
+ TH1D *diffFlowCumulants[2][2][4] = {{{NULL}}};
+ // differential flow estimates from Q-cumulants:
+ TList *diffFlowHistList[2][2] = {{NULL}};
+ TString diffFlowName = "fDiffFlow";
+ diffFlowName += fAnalysisLabel->Data();
+ TH1D *diffFlow[2][2][4] = {{{NULL}}};
+ // differential covariances:
+ TList *diffFlowCovariancesHistList[2][2] = {{NULL}};
+ TString diffFlowCovariancesName = "fDiffFlowCovariances";
+ diffFlowCovariancesName += fAnalysisLabel->Data();
+ TH1D *diffFlowCovariances[2][2][5] = {{{NULL}}};
+ for(Int_t t=0;t<2;t++) // type: RP or POI
+ {
+ for(Int_t pe=0;pe<2;pe++) // pt or eta
+ {
+ // reduced correlations:
+ diffFlowCorrelationsHistList[t][pe] = dynamic_cast<TList*>(diffFlowListResults->FindObject(Form("Correlations (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data())));
+ if(!diffFlowCorrelationsHistList[t][pe])
+ {
+ cout<<"WARNING: diffFlowCorrelationsHistList[t][pe] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ exit(0);
+ }
+ for(Int_t index=0;index<4;index++)
+ {
+ diffFlowCorrelationsHist[t][pe][index] = dynamic_cast<TH1D*>(diffFlowCorrelationsHistList[t][pe]->FindObject(Form("%s, %s, %s, %s",diffFlowCorrelationsHistName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),reducedCorrelationIndex[index].Data())));
+ if(diffFlowCorrelationsHist[t][pe][index])
+ {
+ this->SetDiffFlowCorrelationsHist(diffFlowCorrelationsHist[t][pe][index],t,pe,index);
+ } else
+ {
+ cout<<"WARNING: diffFlowCorrelationsHist[t][pe][index] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"index = "<<index<<endl;
+ exit(0);
+ }
+ } // end of for(Int_t index=0;index<4;index++)
+ // corrections:
+ diffFlowCorrectionsHistList[t][pe] = dynamic_cast<TList*>(diffFlowListResults->FindObject(Form("Histograms with correction terms for NUA (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data())));
+ if(!diffFlowCorrectionsHistList[t][pe])
+ {
+ cout<<"WARNING: diffFlowCorrectionsHistList[t][pe] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ exit(0);
+ }
+ // correction terms for NUA:
+ for(Int_t sc=0;sc<2;sc++) // sin or cos
+ {
+ for(Int_t cti=0;cti<9;cti++) // correction term index
+ {
+ diffFlowCorrectionTermsForNUAHist[t][pe][sc][cti] = dynamic_cast<TH1D*>(diffFlowCorrectionsHistList[t][pe]->FindObject(Form("%s, %s, %s, %s, cti = %d",diffFlowCorrectionTermsForNUAHistName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),sinCosFlag[sc].Data(),cti+1)));
+ if(diffFlowCorrectionTermsForNUAHist[t][pe][sc][cti])
+ {
+ this->SetDiffFlowCorrectionTermsForNUAHist(diffFlowCorrectionTermsForNUAHist[t][pe][sc][cti],t,pe,sc,cti);
+ } else
+ {
+ cout<<"WARNING: diffFlowCorrectionTermsForNUAHist[t][pe][sc][cti] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"sc = "<<sc<<endl;
+ cout<<"cti = "<<cti<<endl;
+ }
+ } // end of for(Int_t cti=0;cti<9;cti++) // correction term index
+ } // end of for(Int_t sc=0;sc<2;sc++) // sin or cos
+ // ...
+ // differential Q-cumulants:
+ diffFlowCumulantsHistList[t][pe] = dynamic_cast<TList*>(diffFlowListResults->FindObject(Form("Differential Q-cumulants (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data())));
+ if(!diffFlowCumulantsHistList[t][pe])
+ {
+ cout<<"WARNING: diffFlowCumulantsHistList[t][pe] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ exit(0);
+ }
+ for(Int_t index=0;index<4;index++)
+ {
+ diffFlowCumulants[t][pe][index] = dynamic_cast<TH1D*>(diffFlowCumulantsHistList[t][pe]->FindObject(Form("%s, %s, %s, %s",diffFlowCumulantsName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),differentialCumulantIndex[index].Data())));
+ if(diffFlowCumulants[t][pe][index])
+ {
+ this->SetDiffFlowCumulants(diffFlowCumulants[t][pe][index],t,pe,index);
+ } else
+ {
+ cout<<"WARNING: diffFlowCumulants[t][pe][index] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"index = "<<index<<endl;
+ exit(0);
+ }
+ } // end of for(Int_t index=0;index<4;index++)
+ // differential flow estimates from Q-cumulants:
+ diffFlowHistList[t][pe] = dynamic_cast<TList*>(diffFlowListResults->FindObject(Form("Differential flow (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data())));
+ if(!diffFlowHistList[t][pe])
+ {
+ cout<<"WARNING: diffFlowHistList[t][pe] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ exit(0);
+ }
+ for(Int_t index=0;index<4;index++)
+ {
+ diffFlow[t][pe][index] = dynamic_cast<TH1D*>(diffFlowHistList[t][pe]->FindObject(Form("%s, %s, %s, %s",diffFlowName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),differentialFlowIndex[index].Data())));
+ if(diffFlow[t][pe][index])
+ {
+ this->SetDiffFlow(diffFlow[t][pe][index],t,pe,index);
+ } else
+ {
+ cout<<"WARNING: diffFlow[t][pe][index] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"index = "<<index<<endl;
+ exit(0);
+ }
+ } // end of for(Int_t index=0;index<4;index++)
+ // differential covariances:
+ diffFlowCovariancesHistList[t][pe] = dynamic_cast<TList*>(diffFlowListResults->FindObject(Form("Covariances of correlations (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data())));
+ if(!diffFlowCovariancesHistList[t][pe])
+ {
+ cout<<"WARNING: diffFlowCovariancesHistList[t][pe] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ exit(0);
+ }
+ for(Int_t covIndex=0;covIndex<5;covIndex++)
+ {
+ diffFlowCovariances[t][pe][covIndex] = dynamic_cast<TH1D*>(diffFlowCovariancesHistList[t][pe]->FindObject(Form("%s, %s, %s, %s",diffFlowCovariancesName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),covarianceName[covIndex].Data())));
+ if(diffFlowCovariances[t][pe][covIndex])
+ {
+ this->SetDiffFlowCovariances(diffFlowCovariances[t][pe][covIndex],t,pe,covIndex);
+ } else
+ {
+ cout<<"WARNING: diffFlowCovariances[t][pe][covIndex] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"covIndex = "<<covIndex<<endl;
+ exit(0);
+ }
+ } // end of for(Int_t covIndex=0;covIndex<5;covIndex++) // covariance index
+ } // end of for(Int_t pe=0;pe<2;pe++) // pt or eta
+ } // end of for(Int_t t=0;t<2;t++) // type: RP or POI
+ // sum of event weights for reduced correlations:
+ TList *diffFlowSumOfEventWeightsHistList[2][2][2] = {{{NULL}}};
+ TString diffFlowSumOfEventWeightsName = "fDiffFlowSumOfEventWeights";
+ diffFlowSumOfEventWeightsName += fAnalysisLabel->Data();
+ TH1D *diffFlowSumOfEventWeights[2][2][2][4] = {{{{NULL}}}};
+ for(Int_t t=0;t<2;t++) // type is RP or POI
+ {
+ for(Int_t pe=0;pe<2;pe++) // pt or eta
+ {
+ for(Int_t p=0;p<2;p++) // power of event weights is either 1 or 2
+ {
+ diffFlowSumOfEventWeightsHistList[t][pe][p] = dynamic_cast<TList*>(diffFlowListResults->FindObject(Form("Sum of %s event weights (%s, %s)",powerFlag[p].Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data())));
+ if(!diffFlowSumOfEventWeightsHistList[t][pe][p])
+ {
+ cout<<"WARNING: diffFlowSumOfEventWeightsHistList[t][pe][p] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"power = "<<p<<endl;
+ exit(0);
+ }
+ for(Int_t ew=0;ew<4;ew++) // index of reduced correlation
+ {
+ diffFlowSumOfEventWeights[t][pe][p][ew] = dynamic_cast<TH1D*>(diffFlowSumOfEventWeightsHistList[t][pe][p]->FindObject(Form("%s, %s, %s, %s, %s",diffFlowSumOfEventWeightsName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),powerFlag[p].Data(),reducedCorrelationIndex[ew].Data())));
+ if(diffFlowSumOfEventWeights[t][pe][p][ew])
+ {
+ this->SetDiffFlowSumOfEventWeights(diffFlowSumOfEventWeights[t][pe][p][ew],t,pe,p,ew);
+ } else
+ {
+ cout<<"WARNING: diffFlowSumOfEventWeights[t][pe][p][ew] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"power = "<<p<<endl;
+ cout<<"ew = "<<ew<<endl;
+ exit(0);
+ }
+ }
+ } // end of for(Int_t p=0;p<2;p++) // power of event weights is either 1 or 2
+ } // end of for(Int_t pe=0;pe<2;pe++) // pt or eta
+ } // end of for(Int_t t=0;t<2;t++) // type is RP or POI
+ //
+ TList *diffFlowSumOfProductOfEventWeightsHistList[2][2] = {{NULL}};
+ TString diffFlowSumOfProductOfEventWeightsName = "fDiffFlowSumOfProductOfEventWeights";
+ diffFlowSumOfProductOfEventWeightsName += fAnalysisLabel->Data();
+ TH1D *diffFlowSumOfProductOfEventWeights[2][2][8][8] = {{{{NULL}}}};
+ for(Int_t t=0;t<2;t++) // type is RP or POI
+ {
+ for(Int_t pe=0;pe<2;pe++) // pt or eta
+ {
+ diffFlowSumOfProductOfEventWeightsHistList[t][pe] = dynamic_cast<TList*>(diffFlowListResults->FindObject(Form("Sum of products of event weights (%s, %s)",typeFlag[t].Data(),ptEtaFlag[pe].Data())));
+ if(!diffFlowSumOfProductOfEventWeightsHistList[t][pe])
+ {
+ cout<<"WARNING: diffFlowSumOfProductOfEventWeightsHistList[t][pe] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ exit(0);
+ }
+ for(Int_t mci1=0;mci1<8;mci1++) // mixed correlation index
+ {
+ for(Int_t mci2=mci1+1;mci2<8;mci2++) // mixed correlation index
+ {
+ diffFlowSumOfProductOfEventWeights[t][pe][mci1][mci2] = dynamic_cast<TH1D*>(diffFlowSumOfProductOfEventWeightsHistList[t][pe]->FindObject(Form("%s, %s, %s, %s, %s",diffFlowSumOfProductOfEventWeightsName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),mixedCorrelationIndex[mci1].Data(),mixedCorrelationIndex[mci2].Data())));
+ if(diffFlowSumOfProductOfEventWeights[t][pe][mci1][mci2])
+ {
+ this->SetDiffFlowSumOfProductOfEventWeights(diffFlowSumOfProductOfEventWeights[t][pe][mci1][mci2],t,pe,mci1,mci2);
+ } else
+ {
+ cout<<"WARNING: diffFlowSumOfProductOfEventWeights[t][pe][mci1][mci2] is NULL in AFAWQC::GPFDFH() !!!!"<<endl;
+ cout<<"t = "<<t<<endl;
+ cout<<"pe = "<<pe<<endl;
+ cout<<"mci1 = "<<mci1<<endl;
+ cout<<"mci2 = "<<mci2<<endl;
+ exit(0);
+ }
+ if(mci1%2 == 0) mci2++; // products which DO NOT include reduced correlations are not stored here
+ } // end of for(Int_t mci2=mci1+1;mci2<8;mci2++) // mixed correlation index
+ } // end of for(Int_t mci1=0;mci1<8;mci1++) // mixed correlation index
+ } // end of for(Int_t pe=0;pe<2;pe++) // pt or eta
+ } // end of for(Int_t t=0;t<2;t++) // type is RP or POI
+
+} // end void AliFlowAnalysisWithQCumulants::GetPointersForDiffFlowHistograms()
+
+
+//================================================================================================================================
+
+
+void AliFlowAnalysisWithQCumulants::BookEverythingForDifferentialFlow()
+{
+ // Book all histograms and profiles needed for differential flow.
+ // a) Define flags locally (to be improved: should I promote flags to data members?);
+ // b) Book profile to hold all flags for differential flow;
+ // c) Book e-b-e quantities;
+ // d) Book profiles;
+ // e) Book histograms holding final results.
+
+ // a) Define flags locally (to be improved: should I promote flags to data members?):
+ TString typeFlag[2] = {"RP","POI"};
+ TString ptEtaFlag[2] = {"p_{T}","#eta"};
+ TString powerFlag[2] = {"linear","quadratic"};
+ TString sinCosFlag[2] = {"sin","cos"};
+ TString differentialCumulantIndex[4] = {"QC{2'}","QC{4'}","QC{6'}","QC{8'}"};
+ TString differentialFlowIndex[4] = {"v'{2}","v'{4}","v'{6}","v'{8}"};
+ TString reducedCorrelationIndex[4] = {"<2'>","<4'>","<6'>","<8'>"};
+ TString mixedCorrelationIndex[8] = {"<2>","<2'>","<4>","<4'>","<6>","<6'>","<8>","<8'>"};
+ TString covarianceName[5] = {"Cov(<2>,<2'>)","Cov(<2>,<4'>)","Cov(<4>,<2'>)","Cov(<4>,<4'>)","Cov(<2'>,<4'>)"};
+ Int_t nBinsPtEta[2] = {fnBinsPt,fnBinsEta};
+ Double_t minPtEta[2] = {fPtMin,fEtaMin};
+ Double_t maxPtEta[2] = {fPtMax,fEtaMax};
+
+ // b) Book profile to hold all flags for differential flow:
+ TString diffFlowFlagsName = "fDiffFlowFlags";
+ diffFlowFlagsName += fAnalysisLabel->Data();
+ fDiffFlowFlags = new TProfile(diffFlowFlagsName.Data(),"Flags for Differential Flow",4,0,4);
+ fDiffFlowFlags->SetTickLength(-0.01,"Y");
+ fDiffFlowFlags->SetMarkerStyle(25);
+ fDiffFlowFlags->SetLabelSize(0.05);
+ fDiffFlowFlags->SetLabelOffset(0.02,"Y");
+ (fDiffFlowFlags->GetXaxis())->SetBinLabel(1,"Particle Weights");
+ (fDiffFlowFlags->GetXaxis())->SetBinLabel(2,"Event Weights");
+ (fDiffFlowFlags->GetXaxis())->SetBinLabel(3,"Corrected for NUA?");
+ (fDiffFlowFlags->GetXaxis())->SetBinLabel(4,"Calculated 2D flow?");
+ fDiffFlowList->Add(fDiffFlowFlags);
+
+ // c) Book e-b-e quantities:
+ // Event-by-event r_{m*n,k}(pt,eta), p_{m*n,k}(pt,eta) and q_{m*n,k}(pt,eta)
+ // Explanantion of notation:
+ // 1.) n is harmonic, m is multiple of harmonic;
+ // 2.) k is power of particle weight;
+ // 3.) r_{m*n,k}(pt,eta) = Q-vector evaluated in harmonic m*n for RPs in particular (pt,eta) bin (i-th RP is weighted with w_i^k);
+ // 4.) p_{m*n,k}(pt,eta) = Q-vector evaluated in harmonic m*n for POIs in particular (pt,eta) bin
+ // (if i-th POI is also RP, than it is weighted with w_i^k);
+ // 5.) q_{m*n,k}(pt,eta) = Q-vector evaluated in harmonic m*n for particles which are both RPs and POIs in particular (pt,eta) bin
+ // (i-th RP&&POI is weighted with w_i^k)
+
+ // 1D:
+ for(Int_t t=0;t<3;t++) // typeFlag (0 = RP, 1 = POI, 2 = RP && POI )
+ {
+ for(Int_t pe=0;pe<2;pe++) // pt or eta
+ {
+ for(Int_t m=0;m<4;m++) // multiple of harmonic
+ {
+ for(Int_t k=0;k<9;k++) // power of particle weight
+ {
+ fReRPQ1dEBE[t][pe][m][k] = new TProfile(Form("TypeFlag%dpteta%dmultiple%dpower%dRe",t,pe,m,k),
+ Form("TypeFlag%dpteta%dmultiple%dpower%dRe",t,pe,m,k),nBinsPtEta[pe],minPtEta[pe],maxPtEta[pe]);
+ fImRPQ1dEBE[t][pe][m][k] = new TProfile(Form("TypeFlag%dpteta%dmultiple%dpower%dIm",t,pe,m,k),
+ Form("TypeFlag%dpteta%dmultiple%dpower%dIm",t,pe,m,k),nBinsPtEta[pe],minPtEta[pe],maxPtEta[pe]);
+ }
+ }
+ }
+ }
+ // to be improved (add explanation of fs1dEBE[t][pe][k]):
+ for(Int_t t=0;t<3;t++) // typeFlag (0 = RP, 1 = POI, 2 = RP&&POI )
+ {
+ for(Int_t pe=0;pe<2;pe++) // pt or eta
+ {
+ for(Int_t k=0;k<9;k++) // power of particle weight
+ {
+ fs1dEBE[t][pe][k] = new TProfile(Form("TypeFlag%dpteta%dmultiple%d",t,pe,k),
+ Form("TypeFlag%dpteta%dmultiple%d",t,pe,k),nBinsPtEta[pe],minPtEta[pe],maxPtEta[pe]);
+ }
+ }
+ }
+ // correction terms for nua:
+ for(Int_t t=0;t<2;t++) // typeFlag (0 = RP, 1 = POI)
+ {
+ for(Int_t pe=0;pe<2;pe++) // pt or eta
+ {
+ for(Int_t sc=0;sc<2;sc++) // sin or cos
+ {
+ for(Int_t cti=0;cti<9;cti++) // correction term index
+ {
+ fDiffFlowCorrectionTermsForNUAEBE[t][pe][sc][cti] = new TH1D(Form("typeFlag%d pteta%d sincos%d cti%d",t,pe,sc,cti),
+ Form("typeFlag%d pteta%d sincos%d cti%d",t,pe,sc,cti),nBinsPtEta[pe],minPtEta[pe],maxPtEta[pe]);
+ }
+ }
+ }
+ }
+ // 2D:
+ TProfile2D styleRe("typeMultiplePowerRe","typeMultiplePowerRe",fnBinsPt,fPtMin,fPtMax,fnBinsEta,fEtaMin,fEtaMax);
+ TProfile2D styleIm("typeMultiplePowerIm","typeMultiplePowerIm",fnBinsPt,fPtMin,fPtMax,fnBinsEta,fEtaMin,fEtaMax);
+ for(Int_t t=0;t<3;t++) // typeFlag (0 = RP, 1 = POI, 2 = RP&&POI )
+ {
+ for(Int_t m=0;m<4;m++)
+ {
+ for(Int_t k=0;k<9;k++)
+ {
+ fReRPQ2dEBE[t][m][k] = (TProfile2D*)styleRe.Clone(Form("typeFlag%dmultiple%dpower%dRe",t,m,k));
+ fImRPQ2dEBE[t][m][k] = (TProfile2D*)styleIm.Clone(Form("typeFlag%dmultiple%dpower%dIm",t,m,k));
+ }
+ }
+ }
+ TProfile2D styleS("typePower","typePower",fnBinsPt,fPtMin,fPtMax,fnBinsEta,fEtaMin,fEtaMax);
+ for(Int_t t=0;t<3;t++) // typeFlag (0 = RP, 1 = POI, 2 = RP&&POI )
+ {
+ for(Int_t k=0;k<9;k++)
+ {
+ fs2dEBE[t][k] = (TProfile2D*)styleS.Clone(Form("typeFlag%dpower%d",t,k));
+ }
+ }
+ // reduced correlations e-b-e:
+ TString diffFlowCorrelationsEBEName = "fDiffFlowCorrelationsEBE";
+ diffFlowCorrelationsEBEName += fAnalysisLabel->Data();
+ for(Int_t t=0;t<2;t++) // type: RP or POI
+ {
+ for(Int_t pe=0;pe<2;pe++) // pt or eta
+ {
+ for(Int_t rci=0;rci<4;rci++) // reduced correlation index
+ {
+ fDiffFlowCorrelationsEBE[t][pe][rci] = new TH1D(Form("%s, %s, %s, %s",diffFlowCorrelationsEBEName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),reducedCorrelationIndex[rci].Data()),Form("%s, %s, %s, %s",diffFlowCorrelationsEBEName.Data(),typeFlag[t].Data(),ptEtaFlag[pe].Data(),reducedCorrelationIndex[rci].Data()),nBinsPtEta[pe],minPtEta[pe],maxPtEta[pe]);
+ } // end of for(Int_t ci=0;ci<4;ci++) // correlation index
+ } // end of for(Int_t pe=0;pe<2;pe++) // pt or eta
+ } // end of for(Int_t t=0;t<2;t++) // type: RP or POI
+ // event weights for reduced correlations e-b-e:
+ TString diffFlowEventWeightsForCorrelationsEBEName = "fDiffFlowEventWeightsForCorrelationsEBE";
+ diffFlowEventWeightsForCorrelationsEBEName += fAnalysisLabel->Data();
+ for(Int_t t=0;t<2;t++) // type: RP or POI
+ {
+ for(Int_t pe=0;pe<2;pe++) // pt or eta