I am trying to calculate the yields of both signal and background using gaussian and pol3. I can directly get the integral and its error of both “signal +background” from my fit function “f1_exp” in my code. Then for the yield of the signal only I defined gaussian function, “fsig_exp” in my code and for the yield of the background only I defined pol3 function, “fbac_exp”. Also I defined the covariance matrix to calculate their errors but I have notice that the yields of the background only is greater than the yields of both signal and background. I can’t figure out the issue.
Here is my code
bool_use_partial_covariances = False # 0 ==false, 1= true
fit_with_fixed_background = False
ffile = ROOT.TFile(sys.argv[1])
h1_exp = ffile.Get("eff_ECAL/h1_expMM_pmmbin5")
#define fit function
f1_exp=ROOT.TF1("f1_exp","gaus(0)+pol3(3)",0.75,1.19)
f1_exp.SetParameters(1,1,0.2,1,1,1,1)
f1_exp.SetLineWidth(1)
# open canvas
can_exp = ROOT.TCanvas('can_exp','can_exp',1100,800)
h1_exp.SetTitle("Exp_MM [Pmm = 1.00 - 1.20] ")
h1_exp.Draw("pe")
axisY = h1_exp.GetYaxis()
axisY.SetTitle("Number of events")
r_exp = h1_exp.Fit(f1_exp,"SR")
mu_exp,sig_exp = f1_exp.GetParameter(1), f1_exp.GetParameter(2)
muerr_exp,sigerr_exp = f1_exp.GetParError(1), f1_exp.GetParError(2)
r1,r2 = mu_exp-6*sig_exp,mu_exp+6*sig_exp
#define fit function for only signal
fsig_exp = ROOT.TF1("fsig_exp","gaus",r1,r2)
#define fit function for only bg
fbac_exp = ROOT.TF1("fbac_exp", "pol3(3)",r1,r2)
for x in range(3):
fsig_exp.SetParameter(x,f1_exp.GetParameter(x))
for x in range(4):
fbac_exp.SetParameter(x+3,f1_exp.GetParameter(x+3))
fsig_exp.SetLineColor(4)
fsig_exp.SetLineWidth(1)
fsig_exp.Draw('same')
fbac_exp.SetLineColor(8)
fbac_exp.SetLineWidth(1)
fbac_exp.Draw('same')
# get cov Matrix to find error
c_exp = r_exp.GetCovarianceMatrix()
# c1.Print()
#++++ gaua signal +++++
c_exp_sig = c_exp.GetSub(0, 2, 0, 2)
# c_exp_sig.Print()
if bool_use_partial_covariances:
c_exp_sig -=c_exp.GetSub(0, 2, 3, 6) *c_exp.GetSub(3, 6, 3, 6).InvertFast() *c_exp.GetSub(3, 6, 0, 2)
#+++++ background signal +++++
c_exp_bg = c_exp.GetSub(3, 6, 3, 6)
# c_exp_bg.Print()
if bool_use_partial_covariances:
c_exp_bg -=c_exp.GetSub(3, 6, 0, 2) *c_exp.GetSub(0, 2, 0, 2).InvertFast() *c_exp.GetSub(0, 2, 3, 6)
N_tot_exp = (abs(f1_exp.Integral(r1,r2)))/h1_exp.GetBinWidth(1)
N_totErr_exp = (abs(f1_exp.IntegralError(r1,r2)))/h1_exp.GetBinWidth(1)
N_totErr_exp_test = (abs(f1_exp.IntegralError(r1,r2, f1_exp.GetParameters(),c_exp.GetMatrixArray())))/h1_exp.GetBinWidth(1)
# fsig.SetParameters(f1_exp.GetParameters())
N_sig_exp = (abs(fsig_exp.Integral(r1,r2)))/h1_exp.GetBinWidth(1)
N_sigErr_exp = (abs(fsig_exp.IntegralError(r1,r2,fsig_exp.GetParameters(),c_exp_sig.GetMatrixArray())))/h1_exp.GetBinWidth(1)
N_bac_exp = (abs(fbac_exp.Integral(r1,r2)))/h1_exp.GetBinWidth(1)
N_bacErr_exp = (abs(fbac_exp.IntegralError(r1,r2, fbac_exp.GetParameters(),c_exp_bg.GetMatrixArray() )))/h1_exp.GetBinWidth(1)
# fbac_exp_bin2.SetParameters(f1_exp.GetParameters() + 3)
print( "N_tot_exp" ,N_tot_exp, " +/-", N_totErr_exp)
print( "N_tot_exp_test" ,N_tot_exp, " +/-", N_totErr_exp_test)
print( "N_sig_exp" ,N_sig_exp, " +/-", N_sigErr_exp)
print( "N_bac_exp" ,N_bac_exp, " +/-", N_bacErr_exp)
Could you please help me to get them correctly.
Many thanks !!