#include <TGraphErrors.h>
#include <TGraphAsymmErrors.h>
-
#include "AliAODRecoDecayHF.h"
#include "AliRDHFCuts.h"
#include "AliRDHFCutsDplustoKpipi.h"
+#include "AliRDHFCutsDStartoKpipi.h"
#include "AliRDHFCutsD0toKpi.h"
#include "AliHFMassFitter.h"
-#include "AliEventPlaneResolution.h"
+#include "AliEventPlaneResolutionHandler.h"
+#include "AliVertexingHFUtils.h"
#endif
// Common variables: to be configured by the user
-TString filname="AnalysisResultsgoodruns.root";
-TString listname="coutputv2D0Std";//"coutputv2D0RAACuts"; //"coutputv2D0Std";
-TString filcutsname="AnalysisResultsgoodruns.root";
+TString filename="AnalysisResults_train63.root";
+TString suffix="v2Dplus3050Cut4upcutPIDTPC";
+TString partname="Dplus";
Int_t minCent=30;
Int_t maxCent=50;
-Int_t mesonPDG=421;
-
-const Int_t nFinalPtBins=3;
-Double_t ptlims[nFinalPtBins+1]={2.,5.,8.,12.};
+//evPlane flag from AliEventPlaneResolutionHandler:
+//kTPCFullEta, kTPCPosEta,kVZERO,kVZEROA,kVZEROC
+Int_t evPlane=AliEventPlaneResolutionHandler::kTPCPosEta;
+//resolution flag fromAliEventPlaneResolutionHandler:
+//kTwoRandSub,kTwoChargeSub,kTwoEtaSub,kThreeSub,kThreeSubTPCGap
+Int_t evPlaneRes=AliEventPlaneResolutionHandler::kTwoEtaSub;
+Bool_t useNcollWeight=kFALSE;
+
+const Int_t nFinalPtBins=4;
+Double_t ptlims[nFinalPtBins+1]={3.,4.,6.,8.,12.};
Double_t sigmaRangeForSig=2.5;
Double_t sigmaRangeForBkg=4.5;
-Int_t rebinHistoSideBands[nFinalPtBins]={2,2,2};
+Int_t rebinHistoSideBands[nFinalPtBins]={2,2,2,2};
Bool_t useConstantvsBkgVsMass=kFALSE;
-Int_t rebinHistov2Mass[nFinalPtBins]={2,2,5};
+Int_t rebinHistov2Mass[nFinalPtBins]={2,2,2,2};
Int_t factor4refl=0;
-Int_t minPtBin[nFinalPtBins]={-1,-1,-1};
-Int_t maxPtBin[nFinalPtBins]={-1,-1,-1};
+Int_t minPtBin[nFinalPtBins]={-1,-1,-1,-1};
+Int_t maxPtBin[nFinalPtBins]={-1,-1,-1,-1};
+Bool_t saveAllCanvas=kFALSE;
+
+// systematic errors for D+, 30-50 TPC eventplane
+Double_t systErrMeth1[nFinalPtBins]={
+ (0.247359-0.152503)/2.,
+ (0.337073-0.290978)/2.,
+ (0.286888-0.222045)/2.,
+ (-0.078888+0.178259)/2.
+};
+
+Double_t systErrMeth2[nFinalPtBins]={
+ (0.151019-0.066624)/2.,
+ (0.272252-0.232044)/2.,
+ (0.246019-0.143002)/2.,
+ (-0.026077+0.292956)/2.
+};
// systematic errors for 2-5, 5-8 and 8-12 (no run-by-run weights)
/*
*/
// systematic errors for 2-5, 5-8 and 8-12 (93 runs with run-by-run weights)
+/*
Double_t systErrMeth1[nFinalPtBins]={
(0.23-0.10)/2.,
(0.152-0.078)/2.,
(0.164-0.097)/2.,
(0.110-0.012)/2.,
(0.131-0.036)/2.
-};
+;
+*/
+
// systematic errors for 2-5, 5-8 and 8-12 (93 runs with run-by-run weights, RAA cuts)
Bool_t DefinePtBins(TDirectoryFile* df);
Double_t v2vsMass(Double_t *x, Double_t *par);
Double_t GetEPResolution(TList* lst, Double_t &rcflow, Double_t &rcfhigh);
+Double_t ComputeEventPlaneResolution(Double_t& error, TH1F* hsubev1, TH1F* hsubev2, TH1F* hsubev3);
Bool_t FitMassSpectrum(TH1F* hMass, TPad* pad);
TH1F* DoSideBands(Int_t iFinalPtBin, TH2F* hMassDphi, TH1F* hMass, Int_t rebin, TCanvas* c1, Int_t optBkg=0);
TF1* DoFitv2vsMass(Int_t iFinalPtBin, TH2F* hMassDphi, TH1F* hMass, Int_t rebin, TCanvas* c2, Int_t optErrors=0, Bool_t useConst=kTRUE);
+//______________________________________________________________________________
void Extractv2from2Dhistos(){
+ // main function: computes v2 with side band and v2(M) methods
gInterpreter->ExecuteMacro("$ALICE_ROOT/PWG3/vertexingHF/macros/LoadLibraries.C");
- TFile* filcuts=new TFile(filcutsname.Data());
- TDirectoryFile* dfcuts=(TDirectoryFile*)filcuts->Get("PWG3_D2H_HFv2");
- Bool_t binOK=DefinePtBins(dfcuts);
+ TString dirname=Form("PWGHF_D2H_HFv2_%s%s",partname.Data(),suffix.Data());
+ TString listname=Form("coutputv2%s%s",partname.Data(),suffix.Data());
+
+ TFile *f = TFile::Open(filename.Data());
+ if(!f){
+ printf("file %s not found, please check file name\n",filename.Data());
+ return;
+ }
+ TDirectoryFile* dir=(TDirectoryFile*)f->Get(dirname.Data());
+ if(!dir){
+ printf("Directory %s not found, please check dir name\n",dirname.Data());
+ return;
+ }
+ Bool_t binOK=DefinePtBins(dir);
if(!binOK){
printf("ERROR: mismatch in pt binning\n");
return;
}
-
- TFile* fil=new TFile(filname.Data());
- TDirectoryFile* df=(TDirectoryFile*)fil->Get("PWG3_D2H_HFv2");
- TList* lst=(TList*)df->Get(listname.Data());
+ TList *lst =(TList*)dir->Get(listname.Data());
if(!lst){
- printf("ERROR: list %s not found in file\n",listname.Data());
+ printf("list %s not found in file, please check list name\n",listname.Data());
return;
}
- Double_t rcfmin,rcfmax;
- Double_t resolFull=GetEPResolution(lst,rcfmin,rcfmax);
- Double_t resolSyst=(rcfmax-rcfmin)/2./resolFull;
- printf("Relative Systematic Error on RCF=%f\n",resolSyst);
+
+ // Double_t rcfmin,rcfmax;
+ // Double_t resolFull=GetEPResolution(lst,rcfmin,rcfmax);
+ // Double_t resolSyst=(rcfmax-rcfmin)/2./resolFull;
+ // printf("Relative Systematic Error on RCF=%f\n",resolSyst);
+ AliEventPlaneResolutionHandler* epres=new AliEventPlaneResolutionHandler();
+ epres->SetEventPlane(evPlane);
+ epres->SetResolutionOption(evPlaneRes);
+ epres->SetUseNcollWeights(useNcollWeight);
+ Double_t resolFull=epres->GetEventPlaneResolution(minCent,maxCent);
+ Double_t rcfmin=epres->GetEventPlaneResolution(minCent,minCent+2.5);
+ Double_t rcfmax=epres->GetEventPlaneResolution(maxCent-2.5,maxCent);
+ Double_t resolSyst=TMath::Abs(rcfmax-rcfmin)/2./resolFull;
+ delete epres;
+ printf("Event plane resolution %f\n",resolFull);
TH2F** hMassDphi=new TH2F*[nFinalPtBins];
for(Int_t iFinalPtBin=0; iFinalPtBin<nFinalPtBins; iFinalPtBin++){
v2M2[iFinalPtBin]=v2obsM2/resolFull;
errv2M2[iFinalPtBin]=errv2obsM2/resolFull;
printf("v2 = %f +- %f\n",v2M2[iFinalPtBin],errv2M2[iFinalPtBin]);
- c1[iFinalPtBin]->SaveAs(Form("cMeth1PtBin%d.root",iFinalPtBin));
- c2[iFinalPtBin]->SaveAs(Form("cMeth2PtBin%d.root",iFinalPtBin));
+ if(saveAllCanvas){
+ c1[iFinalPtBin]->SaveAs(Form("cMeth1PtBin%d.root",iFinalPtBin));
+ c2[iFinalPtBin]->SaveAs(Form("cMeth2PtBin%d.root",iFinalPtBin));
+ c1[iFinalPtBin]->SaveAs(Form("cMeth1PtBin%d.eps",iFinalPtBin));
+ c2[iFinalPtBin]->SaveAs(Form("cMeth2PtBin%d.eps",iFinalPtBin));
+ }
}
printf("\n--------- Summary ------------\n");
ent->SetTextColor(hv2m2->GetMarkerColor());
leg2->Draw();
cv2->Update();
- cv2->SaveAs("Dzero-v2-2Dmethods.gif");
+ cv2->SaveAs(Form("v2%s-2Dmethods-%s.png",partname.Data(),suffix.Data()));
+ cv2->SaveAs(Form("v2%s-2Dmethods-%s.eps",partname.Data(),suffix.Data()));
- TFile* outfil=new TFile("Dzero-v2-2Dmethods.root","recreate");
+ TFile* outfil=new TFile(Form("v2Output2DMeth_%s_%d_%d_%s.root",partname.Data(),minCent,maxCent,suffix.Data()),"recreate");
outfil->cd();
hv2m1->Write();
hv2m2->Write();
}
+//______________________________________________________________________________
Double_t v2vsMass(Double_t *x, Double_t *par){
// Fit function for signal+background
// par[0] = S/B at the mass peak
return par[3]*fracsig+par[4]*fracbkg;
}
+//______________________________________________________________________________
Double_t v2vsMassLin(Double_t *x, Double_t *par){
// Fit function for signal+background
// par[0] = S/B at the mass peak
return par[3]*fracsig+(par[4]+par[5]*x[0])*fracbkg;
}
+//______________________________________________________________________________
Bool_t FitMassSpectrum(TH1F* hMass, TPad* pad){
Int_t nMassBins=hMass->GetNbinsX();
hMaxBin=nMassBins-2;
Double_t hmin=hMass->GetBinLowEdge(hMinBin);
Double_t hmax=hMass->GetBinLowEdge(hMaxBin)+hMass->GetBinWidth(hMaxBin);
- Float_t massD=TDatabasePDG::Instance()->GetParticle(mesonPDG)->Mass();
+ Double_t massD;
+ if(partname.Contains("Dzero")) {
+ massD=TDatabasePDG::Instance()->GetParticle(421)->Mass();
+ }else if(partname.Contains("Dplus")){
+ massD=TDatabasePDG::Instance()->GetParticle(411)->Mass();
+ }else if(partname.Contains("Dstar")) {
+ massD=(TDatabasePDG::Instance()->GetParticle(413)->Mass() - TDatabasePDG::Instance()->GetParticle(421)->Mass());
+ }else{
+ printf("Wrong particle name\n");
+ massD=TDatabasePDG::Instance()->GetParticle(421)->Mass();
+ }
AliHFMassFitter* fitter=new AliHFMassFitter(hMass,hmin,hmax,2,0,0);
fitter->SetReflectionSigmaFactor(factor4refl);
return kTRUE;
}
+//______________________________________________________________________________
TH1F* DoSideBands(Int_t iFinalPtBin,
TH2F* hMassDphi,
TH1F* hMass,
}
+//______________________________________________________________________________
TF1* DoFitv2vsMass(Int_t iFinalPtBin, TH2F* hMassDphi, TH1F* hMass, Int_t rebin, TCanvas* c2, Int_t optErrors, Bool_t useConst){
Int_t npars=fSB->GetNpar();
TH1F* hAverCos2Phi=new TH1F(Form("hAverCos2PhiBin%d",iFinalPtBin),"",nbinsmass,minmass,maxmass);
TH1F* hFractionSig=new TH1F(Form("hFractionSigBin%d",iFinalPtBin),"",nbinsmass,minmass,maxmass);
TH1F* hFractionBkg=new TH1F(Form("hFractionBkgBin%d",iFinalPtBin),"",nbinsmass,minmass,maxmass);
-
for(Int_t iBin=1; iBin<=hMassDphi->GetNbinsY(); iBin++){
TH1F* htemp=(TH1F*)hMassDphi->ProjectionX("htemp",iBin,iBin);
hAverCos2Phi->SetBinContent(iBin,htemp->GetMean());
fv2->FixParameter(1,massSB);
fv2->FixParameter(2,sigmaSB);
+ printf("LOOP %d\n",rebin);
if((hAverCos2Phi->GetNbinsX()%rebin)==0){
hAverCos2Phi->Rebin(rebin);
hAverCos2Phi->Scale(1./(Double_t)rebin);
}
+//______________________________________________________________________________
Double_t GetEPResolution(TList* lst, Double_t &rcflow, Double_t &rcfhigh){
// Event plane resolution syst err (from wide centrality bin
TH1F* hResolSubAB=0x0;
+ TH1F* hResolSubBC=0x0;
+ TH1F* hResolSubAC=0x0;
Double_t xmin=1.;
Double_t xmax=-1.;
TGraphAsymmErrors* grSingle=new TGraphAsymmErrors(0);
TGraphAsymmErrors* grInteg=new TGraphAsymmErrors(0);
Int_t iPt=0;
- for(Int_t iHisC=minCent; iHisC<=maxCent-5; iHisC+=5){
- TString hisnameEP=Form("hEvPlaneResocentr%d_%d",iHisC,iHisC+5);
- TH1F* hResolSubABsing=(TH1F*)lst->FindObject(hisnameEP.Data());
- Double_t resolFull=AliEventPlaneResolution::GetFullEvResol(hResolSubABsing,1);
- Double_t resolFullmin=AliEventPlaneResolution::GetFullEvResolLowLim(hResolSubABsing,1);
- Double_t resolFullmax=AliEventPlaneResolution::GetFullEvResolHighLim(hResolSubABsing,1);
- grSingle->SetPoint(iPt,iHisC+2.5,resolFull);
- grSingle->SetPointEXlow(iPt,2.5);
- grSingle->SetPointEXhigh(iPt,2.5);
+ Int_t minCentTimesTen=minCent*10;
+ Int_t maxCentTimesTen=maxCent*10;
+ Double_t binHalfWid=25./20.;
+ for(Int_t iHisC=minCentTimesTen; iHisC<=maxCentTimesTen-25; iHisC+=25){
+ TString hisnameAB=Form("hEvPlaneResocentr%d_%d",iHisC,iHisC+25);
+ TString hisnameBC=Form("hEvPlaneReso2centr%d_%d",iHisC,iHisC+25);
+ TString hisnameAC=Form("hEvPlaneReso3centr%d_%d",iHisC,iHisC+25);
+ TH1F* hResolSubABsing=(TH1F*)lst->FindObject(hisnameAB.Data());
+ cout<<hResolSubABsing->GetName()<<endl;
+ TH1F* hResolSubBCsing=0x0;
+ TH1F* hResolSubACsing=0x0;
+ if(evPlaneRes==AliEventPlaneResolutionHandler::kThreeSub||
+ evPlaneRes==AliEventPlaneResolutionHandler::kThreeSubTPCGap){
+ hResolSubBCsing=(TH1F*)lst->FindObject(hisnameBC.Data());
+ cout<<hResolSubBCsing->GetName()<<endl;
+ hResolSubACsing=(TH1F*)lst->FindObject(hisnameAC.Data());
+ cout<<hResolSubACsing->GetName()<<endl;
+ }
+ Double_t error;
+ Double_t resolFull=ComputeEventPlaneResolution(error,hResolSubABsing,hResolSubBCsing,hResolSubACsing);
+ Double_t resolFullmin=resolFull-error;
+ Double_t resolFullmax=resolFull+error;
+
+ Double_t binCentr=(Double_t)iHisC/10.+binHalfWid;
+ grSingle->SetPoint(iPt,binCentr,resolFull);
+ grSingle->SetPointEXlow(iPt,binHalfWid);
+ grSingle->SetPointEXhigh(iPt,binHalfWid);
grSingle->SetPointEYlow(iPt,resolFullmax-resolFull);
grSingle->SetPointEYhigh(iPt,resolFull-resolFullmin);
++iPt;
if(resolFullmin<xmin) xmin=resolFullmin;
if(resolFullmax>xmax) xmax=resolFullmax;
- if(iHisC==minCent){
+ if(iHisC==minCentTimesTen){
hResolSubAB=(TH1F*)hResolSubABsing->Clone("hResolSubAB");
+ if(evPlaneRes==AliEventPlaneResolutionHandler::kThreeSub||
+ evPlaneRes==AliEventPlaneResolutionHandler::kThreeSubTPCGap){
+ hResolSubBC=(TH1F*)hResolSubBCsing->Clone("hResolSubAB");
+ hResolSubAC=(TH1F*)hResolSubACsing->Clone("hResolSubAB");
+ }
}else{
hResolSubAB->Add(hResolSubABsing);
+ if(evPlaneRes==AliEventPlaneResolutionHandler::kThreeSub||
+ evPlaneRes==AliEventPlaneResolutionHandler::kThreeSubTPCGap){
+ hResolSubBC->Add(hResolSubBCsing);
+ hResolSubAC->Add(hResolSubACsing);
+ }
}
- printf("Adding histogram %s\n",hisnameEP.Data());
+ printf("Adding histogram %s\n",hResolSubAB->GetName());
}
rcflow=xmin;
rcfhigh=xmax;
- Double_t resolSub=AliEventPlaneResolution::GetSubEvResol(hResolSubAB);
- printf("Resolution on sub event = %.4f\n",resolSub);
- Double_t chisub=AliEventPlaneResolution::FindChi(resolSub,1);
- printf("Chi Subevent = %.4f\n",chisub);
- Double_t chifull=chisub*TMath::Sqrt(2);
- printf("Chi Full Event = %.4f\n",chifull);
- Double_t resolFull=AliEventPlaneResolution::GetFullEvResol(hResolSubAB,1);
- Double_t resolFullmin=AliEventPlaneResolution::GetFullEvResolLowLim(hResolSubAB,1);
- Double_t resolFullmax=AliEventPlaneResolution::GetFullEvResolHighLim(hResolSubAB,1);
-
- AliEventPlaneResolution *resol=new AliEventPlaneResolution(1);
- resol->SetSubEventHisto(hResolSubAB);
- Double_t resolFull2=resol->GetFullEvResol();
- printf("Resolution on full event = %.4f %.4f\n",resolFull,resolFull2);
+
+ Double_t error;
+ Double_t resolFull=ComputeEventPlaneResolution(error,hResolSubAB,hResolSubBC,hResolSubAC);
grInteg->SetPoint(0,40.,resolFull);
grInteg->SetPointEXlow(0,10);
grInteg->SetPointEXhigh(0,10.);
- grInteg->SetPointEYlow(0,resolFullmax-resolFull);
- grInteg->SetPointEYhigh(0,resolFull-resolFullmin);
+ grInteg->SetPointEYlow(0,error);
+ grInteg->SetPointEYhigh(0,error);
-
TCanvas* cEP=new TCanvas("cEP","EvPlaneResol");
cEP->Divide(1,2);
cEP->cd(1);
hResolSubAB->Draw();
+ if(evPlaneRes==AliEventPlaneResolutionHandler::kThreeSub||
+ evPlaneRes==AliEventPlaneResolutionHandler::kThreeSubTPCGap){
+ hResolSubBC->SetLineColor(2);
+ hResolSubBC->Draw("same");
+ hResolSubAC->SetLineColor(4);
+ hResolSubAC->Draw("same");
+ }
TLatex* tres=new TLatex(0.15,0.7,Form("Resolution on full event = %.4f\n",resolFull));
tres->SetNDC();
tres->Draw();
}
+//______________________________________________________________________________
+Double_t ComputeEventPlaneResolution(Double_t& error, TH1F* hsubev1, TH1F* hsubev2, TH1F* hsubev3){
+ Double_t resolFull;
+ if(evPlaneRes==AliEventPlaneResolutionHandler::kTwoRandSub ||
+ evPlaneRes==AliEventPlaneResolutionHandler::kTwoChargeSub){
+ resolFull=AliVertexingHFUtils::GetFullEvResol(hsubev1);
+ error = TMath::Abs(resolFull-AliVertexingHFUtils::GetFullEvResolLowLim(hsubev1));
+ }else if(evPlaneRes==AliEventPlaneResolutionHandler::kTwoEtaSub){
+ if(evPlane==AliEventPlaneResolutionHandler::kTPCFullEta){
+ resolFull=AliVertexingHFUtils::GetFullEvResol(hsubev1);
+ error = TMath::Abs(resolFull-AliVertexingHFUtils::GetFullEvResolLowLim(hsubev1));
+ }else if(evPlane==AliEventPlaneResolutionHandler::kTPCPosEta){
+ resolFull=AliVertexingHFUtils::GetSubEvResol(hsubev1);
+ error = TMath::Abs(resolFull-AliVertexingHFUtils::GetSubEvResolLowLim(hsubev1));
+ }
+ }else{
+ Double_t resolSub[3];
+ Double_t errors[3];
+ TH1F* hevplresos[3];
+ hevplresos[0]=hsubev1;
+ hevplresos[1]=hsubev2;
+ hevplresos[2]=hsubev3;
+ for(Int_t ires=0;ires<3;ires++){
+ resolSub[ires]=hevplresos[ires]->GetMean();
+ errors[ires]=hevplresos[ires]->GetMeanError();
+ }
+ Double_t lowlim[3];for(Int_t ie=0;ie<3;ie++)lowlim[ie]=TMath::Abs(resolSub[ie]-errors[ie]);
+ if(evPlane==AliEventPlaneResolutionHandler::kVZEROC ||
+ evPlane==AliEventPlaneResolutionHandler::kVZERO){
+ resolFull=TMath::Sqrt(resolSub[1]*resolSub[2]/resolSub[0]);
+ error=resolFull-TMath::Sqrt(lowlim[2]*lowlim[1]/lowlim[0]);
+ }
+ else if(evPlane==AliEventPlaneResolutionHandler::kVZEROA){
+ resolFull=TMath::Sqrt(resolSub[0]*resolSub[2]/resolSub[1]);
+ error=resolFull-TMath::Sqrt(lowlim[2]*lowlim[0]/lowlim[1]);
+ }
+ else if(evPlane==AliEventPlaneResolutionHandler::kTPCFullEta ||
+ evPlane==AliEventPlaneResolutionHandler::kTPCPosEta){
+ resolFull=TMath::Sqrt(resolSub[0]*resolSub[1]/resolSub[2]);
+ error=resolFull-TMath::Sqrt(lowlim[0]*lowlim[1]/lowlim[2]);
+ }
+ }
+ return resolFull;
+}
+//______________________________________________________________________________
void LoadMassHistos(TList* lst, TH2F** hMassDphi){
- for(Int_t iHisC=minCent; iHisC<=maxCent-5; iHisC+=5){
+ Int_t minCentTimesTen=minCent*10;
+ Int_t maxCentTimesTen=maxCent*10;
+ for(Int_t iHisC=minCentTimesTen; iHisC<=maxCentTimesTen-25; iHisC+=25){
for(Int_t iFinalPtBin=0; iFinalPtBin<nFinalPtBins; iFinalPtBin++){
for(Int_t iPtBin=minPtBin[iFinalPtBin]; iPtBin<=maxPtBin[iFinalPtBin]; iPtBin++){
- TString hisname=Form("hMc2phi_pt%dcentr%d_%d",iPtBin,iHisC,iHisC+5);
+ TString hisname=Form("hMc2deltaphi_pt%dcentr%d_%d",iPtBin,iHisC,iHisC+25);
TH2F* htmp=(TH2F*)lst->FindObject(hisname.Data());
if(hMassDphi[iFinalPtBin]==0x0){
hMassDphi[iFinalPtBin]=(TH2F*)htmp->Clone(Form("hMassCos2DphiBin%d",iFinalPtBin));
}
}
+//______________________________________________________________________________
Bool_t DefinePtBins(TDirectoryFile* df){
- AliRDHFCuts *d0cuts;//=(AliRDHFCutsD0toKpi*)df->Get("D0toKpiCuts");
- if(mesonPDG==421)d0cuts=(AliRDHFCutsD0toKpi*)df->Get("D0toKpiCuts");
- if(mesonPDG==411)d0cuts=(AliRDHFCutsDplustoKpipi*)df->Get("AnalysisCuts");
+
+ AliRDHFCuts *d0cuts;
+ if(partname.Contains("Dzero")) d0cuts=(AliRDHFCutsD0toKpi*)df->Get(df->GetListOfKeys()->At(2)->GetName());
+ if(partname.Contains("Dplus")) d0cuts=(AliRDHFCutsDplustoKpipi*)df->Get(df->GetListOfKeys()->At(2)->GetName());
+ if(partname.Contains("Dstar")) d0cuts=((AliRDHFCutsDStartoKpipi*)df->Get(df->GetListOfKeys()->At(2)->GetName()));
+
Int_t nPtBinsCuts=d0cuts->GetNPtBins();
printf("Number of pt bins for cut object = %d\n",nPtBinsCuts);
Float_t *ptlimsCuts=d0cuts->GetPtBinLimits();
}
+//______________________________________________________________________________
void SystForSideBands(){
+ // Compute systematic ucnertainty for the side band subtraction method
gInterpreter->ExecuteMacro("$ALICE_ROOT/PWG3/vertexingHF/macros/LoadLibraries.C");
- TFile* filcuts=new TFile(filcutsname.Data());
- TDirectoryFile* dfcuts=(TDirectoryFile*)filcuts->Get("PWG3_D2H_HFv2");
- Bool_t binOK=DefinePtBins(dfcuts);
- if(!binOK) return;
+ TString dirname=Form("PWGHF_D2H_HFv2_%s%s",partname.Data(),suffix.Data());
+ TString listname=Form("coutputv2%s%s",partname.Data(),suffix.Data());
- TFile* fil=new TFile(filname.Data());
- TDirectoryFile* df=(TDirectoryFile*)fil->Get("PWG3_D2H_HFv2");
- TList* lst=(TList*)df->Get(listname.Data());
+ TFile *f = TFile::Open(filename.Data());
+ if(!f){
+ printf("file %s not found, please check file name\n",filename.Data());
+ return;
+ }
+ TDirectoryFile* dir=(TDirectoryFile*)f->Get(dirname.Data());
+ if(!dir){
+ printf("Directory %s not found, please check dir name\n",dirname.Data());
+ return;
+ }
+ Bool_t binOK=DefinePtBins(dir);
+ if(!binOK){
+ printf("ERROR: mismatch in pt binning\n");
+ return;
+ }
+ TList *lst =(TList*)dir->Get(listname.Data());
if(!lst){
- printf("ERROR: list %s not found in file\n",listname.Data());
+ printf("list %s not found in file, please check list name\n",listname.Data());
return;
}
- Double_t rcfmin,rcfmax;
- Double_t resolFull=GetEPResolution(lst,rcfmin,rcfmax);
+ AliEventPlaneResolutionHandler* epres=new AliEventPlaneResolutionHandler();
+ epres->SetEventPlane(evPlane);
+ epres->SetResolutionOption(evPlaneRes);
+ epres->SetUseNcollWeights(useNcollWeight);
+ Double_t resolFull=epres->GetEventPlaneResolution(minCent,maxCent);
+ delete epres;
+ printf("Event plane resolution %f\n",resolFull);
TH2F** hMassDphi=new TH2F*[nFinalPtBins];
for(Int_t iFinalPtBin=0; iFinalPtBin<nFinalPtBins; iFinalPtBin++){
}
}
+//______________________________________________________________________________
void SystForFitv2Mass(){
+ // Compute systematic ucnertainty for the v2(M) method
gInterpreter->ExecuteMacro("$ALICE_ROOT/PWG3/vertexingHF/macros/LoadLibraries.C");
- TFile* filcuts=new TFile(filcutsname.Data());
- TDirectoryFile* dfcuts=(TDirectoryFile*)filcuts->Get("PWG3_D2H_HFv2");
- Bool_t binOK=DefinePtBins(dfcuts);
- if(!binOK) return;
+ TString dirname=Form("PWGHF_D2H_HFv2_%s%s",partname.Data(),suffix.Data());
+ TString listname=Form("coutputv2%s%s",partname.Data(),suffix.Data());
- TFile* fil=new TFile(filname.Data());
- TDirectoryFile* df=(TDirectoryFile*)fil->Get("PWG3_D2H_HFv2");
- TList* lst=(TList*)df->Get(listname.Data());
+ TFile *f = TFile::Open(filename.Data());
+ if(!f){
+ printf("file %s not found, please check file name\n",filename.Data());
+ return;
+ }
+ TDirectoryFile* dir=(TDirectoryFile*)f->Get(dirname.Data());
+ if(!dir){
+ printf("Directory %s not found, please check dir name\n",dirname.Data());
+ return;
+ }
+ Bool_t binOK=DefinePtBins(dir);
+ if(!binOK){
+ printf("ERROR: mismatch in pt binning\n");
+ return;
+ }
+ TList *lst =(TList*)dir->Get(listname.Data());
if(!lst){
- printf("ERROR: list %s not found in file\n",listname.Data());
+ printf("list %s not found in file, please check list name\n",listname.Data());
return;
}
- Double_t rcfmin,rcfmax;
- Double_t resolFull=GetEPResolution(lst,rcfmin,rcfmax);
-
+ AliEventPlaneResolutionHandler* epres=new AliEventPlaneResolutionHandler();
+ epres->SetEventPlane(evPlane);
+ epres->SetResolutionOption(evPlaneRes);
+ epres->SetUseNcollWeights(useNcollWeight);
+ Double_t resolFull=epres->GetEventPlaneResolution(minCent,maxCent);
+ delete epres;
+ printf("Event plane resolution %f\n",resolFull);
+
TH2F** hMassDphi=new TH2F*[nFinalPtBins];
for(Int_t iFinalPtBin=0; iFinalPtBin<nFinalPtBins; iFinalPtBin++){
hMassDphi[iFinalPtBin]=0x0;