From: miweber Date: Wed, 28 Nov 2012 19:50:28 +0000 (+0000) Subject: 1) correcting the normalization for per-trigger-yield for multidimensional analysis... X-Git-Url: http://git.uio.no/git/?a=commitdiff_plain;h=47e040b5e71724fe7be848a003edb5a17fe837d5;p=u%2Fmrichter%2FAliRoot.git 1) correcting the normalization for per-trigger-yield for multidimensional analysis 2) adding the possibility to correct with event mixing for 1D analysis 3) minor changes --- diff --git a/PWGCF/EBYE/BalanceFunctions/AliBalance.cxx b/PWGCF/EBYE/BalanceFunctions/AliBalance.cxx index ec7d21ca58b..064feaa0a99 100644 --- a/PWGCF/EBYE/BalanceFunctions/AliBalance.cxx +++ b/PWGCF/EBYE/BalanceFunctions/AliBalance.cxx @@ -55,7 +55,11 @@ AliBalance::AliBalance() : fAnalyzedEvents(0) , fCentralityId(0) , fCentStart(0.), - fCentStop(0.) + fCentStop(0.), + fHistHBTbefore(NULL), + fHistHBTafter(NULL), + fHistConversionbefore(NULL), + fHistConversionafter(NULL) { // Default constructor @@ -99,13 +103,6 @@ AliBalance::AliBalance() : fHistNN[i] = NULL; } - - //QA histograms - fHistHBTbefore = NULL; - fHistHBTafter = NULL; - fHistConversionbefore = NULL; - fHistConversionafter = NULL; - } @@ -119,7 +116,11 @@ AliBalance::AliBalance(const AliBalance& balance): fAnalyzedEvents(balance.fAnalyzedEvents), fCentralityId(balance.fCentralityId), fCentStart(balance.fCentStart), - fCentStop(balance.fCentStop) { + fCentStop(balance.fCentStop), + fHistHBTbefore(balance.fHistHBTbefore), + fHistHBTafter(balance.fHistHBTafter), + fHistConversionbefore(balance.fHistConversionbefore), + fHistConversionafter(balance.fHistConversionafter) { //copy constructor for(Int_t i = 0; i < ANALYSIS_TYPES; i++){ fNn[i] = balance.fNn[i]; @@ -789,12 +790,12 @@ void AliBalance::PrintResults(Int_t iAnalysisType, TH1D *gHistBalance) { Double_t deltaBalP2 = 0.0, integral = 0.0; Double_t deltaErrorNew = 0.0; - cout<<"=================================================="<GetBinContent(i)<<"\t Error: "<GetBinError(i)<<"\t bin: "<GetBinCenter(i)<GetBinContent(i)<<"\t Error: "<GetBinError(i)<<"\t bin: "<GetBinCenter(i)<GetBinCenter(i); gSumBi += gHistBalance->GetBinContent(i); @@ -814,15 +815,15 @@ void AliBalance::PrintResults(Int_t iAnalysisType, TH1D *gHistBalance) { Double_t delta = gSumBiXi / gSumBi; Double_t deltaError = (gSumBiXi / gSumBi) * TMath::Sqrt(TMath::Power((TMath::Sqrt(gSumXi2DeltaBi2)/gSumBiXi),2) + TMath::Power((gSumDeltaBi2/gSumBi),2) ); - cout<<"Analysis type: "< 0 && correctWithEfficiency){ + if(iAnalysisType == kEta && etaWindow > 0 && correctWithEfficiency && !correctWithMixed){ TH1F* hEffP = NULL; TH1F* hEffN = NULL; @@ -974,7 +975,7 @@ TH1D *AliBalance::GetBalanceFunctionHistogram(Int_t iAnalysisType,Double_t centr // - single particle efficiencies from MC (AliAnalysiTaskEfficiency) // - two particle efficiencies from convolution of data single particle distributions // (normalized to single particle efficiency) - if(iAnalysisType == kPhi && correctWithEfficiency){ + if(iAnalysisType == kPhi && correctWithEfficiency && !correctWithMixed){ TH1F* hEffPhiP = NULL; TH1F* hEffPhiN = NULL; @@ -1020,6 +1021,48 @@ TH1D *AliBalance::GetBalanceFunctionHistogram(Int_t iAnalysisType,Double_t centr } } + // do the correction with the event mixing directly! + if(correctWithMixed){ + + if(hMixed[0] && hMixed[1] && hMixed[2] && hMixed[3]){ + + // scale that EM is 1 at 0 for Deta + // in the region 0-10degree (one 1/2 sector) for Dphi + if(iAnalysisType==6){ + hMixed[0]->Scale(1./(Double_t)hMixed[0]->Integral(hMixed[0]->FindBin(0),hMixed[0]->FindBin(10))*(Double_t)(hMixed[0]->FindBin(10)-hMixed[0]->FindBin(0)+1)); + hMixed[2]->Scale(1./(Double_t)hMixed[2]->Integral(hMixed[2]->FindBin(0),hMixed[2]->FindBin(10))*(Double_t)(hMixed[0]->FindBin(10)-hMixed[0]->FindBin(0)+1)); + hMixed[3]->Scale(1./(Double_t)hMixed[3]->Integral(hMixed[3]->FindBin(0),hMixed[3]->FindBin(10))*(Double_t)(hMixed[0]->FindBin(10)-hMixed[0]->FindBin(0)+1)); + } + else{ + hMixed[0]->Scale(1./(Double_t)hMixed[0]->GetBinContent(1)); + hMixed[2]->Scale(1./(Double_t)hMixed[2]->GetBinContent(1)); + hMixed[3]->Scale(1./(Double_t)hMixed[3]->GetBinContent(1)); + } + + // scale to average efficiency in the pt region (0.3-1.5) and |eta| < 0.8 + // by multiplying the average single particle efficiencies from HIJING + Double_t normPMC = 0.847546; + Double_t normNMC = 0.83827; + hMixed[0]->Scale(normNMC*normPMC); + hMixed[2]->Scale(normNMC*normNMC); + hMixed[3]->Scale(normPMC*normPMC); + + // divide by event mixing + hTemp1->Divide(hMixed[0]); + hTemp2->Divide(hMixed[0]); + hTemp3->Divide(hMixed[2]); + hTemp4->Divide(hMixed[3]); + + // scale also single histograms with average efficiency + hTemp5->Scale(1./normNMC); + hTemp6->Scale(1./normPMC); + + } + else{ + AliError("Correction with EventMixing requested, but not all Histograms there!"); + return NULL; + } + } if((hTemp1)&&(hTemp2)&&(hTemp3)&&(hTemp4)) { @@ -1036,7 +1079,7 @@ TH1D *AliBalance::GetBalanceFunctionHistogram(Int_t iAnalysisType,Double_t centr } // do the acceptance correction (only for Eta and etaWindow > 0) - if(iAnalysisType == kEta && etaWindow > 0 && !correctWithEfficiency){ + if(iAnalysisType == kEta && etaWindow > 0 && !correctWithEfficiency && !correctWithMixed){ for(Int_t iBin = 0; iBin < gHistBalanceFunctionHistogram->GetNbinsX(); iBin++){ Double_t notCorrected = gHistBalanceFunctionHistogram->GetBinContent(iBin+1); diff --git a/PWGCF/EBYE/BalanceFunctions/AliBalance.h b/PWGCF/EBYE/BalanceFunctions/AliBalance.h index fb2e5e29736..acd6e074b8c 100644 --- a/PWGCF/EBYE/BalanceFunctions/AliBalance.h +++ b/PWGCF/EBYE/BalanceFunctions/AliBalance.h @@ -116,7 +116,7 @@ class AliBalance : public TObject { void SetHistNnn(Int_t iAnalysisType, TH2D *gHist) { fHistNN[iAnalysisType] = gHist;} - TH1D *GetBalanceFunctionHistogram(Int_t iAnalysisType,Double_t centrMin, Double_t centrMax, Double_t etaWindow = -1, Bool_t correctWithEfficiency = kFALSE, Bool_t correctWithAcceptanceOnly = kFALSE); + TH1D *GetBalanceFunctionHistogram(Int_t iAnalysisType,Double_t centrMin, Double_t centrMax, Double_t etaWindow = -1, Bool_t correctWithEfficiency = kFALSE, Bool_t correctWithAcceptanceOnly = kFALSE, Bool_t correctWithMixed = kFALSE, TH1D *hMixed[4]=NULL); void PrintResults(Int_t iAnalysisType, TH1D *gHist); private: diff --git a/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.cxx b/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.cxx index 96dd89b10f8..02d2516d2aa 100644 --- a/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.cxx +++ b/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.cxx @@ -929,8 +929,7 @@ TH2D *AliBalancePsi::GetCorrelationFunctionPN(Double_t psiMin, Double_t ptTriggerMin, Double_t ptTriggerMax, Double_t ptAssociatedMin, - Double_t ptAssociatedMax, - Bool_t normToTrig) { + Double_t ptAssociatedMax) { //Returns the 2D correlation function for (+-) pairs // Psi_2: axis 0 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax); @@ -975,7 +974,7 @@ TH2D *AliBalancePsi::GetCorrelationFunctionPN(Double_t psiMin, //c2->cd(); //fHistPN->Project(0,1,2)->DrawCopy("colz"); - if(normToTrig && (Double_t)(fHistP->Project(0,1)->GetEntries())!=0) + if((Double_t)(fHistP->Project(0,1)->GetEntries())!=0) gHist->Scale(1./(Double_t)(fHistP->Project(0,1)->GetEntries())); return gHist; @@ -987,8 +986,7 @@ TH2D *AliBalancePsi::GetCorrelationFunctionNP(Double_t psiMin, Double_t ptTriggerMin, Double_t ptTriggerMax, Double_t ptAssociatedMin, - Double_t ptAssociatedMax, - Bool_t normToTrig) { + Double_t ptAssociatedMax) { //Returns the 2D correlation function for (+-) pairs // Psi_2: axis 0 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax); @@ -1013,7 +1011,7 @@ TH2D *AliBalancePsi::GetCorrelationFunctionNP(Double_t psiMin, //Printf("Entries (1D): %lf",(Double_t)(fHistN->Project(0,2)->GetEntries())); //Printf("Entries (2D): %lf",(Double_t)(fHistNP->Project(0,2,3)->GetEntries())); - if(normToTrig && (Double_t)(fHistN->Project(0,1)->GetEntries())!=0) + if((Double_t)(fHistN->Project(0,1)->GetEntries())!=0) gHist->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries())); return gHist; @@ -1025,8 +1023,7 @@ TH2D *AliBalancePsi::GetCorrelationFunctionPP(Double_t psiMin, Double_t ptTriggerMin, Double_t ptTriggerMax, Double_t ptAssociatedMin, - Double_t ptAssociatedMax, - Bool_t normToTrig) { + Double_t ptAssociatedMax) { //Returns the 2D correlation function for (+-) pairs // Psi_2: axis 0 fHistP->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax); @@ -1051,7 +1048,7 @@ TH2D *AliBalancePsi::GetCorrelationFunctionPP(Double_t psiMin, //Printf("Entries (1D): %lf",(Double_t)(fHistP->Project(0,2)->GetEntries())); //Printf("Entries (2D): %lf",(Double_t)(fHistPP->Project(0,2,3)->GetEntries())); - if(normToTrig && (Double_t)(fHistP->Project(0,1)->GetEntries())!=0) + if((Double_t)(fHistP->Project(0,1)->GetEntries())!=0) gHist->Scale(1./(Double_t)(fHistP->Project(0,1)->GetEntries())); return gHist; @@ -1063,8 +1060,7 @@ TH2D *AliBalancePsi::GetCorrelationFunctionNN(Double_t psiMin, Double_t ptTriggerMin, Double_t ptTriggerMax, Double_t ptAssociatedMin, - Double_t ptAssociatedMax, - Bool_t normToTrig) { + Double_t ptAssociatedMax) { //Returns the 2D correlation function for (+-) pairs // Psi_2: axis 0 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax); @@ -1089,7 +1085,7 @@ TH2D *AliBalancePsi::GetCorrelationFunctionNN(Double_t psiMin, //Printf("Entries (1D): %lf",(Double_t)(fHistN->Project(0,2)->GetEntries())); //Printf("Entries (2D): %lf",(Double_t)(fHistNN->Project(0,2,3)->GetEntries())); - if(normToTrig && (Double_t)(fHistN->Project(0,1)->GetEntries())!=0) + if((Double_t)(fHistN->Project(0,1)->GetEntries())!=0) gHist->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries())); return gHist; @@ -1101,8 +1097,7 @@ TH2D *AliBalancePsi::GetCorrelationFunctionChargeIndependent(Double_t psiMin, Double_t ptTriggerMin, Double_t ptTriggerMax, Double_t ptAssociatedMin, - Double_t ptAssociatedMax, - Bool_t normToTrig) { + Double_t ptAssociatedMax) { //Returns the 2D correlation function for the sum of all charge combination pairs // Psi_2: axis 0 fHistN->GetGrid(0)->GetGrid()->GetAxis(0)->SetRangeUser(psiMin,psiMax); @@ -1156,8 +1151,8 @@ TH2D *AliBalancePsi::GetCorrelationFunctionChargeIndependent(Double_t psiMin, gHistNN->Add(gHistNP); gHistNN->Add(gHistPN); - // if normalization to trigger particle, divide by sum of + and - triggers - if(normToTrig && (Double_t)(fHistN->Project(0,1)->GetEntries())!=0 && (Double_t)(fHistP->Project(0,1)->GetEntries())!=0) + // divide by sum of + and - triggers + if((Double_t)(fHistN->Project(0,1)->GetEntries())!=0 && (Double_t)(fHistP->Project(0,1)->GetEntries())!=0) gHistNN->Scale(1./(Double_t)(fHistN->Project(0,1)->GetEntries() + fHistN->Project(0,1)->GetEntries())); return gHistNN; diff --git a/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.h b/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.h index 370a3c28f93..f206c1893ff 100644 --- a/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.h +++ b/PWGCF/EBYE/BalanceFunctions/AliBalancePsi.h @@ -71,32 +71,27 @@ class AliBalancePsi : public TObject { Double_t ptTriggerMin=-1., Double_t ptTriggerMax=-1., Double_t ptAssociatedMin=-1., - Double_t ptAssociatedMax=-1, - Bool_t normToTrig=kTRUE); + Double_t ptAssociatedMax=-1); TH2D *GetCorrelationFunctionNP(Double_t psiMin, Double_t psiMax, Double_t ptTriggerMin=-1., Double_t ptTriggerMax=-1., Double_t ptAssociatedMin=-1., - Double_t ptAssociatedMax=-1, - Bool_t normToTrig=kTRUE); + Double_t ptAssociatedMax=-1); TH2D *GetCorrelationFunctionPP(Double_t psiMin, Double_t psiMax, Double_t ptTriggerMin=-1., Double_t ptTriggerMax=-1., Double_t ptAssociatedMin=-1., - Double_t ptAssociatedMax=-1, - Bool_t normToTrig=kTRUE); + Double_t ptAssociatedMax=-1); TH2D *GetCorrelationFunctionNN(Double_t psiMin, Double_t psiMax, Double_t ptTriggerMin=-1., Double_t ptTriggerMax=-1., Double_t ptAssociatedMin=-1., - Double_t ptAssociatedMax=-1, - Bool_t normToTrig=kTRUE); + Double_t ptAssociatedMax=-1); TH2D *GetCorrelationFunctionChargeIndependent(Double_t psiMin, Double_t psiMax, Double_t ptTriggerMin=-1., Double_t ptTriggerMax=-1., Double_t ptAssociatedMin=-1., - Double_t ptAssociatedMax=-1, - Bool_t normToTrig=kTRUE); + Double_t ptAssociatedMax=-1); AliTHn *GetHistNp() {return fHistP;} AliTHn *GetHistNn() {return fHistN;} diff --git a/PWGCF/EBYE/macros/drawCorrelationFunctionPsi.C b/PWGCF/EBYE/macros/drawCorrelationFunctionPsi.C index b095a2bb384..c9357198713 100644 --- a/PWGCF/EBYE/macros/drawCorrelationFunctionPsi.C +++ b/PWGCF/EBYE/macros/drawCorrelationFunctionPsi.C @@ -432,7 +432,7 @@ void draw(TList *list, TList *listBFShuffled, TList *listBFMixed, histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})"; // if normalization to trigger then do not divide Event mixing by number of trigger particles - gHistPN[2] = bMixed->GetCorrelationFunctionPN(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,normToTrig); + gHistPN[2] = bMixed->GetCorrelationFunctionPN(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax); if(rebinEta > 1 || rebinPhi > 1) gHistPN[2]->Rebin2D(rebinEta,rebinPhi); // normalization to 1 at (0,0) --> Jan Fietes method @@ -558,7 +558,7 @@ void draw(TList *list, TList *listBFShuffled, TList *listBFMixed, histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})"; // if normalization to trigger then do not divide Event mixing by number of trigger particles - gHistNP[2] = bMixed->GetCorrelationFunctionNP(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,normToTrig); + gHistNP[2] = bMixed->GetCorrelationFunctionNP(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax); if(rebinEta > 1 || rebinPhi > 1) gHistNP[2]->Rebin2D(rebinEta,rebinPhi); // normalization to 1 at (0,0) --> Jan Fietes method @@ -684,7 +684,7 @@ void draw(TList *list, TList *listBFShuffled, TList *listBFMixed, histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})"; // if normalization to trigger then do not divide Event mixing by number of trigger particles - gHistPP[2] = bMixed->GetCorrelationFunctionPP(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,normToTrig); + gHistPP[2] = bMixed->GetCorrelationFunctionPP(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax); if(rebinEta > 1 || rebinPhi > 1) gHistPP[2]->Rebin2D(rebinEta,rebinPhi); // normalization to 1 at (0,0) --> Jan Fietes method @@ -810,7 +810,7 @@ void draw(TList *list, TList *listBFShuffled, TList *listBFMixed, histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})"; // if normalization to trigger then do not divide Event mixing by number of trigger particles - gHistNN[2] = bMixed->GetCorrelationFunctionNN(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,normToTrig); + gHistNN[2] = bMixed->GetCorrelationFunctionNN(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax); if(rebinEta > 1 || rebinPhi > 1) gHistNN[2]->Rebin2D(rebinEta,rebinPhi); // normalization to 1 at (0,0) --> Jan Fietes method diff --git a/PWGCF/EBYE/macros/drawCorrelationFunctionPsiChargeIndependent.C b/PWGCF/EBYE/macros/drawCorrelationFunctionPsiChargeIndependent.C index 034dbc81a13..2d5504f6875 100644 --- a/PWGCF/EBYE/macros/drawCorrelationFunctionPsiChargeIndependent.C +++ b/PWGCF/EBYE/macros/drawCorrelationFunctionPsiChargeIndependent.C @@ -426,7 +426,7 @@ void draw(TList *list, TList *listBFShuffled, TList *listBFMixed, histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})"; // if normalization to trigger then do not divide Event mixing by number of trigger particles - gHist[2] = bMixed->GetCorrelationFunctionChargeIndependent(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,normToTrig); + gHist[2] = bMixed->GetCorrelationFunctionChargeIndependent(psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax); if(rebinEta > 1 || rebinPhi > 1) gHist[2]->Rebin2D(rebinEta,rebinPhi); // normalization to 1 at (0,0) --> Jan Fietes method diff --git a/PWGCF/EBYE/macros/readBalanceFunction.C b/PWGCF/EBYE/macros/readBalanceFunction.C index ab2b026eeaf..a88e70b9222 100644 --- a/PWGCF/EBYE/macros/readBalanceFunction.C +++ b/PWGCF/EBYE/macros/readBalanceFunction.C @@ -6,7 +6,7 @@ const Double_t centralityArray[nrOfCentralities+1] = {0.,5.,10.,20.,30.,40.,50., const Double_t cent[nrOfCentralities] = {382.8,329.7,260.5,186.4,128.9,85.,52.8,30.,15.8}; // hard coded at the moment for centrality percentiles const Double_t centE[nrOfCentralities] = {3.1,4.6,4.4,3.9,3.3,2.6,2.0,1.3,0.6}; // (0-5,5-10,10-20,20-30,...,70-80) -void readBalanceFunction(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root",Int_t fStartBinBFWidth = 3, Int_t fRebin = 2,Int_t fStartBinBFWidthPhi = 2, Int_t fRebinPhi = 2,TString centEst = "V0M",Double_t etaWindow = -1, Int_t etaBins = -1, Bool_t correctWithEfficiency = kFALSE, Bool_t correctWithAcceptanceOnly = kFALSE) { +void readBalanceFunction(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root",Int_t fStartBinBFWidth = 3, Int_t fRebin = 2,Int_t fStartBinBFWidthPhi = 2, Int_t fRebinPhi = 2,TString centEst = "V0M",Double_t etaWindow = -1, Int_t etaBins = -1, Bool_t correctWithEfficiency = kFALSE, Bool_t correctWithAcceptanceOnly = kFALSE, Bool_t correctWithMixed = kFALSE) { // Macro to read the output of the BF analysis: MW: CHANGE THIS!!!! //i) Prints and draws the final BF output //ii) Plots the QA part of the analysis @@ -20,7 +20,7 @@ void readBalanceFunction(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResu gSystem->Load("libPWGCFebye.so"); //Draw BF - drawBF(bHistos,inFile, fStartBinBFWidth, fRebin,fStartBinBFWidthPhi, fRebinPhi,centEst, "", etaWindow,etaBins, correctWithEfficiency,correctWithAcceptanceOnly); + drawBF(bHistos,inFile, fStartBinBFWidth, fRebin,fStartBinBFWidthPhi, fRebinPhi,centEst, "", etaWindow,etaBins, correctWithEfficiency,correctWithAcceptanceOnly,correctWithMixed); //Merge the output @@ -28,7 +28,7 @@ void readBalanceFunction(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResu } //___________________________________________________________// -void drawBF(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root", Int_t fStartBinBFWidth = 1, Int_t fRebin = 1, Int_t fStartBinBFWidthPhi = 1, Int_t fRebinPhi = 1, TString centEst = "V0M",TString extraString = "", Double_t etaWindow = -1, Int_t etaBins = -1, Bool_t correctWithEfficiency = kFALSE, Bool_t correctWithAcceptanceOnly = kFALSE) { +void drawBF(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root", Int_t fStartBinBFWidth = 1, Int_t fRebin = 1, Int_t fStartBinBFWidthPhi = 1, Int_t fRebinPhi = 1, TString centEst = "V0M",TString extraString = "", Double_t etaWindow = -1, Int_t etaBins = -1, Bool_t correctWithEfficiency = kFALSE, Bool_t correctWithAcceptanceOnly = kFALSE, Bool_t correctWithMixed = kFALSE) { //Function to draw the BF objects and write them into the output file Int_t maximumCanvases = 13; @@ -43,6 +43,41 @@ void drawBF(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root", In TCanvas *cBF[13][10]; TCanvas *cBFS[13][10]; + // only one correction method allowed (correctWithEfficiency/Acceptance OR divide by event mixing) + if((correctWithEfficiency||correctWithAcceptanceOnly) && correctWithMixed){ + Printf("only one correction method allowed (correctWithEfficiency/Acceptance OR divide by event mixing)"); + return; + } + + // if divide by event mixing, first get all files for all centralities + // here the event mixing file is more or less hard coded (has to be changed if somebody else wants to use this) + TH1D *hMixedEta[nrOfCentralities][4]; + TH1D *hMixedPhi[nrOfCentralities][4]; + + if(correctWithMixed){ + TFile* fMixed[nrOfCentralities]; + + for (Int_t iCent = 0; iCent < nrOfCentralities ; iCent++){ + fMixed[iCent] = TFile::Open(Form("/Users/physics/ALICE/balance/BF_EM_new/gridOut/Centrality%.0f_%.0f_oldMethod/Histograms_WMstart%d_rebin%d_WMstartPhi%d_rebinPhi%d_EventMixing_%.0f-%.0f_AnalysisResults.root",centralityArray[iCent],centralityArray[iCent+1],fStartBinBFWidth, fRebin,fStartBinBFWidthPhi, fRebinPhi,centralityArray[iCent],centralityArray[iCent+1]),"READ"); + + if(!fMixed[iCent]){ + Printf("Event Mixing file not available!"); + return; + } + + // at the time of production the event mixing had no centrality info (all in most central bin, centarlity is selected before) + hMixedEta[iCent][0] = (TH1D*)fMixed[iCent]->Get(Form("hPN_eta_0")); + hMixedPhi[iCent][0] = (TH1D*)fMixed[iCent]->Get(Form("hPN_phi_0")); + hMixedEta[iCent][1] = (TH1D*)fMixed[iCent]->Get(Form("hNP_eta_0")); + hMixedPhi[iCent][1] = (TH1D*)fMixed[iCent]->Get(Form("hNP_phi_0")); + hMixedEta[iCent][2] = (TH1D*)fMixed[iCent]->Get(Form("hNN_eta_0")); + hMixedPhi[iCent][2] = (TH1D*)fMixed[iCent]->Get(Form("hNN_phi_0")); + hMixedEta[iCent][3] = (TH1D*)fMixed[iCent]->Get(Form("hPP_eta_0")); + hMixedPhi[iCent][3] = (TH1D*)fMixed[iCent]->Get(Form("hPP_phi_0")); + + } + } + // get the file TFile *f = TFile::Open(inFile.Data()); if(!f) { @@ -117,7 +152,7 @@ void drawBF(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root", In list = (TList*)key->ReadObj(); listName = TString(list->GetName()); - cout<<"Processing list "<SetHistNnn(a, fHistNN[iList][a]); for(iCanvas = 0; iCanvas < nrOfCentralities; iCanvas++){ + + if(correctWithMixed && a == 1) gbf[iList][iCanvas][a] = bf[iList][a]->GetBalanceFunctionHistogram(a,centralityArray[iCanvas],centralityArray[iCanvas+1],etaWindow,correctWithEfficiency,correctWithAcceptanceOnly,correctWithMixed,hMixedEta[iCanvas]); + else if(correctWithMixed && a == 6) gbf[iList][iCanvas][a] = bf[iList][a]->GetBalanceFunctionHistogram(a,centralityArray[iCanvas],centralityArray[iCanvas+1],etaWindow,correctWithEfficiency,correctWithAcceptanceOnly,correctWithMixed,hMixedPhi[iCanvas]); + else gbf[iList][iCanvas][a] = bf[iList][a]->GetBalanceFunctionHistogram(a,centralityArray[iCanvas],centralityArray[iCanvas+1],etaWindow,correctWithEfficiency,correctWithAcceptanceOnly); - gbf[iList][iCanvas][a] = bf[iList][a]->GetBalanceFunctionHistogram(a,centralityArray[iCanvas],centralityArray[iCanvas+1],etaWindow,correctWithEfficiency,correctWithAcceptanceOnly); gbf[iList][iCanvas][a]->SetName(Form("BF_%s_Cent_%.0f_%.0f_%d",gBFAnalysisType[a].Data(),centralityArray[iCanvas],centralityArray[iCanvas+1],iList)); cBF[iList][iCanvas]->cd(a+1); @@ -431,7 +469,9 @@ void drawBF(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root", In for(iCanvas = 0; iCanvas < nrOfCentralities; iCanvas++){ - gbfs[iList][iCanvas][a] = bfs[iList][a]->GetBalanceFunctionHistogram(a,centralityArray[iCanvas],centralityArray[iCanvas+1],etaWindow, correctWithEfficiency,correctWithAcceptanceOnly); + if(correctWithMixed && a == 1) gbfs[iList][iCanvas][a] = bfs[iList][a]->GetBalanceFunctionHistogram(a,centralityArray[iCanvas],centralityArray[iCanvas+1],etaWindow, correctWithEfficiency,correctWithAcceptanceOnly,correctWithMixed,hMixedEta[iCanvas]); + else if(correctWithMixed && a == 6) gbfs[iList][iCanvas][a] = bfs[iList][a]->GetBalanceFunctionHistogram(a,centralityArray[iCanvas],centralityArray[iCanvas+1],etaWindow, correctWithEfficiency,correctWithAcceptanceOnly,correctWithMixed,hMixedPhi[iCanvas]); + else gbfs[iList][iCanvas][a] = bfs[iList][a]->GetBalanceFunctionHistogram(a,centralityArray[iCanvas],centralityArray[iCanvas+1],etaWindow, correctWithEfficiency,correctWithAcceptanceOnly); gbfs[iList][iCanvas][a]->SetName(Form("BFS_%s_Cent_%.0f_%.0f_%d",gBFAnalysisType[a].Data(),centralityArray[iCanvas],centralityArray[iCanvas+1],iList)); cBFS[iList][iCanvas]->cd(a+1); @@ -541,6 +581,9 @@ void drawBF(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root", In fOut = TFile::Open(Form("Histograms_AccCorrWithEfficiency_Window%.1f_Bins%d_WMstart%d_rebin%d_WMstartPhi%d_rebinPhi%d_%s", etaWindow, etaBins, fStartBinBFWidth, fRebin,fStartBinBFWidthPhi, fRebinPhi,inFile.Data()),"RECREATE"); } } + else if(correctWithMixed){ + fOut = TFile::Open(Form("Histograms_CorrEM_Window%.1f_Bins%d_WMstart%d_rebin%d_WMstartPhi%d_rebinPhi%d_%s", etaWindow, etaBins, fStartBinBFWidth, fRebin,fStartBinBFWidthPhi, fRebinPhi,inFile.Data()),"RECREATE"); + } else{ fOut = TFile::Open(Form("Histograms_AccCorr_Window%.1f_Bins%d_WMstart%d_rebin%d_WMstartPhi%d_rebinPhi%d_%s", etaWindow, etaBins, fStartBinBFWidth, fRebin,fStartBinBFWidthPhi, fRebinPhi,inFile.Data()),"RECREATE"); } @@ -550,9 +593,11 @@ void drawBF(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root", In } fOut->cd(); for(Int_t i = 0; i < iList+1; i++){ + if(gbf[i][0][1]){ cout<<"PROCESS LIST "<ProjectionY(Form("hPN_%s_%d",gBFAnalysisType[a].Data(),i)))->Write(); @@ -619,6 +664,7 @@ void drawBF(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root", In cout<<"//=========================Centrality "<Write(); } + } } fOut->Close(); f->Close(); @@ -786,13 +833,13 @@ void GetWeightedMean(TH1D *gHistBalance, Int_t fStartBin = 1, Int_t fStopBin = - if(fStopBin > -1) fNumberOfBins = fStopBin; Double_t fP2Step = gHistBalance->GetBinWidth(1); // assume equal binning! - cout<<"=================================================="<GetBinContent(i)<<"\t Error: "<GetBinError(i)<<"\t bin: "<GetBinCenter(i)<GetBinContent(i)<<"\t Error: "<GetBinError(i)<<"\t bin: "<GetBinCenter(i)<GetBinCenter(i); gSumBi += gHistBalance->GetBinContent(i); @@ -813,10 +860,10 @@ void GetWeightedMean(TH1D *gHistBalance, Int_t fStartBin = 1, Int_t fStopBin = - Double_t delta = gSumBiXi / gSumBi; Double_t deltaError = (gSumBiXi / gSumBi) * TMath::Sqrt(TMath::Power((TMath::Sqrt(gSumXi2DeltaBi2)/gSumBiXi),2) + TMath::Power((gSumDeltaBi2/gSumBi),2) ); - cout<<"Width: "<GetNbinsX(); Double_t fP2Step = gHistBalance->GetBinWidth(1); // assume equal binning! - cout<<"=================================================="<GetBinContent(i)<<"\t Error: "<GetBinError(i)<<"\t bin: "<GetBinCenter(i)<GetBinContent(i)<<"\t Error: "<GetBinError(i)<<"\t bin: "<GetBinCenter(i)<GetBinCenter(i); gSumBi += gHistBalance->GetBinContent(i); @@ -866,10 +913,10 @@ void GetIntegral(TH1D *gHistBalance, Double_t &integ, Double_t &integE) { Double_t delta = gSumBiXi / gSumBi; Double_t deltaError = (gSumBiXi / gSumBi) * TMath::Sqrt(TMath::Power((TMath::Sqrt(gSumXi2DeltaBi2)/gSumBiXi),2) + TMath::Power((gSumDeltaBi2/gSumBi),2) ); - cout<<"Width: "<