代码拉取完成,页面将自动刷新
//macro for dN/dphi,v1/v2 vs. rapidity
//write by Chenlu Hu
//Created on :
#include "RooRealVar.h"
#include "RooDataHist.h"
#include "RooGaussian.h"
#include "RooExtendPdf.h"
#include "TCanvas.h"
#include "RooPlot.h"
using namespace RooFit ;
const double pi=TMath::Pi();
void FormatFrame(RooPlot*);
void FormatFrame(TH1F*);
void HistStyle(TGraphErrors*);
void HistAxiasName(TGraphErrors*, string, string);
void HistStyle(TH1F*);
void HistAxiasName(TH1F*, string, string);
const int nrapbins = 3;
float rapbin[nrapbins+1] = {-0.2,-0.4,-0.6,-0.8};
float rapbincent[nrapbins+1] = {-0.3, -0.5, -0.7};
//double the_gausppar(double* vars, double* pars){
//return pars[0]*TMath::Gaus(vars[0],pars[1],pars[2])+
//pars[3]+pars[4]*vars[0]+pars[5]*vars[0]*vars[0]+pars[6]*vars[0]*vars[0]*vars[0];}
double the_gausppar(double* vars, double* pars){
return pars[0]*TMath::Gaus(vars[0],pars[1],pars[2])+
pars[3]+pars[4]*vars[0];}
//const int nrapbins = 2;
//float rapbin[nrapbins+1] = {-0.2,-0.5,-0.8};
//float rapbincent[nrapbins+1] = {-0.35, -0.65};
int opt=0; // without efficiency correction
//int opt=1; // with efficiency correction
int nol = 1; // normalize range left
//int nol = 2; // normalize range right
//int nol = 3; // normalize range left + right
//int nb = 0;
int nb = 1;
//int Rot=0; // bkg deuteron rotation 180 degree
int Rot=1; // bkg deuteron rotation 150-210 degree
//int Rot=2; // bkg deuteron rotation 0-180 degree
int mm=0; // no lambda mass cut
//int mm=1; // with lambda mass cut
//int nbin = 50;
//int nbin = 37;
int nbin = 25;
//int pid=3004;
//int pid=3005;
int pid=103004;
float dL=8.0;
float pidca=1.5;
//float pidca=2.0;
const int ncutphibins=4;
const char* dphi1[ncutphibins+1] = {"0", "#frac{#pi}{4}", "#frac{#pi}{2}", "#frac{3#pi}{4}", "#pi"};
const char* dphi2[ncutphibins+1] = {"0", "#frac{#pi}{8}", "#frac{#pi}{4}", "#frac{3#pi}{8}", "#frac{#pi}{2}"};
void flow_ht3_rpm_fit_3mass()
{
TH1::SetDefaultSumw2();
TH2::SetDefaultSumw2();
gStyle->SetOptFit(111);
//gStyle->SetOptFit(0);
gStyle->SetOptStat(0);
gStyle->SetStatX(.97);
gStyle->SetStatY(0.31);
float Res1;
float Res2;
TString InFile1;
TString InFile2;
TString InFile3;
TString InFile4;
if(pid==103004){
Res1=1.91; // Lambda pt:0.4-2.0
Res2=6.01;
//InFile2=Form("vnm_dNdp_noshift_ptmi1.0_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f.root",nbin,dL,pidca);
//InFile1=Form("vnmb_dNdp_noshift_ptmi1.0_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f.root",nbin,dL,pidca);
//InFile3=Form("vnmb_dNdp_noshift_ptmi1.0_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f_150210.root",nbin,dL,pidca);
//InFile4=Form("vnmb_dNdp_noshift_ptmi1.0_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f_0180.root",nbin,dL,pidca);
if(mm==0){
InFile2=Form("vnm_dNdp_noshift_ptmi1.2_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f.root",nbin,dL,pidca);
InFile1=Form("vnmb_dNdp_noshift_ptmi1.2_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f.root",nbin,dL,pidca);
InFile3=Form("vnmb_dNdp_noshift_ptmi1.2_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f_150210.root",nbin,dL,pidca);
InFile4=Form("vnmb_dNdp_noshift_ptmi1.2_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f_0180.root",nbin,dL,pidca);
}
else{
InFile2=Form("vnm_dNdp_noshift_ptmi1.2_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f_Mppi.root",nbin,dL,pidca);
InFile1=Form("vnmb_dNdp_noshift_ptmi1.2_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f_Mppi.root",nbin,dL,pidca);
InFile3=Form("vnmb_dNdp_noshift_ptmi1.2_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f_150210_Mppi.root",nbin,dL,pidca);
InFile4=Form("vnmb_dNdp_noshift_ptmi1.2_ptma2.5_rapmi-0.8_rapma-0.2_sig_103004_PID2_dE1.5_eff1_nbin%d_dL%3.1f_pidca%3.1f_0180_Mppi.root",nbin,dL,pidca);
}
}
TString Outpdf=InFile2;
if(opt==0){
if(Rot==0) Outpdf.ReplaceAll(".root",Form("_withouteff_3mass_norm%d.pdf",nol));
else if(Rot==1) Outpdf.ReplaceAll(".root",Form("_withouteff_3mass_norm%d_150210.pdf",nol));
else Outpdf.ReplaceAll(".root",Form("_withouteff_3mass_norm%d_0180.pdf",nol));
}
else{
if(Rot==0) Outpdf.ReplaceAll(".root",Form("_witheff_3mass_norm%d.pdf",nol));
else if(Rot==1) Outpdf.ReplaceAll(".root",Form("_witheff_3mass_norm%d_150210.pdf",nol));
else Outpdf.ReplaceAll(".root",Form("_witheff_3mass_norm%d_0180.pdf",nol));
}
TFile *t1 = new TFile(InFile1);
TFile *t2 = new TFile(InFile2);
TFile *t3 = new TFile(InFile3);
TFile *t4 = new TFile(InFile4);
float hyperdN1[nrapbins][ncutphibins];
float dN1[nrapbins][ncutphibins];
float dN2[nrapbins][ncutphibins];
float dN3[nrapbins][ncutphibins];
float dN1err[nrapbins][ncutphibins];
float dN2err[nrapbins][ncutphibins];
float dN3err[nrapbins][ncutphibins];
float hyperdN1err[nrapbins][ncutphibins];
float hyperdN2[nrapbins][ncutphibins];
float hyperdN2err[nrapbins][ncutphibins];
float phibincent1[ncutphibins];
float phibincent2[ncutphibins];
float phibindN[ncutphibins];
float phistep1 = pi / ncutphibins; // for v1
float phistep2 = pi / (2*ncutphibins); // for v2
for(int i=0;i<ncutphibins;i++){
phibincent1[i]=(i+0.5)*phistep1;
phibincent2[i]=(i+0.5)*phistep2;
}
TLine *li = new TLine();
li->SetLineColor(4);
li->SetLineStyle(4);
TLatex* t=new TLatex();
t->SetTextSize(0.04);
TLegend* tl = new TLegend(0.2, 0.58, 0.45, 0.9);
TLegend* leg = new TLegend(0.18, 0.55, 0.5, 0.85);
//TLegend* leg = new TLegend(0.6, 0.1, 0.95, 0.5);
leg->SetTextSize(0.04);
Float_t fxmi=-1.0;
Float_t fxma=1.0;
//Centrality 10-40%
/*
TF1* fproton = new TF1("fproton","1.0*0.351*x",fxmi,fxma);
TF1* fdeuteron = new TF1("fdeuteron","2*0.351*x",fxmi,fxma);
TF1* ftriton = new TF1("ftriton","3*0.339*x",fxmi,fxma);
TF1* fHe3 = new TF1("fHe3","3*0.323*x",fxmi,fxma);
TF1* fHe4 = new TF1("fHe4","4*0.313*x",fxmi,fxma);
*/
//Centrality 0-60%
/* 2020.7.1 from Xionghong
TF1* fproton = new TF1("fproton","0.338111*x",fxmi,fxma);
TF1* fdeuteron = new TF1("fdeuteron","0.603586*x",fxmi,fxma);
TF1* ftriton = new TF1("ftriton","0.858492*x",fxmi,fxma);
TF1* fHe3 = new TF1("fHe3","0.81407*x",fxmi,fxma);
TF1* fHe4 = new TF1("fHe4","1.05194*x",fxmi,fxma);
*/
//Centrality 0-60%
// 2020.7.9 from Xionghong
TF1* fproton = new TF1("fproton","0.331*x",fxmi,fxma);
TF1* fdeuteron = new TF1("fdeuteron","0.61*x",fxmi,fxma);
TF1* ftriton = new TF1("ftriton","0.914*x",fxmi,fxma);
TF1* fHe3 = new TF1("fHe3","0.828*x",fxmi,fxma);
TF1* fHe4 = new TF1("fHe4","1.151*x",fxmi,fxma);
TH1F* sHtm_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* bHtm_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* Htm_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* bHtm_150210_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* bHtm_0180_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TCanvas* c1 = new TCanvas("c1","c1",1000,200*nrapbins);
TCanvas* c1s = new TCanvas("c1s","c1s",1000,200*nrapbins);
c1->Print((Outpdf+"[").Data());
c1->Divide(ncutphibins,nrapbins,0.12,0.0);
c1s->Divide(ncutphibins,nrapbins,0.12,0.0);
float dNN1[nrapbins][ncutphibins];
float sN1[nrapbins][ncutphibins];
float sN1e[nrapbins][ncutphibins];
float dnn1[nrapbins][ncutphibins];
TF1 gaussian("gaussian","[0]*TMath::Gaus(x,[1],[2])",2.9925,2.9935);
TF1 gausppar("gausppar",the_gausppar,2.97,20,7);
double a=0; double b=0.5; double c=-0.2; double d=0.0;
double norm = 100; double mmean = 2.9927; double ssigma = 0.0019;
gausppar.SetParameters(norm,mmean,ssigma,a,b,c,d);
gausppar.SetParNames("Norm","Mean","Sigma","a","b","c","d");
float nx1=2.97;
float nx2=2.985;
//float nx2=2.9873;
//float nx3=2.9981;
float nx3=3.00;
float nx4=3.02;
Int_t npad1=1;
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
bHtm_150210_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t3->Get(Form("sHtm_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
bHtm_0180_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t4->Get(Form("sHtm_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t1->Get(Form("sHtm_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t2->Get(Form("sHtm_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
//bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Rebin(2);
//sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Rebin(2);
int nbin1 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx1);
int nbin2 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx2);
int nbin3 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx3);
int nbin4 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx4);
int sNb = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int bNb = bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int bNb1 = bHtm_150210_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int bNb2 = bHtm_0180_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int sNu = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin3,nbin4);
int bNu = bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin3,nbin4);
float eta;
float eta1;
float eta2;
if(nol == 1) eta = 1.0*sNb/bNb;
else if(nol == 2) eta = 1.0*sNu/bNu;
else eta = 1.0*(sNb+sNu)/(bNb+bNu);
cout<<"eta = "<<eta<<endl;
eta1 = 1.0*sNb/bNb1;
eta2 = 1.0*sNb/bNb2;
bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Scale(eta);
bHtm_150210_phi1_rap_piDca2_DecayL3[iphi][irap]->Scale(eta1);
cout<<"eta1 = "<<eta1<<endl;
bHtm_0180_phi1_rap_piDca2_DecayL3[iphi][irap]->Scale(eta2);
cout<<"eta2 = "<<eta2<<endl;
Htm_phi1_rap_piDca2_DecayL3[iphi][irap] = (TH1F*)sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Clone(Form("Htm_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->Reset();
Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->Sumw2();
if(Rot==0) Htm_phi1_rap_piDca2_DecayL3[iphi][irap] -> Add(sHtm_phi1_rap_piDca2_DecayL3[iphi][irap],bHtm_phi1_rap_piDca2_DecayL3[iphi][irap],1.0,-1.0);
else if(Rot==1) Htm_phi1_rap_piDca2_DecayL3[iphi][irap] -> Add(sHtm_phi1_rap_piDca2_DecayL3[iphi][irap],bHtm_150210_phi1_rap_piDca2_DecayL3[iphi][irap],1.0,-1.0);
else Htm_phi1_rap_piDca2_DecayL3[iphi][irap] -> Add(sHtm_phi1_rap_piDca2_DecayL3[iphi][irap],bHtm_0180_phi1_rap_piDca2_DecayL3[iphi][irap],1.0,-1.0);
int nbin = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetNbinsX();
for(int ii=1; ii<nbin+1; ii++){
float Nsbin=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinContent(ii);
float Nbbin=bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinContent(ii);
Nbbin = Nbbin/10;
Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetBinError(ii,sqrt(Nsbin+Nbbin));
}
float nmax = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetMaximum();
float xmin, xmax;
float m1, m2, m3;
float s1, s2, s3;
if(pid==103004){
//xmin=2.975; xmax=3.015; // for H3L
//xmin=2.9873; xmax=2.9981; // for H3L
xmin=2.97; xmax=3.02; // for H3L
//xmin=2.975; xmax=3.01; // for H3L
//xmin=2.97; xmax=3.01; // for H3L
m1=2.9924; m2=2.992; m3=2.9927; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
s1=0.0016; s2=0.0014; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
if(irap==0){
m1=2.9924; m2=2.9920; m3=2.9927;
s1=0.0016; s2=0.0014; s3=0.0018;
}
if(irap==1){
m1=2.9927; m2=2.9927; m3=2.9927;
s1=0.0012; s2=0.0010; s3=0.0014;
}
if(irap==2){
m1=2.9920; m2=2.9920; m3=2.9920;
s1=0.0016; s2=0.0015; s3=0.0019;
//s1=0.0016; s2=0.0015; s3=0.0020;
//s1=0.0015; s2=0.0012; s3=0.0022;
//s1=0.0014; s2=0.0014; s3=0.0014;
}
//m1=2.9925, m2=2.989, m3=2.996;
//m1=2.9927; m2=2.9926; m3=2.9928; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//m1=2.9918; m2=2.9915; m3=2.9919; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 2ybin
//m1=2.9924; m2=2.9924; m3=2.9924;
//m1=2.9922; m2=2.9922; m3=2.9922;
//m1=2.9928; m2=2.9920; m3=2.9930;
//s1=0.0021; s2=0.0016; s3=0.0024;
//s1=0.0019; s2=0.0015; s3=0.0024;
//s1=0.002; s2=0.001; s3=0.0030;
//s1=0.0018; s2=0.0015; s3=0.002; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0019; s2=0.0014; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
}
RooRealVar x("x", "x", xmin, xmax); // for H3L
RooRealVar mean("mean", "mean", m1, m2, m3);
RooRealVar sigma("sigma", "sigma", s1, s2, s3);//H4L
RooGaussian sig("sig", "sig", x, mean, sigma);
RooRealVar par0("par0", "par0", 0.5, -1, 1);
RooRealVar par1("par1", "par1", -0.2, -1, 1);
RooRealVar par2("par2", "par2", 0, -1, 1);
RooRealVar par3("par3", "par3", 0, -1, 1);
RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1, par2));
//RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1));
//RooChebychev bkg("bkg", "bkg", x, RooArgList(par0));
RooRealVar nsig("nsig","number of signal events",100, 30, 8000);
//RooRealVar nbkg("nbkg","number of bkg events", 100, 0, 10e4);
RooRealVar nbkg("nbkg","number of bkg events", 100, -100, 10e4);
RooExtendPdf esig("esig", "esig", sig, nsig);
RooExtendPdf ebkg("ebkg", "ebkg", bkg, nbkg);
RooAddPdf model("model", "model", RooArgList(esig, ebkg));
//RooAddPdf model("model", "model", RooArgList(esig));
sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(2);
sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(1);
bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
float nmax1=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetMaximum();
/*
RooDataHist dh("dh", "dh", x, Import(*(Htm_phi1_rap_piDca2_DecayL3[iphi][irap])));
//RooDataHist dh("dh", "dh", x, Import(*(bHtm_phi1_rap_piDca2_DecayL3[iphi][irap])));
model.fitTo(dh, Extended());
RooPlot* frame = x.frame(" ");
dh.plotOn(frame, Name("dh"), MarkerSize(0.6), XErrorSize(0), DrawOption("PZ"));
model.plotOn(frame, Components(esig), LineStyle(kSolid), LineColor(2), LineWidth(2));
model.plotOn(frame, Components(ebkg), LineStyle(kDashed), LineColor(4), LineWidth(2));
model.plotOn(frame, Name("model"), LineStyle(kDashed), LineColor(2), LineWidth(2));
*/
//Float_t fmean = 2.9927;
//Float_t fsig = 0.0018;
//Float_t ffac=3.;
//Float_t fsmi=fmean-ffac*fsig;
//Float_t fsma=fmean+ffac*fsig;
//TF1 f2("f2","gaus",fsmi,fsma);
//f2.SetParameters(1,2.9927,0.0018);
//Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->Fit("f2","S","",fsmi,fsma);
float gmean=mean.getVal();
float gsig=sigma.getVal();
//float tt=3.;
//float tt=2.;
float tt=1.5;
float a1=gmean-tt*gsig;
float a2=gmean+tt*gsig;
int binx1 = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(gmean);
int binx2 = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a1);
int binx3 = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a2);
float Ntt = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx1-1,binx1+1);
float xx1 = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1-1);
float xx2 = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1+2);
float xx3 = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx2);
float xx4 = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx3+1);
float Ntt1 = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx2,binx3);
float NT = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx2,binx3);
float NB = bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx2,binx3);
float NT1 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx1-1,binx1+1);
float NB1 = bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx1-1,binx1+1);
x.setRange("bkg_range",xx1,xx2);
RooAbsReal* ib = bkg.createIntegral(x,NormSet(x),Range("bkg_range"));
float Nbb = ib->getVal()*nbkg.getVal();
float Nss = Ntt-Nbb;
x.setRange("bkg1_range",xx3,xx4);
RooAbsReal* ib1 = bkg.createIntegral(x,NormSet(x),Range("bkg1_range"));
float Nbb1 = ib1->getVal()*nbkg.getVal();
if(Nbb1<0) Nbb1=0;
float Nss1 = Ntt1-Nbb1;
float eNS = 1.0*sqrt(NT+NB+Nbb1);
//float eNS = 1.0*sqrt(NT+NB);
float eNS1 = 1.0*sqrt(NT1+NB1);
//dN1[irap][iphi] = Nss;
dN1[irap][iphi] = nsig.getVal();
dN1err[irap][iphi] = nsig.getError();
//dN1err[irap][iphi] = sqrt(Ntt+Nbb);
dNN1[irap][iphi] = dN1err[irap][iphi]/dN1[irap][iphi];
if(nb==0) sN1[irap][iphi] = Ntt1;
else sN1[irap][iphi] = Nss1;
//sN1[irap][iphi] = Nss1;
//sN1[irap][iphi] = Ntt1;
sN1e[irap][iphi] = eNS;
//sN1e[irap][iphi] = sqrt(Ntt1 + Nbb1);
//sN1[irap][iphi] = Nss;
//sN1e[irap][iphi] = eNS1;
//sN1e[irap][iphi] = sqrt(Ntt + Nbb);
dnn1[irap][iphi] = sN1e[irap][iphi]/sN1[irap][iphi];
int binmi = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmin);
int binma = Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmax);
int i=0;
/*
TGraphErrors* tra1=new TGraphErrors();
for(int ibin=binmi; ibin<binma+1; ibin++){
float Nt=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinContent(ibin);
float bc=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinCenter(ibin);
float bw=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinWidth(ibin);
float xmi=bc-0.5*bw;
float xma=bc+0.5*bw;
x.setRange("bg_range",xmi,xma);
RooAbsReal* ibg = bkg.createIntegral(x,NormSet(x),Range("bg_range"));
float Nb = ibg->getVal()*nbkg.getVal();
float Ns=Nt-Nb;
float Ne=1.0*sqrt(Nt+Nb);
float errbin=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinError(ibin);
tra1->SetPoint(i,bc,Ns);
//tra1->SetPointError(i,0,Ne);
tra1->SetPointError(i,0,errbin);
i++;
}
*/
//frame->GetYaxis()->SetRangeUser(-0.08*nmax,1.1*nmax);
gausppar.SetParLimits(1,2.9922,2.9929);
gausppar.SetParLimits(2,0.0014,0.0019);
t->SetTextSize(0.05);
c1->cd(npad1);
//gStyle->SetOptFit(1110);
gPad->SetMargin(0.135, 0.01, 0., 0.);
if(irap==nrapbins-1) gPad->SetMargin(0.135, 0.01, 0.1, 0.);
//FormatFrame(frame);
//frame->Draw();
FormatFrame(Htm_phi1_rap_piDca2_DecayL3[iphi][irap]);
Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerColor(1);
Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerSize(0.8);
Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetTitle("");
Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("ep");
gausppar.SetParLimits(1,2.9922,2.9929);
gausppar.SetParLimits(2,0.0012,0.0019);
Htm_phi1_rap_piDca2_DecayL3[iphi][irap]->Fit(&gausppar, "S");
double aa1 = gausppar.GetParameter(3);
double bb1 = gausppar.GetParameter(4);
//double cc1 = gausppar.GetParameter(5);
//double dd1 = gausppar.GetParameter(6);
double mmm = gausppar.GetParameter(1);
double sss = gausppar.GetParameter(2);
double xb = mmm - 2.0*sss;
double xu = mmm + 2.0*sss;
//TF1* poly = new TF1("poly","[0]+[1]*x+[2]*x*x+[3]*x*x*x",2.97,23.02);
TF1* poly = new TF1("poly","[0]+[1]*x",2.97,23.02);
//poly->SetParameters(aa1,bb1,cc1,dd1);
poly->SetParameters(aa1,bb1);
//poly->SetParameters(1,20,0,0);
poly->SetLineStyle(2);
poly->SetLineColor(4);
poly->Draw("same");
li->SetLineStyle(1);
li->SetLineColor(1);
li->DrawLine(2.97,0,3.02,0);
li->SetLineStyle(2);
li->SetLineColor(4);
double NNt = gausppar.Integral(xb,xu);
NNt = NNt/0.002;
double NNb = poly->Integral(xb,xu);
NNb = NNb/0.002;
double NNs = NNt - NNb;
t->DrawLatexNDC(0.2,0.7,Form("Nsig = %3.1f",NNs));
int binb = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xb);
int binu = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xu);
float edg1 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binb);
float edg2 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binu+1);
NT = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binb,binu);
NB = bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binb,binu);
//eNS = 1.0*sqrt(NT + NB);
eNS = 1.0*sqrt(NT + NB/10.0);
li->DrawLine(edg1, -0.08*nmax, edg1, 1.1*nmax);
li->DrawLine(edg2, -0.08*nmax, edg2, 1.1*nmax);
sN1[irap][iphi] = NNs;
sN1e[irap][iphi] = eNS;
dnn1[irap][iphi] = sN1e[irap][iphi]/sN1[irap][iphi];
//li->DrawLine(xx3, -0.08*nmax, xx3, 1.1*nmax);
//li->DrawLine(xx4, -0.08*nmax, xx4, 1.1*nmax);
//li->DrawLine(nx1, -0.08*nmax, nx1, 1.1*nmax);
//li->DrawLine(nx2, -0.08*nmax, nx2, 1.1*nmax);
//li->DrawLine(nx3, -0.08*nmax, nx3, 1.1*nmax);
//li->DrawLine(nx4, -0.08*nmax, nx4, 1.1*nmax);
//tra1->SetMarkerStyle(20);
//tra1->SetMarkerSize(0.6);
//tra1->SetMarkerColor(1);
//tra1->GetXaxis()->SetLimits(xmin,xmax);
//tra1->Draw("psame");
//t->DrawLatexNDC(0.21,0.55,Form("Ns = %4.1f, eNs = %4.1f, dN/N = %5.4f",dN1[irap][iphi],dN1err[irap][iphi],dNN1[irap][iphi]));
//t->DrawLatexNDC(0.21,0.45,Form("3Ns = %4.1f, 3eNs = %4.1f, 3dN/N = %5.4f",sN1[irap][iphi],sN1e[irap][iphi],dnn1[irap][iphi]));
//t->DrawLatexNDC(0.22,0.45,Form("dN/N = %5.4f",dNN1[irap][iphi]));
TPaveText *label_rap = new TPaveText(0.22,0.85,0.52,1.0,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.07);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
//TPaveText *label_phi = new TPaveText(0.7,0.3,0.9,0.4,"blNDC");
TPaveText *label_phi = new TPaveText(0.7,0.85,0.9,1.0,"blNDC");
label_phi->SetBorderSize(0);
label_phi->SetFillStyle(0);
label_phi->SetTextFont(62);
label_phi->SetTextSize(0.07);
// label_phi->SetTextColor(2);
label_phi->AddText(Form("%s<#phi-#Psi_{1}<%s", dphi1[iphi], dphi1[iphi+1]));
label_phi->Draw();
float chi2;
//chi2=frame->chiSquare("model", "dh");
t->SetTextSize(0.05);
//t->DrawLatexNDC(0.72,0.4,Form("#chi^{2} = %3.4f",chi2));
//t->DrawLatexNDC(0.22,0.75,Form("#mu = %3.4f",gmean));
//t->DrawLatexNDC(0.22,0.65,Form("#sigma = %3.4f",gsig));
//cout << frame->chiSquare("model", "dh") << endl;
cout <<"par0 = "<< par0.getVal() << endl;
//label->AddText(Form("#chi^{2}/dof = ", frame->chiSquare("model", "dh")));
//label->Draw();
c1s->cd(npad1);npad1++;
gPad->SetMargin(0.135, 0.01, 0., 0.);
if(irap==nrapbins-1) gPad->SetMargin(0.135, 0.01, 0.1, 0.);
FormatFrame(sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]);
FormatFrame(bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]);
sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetXaxis()->SetRangeUser(2.97,3.02);
sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerColor(2);
//sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetXaxis()->SetTitle("Inv.mass(#pi^{-}-d-p) [GeV/c^{2}]");
//sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetYaxis()->SetTitle("Counts");
sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerSize(0.8);
sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("e");
bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(1);
//bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("histsame");
bHtm_150210_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(4);
bHtm_150210_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("histsame");
bHtm_0180_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(6);
bHtm_0180_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("histsame");
li->DrawLine(xx3, 0, xx3, 1.1*nmax1);
li->DrawLine(xx4, 0, xx4, 1.1*nmax1);
if(irap==0&&iphi==0)
{
leg->AddEntry(sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]," N_{tot}","PE");
leg->AddEntry(bHtm_phi1_rap_piDca2_DecayL3[iphi][irap]," N_{bkg}_Rot(180)","l");
leg->AddEntry(bHtm_150210_phi1_rap_piDca2_DecayL3[iphi][irap]," N_{bkg}_Rot(150,210)","l");
leg->AddEntry(bHtm_0180_phi1_rap_piDca2_DecayL3[iphi][irap]," N_{bkg}_Rot(0,180)","l");
leg->Draw("same");
}
label_rap->Draw();
label_phi->Draw();
}
}
c1s->Print(Outpdf);
c1->Print(Outpdf);
TH1F* sHtmass_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* bHtmass_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* Htmass_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* bHtmass_150210_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* bHtmass_0180_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
float dNN2[nrapbins][ncutphibins];
float sN2[nrapbins][ncutphibins];
float sN2e[nrapbins][ncutphibins];
float dnn2[nrapbins][ncutphibins];
TCanvas* c2 = new TCanvas("c2","c2",1000,200*nrapbins);
TCanvas* c2s = new TCanvas("c2s","c2s",1000,200*nrapbins);
c2->Divide(ncutphibins,nrapbins,0.12,0.0);
c2s->Divide(ncutphibins,nrapbins,0.12,0.0);
Int_t npad2=1;
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t2->Get(Form("sHtmass_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
bHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t1->Get(Form("sHtmass_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
bHtmass_150210_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t3->Get(Form("sHtmass_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
bHtmass_0180_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t4->Get(Form("sHtmass_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
//sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Rebin(2);
//bHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Rebin(2);
int nbin1 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx1);
int nbin2 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx2);
int nbin3 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx3);
int nbin4 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx4);
int sNb = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int bNb = bHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int bNb1 = bHtmass_150210_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int bNb2 = bHtmass_0180_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int sNu = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin3,nbin4);
int bNu = bHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin3,nbin4);
float eta;
float eta1;
float eta2;
if(nol == 1) eta = 1.0*sNb/bNb;
else if(nol == 2)eta = 1.0*sNu/bNu;
else eta = 1.0*(sNb+sNu)/(bNb+bNu);
cout<<"eta = "<<eta<<endl;
eta1 = 1.0*sNb/bNb1;
eta2 = 1.0*sNb/bNb2;
bHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Scale(eta);
bHtmass_150210_phi1_rap_piDca2_DecayL3[iphi][irap]->Scale(eta1);
bHtmass_0180_phi1_rap_piDca2_DecayL3[iphi][irap]->Scale(eta2);
Htmass_phi1_rap_piDca2_DecayL3[iphi][irap] = (TH1F*)sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Clone(Form("Htmass_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Reset();
Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Sumw2();
if(Rot==0) Htmass_phi1_rap_piDca2_DecayL3[iphi][irap] -> Add(sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap],bHtmass_phi1_rap_piDca2_DecayL3[iphi][irap],1.0,-1.0);
else if(Rot==1) Htmass_phi1_rap_piDca2_DecayL3[iphi][irap] -> Add(sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap],bHtmass_150210_phi1_rap_piDca2_DecayL3[iphi][irap],1.0,-1.0);
else Htmass_phi1_rap_piDca2_DecayL3[iphi][irap] -> Add(sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap],bHtmass_0180_phi1_rap_piDca2_DecayL3[iphi][irap],1.0,-1.0);
float nmax = Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetMaximum();
float nmax1 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetMaximum();
float xmin, xmax;
float m1, m2, m3;
float s1, s2, s3;
if(pid==103004){
xmin=2.97; xmax=3.02; // for H3L
//xmin=2.975; xmax=3.01; // for H3L
//xmin=2.975; xmax=3.015; // for H3L
m1=2.9924; m2=2.992; m3=2.9927; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
s1=0.0016; s2=0.0014; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
if(irap==0){
m1=2.9924; m2=2.9920; m3=2.9927;
s1=0.0016; s2=0.0014; s3=0.0018;
}
if(irap==1){
m1=2.9927; m2=2.9927; m3=2.9927;
s1=0.0014; s2=0.0014; s3=0.0014;
}
if(irap==2){
m1=2.9924; m2=2.9924; m3=2.9924;
s1=0.0014; s2=0.0014; s3=0.0014;
//s1=0.0015; s2=0.0014; s3=0.0022;
}
//m1=2.9925, m2=2.989, m3=2.996;
//m1=2.9927; m2=2.9926; m3=2.9928; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//m1=2.9927; m2=2.9927; m3=2.9927; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//m1=2.9918; m2=2.9915; m3=2.9919; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 2ybin
//m1=2.9924; m2=2.9924; m3=2.9924;
//m1=2.9922; m2=2.9922; m3=2.9922;
//m1=2.9928; m2=2.9920; m3=2.9930;
//s1=0.0021; s2=0.0016; s3=0.0024;
//s1=0.0019; s2=0.0015; s3=0.0024;
//s1=0.002; s2=0.001; s3=0.0030;
//s1=0.0016; s2=0.0014; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0018; s2=0.0015; s3=0.002; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0019; s2=0.0014; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
}
RooRealVar x("x", "x", xmin, xmax); // for H3L
RooRealVar mean("mean", "mean", m1, m2, m3);
RooRealVar sigma("sigma", "sigma", s1, s2, s3);//H4L
RooGaussian sig("sig", "sig", x, mean, sigma);
RooRealVar par0("par0", "par0", 0.5, -1, 1);
RooRealVar par1("par1", "par1", -0.2, -1, 1);
RooRealVar par2("par2", "par2", 0, -1, 1);
RooRealVar par3("par3", "par3", 0, -1, 1);
RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1, par2));
//RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1));
//RooChebychev bkg("bkg", "bkg", x, RooArgList(par0));
RooRealVar nsig("nsig","number of signal events",100, 0, 8000);
RooRealVar nbkg("nbkg","number of bkg events", 100, 0, 10e4);
RooExtendPdf esig("esig", "esig", sig, nsig);
RooExtendPdf ebkg("ebkg", "ebkg", bkg, nbkg);
RooAddPdf model("model", "model", RooArgList(esig, ebkg));
sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(2);
sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
/*
RooDataHist dh("dh", "dh", x, Import(*(Htmass_phi1_rap_piDca2_DecayL3[iphi][irap])));
model.fitTo(dh, Extended());
RooPlot* frame = x.frame(" ");
dh.plotOn(frame, Name("dh"), MarkerSize(0.5), XErrorSize(0), DrawOption("PZ"));
model.plotOn(frame, Components(esig), LineStyle(kSolid), LineColor(2), LineWidth(2));
model.plotOn(frame, Components(ebkg), LineStyle(kDashed), LineColor(4), LineWidth(2));
model.plotOn(frame, Name("model"), LineStyle(kDashed), LineColor(2), LineWidth(2));
*/
float gmean=mean.getVal();
float gsig=sigma.getVal();
//float tt=3.;
//float tt=2.;
float tt=1.5;
float a1=gmean-tt*gsig;
float a2=gmean+tt*gsig;
int binx1 = Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(gmean);
int binx2 = Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a1);
int binx3 = Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a2);
float Ntt = Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx1-1,binx1+1);
float xx1 = Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1-1);
float xx2 = Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1+2);
float xx3 = Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx2);
float xx4 = Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx3+1);
float Ntt1=0;
for(int ibin=binx2; ibin<binx3+1; ibin++){
Ntt1 += Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinContent(ibin);
}
//float Ntt1 = Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx2,binx3);
x.setRange("bkg_range",xx1,xx2);
RooAbsReal* ib = bkg.createIntegral(x,NormSet(x),Range("bkg_range"));
float Nbb = ib->getVal()*nbkg.getVal();
float Nss = Ntt-Nbb;
x.setRange("bkg1_range",xx3,xx4);
RooAbsReal* ib1 = bkg.createIntegral(x,NormSet(x),Range("bkg1_range"));
float Nbb1 = ib1->getVal()*nbkg.getVal();
if(Nbb<0) Nbb1=0;
float Nss1 = Ntt1-Nbb1;
//dN2[irap][iphi] = Nss;
dN2[irap][iphi] = nsig.getVal();
dN2err[irap][iphi] = dN2[irap][iphi]*dnn1[irap][iphi];
//dN2err[irap][iphi] = nsig.getError();
//dN2err[irap][iphi] = sqrt(Ntt+Nbb);
//dN2err[irap][iphi] = dN2[irap][iphi]*dNN1[irap][iphi];
dNN2[irap][iphi]=dN2err[irap][iphi]/dN2[irap][iphi];
if(nb==0) sN2[irap][iphi] = Ntt1;
else sN2[irap][iphi] = Nss1;
//sN2[irap][iphi] = Ntt1;
//sN2[irap][iphi] = Nss1;
//sN2[irap][iphi] = Nss;
//sN2[irap][iphi] = Nss1;
//dN1[irap][iphi] = nsig.getVal();
//dN1err[irap][iphi] = nsig.getError();
//sN2e[irap][iphi] = sqrt(Ntt1+Nbb1);
sN2e[irap][iphi] = sN2[irap][iphi]*dnn1[irap][iphi];
dnn2[irap][iphi] = sN2e[irap][iphi]/sN2[irap][iphi];
//t->DrawLatexNDC(0.22,0.65,Form("Ns = %f, eNs = %f",dN1[irap][iphi],dN1err[irap][iphi]));
int binmi=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmin);
int binma=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmax);
/*
TGraphErrors* tra2=new TGraphErrors();
int i=1;
for(int ibin=binmi; ibin<binma+1; ibin++){
float Nt=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinContent(ibin);
float bc=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinCenter(ibin);
float bw=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinWidth(ibin);
float xmi=bc-0.5*bw;
float xma=bc+0.5*bw;
x.setRange("bg_range",xmi,xma);
RooAbsReal* ibg = bkg.createIntegral(x,NormSet(x),Range("bg_range"));
float Nb = ibg->getVal()*nbkg.getVal();
float Ns=Nt-Nb;
float Ne=1.0*sqrt(Nt+Nb);
float errbin=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinError(ibin);
tra2->SetPoint(i,bc,Ns);
//tra2->SetPointError(i,0,Ne);
tra2->SetPointError(i,0,errbin);
i++;
}
*/
//frame->GetYaxis()->SetRangeUser(-0.08*nmax,1.1*nmax);
gausppar.SetParLimits(1,2.9922,2.9929);
gausppar.SetParLimits(2,0.0014,0.0019);
t->SetTextSize(0.05);
c2->cd(npad2);
gPad->SetMargin(0.135, 0.01, 0., 0.);
if(irap==nrapbins-1) gPad->SetMargin(0.135, 0.01, 0.1, 0.);
//FormatFrame(frame);
//frame->Draw();
FormatFrame(Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]);
Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerColor(1);
Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerSize(0.8);
Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetTitle("");
Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("e");
gausppar.SetParLimits(1,2.9922,2.9929);
gausppar.SetParLimits(2,0.0012,0.0019);
Htmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Fit(&gausppar, "S");
double aa1 = gausppar.GetParameter(3);
double bb1 = gausppar.GetParameter(4);
//double cc1 = gausppar.GetParameter(5);
//double dd1 = gausppar.GetParameter(6);
double mmm = gausppar.GetParameter(1);
double sss = gausppar.GetParameter(2);
double xb = mmm - 2.0*sss;
double xu = mmm + 2.0*sss;
//TF1* poly = new TF1("poly","[0]+[1]*x+[2]*x*x+[3]*x*x*x",2.97,23.02);
TF1* poly = new TF1("poly","[0]+[1]*x",2.97,23.02);
//poly->SetParameters(aa1,bb1,cc1,dd1);
poly->SetParameters(aa1,bb1);
//poly->SetParameters(1,20,0,0);
poly->SetLineStyle(2);
poly->SetLineColor(4);
poly->Draw("same");
li->SetLineStyle(1);
li->SetLineColor(1);
li->DrawLine(2.97,0,3.02,0);
li->SetLineStyle(2);
li->SetLineColor(4);
double NNt = gausppar.Integral(xb,xu);
NNt = NNt/0.002;
double NNb = poly->Integral(xb,xu);
NNb = NNb/0.002;
double NNs = NNt - NNb;
t->DrawLatexNDC(0.2,0.7,Form("Nsig = %3.1f",NNs));
int binb = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xb);
int binu = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xu);
float edg1 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binb);
float edg2 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binu+1);
float NT = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binb,binu);
float NB = bHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binb,binu);
//float eNS = 1.0*sqrt(NT + NB);
float eNS = 1.0*sqrt(NT + NB/10.0);
li->DrawLine(edg1, -0.08*nmax, edg1, 1.1*nmax);
li->DrawLine(edg2, -0.08*nmax, edg2, 1.1*nmax);
sN2[irap][iphi] = NNs;
sN2e[irap][iphi] = sN2[irap][iphi]*dnn1[irap][iphi];
dnn2[irap][iphi] = sN1e[irap][iphi]/sN1[irap][iphi];
li->SetLineColor(4);
//li->DrawLine(xx3, -0.08*nmax, xx3, 1.1*nmax);
//li->DrawLine(xx4, -0.08*nmax, xx4, 1.1*nmax);
//li->DrawLine(nx1, -0.08*nmax, nx1, 1.1*nmax);
//li->DrawLine(nx2, -0.08*nmax, nx2, 1.1*nmax);
//li->DrawLine(nx3, -0.08*nmax, nx3, 1.1*nmax);
//li->DrawLine(nx4, -0.08*nmax, nx4, 1.1*nmax);
/*
tra2->SetMarkerStyle(20);
tra2->SetMarkerSize(0.6);
tra2->SetMarkerColor(1);
tra2->GetXaxis()->SetLimits(xmin,xmax);
tra2->Draw("psame");
*/
//t->DrawLatexNDC(0.21,0.55,Form("Ns = %4.1f, eNs = %4.1f, dN/N = %5.4f",dN2[irap][iphi],dN2err[irap][iphi],dNN2[irap][iphi]));
//t->DrawLatexNDC(0.21,0.45,Form("3Ns = %4.1f, 3eNs = %4.1f, 3dN/N = %5.4f",sN2[irap][iphi],sN2e[irap][iphi],dnn2[irap][iphi]));
//t->DrawLatexNDC(0.21,0.35,Form("Ntt1 = %5.1f, Nbb1 = %5.1f, Nss1 = %5.1f",Ntt1,Nbb1,Nss1));
//t->DrawLatexNDC(0.21,0.25,Form("Nball = %5.1f",nbkg.getVal()));
//t->DrawLatexNDC(0.22,0.45,Form("dN/N = %5.4f",dNN2[irap][iphi]));
TPaveText *label_rap = new TPaveText(0.22,0.85,0.52,1.0,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.07);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
//TPaveText *label_phi = new TPaveText(0.7,0.3,0.9,0.4,"blNDC");
TPaveText *label_phi = new TPaveText(0.7,0.85,0.9,1.0,"blNDC");
label_phi->SetBorderSize(0);
label_phi->SetFillStyle(0);
label_phi->SetTextFont(62);
label_phi->SetTextSize(0.07);
// label_phi->SetTextColor(2);
label_phi->AddText(Form("%s<#phi-#Psi_{1}<%s", dphi1[iphi], dphi1[iphi+1]));
label_phi->Draw();
float chi2;
//chi2=frame->chiSquare("model", "dh");
t->SetTextSize(0.05);
//t->DrawLatexNDC(0.72,0.4,Form("#chi^{2} = %3.4f",chi2));
//t->DrawLatexNDC(0.22,0.75,Form("#mu = %3.4f",gmean));
//t->DrawLatexNDC(0.22,0.65,Form("#sigma = %3.4f",gsig));
//cout << frame->chiSquare("model", "dh") << endl;
//cout <<"par0 = "<< par0.getVal() << endl;
c2s->cd(npad2);npad2++;
gPad->SetMargin(0.135, 0.01, 0., 0.);
if(irap==nrapbins-1) gPad->SetMargin(0.135, 0.01, 0.1, 0.);
FormatFrame(sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]);
FormatFrame(bHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]);
sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetXaxis()->SetRangeUser(2.97,3.02);
sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerColor(2);
//sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetXaxis()->SetTitle("Inv.mass(#pi^{-}-d-p) [GeV/c^{2}]");
//sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetYaxis()->SetTitle("Counts");
sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerSize(0.8);
sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("e");
bHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(1);
//bHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("histsame");
bHtmass_150210_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(4);
bHtmass_150210_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("histsame");
bHtmass_0180_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(6);
bHtmass_0180_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("histsame");
label_rap->Draw();
label_phi->Draw();
li->DrawLine(xx3, 0, xx3, 1.1*nmax1);
li->DrawLine(xx4, 0, xx4, 1.1*nmax1);
if(irap==0&&iphi==0)
leg->Draw("same");
}
}
c2s->Print(Outpdf);
c2->Print(Outpdf);
TH1F* sHtmass_eff_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* bHtmass_eff_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* Htmass_eff_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* bHtmass_150210_eff_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
TH1F* bHtmass_0180_eff_phi1_rap_piDca2_DecayL3[ncutphibins][nrapbins];
float dNN3[nrapbins][ncutphibins];
float sN3[nrapbins][ncutphibins];
float sN3e[nrapbins][ncutphibins];
float dnn3[nrapbins][ncutphibins];
TCanvas* c3 = new TCanvas("c3","c3",1000,200*nrapbins);
TCanvas* c3s = new TCanvas("c3s","c3s",1000,200*nrapbins);
c3->Divide(ncutphibins,nrapbins,0.12,0.0);
c3s->Divide(ncutphibins,nrapbins,0.12,0.0);
Int_t npad3=1;
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t2->Get(Form("sHtmass_eff_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
bHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t1->Get(Form("sHtmass_eff_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
bHtmass_150210_eff_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t3->Get(Form("sHtmass_eff_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
bHtmass_0180_eff_phi1_rap_piDca2_DecayL3[iphi][irap]=(TH1F*)t4->Get(Form("sHtmass_eff_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
int nbin1 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx1);
int nbin2 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx2);
int nbin3 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx3);
int nbin4 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(nx4);
int sNb = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int bNb = bHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int bNb1 = bHtmass_150210_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int bNb2 = bHtmass_0180_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin1,nbin2);
int sNu = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin3,nbin4);
int bNu = bHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(nbin3,nbin4);
float eta,eta1,eta2;
if(nol == 1) eta = 1.0*sNb/bNb;
else if(nol == 2)eta = 1.0*sNu/bNu;
else eta = 1.0*(sNb+sNu)/(bNb+bNu);
eta1 = 1.0*sNb/bNb1;
eta2 = 1.0*sNb/bNb2;
bHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Scale(eta);
bHtmass_150210_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Scale(eta1);
bHtmass_0180_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Scale(eta2);
Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap] = (TH1F*)sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Clone(Form("Htmass_eff_phi1_rap_piDca2_DecayL3_phi1%d_rap%d",iphi,irap));
Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Reset();
Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Sumw2();
if(Rot==0) Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap] -> Add(sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap],bHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap],1.0,-1.0);
else if(Rot==1) Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap] -> Add(sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap],bHtmass_150210_eff_phi1_rap_piDca2_DecayL3[iphi][irap],1.0,-1.0);
else Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap] -> Add(sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap],bHtmass_0180_eff_phi1_rap_piDca2_DecayL3[iphi][irap],1.0,-1.0);
float nmax = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetMaximum();
float nmax1 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetMaximum();
float xmin, xmax;
float m1, m2, m3;
float s1, s2, s3;
if(pid==103004){
//xmin=2.975; xmax=3.01; // for H3L
xmin=2.97; xmax=3.02; // for H3L
//m1=2.9925, m2=2.989, m3=2.996;
//m1=2.9927; m2=2.9926; m3=2.9928; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//m1=2.993; m2=2.9928; m3=2.9940; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
m1=2.9927; m2=2.9927; m3=2.9927; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//m1=2.9918; m2=2.9915; m3=2.9919; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 2ybin
//m1=2.9924; m2=2.9924; m3=2.9924;
//m1=2.9922; m2=2.9922; m3=2.9922;
//m1=2.9928; m2=2.9920; m3=2.9930;
//s1=0.0021; s2=0.0016; s3=0.0024;
//s1=0.0019; s2=0.0015; s3=0.0024;
//s1=0.002; s2=0.001; s3=0.0030;
s1=0.0018; s2=0.0017; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0018; s2=0.0015; s3=0.002; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0019; s2=0.0014; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
}
RooRealVar x("x", "x", xmin, xmax); // for H3L
RooRealVar mean("mean", "mean", m1, m2, m3);
RooRealVar sigma("sigma", "sigma", s1, s2, s3);//H4L
RooGaussian sig("sig", "sig", x, mean, sigma);
RooRealVar par0("par0", "par0", 0.5, -1, 1);
RooRealVar par1("par1", "par1", -0.2, -1, 1);
RooRealVar par2("par2", "par2", 0, -1, 1);
RooRealVar par3("par3", "par3", 0, -1, 1);
RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1, par2));
//RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1));
RooRealVar nsig("nsig","number of signal events",200, 0, 8000000);
RooRealVar nbkg("nbkg","number of bkg events", 10e4,100,10e8);
RooExtendPdf esig("esig", "esig", sig, nsig);
RooExtendPdf ebkg("ebkg", "ebkg", bkg, nbkg);
RooAddPdf model("model", "model", RooArgList(esig, ebkg));
Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(2);
Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
RooPlot* frame = x.frame(" ");
//RooDataHist dh1("dh1", "dh1", x, Import(*(sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap])));
RooDataHist dh1("dh1", "dh1", x, Import(*(Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap])));
model.fitTo(dh1, Extended());
dh1.plotOn(frame, Name("dh1"), MarkerSize(0.5), XErrorSize(0), DrawOption("PZ"));
model.plotOn(frame, Components(esig), LineStyle(kSolid), LineColor(2), LineWidth(2));
model.plotOn(frame, Components(ebkg), LineStyle(kDashed), LineColor(4), LineWidth(2));
model.plotOn(frame, Name("model"), LineStyle(kDashed), LineColor(2), LineWidth(2));
float gmean=mean.getVal();
float gsig=sigma.getVal();
//float tt=3.;
//float tt=2.;
float tt=1.5;
float a1=gmean-tt*gsig;
float a2=gmean+tt*gsig;
int binx1 = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(gmean);
int binx2 = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a1);
int binx3 = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a2);
float Ntt = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx1-1,binx1+1);
float xx1 = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1-1);
float xx2 = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1+2);
float xx3 = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx2);
float xx4 = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx3+1);
float Ntt1 = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx2,binx3);
float NT = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx2,binx3);
float NB = bHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx2,binx3);
float NT1 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx1-1,binx1+1);
float NB1 = bHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx1-1,binx1+1);
x.setRange("bkg_range",xx1,xx2);
RooAbsReal* ib = bkg.createIntegral(x,NormSet(x),Range("bkg_range"));
float Nbb = ib->getVal()*nbkg.getVal();
float Nss = Ntt-Nbb;
x.setRange("bkg1_range",xx3,xx4);
RooAbsReal* ib1 = bkg.createIntegral(x,NormSet(x),Range("bkg1_range"));
float Nbb1 = ib1->getVal()*nbkg.getVal();
if(Nbb1<0) Nbb1=0;
float Nss1 = Ntt1-Nbb1;
float eNS = 1.0*sqrt(NT+NB+Nbb1);
//float eNS = 1.0*sqrt(NT+NB);
float eNS1 = 1.0*sqrt(NT1+NB1);
sN3[0][0] = sN2[0][0]/0.0297078;
sN3[0][1] = sN2[0][1]/0.0324018;
sN3[0][2] = sN2[0][2]/0.0323103;
sN3[0][3] = sN2[0][3]/0.0324594;
sN3[1][0] = sN2[1][0]/0.0460016;
sN3[1][1] = sN2[1][1]/0.0475578;
sN3[1][2] = sN2[1][2]/0.0414245;
sN3[1][3] = sN2[1][3]/0.0468818;
sN3[2][0] = sN2[2][0]/0.047831;
sN3[2][1] = sN2[2][1]/0.0490432;
sN3[2][2] = sN2[2][2]/0.0480522;
sN3[2][3] = sN2[2][3]/0.0532116;
sN3e[0][0] = sN2e[0][0]/0.0297078;
sN3e[0][1] = sN2e[0][1]/0.0324018;
sN3e[0][2] = sN2e[0][2]/0.0323103;
sN3e[0][3] = sN2e[0][3]/0.0324594;
sN3e[1][0] = sN2e[1][0]/0.0460016;
sN3e[1][1] = sN2e[1][1]/0.0475578;
sN3e[1][2] = sN2e[1][2]/0.0414245;
sN3e[1][3] = sN2e[1][3]/0.0468818;
sN3e[2][0] = sN2e[2][0]/0.047831;
sN3e[2][1] = sN2e[2][1]/0.0490432;
sN3e[2][2] = sN2e[2][2]/0.0480522;
sN3e[2][3] = sN2e[2][3]/0.0532116;
if(nb==0) dN3[irap][iphi] = Ntt1;
else dN3[irap][iphi] = Nss1;
//dN3[irap][iphi] = Ntt1;
//dN3[irap][iphi] = Nss1;
dN3err[irap][iphi] = dN3[irap][iphi]*dnn1[irap][iphi];
dNN3[irap][iphi] = dN3err[irap][iphi]/dN3[irap][iphi];
//dN3[irap][iphi] = Nss;
//dN3err[irap][iphi] = nsig.getError();
//dN3err[irap][iphi] = sqrt(Ntt+Nbb);
//sN3[irap][iphi] = Nss1;
//sN3e[irap][iphi] = sN3[irap][iphi]*dnn1[irap][iphi];
//dnn3[irap][iphi] = sN3e[irap][iphi]/sN3[irap][iphi];
//sN3e[irap][iphi] = sqrt(Ntt3 + Nbb3);
//sN3[irap][iphi] = Nss;
//sN3e[irap][iphi] = eNS3;
//sN3e[irap][iphi] = sqrt(Ntt + Nbb);
int binmi = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmin);
int binma = Htmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmax);
frame->GetYaxis()->SetRangeUser(-0.08*nmax,1.1*nmax);
t->SetTextSize(0.05);
c3->cd(npad3);
gPad->SetMargin(0.135, 0.01, 0., 0.);
if(irap==nrapbins-1) gPad->SetMargin(0.135, 0.01, 0.1, 0.);
FormatFrame(frame);
frame->Draw();
li->DrawLine(xx3, -0.08*nmax, xx3, 1.1*nmax);
li->DrawLine(xx4, -0.08*nmax, xx4, 1.1*nmax);
/*
tra3->SetMarkerStyle(20);
tra3->SetMarkerSize(0.6);
tra3->SetMarkerColor(1);
tra3->GetXaxis()->SetLimits(xmin,xmax);
tra3->Draw("psame");
t->DrawLatexNDC(0.21,0.55,Form("Ns = %4.1f, eNs = %4.1f, dN/N = %5.4f",dN3[irap][iphi],dN3err[irap][iphi],dNN3[irap][iphi]));
t->DrawLatexNDC(0.21,0.45,Form("3Ns = %4.1f, 3eNs = %4.1f, 3dN/N = %5.4f",sN3[irap][iphi],sN3e[irap][iphi],dnn3[irap][iphi]));
//t->DrawLatexNDC(0.22,0.45,Form("dN/N = %5.4f",dNN3[irap][iphi]));
*/
TPaveText *label_rap = new TPaveText(0.22,0.85,0.52,1.0,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.07);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
TPaveText *label_phi = new TPaveText(0.7,0.3,0.9,0.4,"blNDC");
label_phi->SetBorderSize(0);
label_phi->SetFillStyle(0);
label_phi->SetTextFont(62);
label_phi->SetTextSize(0.07);
// label_phi->SetTextColor(2);
label_phi->AddText(Form("%s<#phi-#Psi_{1}<%s", dphi1[iphi], dphi1[iphi+1]));
label_phi->Draw();
float chi2;
chi2=frame->chiSquare("model", "dh");
t->SetTextSize(0.05);
//t->DrawLatexNDC(0.72,0.4,Form("#chi^{2} = %3.4f",chi2));
t->DrawLatexNDC(0.22,0.75,Form("#mu = %3.4f",gmean));
t->DrawLatexNDC(0.22,0.65,Form("#sigma = %3.4f",gsig));
cout << frame->chiSquare("model", "dh") << endl;
cout <<"par0 = "<< par0.getVal() << endl;
c3s->cd(npad3);npad3++;
gPad->SetMargin(0.135, 0.01, 0., 0.);
if(irap==nrapbins-1) gPad->SetMargin(0.135, 0.01, 0.1, 0.);
FormatFrame(sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]);
FormatFrame(bHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]);
sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetXaxis()->SetRangeUser(2.97,3.02);
sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerColor(2);
//sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetXaxis()->SetTitle("Inv.mass(#pi^{-}-d-p) [GeV/c^{2}]");
//sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetYaxis()->SetTitle("Counts");
sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerSize(0.8);
sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("e");
bHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(1);
//bHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("histsame");
bHtmass_150210_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(4);
bHtmass_150210_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("histsame");
bHtmass_0180_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(6);
bHtmass_0180_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Draw("histsame");
li->DrawLine(xx3, 0, xx3, 1.1*nmax1);
li->DrawLine(xx4, 0, xx4, 1.1*nmax1);
if(irap==0&&iphi==0)
leg->Draw("same");
label_rap->Draw();
label_phi->Draw();
}
}
c3s->Print(Outpdf);
c3->Print(Outpdf);
float sv1[nrapbins];
float sv1_err[nrapbins];
float sv1c[nrapbins];
float sv1c_err[nrapbins];
float sv2[nrapbins];
float sv2_err[nrapbins];
float sv2c[nrapbins];
float sv2c_err[nrapbins];
float sv3[nrapbins];
float sv3_err[nrapbins];
float sv3c[nrapbins];
float sv3c_err[nrapbins];
TF1 *f1[nrapbins];
TF1 *f1c[nrapbins];
TF1 *f2[nrapbins];
TF1 *f2c[nrapbins];
TF1 *f3[nrapbins];
TF1 *f3c[nrapbins];
TGraphErrors* v1y = new TGraphErrors();
TGraphErrors* v1ym = new TGraphErrors();
TGraphErrors* v1cy = new TGraphErrors();
TGraphErrors* v1cym = new TGraphErrors();
TGraphErrors* v2y = new TGraphErrors();
TGraphErrors* v2ym = new TGraphErrors();
TGraphErrors* v2cy = new TGraphErrors();
TGraphErrors* v2cym = new TGraphErrors();
TGraphErrors* v3y = new TGraphErrors();
TGraphErrors* v3ym = new TGraphErrors();
TGraphErrors* v3cy = new TGraphErrors();
TGraphErrors* v3cym = new TGraphErrors();
TGraphErrors *dN1rap[nrapbins];
TGraphErrors *dN1crap[nrapbins];
TGraphErrors *dN2rap[nrapbins];
TGraphErrors *dN2crap[nrapbins];
TGraphErrors *dN3rap[nrapbins];
TGraphErrors *dN3crap[nrapbins];
for(int i=0;i<nrapbins;i++){
dN1rap[i]=new TGraphErrors();
dN1rap[i]->SetName(Form("dN1rap_%d",i));
dN1crap[i]=new TGraphErrors();
dN1crap[i]->SetName(Form("dN1crap_%d",i));
dN2rap[i]=new TGraphErrors();
dN2rap[i]->SetName(Form("dN2rap_%d",i));
dN2crap[i]=new TGraphErrors();
dN2crap[i]->SetName(Form("dN2crap_%d",i));
dN3rap[i]=new TGraphErrors();
dN3rap[i]->SetName(Form("dN3rap_%d",i));
dN3crap[i]=new TGraphErrors();
dN3crap[i]->SetName(Form("dN3crap_%d",i));
}
TCanvas* c4= new TCanvas("c4","c4",400*nrapbins,800);
c4->Divide(nrapbins,2);
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
dN1rap[irap]->SetPoint(iphi,phibincent1[iphi],dN1[irap][iphi]);
dN1rap[irap]->SetPointError(iphi,0,dN1err[irap][iphi]);
dN1crap[irap]->SetPoint(iphi,phibincent1[iphi],sN1[irap][iphi]);
dN1crap[irap]->SetPointError(iphi,0,sN1e[irap][iphi]);
}
c4->cd(irap+1);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN1rap[irap]->SetMarkerStyle(20);
dN1rap[irap]->SetMarkerColor(1);
dN1rap[irap]->SetMarkerSize(1.);
HistStyle(dN1rap[irap]);
HistAxiasName(dN1rap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f1[irap] = new TF1(Form("f1_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f1[irap] ->SetParameters(200, 0.04);
dN1rap[irap]->Fit(Form("f1_%d", irap),"i");
sv1[irap] = f1[irap]->GetParameter(1);
sv1_err[irap] = f1[irap]->GetParError(1);
v1y->SetPoint(irap, rapbincent[irap], sv1[irap]*Res1);
v1y->SetPointError(irap, 0, sv1_err[irap]*Res1);
v1ym->SetPoint(irap, -rapbincent[irap], -sv1[irap]*Res1);
v1ym->SetPointError(irap, 0, sv1_err[irap]*Res1);
dN1rap[irap]->Draw("ap");
TPaveText *label_rap = new TPaveText(0.4,0.8,0.75,0.9,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.06);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
c4->cd(irap+4);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN1crap[irap]->SetMarkerStyle(20);
dN1crap[irap]->SetMarkerColor(1);
dN1crap[irap]->SetMarkerSize(1.);
HistStyle(dN1crap[irap]);
HistAxiasName(dN1crap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f1c[irap] = new TF1(Form("f1c_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f1c[irap] ->SetParameters(200, 0.04);
dN1crap[irap]->Fit(Form("f1c_%d", irap),"i");
sv1c[irap] = f1c[irap]->GetParameter(1);
sv1c_err[irap] = f1c[irap]->GetParError(1);
v1cy->SetPoint(irap, rapbincent[irap], sv1c[irap]*Res1);
v1cy->SetPointError(irap, 0, sv1c_err[irap]*Res1);
v1cym->SetPoint(irap, -rapbincent[irap], -sv1c[irap]*Res1);
v1cym->SetPointError(irap, 0, sv1c_err[irap]*Res1);
dN1crap[irap]->Draw("ap");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.06);
// label_rap->SetTextColor(2);
label_rap->Draw();
}
c4->Print(Outpdf);
fproton->SetLineStyle(2);
fproton->SetLineColor(1);
fdeuteron->SetLineStyle(2);
fdeuteron->SetLineColor(2);
fHe3->SetLineStyle(2);
fHe3->SetLineColor(3);
fHe4->SetLineStyle(2);
fHe4->SetLineColor(6);
ftriton->SetLineStyle(2);
ftriton->SetLineColor(4);
tl->SetTextSize(0.035);
TCanvas* c5= new TCanvas("c5","c5",900,400);
c5->Divide(2,1,0.02,0.02);
gPad->SetMargin(0.17, 0.0, 0.2, 0.03);
TH2D* htest = new TH2D("htest", "", 9, -0.9, 0.9,100,-1.5,1.2);
v1y->SetMarkerStyle(20);
v1y->SetMarkerSize(1.0);
v1y->SetMarkerColor(2);
v1ym->SetMarkerStyle(24);
v1ym->SetMarkerSize(1.0);
v1ym->SetMarkerColor(2);
v1cy->SetMarkerStyle(20);
v1cy->SetMarkerSize(1.0);
v1cy->SetMarkerColor(2);
v1cym->SetMarkerStyle(24);
v1cym->SetMarkerSize(1.0);
v1cym->SetMarkerColor(2);
htest->GetYaxis()->SetRangeUser(-1.2, 1.2);
htest->GetXaxis()->CenterTitle();
htest->GetYaxis()->CenterTitle();
htest->GetYaxis()->SetTitle("v_{1}");
htest->GetXaxis()->SetTitle("Rapidity (y_{c.m.})");
htest->GetXaxis()->SetTitleSize(0.05);
htest->GetYaxis()->SetTitleSize(0.05);
htest->GetXaxis()->SetTitleOffset(1.1);
htest->GetYaxis()->SetTitleOffset(0.8);
c5->cd(1);
htest->Draw();
v1y->Draw("pesame");
v1ym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
if(opt==0) t->DrawLatex(-0.4,1.25,"Without efficieny correction");
else t->DrawLatex(-0.4,1.25,"With efficieny correction");
//TF1* fitv1y = new TF1("fitv1y", "[0]*x+[1]*x*x*x", -0.9, 0.9);
TF1* fitv1y = new TF1("fitv1y", "[0]*x", -0.9, 0.9);
//TF1* fitv1y = new TF1("fitv1y", "[0]*x", -0.5, 0.);
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
if(pid==3005)
tl->AddEntry(v1y,"{}^{4}_{#Lambda}H_dNdphi","PE");
if(pid==3004)
tl->AddEntry(v1y,"{}^{3}_{#Lambda}H","PE");
if(pid==103004)
tl->AddEntry(v1y,"{}^{3}_{#Lambda}H(3-body)","PE");
tl->AddEntry(ftriton,"Triton","l");
tl->AddEntry(fproton,"Proton","l");
tl->AddEntry(fdeuteron,"Deuteron","l");
tl->AddEntry(fHe3,"Helium3","l");
tl->AddEntry(fHe4,"Helium4","l");
tl->Draw();
v1y->Fit("fitv1y", "r");
c5->cd(2);
htest->Draw();
v1cy->Draw("pesame");
v1cym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
if(opt==0) t->DrawLatex(-0.4,1.25,"Without efficieny correction");
else t->DrawLatex(-0.4,1.25,"With efficieny correction");
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
tl->Draw();
v1cy->Fit("fitv1y", "r");
c5->Print(Outpdf);
TCanvas* c6= new TCanvas("c6","c6",400*nrapbins,800);
c6->Divide(nrapbins,2);
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
dN2rap[irap]->SetPoint(iphi,phibincent1[iphi],dN2[irap][iphi]);
dN2rap[irap]->SetPointError(iphi,0,dN2err[irap][iphi]);
dN2crap[irap]->SetPoint(iphi,phibincent1[iphi],sN2[irap][iphi]);
dN2crap[irap]->SetPointError(iphi,0,sN2e[irap][iphi]);
}
c6->cd(irap+1);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN2rap[irap]->SetMarkerStyle(20);
dN2rap[irap]->SetMarkerColor(1);
dN2rap[irap]->SetMarkerSize(1.);
HistStyle(dN2rap[irap]);
HistAxiasName(dN2rap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f2[irap] = new TF1(Form("f2_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f2[irap] ->SetParameters(200, 0.04);
dN2rap[irap]->Fit(Form("f2_%d", irap),"i");
sv2[irap] = f2[irap]->GetParameter(1);
sv2_err[irap] = f2[irap]->GetParError(1);
v2y->SetPoint(irap, rapbincent[irap], sv2[irap]*Res1);
v2y->SetPointError(irap, 0, sv2_err[irap]*Res1);
v2ym->SetPoint(irap, -rapbincent[irap], -sv2[irap]*Res1);
v2ym->SetPointError(irap, 0, sv2_err[irap]*Res1);
dN2rap[irap]->Draw("ap");
TPaveText *label_rap = new TPaveText(0.4,0.8,0.75,0.9,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.06);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
c6->cd(irap+4);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN2crap[irap]->SetMarkerStyle(20);
dN2crap[irap]->SetMarkerColor(1);
dN2crap[irap]->SetMarkerSize(1.);
HistStyle(dN2crap[irap]);
HistAxiasName(dN2crap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f2c[irap] = new TF1(Form("f2c_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f2c[irap] ->SetParameters(200, 0.04);
dN2crap[irap]->Fit(Form("f2c_%d", irap),"i");
sv2c[irap] = f2c[irap]->GetParameter(1);
sv2c_err[irap] = f2c[irap]->GetParError(1);
v2cy->SetPoint(irap, rapbincent[irap], sv2c[irap]*Res1);
v2cy->SetPointError(irap, 0, sv2c_err[irap]*Res1);
v2cym->SetPoint(irap, -rapbincent[irap], -sv2c[irap]*Res1);
v2cym->SetPointError(irap, 0, sv2c_err[irap]*Res1);
dN2crap[irap]->Draw("ap");
label_rap->Draw();
}
c6->Print(Outpdf);
TCanvas* c7= new TCanvas("c7","c7",900,400);
c7->Divide(2,1);
gPad->SetMargin(0.17, 0.0, 0.2, 0.03);
v2y->SetMarkerStyle(20);
v2y->SetMarkerSize(1.0);
v2y->SetMarkerColor(2);
v2ym->SetMarkerStyle(24);
v2ym->SetMarkerSize(1.0);
v2ym->SetMarkerColor(2);
v2cy->SetMarkerStyle(20);
v2cy->SetMarkerSize(1.0);
v2cy->SetMarkerColor(2);
v2cym->SetMarkerStyle(24);
v2cym->SetMarkerSize(1.0);
v2cym->SetMarkerColor(2);
c7->cd(1);
htest->Draw();
v2y->Draw("pesame");
v2ym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
if(opt==0) t->DrawLatex(-0.4,1.25,"Without efficieny correction");
else t->DrawLatex(-0.4,1.25,"With efficieny correction");
//TF1* fitv1y = new TF1("fitv1y", "[0]*x", -0.9, 0.9);
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
tl->Draw();
v2y->Fit("fitv1y", "r");
c7->cd(2);
htest->Draw();
v2cy->Draw("pesame");
v2cym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
if(opt==0) t->DrawLatex(-0.4,1.25,"Without efficieny correction");
else t->DrawLatex(-0.4,1.25,"With efficieny correction");
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
tl->Draw();
v2cy->Fit("fitv1y", "r");
c7->Print(Outpdf);
TCanvas* c8= new TCanvas("c8","c8",400*nrapbins,800);
c8->Divide(nrapbins,2);
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
dN3rap[irap]->SetPoint(iphi,phibincent1[iphi],dN3[irap][iphi]);
dN3rap[irap]->SetPointError(iphi,0,dN3err[irap][iphi]);
dN3crap[irap]->SetPoint(iphi,phibincent1[iphi],sN3[irap][iphi]);
dN3crap[irap]->SetPointError(iphi,0,sN3e[irap][iphi]);
}
c8->cd(irap+1);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN3rap[irap]->SetMarkerStyle(20);
dN3rap[irap]->SetMarkerColor(1);
dN3rap[irap]->SetMarkerSize(1.);
HistStyle(dN3rap[irap]);
HistAxiasName(dN1rap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f3[irap] = new TF1(Form("f3_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f3[irap] ->SetParameters(200, 0.04);
dN3rap[irap]->Fit(Form("f3_%d", irap),"i");
sv3[irap] = f3[irap]->GetParameter(1);
sv3_err[irap] = f3[irap]->GetParError(1);
v3y->SetPoint(irap, rapbincent[irap], sv3[irap]*Res1);
v3y->SetPointError(irap, 0, sv3_err[irap]*Res1);
v3ym->SetPoint(irap, -rapbincent[irap], -sv3[irap]*Res1);
v3ym->SetPointError(irap, 0, sv3_err[irap]*Res1);
dN3rap[irap]->Draw("ap");
TPaveText *label_rap = new TPaveText(0.4,0.8,0.75,0.9,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.06);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
c8->cd(irap+4);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN3crap[irap]->SetMarkerStyle(20);
dN3crap[irap]->SetMarkerColor(1);
dN3crap[irap]->SetMarkerSize(1.);
HistStyle(dN3crap[irap]);
HistAxiasName(dN3crap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f3c[irap] = new TF1(Form("f3c_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f3c[irap] ->SetParameters(200, 0.04);
dN3crap[irap]->Fit(Form("f3c_%d", irap),"i");
sv3c[irap] = f3c[irap]->GetParameter(1);
sv3c_err[irap] = f3c[irap]->GetParError(1);
v3cy->SetPoint(irap, rapbincent[irap], sv3c[irap]*Res1);
v3cy->SetPointError(irap, 0, sv3c_err[irap]*Res1);
v3cym->SetPoint(irap, -rapbincent[irap], -sv3c[irap]*Res1);
v3cym->SetPointError(irap, 0, sv3c_err[irap]*Res1);
dN3crap[irap]->Draw("ap");
label_rap->Draw();
}
c8->Print(Outpdf);
TCanvas* c9= new TCanvas("c9","c9",1000,500);
c9->Divide(2,1);
gPad->SetMargin(0.17, 0.0, 0.2, 0.03);
v3y->SetMarkerStyle(20);
v3y->SetMarkerSize(1.0);
v3y->SetMarkerColor(2);
v3ym->SetMarkerStyle(24);
v3ym->SetMarkerSize(1.0);
v3ym->SetMarkerColor(2);
v3cy->SetMarkerStyle(20);
v3cy->SetMarkerSize(1.0);
v3cy->SetMarkerColor(2);
v3cym->SetMarkerStyle(24);
v3cym->SetMarkerSize(1.0);
v3cym->SetMarkerColor(2);
c9->cd(1);
htest->Draw();
v3y->Draw("pesame");
v3ym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
t->DrawLatex(-0.4,1.25,"With efficieny correction");
//TF1* fitv3y = new TF1("fitv3y", "[0]*x+[1]*x*x*x", -0.9, 0.9);
TF1* fitv3y = new TF1("fitv3y", "[0]*x", -0.9, 0.9);
//TF1* fitv3y = new TF1("fitv3y", "[0]*x", -0.5, 0.);
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
tl->Draw();
v3y->Fit("fitv3y", "r");
c9->cd(2);
htest->Draw();
v3cy->Draw("pesame");
v3cym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
t->DrawLatex(-0.4,1.25,"With efficieny correction");
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
tl->Draw();
v3cy->Fit("fitv3y", "r");
c9->Print(Outpdf);
float dN4[nrapbins][ncutphibins];
float dN4err[nrapbins][ncutphibins];
float dNN4[nrapbins][ncutphibins];
float sN4[nrapbins][ncutphibins];
float sN4e[nrapbins][ncutphibins];
float dnn4[nrapbins][ncutphibins];
TCanvas* c10 = new TCanvas("c10","c10",1000,200*nrapbins);
c10->Divide(ncutphibins,nrapbins,0.12,0.0);
npad3=1;
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
float xmin, xmax;
float m1, m2, m3;
float s1, s2, s3;
if(pid==103004){
xmin=2.97; xmax=3.02; // for H3L
//m1=2.9925, m2=2.989, m3=2.996;
//m1=2.9927; m2=2.9926; m3=2.9928; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//m1=2.993; m2=2.9928; m3=2.9940; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
m1=2.9927; m2=2.9927; m3=2.9927; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//m1=2.9918; m2=2.9915; m3=2.9919; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 2ybin
//m1=2.9924; m2=2.9924; m3=2.9924;
//m1=2.9922; m2=2.9922; m3=2.9922;
//m1=2.9928; m2=2.9920; m3=2.9930;
//s1=0.0021; s2=0.0016; s3=0.0024;
//s1=0.0019; s2=0.0015; s3=0.0024;
//s1=0.002; s2=0.001; s3=0.0030;
s1=0.0018; s2=0.0017; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0018; s2=0.0015; s3=0.002; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0019; s2=0.0014; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
}
RooRealVar x("x", "x", xmin, xmax); // for H3L
RooRealVar mean("mean", "mean", m1, m2, m3);
RooRealVar sigma("sigma", "sigma", s1, s2, s3);//H4L
RooGaussian sig("sig", "sig", x, mean, sigma);
RooRealVar par0("par0", "par0", 0.5, -1, 1);
RooRealVar par1("par1", "par1", -0.2, -1, 1);
RooRealVar par2("par2", "par2", 0, -1, 1);
RooRealVar par3("par3", "par3", 0, -1, 1);
RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1, par2));
//RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1));
RooRealVar nsig("nsig","number of signal events",100, 10, 8000);
RooRealVar nbkg("nbkg","number of bkg events", 100, 10, 10e6);
RooExtendPdf esig("esig", "esig", sig, nsig);
RooExtendPdf ebkg("ebkg", "ebkg", bkg, nbkg);
RooAddPdf model("model", "model", RooArgList(esig, ebkg));
sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(2);
sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
float nmax=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetMaximum();
RooPlot* frame = x.frame(" ");
RooDataHist dh1("dh1", "dh1", x, Import(*(sHtm_phi1_rap_piDca2_DecayL3[iphi][irap])));
model.fitTo(dh1, Extended());
dh1.plotOn(frame, Name("dh1"), MarkerSize(0.5), XErrorSize(0), DrawOption("PZ"));
model.plotOn(frame, Components(esig), LineStyle(kSolid), LineColor(2), LineWidth(2));
model.plotOn(frame, Components(ebkg), LineStyle(kDashed), LineColor(4), LineWidth(2));
model.plotOn(frame, Name("model"), LineStyle(kDashed), LineColor(2), LineWidth(2));
int binmi=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmin);
int binma=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmax);
int j=0;
TGraphErrors* tra4=new TGraphErrors();
for(int ibin=binmi; ibin<binma+1; ibin++){
int xbin=ibin-11;
float Nt=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinContent(ibin);
float bc=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinCenter(ibin);
float bw=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinWidth(ibin);
float xmi=bc-0.5*bw;
float xma=bc+0.5*bw;
x.setRange("bg_range",xmi,xma);
RooAbsReal* ibg = bkg.createIntegral(x,NormSet(x),Range("bg_range"));
float Nb = ibg->getVal()*nbkg.getVal();
float Ns=Nt-Nb;
float Ne=1.0*sqrt(Nt+Nb);
float errbin=sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinError(ibin);
if(j==16)
//t->DrawLatexNDC(0.22,0.55,Form("ibin = %d, errbin = %f",ibin,errbin));
cout<<"ibin = "<<ibin<<", Nt = "<<Nt<<", errbin = "<<errbin<<", Ne = "<<Ne<<endl;
tra4->SetPoint(j,bc,Ns);
//tra4->SetPointError(j,0,Ne);
tra4->SetPointError(j,0,errbin);
j++;
}
float gmean = mean.getVal();
float gsig = sigma.getVal();
//float tt=3.;
//float tt=2.;
float tt=1.5;
float a1=gmean-tt*gsig;
float a2=gmean+tt*gsig;
int binx1 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(gmean);
int binx2 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a1);
int binx3 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a2);
float Ntt = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx1-1,binx1+1);
float xx1 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1-1);
float xx2 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1+2);
float xx3 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx2);
float xx4 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx3+1);
float Ntt1 = sHtm_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx2,binx3);
x.setRange("bkge_range",xx1,xx2);
RooAbsReal* ibb = bkg.createIntegral(x,NormSet(x),Range("bkge_range"));
float Nbb = ibb->getVal()*nbkg.getVal();
float Nss = Ntt-Nbb;
x.setRange("bkg1_range",xx3,xx4);
RooAbsReal* ib1 = bkg.createIntegral(x,NormSet(x),Range("bkg1_range"));
float Nbb1 = ib1->getVal()*nbkg.getVal();
float Nss1 = Ntt1-Nbb1;
float eN1 = 1.0*sqrt(Ntt1+Nbb1);
//dN3[irap][iphi] = Nss;
dN4[irap][iphi] = nsig.getVal();
//dN4err[irap][iphi] = dN4[irap][iphi]*dnn1[irap][iphi];
dN4err[irap][iphi] = nsig.getError();
//dN4err[irap][iphi] = sqrt(Ntt+Nbb);
dNN4[irap][iphi] = dN4err[irap][iphi]/dN4[irap][iphi];
sN4[irap][iphi] = Nss1;
//sN4[irap][iphi] = Nss1;
//dN1[irap][iphi] = nsig.getVal();
//dN1err[irap][iphi] = nsig.getError();
//sN4e[irap][iphi] = sqrt(Ntt1+Nbb1);
sN4e[irap][iphi] = eN1;
dnn4[irap][iphi] = sN4e[irap][iphi]/sN4[irap][iphi];
frame->GetYaxis()->SetRangeUser(-0.08*nmax,1.1*nmax);
t->SetTextSize(0.05);
c10->cd(npad3); npad3++;
gPad->SetMargin(0.135, 0.01, 0., 0.);
if(irap==nrapbins-1) gPad->SetMargin(0.135, 0.01, 0.1, 0.);
FormatFrame(frame);
frame->Draw();
li->SetLineColor(4);
li->DrawLine(xx3, -0.08*nmax, xx3, 1.1*nmax);
li->DrawLine(xx4, -0.08*nmax, xx4, 1.1*nmax);
tra4->SetMarkerStyle(20);
tra4->SetMarkerSize(0.6);
tra4->SetMarkerColor(1);
tra4->GetXaxis()->SetLimits(xmin,xmax);
tra4->Draw("psame");
//t->DrawLatexNDC(0.21,0.55,Form("Ns = %4.1f, eNs = %4.1f, dN/N = %5.4f",dN4[irap][iphi],dN4err[irap][iphi],dNN4[irap][iphi]));
//t->DrawLatexNDC(0.21,0.45,Form("3Ns = %4.1f, 3eNs = %4.1f, 3dN/N = %5.4f",sN4[irap][iphi],sN4e[irap][iphi],dnn4[irap][iphi]));
//t->DrawLatexNDC(0.22,0.45,Form("dN/N = %5.4f",dNN3[irap][iphi]));
TPaveText *label_rap = new TPaveText(0.22,0.85,0.52,1.0,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.07);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
TPaveText *label_phi = new TPaveText(0.7,0.3,0.9,0.4,"blNDC");
label_phi->SetBorderSize(0);
label_phi->SetFillStyle(0);
label_phi->SetTextFont(62);
label_phi->SetTextSize(0.07);
// label_phi->SetTextColor(2);
label_phi->AddText(Form("%s<#phi-#Psi_{1}<%s", dphi1[iphi], dphi1[iphi+1]));
label_phi->Draw();
float chi2;
chi2=frame->chiSquare("model", "dh");
t->SetTextSize(0.05);
//t->DrawLatexNDC(0.72,0.4,Form("#chi^{2} = %3.4f",chi2));
//t->DrawLatexNDC(0.22,0.75,Form("#mu = %3.4f",gmean));
//t->DrawLatexNDC(0.22,0.65,Form("#sigma = %3.4f",gsig));
cout << frame->chiSquare("model", "dh") << endl;
cout <<"par0 = "<< par0.getVal() << endl;
}
}
c10->Print(Outpdf);
TGraphErrors *dN4rap[nrapbins];
TGraphErrors *dN4crap[nrapbins];
TGraphErrors* v4y = new TGraphErrors();
TGraphErrors* v4ym = new TGraphErrors();
TGraphErrors* v4cy = new TGraphErrors();
TGraphErrors* v4cym = new TGraphErrors();
float sv4[nrapbins];
float sv4_err[nrapbins];
float sv4c[nrapbins];
float sv4c_err[nrapbins];
for(int i=0;i<nrapbins;i++){
dN4rap[i]=new TGraphErrors();
dN4rap[i]->SetName(Form("dN1rap_%d",i));
dN4crap[i]=new TGraphErrors();
dN4crap[i]->SetName(Form("dN1crap_%d",i));
}
TCanvas* c11= new TCanvas("c11","c11",400*nrapbins,800);
c11->Divide(nrapbins,2);
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
dN4rap[irap]->SetPoint(iphi,phibincent1[iphi],dN4[irap][iphi]);
dN4rap[irap]->SetPointError(iphi,0,dN4err[irap][iphi]);
dN4crap[irap]->SetPoint(iphi,phibincent1[iphi],sN4[irap][iphi]);
dN4crap[irap]->SetPointError(iphi,0,sN4e[irap][iphi]);
}
c11->cd(irap+1);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN4rap[irap]->SetMarkerStyle(20);
dN4rap[irap]->SetMarkerColor(1);
dN4rap[irap]->SetMarkerSize(1.);
HistStyle(dN4rap[irap]);
HistAxiasName(dN4rap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f1[irap] = new TF1(Form("f1_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f1[irap] ->SetParameters(200, 0.04);
dN4rap[irap]->Fit(Form("f1_%d", irap),"i");
sv4[irap] = f1[irap]->GetParameter(1);
sv4_err[irap] = f1[irap]->GetParError(1);
v4y->SetPoint(irap, rapbincent[irap], sv4[irap]*Res1);
v4y->SetPointError(irap, 0, sv4_err[irap]*Res1);
v4ym->SetPoint(irap, -rapbincent[irap], -sv4[irap]*Res1);
v4ym->SetPointError(irap, 0, sv4_err[irap]*Res1);
dN4rap[irap]->Draw("ap");
TPaveText *label_rap = new TPaveText(0.4,0.8,0.75,0.9,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.06);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
c11->cd(irap+4);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN4crap[irap]->SetMarkerStyle(20);
dN4crap[irap]->SetMarkerColor(1);
dN4crap[irap]->SetMarkerSize(1.);
HistStyle(dN4crap[irap]);
HistAxiasName(dN4crap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f1c[irap] = new TF1(Form("f1c_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f1c[irap] ->SetParameters(200, 0.04);
dN4crap[irap]->Fit(Form("f1c_%d", irap),"i");
sv4c[irap] = f1c[irap]->GetParameter(1);
sv4c_err[irap] = f1c[irap]->GetParError(1);
v4cy->SetPoint(irap, rapbincent[irap], sv4c[irap]*Res1);
v4cy->SetPointError(irap, 0, sv4c_err[irap]*Res1);
v4cym->SetPoint(irap, -rapbincent[irap], -sv4c[irap]*Res1);
v4cym->SetPointError(irap, 0, sv4c_err[irap]*Res1);
dN4crap[irap]->Draw("ap");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.06);
// label_rap->SetTextColor(2);
label_rap->Draw();
}
c11->Print(Outpdf);
TCanvas* c12= new TCanvas("c12","c12",900,400);
c12->Divide(2,1,0.02,0.02);
gPad->SetMargin(0.17, 0.0, 0.2, 0.03);
v4y->SetMarkerStyle(20);
v4y->SetMarkerSize(1.0);
v4y->SetMarkerColor(2);
v4ym->SetMarkerStyle(24);
v4ym->SetMarkerSize(1.0);
v4ym->SetMarkerColor(2);
v4cy->SetMarkerStyle(20);
v4cy->SetMarkerSize(1.0);
v4cy->SetMarkerColor(2);
v4cym->SetMarkerStyle(24);
v4cym->SetMarkerSize(1.0);
v4cym->SetMarkerColor(2);
htest->GetYaxis()->SetRangeUser(-1.2, 1.2);
htest->GetXaxis()->CenterTitle();
htest->GetYaxis()->CenterTitle();
htest->GetYaxis()->SetTitle("v_{1}");
htest->GetXaxis()->SetTitle("Rapidity (y_{c.m.})");
htest->GetXaxis()->SetTitleSize(0.05);
htest->GetYaxis()->SetTitleSize(0.05);
htest->GetXaxis()->SetTitleOffset(1.1);
htest->GetYaxis()->SetTitleOffset(0.8);
c12->cd(1);
htest->Draw();
v4y->Draw("pesame");
v4ym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
if(opt==0) t->DrawLatex(-0.4,1.25,"Without efficieny correction");
else t->DrawLatex(-0.4,1.25,"With efficieny correction");
//TF1* fitv1y = new TF1("fitv1y", "[0]*x+[1]*x*x*x", -0.9, 0.9);
fitv1y = new TF1("fitv1y", "[0]*x", -0.9, 0.9);
//TF1* fitv1y = new TF1("fitv1y", "[0]*x", -0.5, 0.);
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
tl->Draw();
v4y->Fit("fitv1y", "r");
c12->cd(2);
htest->Draw();
v4cy->Draw("pesame");
v4cym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
if(opt==0) t->DrawLatex(-0.4,1.25,"Without efficieny correction");
else t->DrawLatex(-0.4,1.25,"With efficieny correction");
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
tl->Draw();
v4cy->Fit("fitv1y", "r");
c12->Print(Outpdf);
float dN5[nrapbins][ncutphibins];
float dN5err[nrapbins][ncutphibins];
float dNN5[nrapbins][ncutphibins];
float sN5[nrapbins][ncutphibins];
float sN5e[nrapbins][ncutphibins];
float dnn5[nrapbins][ncutphibins];
TCanvas* c13 = new TCanvas("c13","c13",1000,200*nrapbins);
c13->Divide(ncutphibins,nrapbins,0.12,0.0);
npad3=1;
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
float xmin, xmax;
float m1, m2, m3;
float s1, s2, s3;
if(pid==103004){
xmin=2.97; xmax=3.02; // for H3L
//m1=2.9925, m2=2.989, m3=2.996;
//m1=2.9927; m2=2.9926; m3=2.9928; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//m1=2.993; m2=2.9928; m3=2.9940; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//m1=2.9927; m2=2.9927; m3=2.9927; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
m1=2.9927; m2=2.9925; m3=2.993; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//m1=2.9918; m2=2.9915; m3=2.9919; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 2ybin
//m1=2.9924; m2=2.9924; m3=2.9924;
//m1=2.9922; m2=2.9922; m3=2.9922;
//m1=2.9928; m2=2.9920; m3=2.9930;
//s1=0.0021; s2=0.0016; s3=0.0024;
//s1=0.0019; s2=0.0015; s3=0.0024;
//s1=0.002; s2=0.001; s3=0.0030;
//s1=0.0017; s2=0.0016; s3=0.002; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
s1=0.0018; s2=0.0018; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0018; s2=0.0015; s3=0.002; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0019; s2=0.0014; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
}
RooRealVar x("x", "x", xmin, xmax); // for H3L
RooRealVar mean("mean", "mean", m1, m2, m3);
RooRealVar sigma("sigma", "sigma", s1, s2, s3);//H4L
RooGaussian sig("sig", "sig", x, mean, sigma);
RooRealVar par0("par0", "par0", 0.5, -1, 1);
RooRealVar par1("par1", "par1", -0.2, -1, 1);
RooRealVar par2("par2", "par2", 0, -1, 1);
RooRealVar par3("par3", "par3", 0, -1, 1);
RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1, par2));
//RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1));
RooRealVar nsig("nsig","number of signal events",100, 10, 8000);
RooRealVar nbkg("nbkg","number of bkg events", 100, 10, 10e6);
RooExtendPdf esig("esig", "esig", sig, nsig);
RooExtendPdf ebkg("ebkg", "ebkg", bkg, nbkg);
RooAddPdf model("model", "model", RooArgList(esig, ebkg));
sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(2);
sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
float nmax=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetMaximum();
RooPlot* frame = x.frame(" ");
RooDataHist dh1("dh1", "dh1", x, Import(*(sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap])));
model.fitTo(dh1, Extended());
dh1.plotOn(frame, Name("dh1"), MarkerSize(0.5), XErrorSize(0), DrawOption("PZ"));
model.plotOn(frame, Components(esig), LineStyle(kSolid), LineColor(2), LineWidth(2));
model.plotOn(frame, Components(ebkg), LineStyle(kDashed), LineColor(4), LineWidth(2));
model.plotOn(frame, Name("model"), LineStyle(kDashed), LineColor(2), LineWidth(2));
int binmi=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmin);
int binma=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmax);
int j=0;
TGraphErrors* tra5=new TGraphErrors();
for(int ibin=binmi; ibin<binma+1; ibin++){
int xbin=ibin-11;
float Nt=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinContent(ibin);
float bc=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinCenter(ibin);
float bw=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinWidth(ibin);
float xmi=bc-0.5*bw;
float xma=bc+0.5*bw;
x.setRange("bg_range",xmi,xma);
RooAbsReal* ibg = bkg.createIntegral(x,NormSet(x),Range("bg_range"));
float Nb = ibg->getVal()*nbkg.getVal();
float Ns=Nt-Nb;
float Ne=1.0*sqrt(Nt+Nb);
float errbin=sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinError(ibin);
if(j==16)
//t->DrawLatexNDC(0.22,0.55,Form("ibin = %d, errbin = %f",ibin,errbin));
cout<<"ibin = "<<ibin<<", Nt = "<<Nt<<", errbin = "<<errbin<<", Ne = "<<Ne<<endl;
tra5->SetPoint(j,bc,Ns);
//tra4->SetPointError(j,0,Ne);
tra5->SetPointError(j,0,errbin);
j++;
}
float gmean = mean.getVal();
float gsig = sigma.getVal();
//float tt=3.;
//float tt=2.;
float tt=1.5;
float a1=gmean-tt*gsig;
float a2=gmean+tt*gsig;
int binx1 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(gmean);
int binx2 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a1);
int binx3 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a2);
float Ntt = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx1-1,binx1+1);
float xx1 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1-1);
float xx2 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1+2);
float xx3 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx2);
float xx4 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx3+1);
float Ntt1 = sHtmass_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx2,binx3);
x.setRange("bkge_range",xx1,xx2);
RooAbsReal* ibb = bkg.createIntegral(x,NormSet(x),Range("bkge_range"));
float Nbb = ibb->getVal()*nbkg.getVal();
float Nss = Ntt-Nbb;
x.setRange("bkg1_range",xx3,xx4);
RooAbsReal* ib1 = bkg.createIntegral(x,NormSet(x),Range("bkg1_range"));
float Nbb1 = ib1->getVal()*nbkg.getVal();
float Nss1 = Ntt1-Nbb1;
float eN1 = 1.0*sqrt(Ntt1+Nbb1);
//dN5[irap][iphi] = Nss;
dN5[irap][iphi] = nsig.getVal();
//dN5err[irap][iphi] = dN5[irap][iphi]*dnn1[irap][iphi];
dN5err[irap][iphi] = nsig.getError();
//dN5err[irap][iphi] = sqrt(Ntt+Nbb);
dNN5[irap][iphi] = dN5err[irap][iphi]/dN5[irap][iphi];
sN5[irap][iphi] = Nss1;
//sN5[irap][iphi] = Nss1;
//dN1[irap][iphi] = nsig.getVal();
//dN1err[irap][iphi] = nsig.getError();
//sN5e[irap][iphi] = sqrt(Ntt1+Nbb1);
sN5e[irap][iphi] = sN5[irap][iphi]*dnn4[irap][iphi];
dnn5[irap][iphi] = sN5e[irap][iphi]/sN5[irap][iphi];
frame->GetYaxis()->SetRangeUser(-0.08*nmax,1.1*nmax);
t->SetTextSize(0.05);
c13->cd(npad3); npad3++;
gPad->SetMargin(0.135, 0.01, 0., 0.);
if(irap==nrapbins-1) gPad->SetMargin(0.135, 0.01, 0.1, 0.);
FormatFrame(frame);
frame->Draw();
li->SetLineColor(4);
li->DrawLine(xx3, -0.08*nmax, xx3, 1.1*nmax);
li->DrawLine(xx4, -0.08*nmax, xx4, 1.1*nmax);
tra5->SetMarkerStyle(20);
tra5->SetMarkerSize(0.6);
tra5->SetMarkerColor(1);
tra5->GetXaxis()->SetLimits(xmin,xmax);
tra5->Draw("psame");
//t->DrawLatexNDC(0.21,0.55,Form("Ns = %4.1f, eNs = %4.1f, dN/N = %5.4f",dN5[irap][iphi],dN5err[irap][iphi],dNN5[irap][iphi]));
//t->DrawLatexNDC(0.21,0.45,Form("3Ns = %4.1f, 3eNs = %4.1f, 3dN/N = %5.4f",sN5[irap][iphi],sN5e[irap][iphi],dnn5[irap][iphi]));
//t->DrawLatexNDC(0.22,0.45,Form("dN/N = %5.4f",dNN3[irap][iphi]));
TPaveText *label_rap = new TPaveText(0.22,0.85,0.52,1.0,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.07);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
TPaveText *label_phi = new TPaveText(0.7,0.3,0.9,0.4,"blNDC");
label_phi->SetBorderSize(0);
label_phi->SetFillStyle(0);
label_phi->SetTextFont(62);
label_phi->SetTextSize(0.07);
// label_phi->SetTextColor(2);
label_phi->AddText(Form("%s<#phi-#Psi_{1}<%s", dphi1[iphi], dphi1[iphi+1]));
label_phi->Draw();
float chi2;
chi2=frame->chiSquare("model", "dh");
t->SetTextSize(0.05);
//t->DrawLatexNDC(0.72,0.4,Form("#chi^{2} = %3.4f",chi2));
//t->DrawLatexNDC(0.22,0.75,Form("#mu = %3.4f",gmean));
//t->DrawLatexNDC(0.22,0.65,Form("#sigma = %3.4f",gsig));
cout << frame->chiSquare("model", "dh") << endl;
cout <<"par0 = "<< par0.getVal() << endl;
}
}
c13->Print(Outpdf);
TGraphErrors *dN5rap[nrapbins];
TGraphErrors *dN5crap[nrapbins];
TGraphErrors* v5y = new TGraphErrors();
TGraphErrors* v5ym = new TGraphErrors();
TGraphErrors* v5cy = new TGraphErrors();
TGraphErrors* v5cym = new TGraphErrors();
float sv5[nrapbins];
float sv5_err[nrapbins];
float sv5c[nrapbins];
float sv5c_err[nrapbins];
for(int i=0;i<nrapbins;i++){
dN5rap[i]=new TGraphErrors();
dN5rap[i]->SetName(Form("dN1rap_%d",i));
dN5crap[i]=new TGraphErrors();
dN5crap[i]->SetName(Form("dN1crap_%d",i));
}
TCanvas* c14= new TCanvas("c14","c14",400*nrapbins,800);
c14->Divide(nrapbins,2);
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
dN5rap[irap]->SetPoint(iphi,phibincent1[iphi],dN5[irap][iphi]);
dN5rap[irap]->SetPointError(iphi,0,dN5err[irap][iphi]);
dN5crap[irap]->SetPoint(iphi,phibincent1[iphi],sN5[irap][iphi]);
dN5crap[irap]->SetPointError(iphi,0,sN5e[irap][iphi]);
}
c14->cd(irap+1);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN5rap[irap]->SetMarkerStyle(20);
dN5rap[irap]->SetMarkerColor(1);
dN5rap[irap]->SetMarkerSize(1.);
HistStyle(dN5rap[irap]);
HistAxiasName(dN5rap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f1[irap] = new TF1(Form("f1_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f1[irap] ->SetParameters(200, 0.04);
dN5rap[irap]->Fit(Form("f1_%d", irap),"i");
sv5[irap] = f1[irap]->GetParameter(1);
sv5_err[irap] = f1[irap]->GetParError(1);
v5y->SetPoint(irap, rapbincent[irap], sv5[irap]*Res1);
v5y->SetPointError(irap, 0, sv5_err[irap]*Res1);
v5ym->SetPoint(irap, -rapbincent[irap], -sv5[irap]*Res1);
v5ym->SetPointError(irap, 0, sv5_err[irap]*Res1);
dN5rap[irap]->Draw("ap");
TPaveText *label_rap = new TPaveText(0.4,0.8,0.75,0.9,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.06);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
c14->cd(irap+4);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN5crap[irap]->SetMarkerStyle(20);
dN5crap[irap]->SetMarkerColor(1);
dN5crap[irap]->SetMarkerSize(1.);
HistStyle(dN5crap[irap]);
HistAxiasName(dN5crap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f1c[irap] = new TF1(Form("f1c_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f1c[irap] ->SetParameters(200, 0.04);
dN5crap[irap]->Fit(Form("f1c_%d", irap),"i");
sv5c[irap] = f1c[irap]->GetParameter(1);
sv5c_err[irap] = f1c[irap]->GetParError(1);
v5cy->SetPoint(irap, rapbincent[irap], sv5c[irap]*Res1);
v5cy->SetPointError(irap, 0, sv5c_err[irap]*Res1);
v5cym->SetPoint(irap, -rapbincent[irap], -sv5c[irap]*Res1);
v5cym->SetPointError(irap, 0, sv5c_err[irap]*Res1);
dN5crap[irap]->Draw("ap");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.06);
// label_rap->SetTextColor(2);
label_rap->Draw();
}
c14->Print(Outpdf);
TCanvas* c15= new TCanvas("c15","c15",900,400);
c15->Divide(2,1,0.02,0.02);
gPad->SetMargin(0.17, 0.0, 0.2, 0.03);
v5y->SetMarkerStyle(20);
v5y->SetMarkerSize(1.0);
v5y->SetMarkerColor(2);
v5ym->SetMarkerStyle(24);
v5ym->SetMarkerSize(1.0);
v5ym->SetMarkerColor(2);
v5cy->SetMarkerStyle(20);
v5cy->SetMarkerSize(1.0);
v5cy->SetMarkerColor(2);
v5cym->SetMarkerStyle(24);
v5cym->SetMarkerSize(1.0);
v5cym->SetMarkerColor(2);
htest->GetYaxis()->SetRangeUser(-1.2, 1.2);
htest->GetXaxis()->CenterTitle();
htest->GetYaxis()->CenterTitle();
htest->GetYaxis()->SetTitle("v_{1}");
htest->GetXaxis()->SetTitle("Rapidity (y_{c.m.})");
htest->GetXaxis()->SetTitleSize(0.05);
htest->GetYaxis()->SetTitleSize(0.05);
htest->GetXaxis()->SetTitleOffset(1.1);
htest->GetYaxis()->SetTitleOffset(0.8);
c15->cd(1);
htest->Draw();
v5y->Draw("pesame");
v5ym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
if(opt==0) t->DrawLatex(-0.4,1.25,"Without efficieny correction");
else t->DrawLatex(-0.4,1.25,"With efficieny correction");
//TF1* fitv1y = new TF1("fitv1y", "[0]*x+[1]*x*x*x", -0.9, 0.9);
fitv1y = new TF1("fitv1y", "[0]*x", -0.9, 0.9);
//TF1* fitv1y = new TF1("fitv1y", "[0]*x", -0.5, 0.);
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
tl->Draw();
v5y->Fit("fitv1y", "r");
c15->cd(2);
htest->Draw();
v5cy->Draw("pesame");
v5cym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
if(opt==0) t->DrawLatex(-0.4,1.25,"Without efficieny correction");
else t->DrawLatex(-0.4,1.25,"With efficieny correction");
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
tl->Draw();
v5cy->Fit("fitv1y", "r");
c15->Print(Outpdf);
float dN6[nrapbins][ncutphibins];
float dN6err[nrapbins][ncutphibins];
float dNN6[nrapbins][ncutphibins];
float sN6[nrapbins][ncutphibins];
float sN6e[nrapbins][ncutphibins];
float dnn6[nrapbins][ncutphibins];
TCanvas* c16 = new TCanvas("c16","c16",1000,200*nrapbins);
c16->Divide(ncutphibins,nrapbins,0.12,0.0);
npad3=1;
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
float xmin, xmax;
float m1, m2, m3;
float s1, s2, s3;
if(pid==103004){
xmin=2.97; xmax=3.02; // for H3L
//m1=2.993; m2=2.9928; m3=2.9940; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//m1=2.9927; m2=2.9927; m3=2.9927; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
m1=2.9927; m2=2.9927; m3=2.9927; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0021; s2=0.0016; s3=0.0024;
//s1=0.0019; s2=0.0015; s3=0.0024;
//s1=0.002; s2=0.001; s3=0.0030;
//s1=0.0017; s2=0.0016; s3=0.002; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
s1=0.0018; s2=0.0018; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0018; s2=0.0015; s3=0.002; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
//s1=0.0019; s2=0.0014; s3=0.0018; // H3L: pt(1.0, 2.5) y(-0.8,-0.2) 3ybin
}
RooRealVar x("x", "x", xmin, xmax); // for H3L
RooRealVar mean("mean", "mean", m1, m2, m3);
RooRealVar sigma("sigma", "sigma", s1, s2, s3);//H4L
RooGaussian sig("sig", "sig", x, mean, sigma);
RooRealVar par0("par0", "par0", 0.5, -1, 1);
RooRealVar par1("par1", "par1", -0.2, -1, 1);
RooRealVar par2("par2", "par2", 0, -1, 1);
RooRealVar par3("par3", "par3", 0, -1, 1);
RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1, par2));
//RooChebychev bkg("bkg", "bkg", x, RooArgList(par0, par1));
RooRealVar nsig("nsig","number of signal events",200, 100, 8000000);
RooRealVar nbkg("nbkg","number of bkg events", 10e4, 10e4, 10e8);
RooExtendPdf esig("esig", "esig", sig, nsig);
RooExtendPdf ebkg("ebkg", "ebkg", bkg, nbkg);
RooAddPdf model("model", "model", RooArgList(esig, ebkg));
sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->SetLineColor(2);
sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->SetMarkerStyle(20);
float nmax=sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetMaximum();
RooPlot* frame = x.frame(" ");
RooDataHist dh1("dh1", "dh1", x, Import(*(sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap])));
model.fitTo(dh1, Extended());
dh1.plotOn(frame, Name("dh1"), MarkerSize(0.5), XErrorSize(0), DrawOption("PZ"));
model.plotOn(frame, Components(esig), LineStyle(kSolid), LineColor(2), LineWidth(2));
model.plotOn(frame, Components(ebkg), LineStyle(kDashed), LineColor(4), LineWidth(2));
model.plotOn(frame, Name("model"), LineStyle(kDashed), LineColor(2), LineWidth(2));
int binmi=sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmin);
int binma=sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(xmax);
int j=0;
TGraphErrors* tra6 = new TGraphErrors();
for(int ibin=binmi; ibin<binma+1; ibin++){
float Nt=sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinContent(ibin);
float bc=sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinCenter(ibin);
float bw=sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinWidth(ibin);
float xmi=bc-0.5*bw;
float xma=bc+0.5*bw;
x.setRange("bg_range",xmi,xma);
RooAbsReal* ibg = bkg.createIntegral(x,NormSet(x),Range("bg_range"));
float Nb = ibg->getVal()*nbkg.getVal();
float Ns=Nt-Nb;
float Ne=1.0*sqrt(Nt+Nb);
float errbin=sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinError(ibin);
if(j==16)
tra6->SetPoint(j,bc,Ns);
tra6->SetPointError(j,0,errbin);
j++;
}
float gmean = mean.getVal();
float gsig = sigma.getVal();
//float tt=3.;
//float tt=2.;
float tt=1.5;
float a1=gmean-tt*gsig;
float a2=gmean+tt*gsig;
int binx1 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(gmean);
int binx2 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a1);
int binx3 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->FindBin(a2);
float Ntt = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx1-1,binx1+1);
float xx1 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1-1);
float xx2 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx1+2);
float xx3 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx2);
float xx4 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->GetBinLowEdge(binx3+1);
float Ntt1 = sHtmass_eff_phi1_rap_piDca2_DecayL3[iphi][irap]->Integral(binx2,binx3);
x.setRange("bkge_range",xx1,xx2);
RooAbsReal* ibb = bkg.createIntegral(x,NormSet(x),Range("bkge_range"));
float Nbb = ibb->getVal()*nbkg.getVal();
float Nss = Ntt-Nbb;
x.setRange("bkg1_range",xx3,xx4);
RooAbsReal* ib1 = bkg.createIntegral(x,NormSet(x),Range("bkg1_range"));
float Nbb1 = ib1->getVal()*nbkg.getVal();
float Nss1 = Ntt1-Nbb1;
float eN1 = 1.0*sqrt(Ntt1+Nbb1);
//dN6[irap][iphi] = Nss;
dN6[irap][iphi] = nsig.getVal();
//dN6err[irap][iphi] = dN6[irap][iphi]*dnn1[irap][iphi];
dN6err[irap][iphi] = nsig.getError();
//dN6err[irap][iphi] = sqrt(Ntt+Nbb);
dNN6[irap][iphi] = dN6err[irap][iphi]/dN6[irap][iphi];
sN6[irap][iphi] = Nss1;
//sN6[irap][iphi] = Nss1;
//dN1[irap][iphi] = nsig.getVal();
//dN1err[irap][iphi] = nsig.getError();
//sN6e[irap][iphi] = sqrt(Ntt1+Nbb1);
sN6e[irap][iphi] = sN6[irap][iphi]*dnn4[irap][iphi];
dnn6[irap][iphi] = sN6e[irap][iphi]/sN6[irap][iphi];
frame->GetYaxis()->SetRangeUser(-0.08*nmax,1.1*nmax);
t->SetTextSize(0.05);
c16->cd(npad3); npad3++;
gPad->SetMargin(0.135, 0.01, 0., 0.);
if(irap==nrapbins-1) gPad->SetMargin(0.135, 0.01, 0.1, 0.);
FormatFrame(frame);
frame->Draw();
li->SetLineColor(4);
li->DrawLine(xx3, -0.08*nmax, xx3, 1.1*nmax);
li->DrawLine(xx4, -0.08*nmax, xx4, 1.1*nmax);
tra6->SetMarkerStyle(20);
tra6->SetMarkerSize(0.6);
tra6->SetMarkerColor(1);
tra6->GetXaxis()->SetLimits(xmin,xmax);
tra6->Draw("psame");
//t->DrawLatexNDC(0.21,0.55,Form("Ns = %4.1f, eNs = %4.1f, dN/N = %5.4f",dN5[irap][iphi],dN5err[irap][iphi],dNN5[irap][iphi]));
//t->DrawLatexNDC(0.21,0.45,Form("3Ns = %4.1f, 3eNs = %4.1f, 3dN/N = %5.4f",sN5[irap][iphi],sN5e[irap][iphi],dnn5[irap][iphi]));
//t->DrawLatexNDC(0.22,0.45,Form("dN/N = %5.4f",dNN3[irap][iphi]));
TPaveText *label_rap = new TPaveText(0.22,0.85,0.52,1.0,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.07);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
TPaveText *label_phi = new TPaveText(0.7,0.3,0.9,0.4,"blNDC");
label_phi->SetBorderSize(0);
label_phi->SetFillStyle(0);
label_phi->SetTextFont(62);
label_phi->SetTextSize(0.07);
// label_phi->SetTextColor(2);
label_phi->AddText(Form("%s<#phi-#Psi_{1}<%s", dphi1[iphi], dphi1[iphi+1]));
label_phi->Draw();
float chi2;
chi2=frame->chiSquare("model", "dh");
t->SetTextSize(0.05);
//t->DrawLatexNDC(0.72,0.4,Form("#chi^{2} = %3.4f",chi2));
//t->DrawLatexNDC(0.22,0.75,Form("#mu = %3.4f",gmean));
//t->DrawLatexNDC(0.22,0.65,Form("#sigma = %3.4f",gsig));
cout << frame->chiSquare("model", "dh") << endl;
cout <<"par0 = "<< par0.getVal() << endl;
}
}
c16->Print(Outpdf);
TGraphErrors *dN6rap[nrapbins];
TGraphErrors *dN6crap[nrapbins];
TGraphErrors* v6y = new TGraphErrors();
TGraphErrors* v6ym = new TGraphErrors();
TGraphErrors* v6cy = new TGraphErrors();
TGraphErrors* v6cym = new TGraphErrors();
float sv6[nrapbins];
float sv6_err[nrapbins];
float sv6c[nrapbins];
float sv6c_err[nrapbins];
for(int i=0;i<nrapbins;i++){
dN6rap[i]=new TGraphErrors();
dN6rap[i]->SetName(Form("dN5rap_%d",i));
dN6crap[i]=new TGraphErrors();
dN6crap[i]->SetName(Form("dN5crap_%d",i));
}
TCanvas* c17= new TCanvas("c17","c17",400*nrapbins,800);
c17->Divide(nrapbins,2);
for(int irap=0; irap<nrapbins; irap++){
for(int iphi=0; iphi<ncutphibins; iphi++){
dN6rap[irap]->SetPoint(iphi,phibincent1[iphi],dN6[irap][iphi]);
dN6rap[irap]->SetPointError(iphi,0,dN6err[irap][iphi]);
dN6crap[irap]->SetPoint(iphi,phibincent1[iphi],sN6[irap][iphi]);
dN6crap[irap]->SetPointError(iphi,0,sN6e[irap][iphi]);
}
c17->cd(irap+1);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN6rap[irap]->SetMarkerStyle(20);
dN6rap[irap]->SetMarkerColor(1);
dN6rap[irap]->SetMarkerSize(1.);
HistStyle(dN6rap[irap]);
HistAxiasName(dN6rap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f1[irap] = new TF1(Form("f1_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f1[irap] ->SetParameters(200, 0.04);
dN6rap[irap]->Fit(Form("f1_%d", irap),"i");
sv6[irap] = f1[irap]->GetParameter(1);
sv6_err[irap] = f1[irap]->GetParError(1);
v6y->SetPoint(irap, rapbincent[irap], sv6[irap]*Res1);
v6y->SetPointError(irap, 0, sv6_err[irap]*Res1);
v6ym->SetPoint(irap, -rapbincent[irap], -sv6[irap]*Res1);
v6ym->SetPointError(irap, 0, sv6_err[irap]*Res1);
dN6rap[irap]->Draw("ap");
TPaveText *label_rap = new TPaveText(0.4,0.8,0.75,0.9,"blNDC");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.06);
// label_rap->SetTextColor(2);
label_rap->AddText(Form("%1.1f < y < %1.1f ", rapbin[irap+1], rapbin[irap]));
label_rap->Draw();
c17->cd(irap+4);
gPad->SetMargin(0.17, 0.03, 0.15, 0.07);
dN6crap[irap]->SetMarkerStyle(20);
dN6crap[irap]->SetMarkerColor(1);
dN6crap[irap]->SetMarkerSize(1.);
HistStyle(dN6crap[irap]);
HistAxiasName(dN6crap[irap], "#phi-#Psi_{1}", "dN/d(#phi-#Psi_{1})");
f1c[irap] = new TF1(Form("f1c_%d", irap), "[0]*(1+[1]*2.*cos(x))", 0, 3.2);
f1c[irap] ->SetParameters(200, 0.04);
dN6crap[irap]->Fit(Form("f1c_%d", irap),"i");
sv6c[irap] = f1c[irap]->GetParameter(1);
sv6c_err[irap] = f1c[irap]->GetParError(1);
v6cy->SetPoint(irap, rapbincent[irap], sv6c[irap]*Res1);
v6cy->SetPointError(irap, 0, sv6c_err[irap]*Res1);
v6cym->SetPoint(irap, -rapbincent[irap], -sv6c[irap]*Res1);
v6cym->SetPointError(irap, 0, sv6c_err[irap]*Res1);
dN6crap[irap]->Draw("ap");
label_rap->SetBorderSize(0);
label_rap->SetFillStyle(0);
label_rap->SetTextFont(62);
label_rap->SetTextSize(0.06);
// label_rap->SetTextColor(2);
label_rap->Draw();
}
c17->Print(Outpdf);
TCanvas* c18= new TCanvas("c18","c18",900,400);
c18->Divide(2,1,0.02,0.02);
gPad->SetMargin(0.17, 0.0, 0.2, 0.03);
v6y->SetMarkerStyle(20);
v6y->SetMarkerSize(1.0);
v6y->SetMarkerColor(2);
v6ym->SetMarkerStyle(24);
v6ym->SetMarkerSize(1.0);
v6ym->SetMarkerColor(2);
v6cy->SetMarkerStyle(20);
v6cy->SetMarkerSize(1.0);
v6cy->SetMarkerColor(2);
v6cym->SetMarkerStyle(24);
v6cym->SetMarkerSize(1.0);
v6cym->SetMarkerColor(2);
htest->GetYaxis()->SetRangeUser(-1.2, 1.2);
htest->GetXaxis()->CenterTitle();
htest->GetYaxis()->CenterTitle();
htest->GetYaxis()->SetTitle("v_{1}");
htest->GetXaxis()->SetTitle("Rapidity (y_{c.m.})");
htest->GetXaxis()->SetTitleSize(0.05);
htest->GetYaxis()->SetTitleSize(0.05);
htest->GetXaxis()->SetTitleOffset(1.1);
htest->GetYaxis()->SetTitleOffset(0.8);
c18->cd(1);
htest->Draw();
v6y->Draw("pesame");
v6ym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
t->DrawLatex(-0.4,1.25,"With efficieny correction");
//TF1* fitv1y = new TF1("fitv1y", "[0]*x+[1]*x*x*x", -0.9, 0.9);
fitv1y = new TF1("fitv1y", "[0]*x", -0.9, 0.9);
//TF1* fitv1y = new TF1("fitv1y", "[0]*x", -0.5, 0.);
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
tl->Draw();
v6y->Fit("fitv1y", "r");
c18->cd(2);
htest->Draw();
v6cy->Draw("pesame");
v6ym->Draw("pesame");
li->SetLineColor(4);
li->SetLineStyle(3);
li->DrawLine(0,-1.2,0,1.2);
li->DrawLine(-0.9,0,0.9,0);
li->DrawLine(-0.5,-1.2,-0.5,1.2);
li->DrawLine(0.5,-1.2,0.5,1.2);
t->SetTextSize(0.04);
t->DrawLatexNDC(0.6,0.8,"0-60%");
t->DrawLatex(-0.4,1.25,"With efficieny correction");
fproton->Draw("same");
fdeuteron->Draw("same");
fHe3->Draw("same");
fHe4->Draw("same");
ftriton->Draw("same");
tl->Draw();
v6cy->Fit("fitv1y", "r");
c18->Print(Outpdf);
c18->Print((Outpdf+"]").Data());
}
void HistStyle(TGraphErrors* hist)
{
hist->GetXaxis()->CenterTitle();
hist->GetYaxis()->CenterTitle();
hist->GetXaxis()->SetTitleFont(62);
hist->GetYaxis()->SetTitleFont(62);
hist->GetXaxis()->SetTitleSize(0.06);
hist->GetYaxis()->SetTitleSize(0.06);
hist->GetXaxis()->SetLabelFont(62);
hist->GetYaxis()->SetLabelFont(62);
hist->GetXaxis()->SetLabelSize(0.05);
hist->GetYaxis()->SetLabelSize(0.05);
hist->GetXaxis()->SetLabelOffset(0.01);
hist->GetYaxis()->SetLabelOffset(0.01);
hist->GetXaxis()->SetTitleOffset(1);
hist->GetYaxis()->SetTitleOffset(1.3);
hist->GetXaxis()->SetNdivisions(510);
hist->GetYaxis()->SetNdivisions(505);
hist->SetMarkerStyle(24);
hist->SetMarkerSize(2);
hist->SetMarkerColor(2);
hist->SetLineWidth(1);
}
void HistAxiasName(TGraphErrors* hist, string xname, string yname)
{
hist->GetXaxis()->SetTitle(xname.c_str());
hist->GetYaxis()->SetTitle(yname.c_str());
}
void FormatFrame(RooPlot* hist)
{
hist->SetTitle("");
hist->GetXaxis()->CenterTitle();
hist->GetYaxis()->CenterTitle();
hist->GetXaxis()->SetTitleFont(62);
hist->GetYaxis()->SetTitleFont(62);
hist->GetXaxis()->SetTitleSize(0.07);
hist->GetYaxis()->SetTitleSize(0.07);
hist->GetXaxis()->SetLabelFont(62);
hist->GetYaxis()->SetLabelFont(62);
hist->GetXaxis()->SetLabelSize(0.06);
hist->GetYaxis()->SetLabelSize(0.06);
hist->GetXaxis()->SetLabelOffset(0.01);
hist->GetYaxis()->SetLabelOffset(0.01);
hist->GetXaxis()->SetTitleOffset(1.1);
hist->GetYaxis()->SetTitleOffset(1.3);
hist->GetXaxis()->SetNdivisions(505);
hist->GetYaxis()->SetNdivisions(506);
//hist->GetYaxis()->SetTitle("Counts");
hist->GetYaxis()->SetTitle("");
if(pid==3004)
hist->GetXaxis()->SetTitle("M(^{3}He#pi^{-})");
if(pid==3005)
hist->GetXaxis()->SetTitle("M(^{4}He#pi^{-})");
}
void FormatFrame(TH1F* hist)
{
hist->SetTitle("");
hist->GetXaxis()->CenterTitle();
hist->GetYaxis()->CenterTitle();
hist->GetXaxis()->SetTitleFont(62);
hist->GetYaxis()->SetTitleFont(62);
hist->GetXaxis()->SetTitleSize(0.07);
hist->GetYaxis()->SetTitleSize(0.07);
hist->GetXaxis()->SetLabelFont(62);
hist->GetYaxis()->SetLabelFont(62);
hist->GetXaxis()->SetLabelSize(0.06);
hist->GetYaxis()->SetLabelSize(0.06);
hist->GetXaxis()->SetLabelOffset(0.01);
hist->GetYaxis()->SetLabelOffset(0.01);
hist->GetXaxis()->SetTitleOffset(1.1);
hist->GetYaxis()->SetTitleOffset(1.3);
hist->GetXaxis()->SetNdivisions(505);
hist->GetYaxis()->SetNdivisions(506);
//hist->GetYaxis()->SetTitle("Counts");
hist->GetYaxis()->SetTitle("");
if(pid==3004)
hist->GetXaxis()->SetTitle("M(^{3}He#pi^{-})");
if(pid==3005)
hist->GetXaxis()->SetTitle("M(^{4}He#pi^{-})");
}
void HistAxiasName(TH1F* hist, string xname, string yname)
{
hist->GetXaxis()->SetTitle(xname.c_str());
hist->GetYaxis()->SetTitle(yname.c_str());
}
void HistStyle(TH1F* hist)
{
hist->SetTitle("");
hist->GetXaxis()->CenterTitle();
hist->GetYaxis()->CenterTitle();
hist->GetXaxis()->SetTitleFont(62);
hist->GetYaxis()->SetTitleFont(62);
hist->GetXaxis()->SetTitleSize(0.06);
hist->GetYaxis()->SetTitleSize(0.06);
hist->GetXaxis()->SetLabelFont(62);
hist->GetYaxis()->SetLabelFont(62);
hist->GetXaxis()->SetLabelSize(0.05);
hist->GetYaxis()->SetLabelSize(0.05);
hist->GetXaxis()->SetLabelOffset(0.01);
hist->GetYaxis()->SetLabelOffset(0.01);
hist->GetXaxis()->SetTitleOffset(1);
hist->GetYaxis()->SetTitleOffset(1.3);
hist->GetXaxis()->SetNdivisions(505);
hist->GetYaxis()->SetNdivisions(505);
hist->SetMarkerStyle(24);
hist->SetMarkerSize(2);
hist->SetMarkerColor(2);
hist->SetLineWidth(1);
}
此处可能存在不合适展示的内容,页面不予展示。您可通过相关编辑功能自查并修改。
如您确认内容无涉及 不当用语 / 纯广告导流 / 暴力 / 低俗色情 / 侵权 / 盗版 / 虚假 / 无价值内容或违法国家有关法律法规的内容,可点击提交进行申诉,我们将尽快为您处理。