23 import matplotlib
as mpl
25 mpl.rcParams[
'text.usetex'] =
True
26 from defaultEvaluationParameters
import categories, rbins, r_subsample, r_size
27 import matplotlib.pyplot
as plt
28 from array
import array
38 workingFilesMC = str(sys.argv[1])
39 workingFilesData = str(sys.argv[2])
41 workingFilesMC = glob.glob(str(workingFilesMC))
42 workingFilesData = glob.glob(str(workingFilesData))
48 treenames = [str(sys.argv[3])]
50 mc_total_notTagged = 0
52 data_total_notTagged = 0
56 ROOT.RooMsgService.instance().setGlobalKillBelow(ROOT.RooFit.WARNING)
59 def efficiencyCalculator(mc, data, method):
63 mc_dataB0 = mc[np.where(mc[:, 1] > 0)]
64 mc_dataB0bar = mc[np.where(mc[:, 1] < 0)]
66 mc_totaldataB0 = np.absolute(mc_dataB0[:, 0])
67 mc_totaldataB0bar = np.absolute(mc_dataB0bar[:, 0])
69 mc_rvalueB0 = (np.histogram(mc_totaldataB0, rbins, weights=mc_totaldataB0)[0] / np.histogram(mc_totaldataB0, rbins)[0])
74 weights=mc_totaldataB0 *
79 mc_rvalueStdB0 = np.sqrt(abs(mc_rvalueMSB0 - mc_rvalueB0 * mc_rvalueB0))
80 mc_entriesB0 = np.histogram(mc_totaldataB0, rbins)[0]
86 weights=mc_totaldataB0bar)[0] /
94 weights=mc_totaldataB0bar *
95 mc_totaldataB0bar)[0] /
99 mc_rvalueStdB0bar = np.sqrt(abs(mc_rvalueMSB0bar - mc_rvalueB0bar * mc_rvalueB0bar))
100 mc_entriesB0bar = np.histogram(mc_totaldataB0bar, rbins)[0]
102 mc_total_entriesB0 = mc_totaldataB0.shape[0] + mc_total_notTagged / 2
103 mc_total_tagged_B0 = mc_totaldataB0.shape[0]
104 mc_event_fractionB0 = mc_entriesB0.astype(float) / mc_total_entriesB0
106 mc_total_entriesB0bar = mc_totaldataB0bar.shape[0] + mc_total_notTagged / 2
107 mc_total_tagged_B0bar = mc_totaldataB0bar.shape[0]
108 mc_event_fractionB0bar = mc_entriesB0bar.astype(float) / mc_total_entriesB0bar
110 mc_total_entries = mc_total_entriesB0 + mc_total_entriesB0bar
111 mc_total_tagged = mc_totaldataB0.shape[0] + mc_totaldataB0bar.shape[0]
112 mc_event_fractionTotal = (mc_entriesB0.astype(float) + mc_entriesB0bar.astype(float)) / mc_total_entries
114 mc_tagging_eff = mc_total_tagged / (mc_total_tagged + mc_total_notTagged)
115 mc_DeltaTagging_eff = math.sqrt(mc_total_tagged * mc_total_notTagged**2 +
116 mc_total_notTagged * mc_total_tagged**2) / (mc_total_entries**2)
118 mc_tot_eff_effB0bar = 0
119 mc_uncertainty_eff_effB0 = 0
120 mc_uncertainty_eff_effB0bar = 0
121 mc_uncertainty_eff_effAverage = 0
122 mc_diff_eff_Uncertainty = 0
123 mc_event_fractionDiff = array(
'f', [0] * r_size)
124 mc_rvalueB0Average = array(
'f', [0] * r_size)
125 mc_rvalueStdB0Average = array(
'f', [0] * r_size)
126 mc_wvalue = array(
'f', [0] * r_size)
127 mc_wvalueUncertainty = array(
'f', [0] * r_size)
128 mc_wvalueB0 = array(
'f', [0] * r_size)
129 mc_wvalueB0bar = array(
'f', [0] * r_size)
130 mc_wvalueDiff = array(
'f', [0] * r_size)
131 mc_wvalueDiffUncertainty = array(
'f', [0] * r_size)
132 mc_entries = array(
'f', [0] * r_size)
133 mc_iEffEfficiency = array(
'f', [0] * r_size)
134 mc_iEffEfficiencyUncertainty = array(
'f', [0] * r_size)
135 mc_iEffEfficiencyB0Uncertainty = array(
'f', [0] * r_size)
136 mc_iEffEfficiencyB0barUncertainty = array(
'f', [0] * r_size)
137 mc_iDeltaEffEfficiency = array(
'f', [0] * r_size)
138 mc_iDeltaEffEfficiencyUncertainty = array(
'f', [0] * r_size)
140 print(
'****************** MEASURED EFFECTIVE EFFICIENCY FOR COMBINER USING ' + method +
' ***************************')
142 print(
'* --> DETERMINATION BASED ON MONTE CARLO *')
144 print(
'* ------------------------------------------------------------------------------------------------ *')
145 print(
'* r-interval <r> Efficiency Delta_Effcy w Delta_w *')
146 print(
'* ------------------------------------------------------------------------------------------------ *')
148 for i
in range(0, r_size - 1):
149 mc_rvalueB0Average[i] = (mc_rvalueB0[i] + mc_rvalueB0bar[i]) / 2
150 mc_rvalueStdB0Average[i] = math.sqrt(mc_rvalueStdB0[i]**2 / (mc_entriesB0[0] - 1) +
151 mc_rvalueStdB0bar[i]**2 / (mc_entriesB0bar[0] - 1)) / 2
152 mc_wvalueB0[i] = (1 - mc_rvalueB0[i]) / 2
153 mc_wvalueB0bar[i] = (1 - mc_rvalueB0bar[i]) / 2
154 mc_wvalueDiff[i] = mc_wvalueB0[i] - mc_wvalueB0bar[i]
155 mc_wvalueDiffUncertainty[i] = math.sqrt((mc_rvalueStdB0[i] / 2)**2 / (mc_entriesB0[0] - 1) +
156 (mc_rvalueStdB0bar[i] / 2)**2 / (mc_entriesB0bar[0] - 1))
157 mc_wvalue[i] = (mc_wvalueB0[i] + mc_wvalueB0bar[i]) / 2
158 mc_wvalueUncertainty[i] = mc_wvalueDiffUncertainty[i] / 2
161 mc_event_fractionDiff[i] = (mc_entriesB0[i] - mc_entriesB0bar[i]) / mc_total_entries
163 mc_iEffEfficiency[i] = mc_tagging_eff * (mc_event_fractionB0[i] * mc_rvalueB0[i] * mc_rvalueB0[i] +
164 mc_event_fractionB0bar[i] * mc_rvalueB0bar[i] * mc_rvalueB0bar[i]) / 2
166 mc_iEffEfficiencyB0Uncertainty[i] = mc_rvalueB0[i] * \
167 math.sqrt((2 * mc_total_entriesB0 * mc_entriesB0[i] * mc_rvalueStdB0[i])**2 / (mc_entriesB0[0] - 1) +
168 mc_rvalueB0[i]**2 * mc_entriesB0[i] *
169 (mc_total_entriesB0 * (mc_total_entriesB0 - mc_entriesB0[i]) +
170 mc_entriesB0[i] * mc_total_notTagged)) / (mc_total_entriesB0**2)
171 mc_iEffEfficiencyB0barUncertainty[i] = mc_rvalueB0bar[i] * \
172 math.sqrt((2 * mc_total_entriesB0bar * mc_entriesB0bar[i] * mc_rvalueStdB0bar[i])**2 / (mc_entriesB0bar[0] - 1) +
173 mc_rvalueB0bar[i]**2 * mc_entriesB0bar[i] *
174 (mc_total_entriesB0bar * (mc_total_entriesB0bar - mc_entriesB0bar[i]) +
175 mc_entriesB0bar[i] * mc_total_notTagged)) / (mc_total_entriesB0bar**2)
177 mc_iDeltaEffEfficiency[i] = mc_event_fractionB0[i] * mc_rvalueB0[i] * \
178 mc_rvalueB0[i] - mc_event_fractionB0bar[i] * mc_rvalueB0bar[i] * mc_rvalueB0bar[i]
180 mc_iDeltaEffEfficiencyUncertainty[i] = math.sqrt(
181 mc_iEffEfficiencyB0Uncertainty[i]**2 +
182 mc_iEffEfficiencyB0barUncertainty[i]**2)
184 mc_iEffEfficiencyUncertainty[i] = mc_iDeltaEffEfficiencyUncertainty[i] / 2
186 mc_diff_eff_Uncertainty = mc_diff_eff_Uncertainty + mc_iDeltaEffEfficiencyUncertainty[i]**2
189 mc_tot_eff_effB0 = mc_tot_eff_effB0 + mc_event_fractionB0[i] * mc_rvalueB0[i] * mc_rvalueB0[i]
190 mc_tot_eff_effB0bar = mc_tot_eff_effB0bar + mc_event_fractionB0bar[i] * mc_rvalueB0bar[i] * mc_rvalueB0bar[i]
191 mc_uncertainty_eff_effAverage = mc_uncertainty_eff_effAverage + mc_iEffEfficiencyUncertainty[i]**2
192 mc_uncertainty_eff_effB0 = mc_uncertainty_eff_effB0 + mc_iEffEfficiencyB0Uncertainty[i]**2
193 mc_uncertainty_eff_effB0bar = mc_uncertainty_eff_effB0bar + mc_iEffEfficiencyB0barUncertainty[i]**2
196 print(
'* ' +
'{:.3f}'.format(r_subsample[i]) +
' - ' +
'{:.3f}'.format(r_subsample[i + 1]) +
' ' +
197 '{:.3f}'.format(mc_rvalueB0Average[i]) +
' +- ' +
'{:.4f}'.format(mc_rvalueStdB0Average[i]) +
' ' +
198 '{:.4f}'.format(mc_event_fractionTotal[i]) +
' ' +
199 '{: .4f}'.format(mc_event_fractionDiff[i]) +
' ' +
200 '{:.4f}'.format(mc_wvalue[i]) +
' +- ' +
'{:.4f}'.format(mc_wvalueUncertainty[i]) +
' ' +
201 '{: .4f}'.format(mc_wvalueDiff[i]) +
' +- ' +
'{:.4f}'.format(mc_wvalueDiffUncertainty[i]) +
' *')
203 mc_average_eff_eff = (mc_tot_eff_effB0 + mc_tot_eff_effB0bar) / 2
204 mc_uncertainty_eff_effAverage = math.sqrt(mc_uncertainty_eff_effAverage)
205 mc_uncertainty_eff_effB0 = math.sqrt(mc_uncertainty_eff_effB0)
206 mc_uncertainty_eff_effB0bar = math.sqrt(mc_uncertainty_eff_effB0bar)
207 mc_diff_eff = mc_tot_eff_effB0 - mc_tot_eff_effB0bar
208 mc_diff_eff_Uncertainty = math.sqrt(mc_diff_eff_Uncertainty)
210 print(
'* ------------------------------------------------------------------------------------------------ *')
212 print(
'* __________________________________________________________________________________________ *')
214 print(
'* | TOTAL NUMBER OF TAGGED EVENTS = ' +
215 '{:<24}'.format(
"%.0f" % mc_total_tagged) +
'{:>36}'.format(
'| *'))
218 '* | TOTAL AVERAGE EFFICIENCY (q=+-1)= ' +
224 mc_DeltaTagging_eff *
229 '* | TOTAL AVERAGE EFFECTIVE EFFICIENCY (q=+-1)= ' +
235 mc_uncertainty_eff_effAverage *
240 '* | TOTAL AVERAGE EFFECTIVE EFFICIENCY ASYMMETRY (q=+-1)= ' +
246 mc_diff_eff_Uncertainty *
250 print(
'* | B0-TAGGER TOTAL EFFECTIVE EFFICIENCIES: ' +
251 '{:.2f}'.format(mc_tot_eff_effB0 * 100) +
" +-" +
'{: 4.2f}'.format(mc_uncertainty_eff_effB0 * 100) +
253 '{:.2f}'.format(mc_tot_eff_effB0bar * 100) +
" +-" +
'{: 4.2f}'.format(mc_uncertainty_eff_effB0bar * 100) +
254 ' % (q=-1) ' +
' | *')
256 print(
'* | FLAVOR PERCENTAGE (MC): ' +
257 '{:.2f}'.format(mc_total_tagged_B0 / mc_total_tagged * 100) +
' % (q=+1) ' +
258 '{:.2f}'.format(mc_total_tagged_B0bar / mc_total_tagged * 100) +
' % (q=-1) Diff=' +
259 '{:^5.2f}'.format((mc_total_tagged_B0 - mc_total_tagged_B0bar) / mc_total_tagged * 100) +
' % | *')
260 print(
'* |__________________________________________________________________________________________| *')
262 print(
'****************************************************************************************************')
264 print(
r'\begin{tabular}{ l r r r r r r r }')
266 print(
r'$r$- Interval & $\varepsilon_i\ $ & $w_i \pm \delta w_i\enskip\, $ ' +
267 r' & $\Delta w_i \pm \delta\Delta w_i $& $\varepsilon_{\text{eff}, i} \pm \delta\varepsilon_{\text{eff}, i}\enskip$ ' +
268 r' & & & $\Delta \varepsilon_{\text{eff}, i} \pm \delta\Delta \varepsilon_{\text{eff}, i} $\\ \hline\hline')
269 for i
in range(0, r_size - 1):
270 print(
'$ ' +
'{:.3f}'.format(r_subsample[i]) +
' - ' +
'{:.3f}'.format(r_subsample[i + 1]) +
'$ & $'
271 '{: 6.1f}'.format(mc_event_fractionTotal[i] * 100) +
'$ & $' +
272 '{: 7.3f}'.format(mc_wvalue[i] * 100) +
r" \pm " +
'{:2.3f}'.format(mc_wvalueUncertainty[i] * 100) +
r' $ & $' +
273 '{: 6.1f}'.format(mc_wvalueDiff[i] * 100) +
r" \pm " +
274 '{:2.3f}'.format(mc_wvalueDiffUncertainty[i] * 100) +
r'\, $ & $' +
275 '{: 8.4f}'.format(mc_iEffEfficiency[i] * 100) +
276 r" \pm " +
'{:2.4f}'.format(mc_iEffEfficiencyUncertainty[i] * 100) +
r'\, $ & & & $' +
277 '{: 6.4f}'.format(mc_iDeltaEffEfficiency[i] * 100) +
278 r" \pm " +
'{:2.4f}'.format(mc_iDeltaEffEfficiencyUncertainty[i] * 100) +
280 print(
r'\hline\hline')
282 r'\multicolumn{1}{r}{Total} & & \multicolumn{3}{r}{ $\varepsilon_\text{eff} = ' +
283 r'\sum_i \varepsilon_i \cdot \langle 1-2w_i\rangle^2 = ' +
289 mc_uncertainty_eff_effAverage *
292 print(
r'& & \multicolumn{2}{r}{ $\Delta \varepsilon_\text{eff} = ' +
293 '{: 6.2f}'.format(mc_diff_eff * 100) +
r" \pm " +
'{: 6.2f}'.format(mc_diff_eff_Uncertainty * 100) +
r'\ \quad\, $ }' +
296 print(
r'\end{tabular}')
300 data_totaldata = np.absolute(data[:, 0])
301 data_splotWeights = data[:, 1]
303 data_entries = np.histogram(data_totaldata, rbins, weights=data_splotWeights)[0]
304 data_rvalueB0Average = (np.histogram(data_totaldata, rbins, weights=data_totaldata * data_splotWeights)[0] / data_entries)
305 data_rvalueMSB0Average = (
309 weights=(data_totaldata)**2 *
310 data_splotWeights)[0] /
312 data_rvalueStdB0Average = np.sqrt(abs(data_rvalueMSB0Average - data_rvalueB0Average * data_rvalueB0Average))
314 data_total_tagged = 0
315 for iEntries
in data_entries:
316 data_total_tagged += iEntries
318 data_total_notTaggedWeighted = data_total_notTagged * data_total_tagged / data_totaldata.shape[0]
320 data_total_entries = data_total_tagged + data_total_notTaggedWeighted
321 data_event_fractionTotal = data_entries.astype(float) / data_total_entries
323 data_tagging_eff = data_total_tagged / data_total_entries
324 data_DeltaTagging_eff = math.sqrt(data_total_tagged * data_total_notTaggedWeighted**2 +
325 data_total_notTaggedWeighted * data_total_tagged**2) / (data_total_entries**2)
326 data_average_eff_eff = 0
327 data_uncertainty_eff_effAverage = 0
328 data_wvalue = array(
'f', [0] * r_size)
329 data_wvalueUncertainty = array(
'f', [0] * r_size)
330 data_iEffEfficiency = array(
'f', [0] * r_size)
331 data_iEffEfficiencyUncertainty = array(
'f', [0] * r_size)
333 print(
'****************** MEASURED EFFECTIVE EFFICIENCY FOR COMBINER USING ' + method +
' ***************************')
335 print(
'* --> DETERMINATION BASED ON DATA *')
337 print(
'* ------------------------------------------------------------------------------------------------ *')
338 print(
'* r-interval <r> Efficiency w *')
339 print(
'* ------------------------------------------------------------------------------------------------ *')
341 for i
in range(0, r_size - 1):
342 data_wvalue[i] = (1 - data_rvalueB0Average[i]) / 2
343 data_wvalueUncertainty[i] = (data_rvalueStdB0Average[i] / math.sqrt(data_entries[i] - 1)) / 2
346 data_iEffEfficiency[i] = data_tagging_eff * (data_event_fractionTotal[i] *
347 data_rvalueB0Average[i] * data_rvalueB0Average[i])
349 data_iEffEfficiencyUncertainty[i] = data_rvalueB0Average[i] * \
350 math.sqrt((2 * data_total_entries * data_entries[i] *
351 (data_rvalueStdB0Average[i] / math.sqrt(data_entries[i] - 1)))**2 +
352 data_rvalueB0Average[i]**2 * data_entries[i] *
353 (data_total_entries * (data_total_entries - data_entries[i]) +
354 data_entries[i] * data_total_notTaggedWeighted)) / (data_total_entries**2)
357 data_average_eff_eff = data_average_eff_eff + \
358 data_event_fractionTotal[i] * data_rvalueB0Average[i] * data_rvalueB0Average[i]
359 data_uncertainty_eff_effAverage = data_uncertainty_eff_effAverage + data_iEffEfficiencyUncertainty[i]**2
362 print(
'* ' +
'{:.3f}'.format(r_subsample[i]) +
' - ' +
'{:.3f}'.format(r_subsample[i + 1]) +
' ' +
363 '{:.3f}'.format(data_rvalueB0Average[i]) +
' +- ' +
'{:.4f}'.format(data_rvalueStdB0Average[i]) +
' ' +
364 '{:.4f}'.format(data_event_fractionTotal[i]) +
' ' +
365 '{:.4f}'.format(data_wvalue[i]) +
' +- ' +
'{:.4f}'.format(data_wvalueUncertainty[i]) +
' ' +
' *')
367 data_uncertainty_eff_effAverage = math.sqrt(data_uncertainty_eff_effAverage)
369 print(
'* ------------------------------------------------------------------------------------------------ *')
371 print(
'* __________________________________________________________________________________________ *')
373 print(
'* | TOTAL NUMBER OF TAGGED EVENTS = ' +
374 '{:<24}'.format(
"%.0f" % data_total_tagged) +
'{:>36}'.format(
'| *'))
377 '* | TOTAL AVERAGE EFFICIENCY (q=+-1)= ' +
383 data_DeltaTagging_eff *
388 '* | TOTAL AVERAGE EFFECTIVE EFFICIENCY (q=+-1)= ' +
390 data_average_eff_eff *
394 data_uncertainty_eff_effAverage *
397 print(
'* |__________________________________________________________________________________________| *')
399 print(
'****************************************************************************************************')
401 print(
r'\begin{tabular}{ l r r r r}')
403 print(
r'$r$- Interval & $\varepsilon_i\ $ & $w_i \pm \delta w_i\ $ ' +
404 r' & $\varepsilon_{\text{eff}, i} \pm \delta\varepsilon_{\text{eff}, i}$ ' +
406 for i
in range(0, len(rbins) - 1):
407 print(
'$ ' +
'{:.3f}'.format(r_subsample[i]) +
' - ' +
'{:.3f}'.format(r_subsample[i + 1]) +
'$ & $'
408 '{: 6.2f}'.format(data_event_fractionTotal[i] * 100) +
'$ & $' +
409 '{: 6.2f}'.format(data_wvalue[i] * 100) +
r" \pm " +
'{:2.2f}'.format(data_wvalueUncertainty[i] * 100) +
r'\, $ & $' +
410 '{: 7.3f}'.format(data_iEffEfficiency[i] * 100) +
411 r" \pm " +
'{:4.3f}'.format(data_iEffEfficiencyUncertainty[i] * 100) +
r'\, $ \\ ')
412 print(
r'\hline\hline')
414 r'\multicolumn{1}{r}{Total} & \multicolumn{3}{r}{ $\varepsilon_\text{eff} = ' +
415 r'\sum_i \varepsilon_i \cdot \langle 1-2w_i\rangle^2 = ' +
417 data_average_eff_eff *
421 data_uncertainty_eff_effAverage *
426 print(
r'\end{tabular}')
431 bins = list(range(-25, 26, 1))
432 for i
in range(0, len(bins)):
433 bins[i] = float(bins[i]) / 25
436 fig1 = plt.figure(1, figsize=(11, 8))
437 ax1 = plt.axes([0.21, 0.15, 0.76, 0.8])
439 qrDataPoints, qrDataBins = np.histogram(data[:, 0], bins=bins, weights=data_splotWeights)
440 qrDataPointsStd = np.sqrt(qrDataPoints)
442 qrDataPointsNorm, qrDataBins = np.histogram(data[:, 0], bins=bins, weights=data_splotWeights, density=1)
443 qrDataPointsStdNorm = qrDataPointsStd * qrDataPointsNorm / qrDataPoints
445 qrBincenters = 0.5 * (qrDataBins[1:] + qrDataBins[:-1])
447 p3 = ax1.errorbar(qrBincenters, qrDataPointsNorm, xerr=0.02, yerr=qrDataPointsStdNorm, elinewidth=2, mew=2, ecolor=
'k',
448 fmt=
'o', mfc=
'k', markersize=6, label=
r'${\rm Data}$')
450 ax1.hist(mc[:, 0], bins=bins, density=1, histtype=
'step',
451 edgecolor=
'b', linewidth=4, linestyle=
'dashed', label=
r'${\rm MC}$')
453 p2, = ax1.plot([], label=
r'${\rm MC}$', linewidth=4, linestyle=
'dashed', c=
'b')
459 ax1.set_ylabel(
r'${\rm Fraction\ \ of\ \ Events\ /\ (\, 0.04\, )}$', fontsize=35)
460 ax1.set_xlabel(
r'$(q\cdot r)_{\rm ' + methodLabel +
'}$', fontsize=35)
461 plt.xticks([-1, -0.75, -0.5, -0.25, 0, 0.25, 0.5, 0.75, 1],
462 [
r'$-1$',
r'',
r'$-0.5$',
r'',
r'$0$',
r'',
r'$0.5$',
r'',
r'$1$'], rotation=0, size=40)
463 ax1.tick_params(axis=
'y', labelsize=35)
464 ax1.legend([p3, p2], [
r'${\rm Data}$',
r'${\rm MC}$'], prop={
'size': 35}, loc=0, numpoints=1)
467 if qrDataPointsNorm[np.argmax(qrDataPointsNorm)] < 1.4:
470 ax1.set_xlim(-1.005, 1.005)
471 plt.savefig(PATH +
'/QR_B0_B0bar' + methodLabel +
'.pdf')
476 fig2 = plt.figure(2, figsize=(11, 8))
477 ax1 = plt.axes([0.21, 0.15, 0.76, 0.8])
479 qrDataPointsCont, qrDataBins = np.histogram(data[:, 0], bins=bins, weights=data[:, 2])
480 qrDataPointsContStd = np.sqrt(qrDataPointsCont)
482 qrDataPointsContNorm, qrDataBins = np.histogram(data[:, 0], bins=bins, weights=data[:, 2], density=1)
483 qrDataPointsContStdNorm = qrDataPointsContStd * qrDataPointsContNorm / qrDataPointsCont
485 qrBincenters = 0.5 * (qrDataBins[1:] + qrDataBins[:-1])
487 p3 = ax1.errorbar(qrBincenters, qrDataPointsContNorm, xerr=0.02, yerr=qrDataPointsContStdNorm, elinewidth=2, mew=2, ecolor=
'k',
488 fmt=
'o', mfc=
'k', mec=
'#990000', markersize=6, label=
r'${\rm Data}$')
490 ax1.set_ylabel(
r'${\rm Fraction\ \ of\ \ Events\ /\ (\, 0.04\, )}$', fontsize=35)
491 ax1.set_xlabel(
r'$(q\cdot r)_{\rm ' + methodLabel +
'}$', fontsize=35)
492 plt.xticks([-1, -0.75, -0.5, -0.25, 0, 0.25, 0.5, 0.75, 1],
493 [
r'$-1$',
r'',
r'$-0.5$',
r'',
r'$0$',
r'',
r'$0.5$',
r'',
r'$1$'], rotation=0, size=40)
494 ax1.tick_params(axis=
'y', labelsize=35)
495 ax1.legend([p3], [
r'${\rm Data}$'], prop={
'size': 35}, loc=1, numpoints=1)
498 if qrDataPointsContNorm[np.argmax(qrDataPointsContNorm)] < 1.4:
501 ax1.set_xlim(-1.005, 1.005)
502 plt.savefig(PATH +
'/QR_B0_B0bar' + methodLabel +
'Continuum.pdf')
505 data_total_tagged_Continuum = 0
506 for iEntries
in qrDataPointsCont:
507 data_total_tagged_Continuum += iEntries
509 print(
'****************************************************************************************************')
510 print(
'* TOTAL NUMBER OF TAGGED CONTINUUM EVENTS = ' +
511 '{:<24}'.format(
"%.0f" % data_total_tagged_Continuum) +
'{:>28}'.format(
' *'))
515 mcHistW = ROOT.TH1F(
"mcHistW",
"mcHistW", int(r_size - 2), r_subsample)
516 dataHistW = ROOT.TH1F(
"dataHistW",
"dataHistW", int(r_size - 2), r_subsample)
518 mcHistEff = ROOT.TH1F(
"mcHistEff",
"mcHistEff", int(r_size - 2), r_subsample)
519 dataHistEff = ROOT.TH1F(
"dataHistEff",
"dataHistEff", int(r_size - 2), r_subsample)
521 print(
'****************************************************************************************************')
523 print(
'* --> TEST OF COMPATIBILITY BETWEEN DATA AND MC *')
525 for i
in range(0, len(rbins) - 1):
527 mcHistW.SetBinContent(i + 1, mc_wvalue[i])
528 mcHistW.SetBinError(i + 1, mc_wvalueUncertainty[i])
530 dataHistW.SetBinContent(i + 1, data_wvalue[i])
531 dataHistW.SetBinError(i + 1, data_wvalueUncertainty[i])
533 mcHistEff.SetBinContent(i + 1, mc_iEffEfficiency[i])
534 mcHistEff.SetBinError(i + 1, mc_iEffEfficiencyUncertainty[i])
536 dataHistEff.SetBinContent(i + 1, data_iEffEfficiency[i])
537 dataHistEff.SetBinError(i + 1, data_iEffEfficiencyUncertainty[i])
539 residualsW = array(
'd', [0] * r_size)
540 residualsEff = array(
'd', [0] * r_size)
542 print(
'* TEST ON w *')
545 chi2W = dataHistW.Chi2Test(mcHistW,
"WW CHI2 P", residualsW)
546 print(
"1-CL for NDF = 7 is = ", ROOT.TMath.Prob(chi2W, 7))
548 print(
'* TEST ON EFFECTIVE EFFICIENCY *')
550 chi2Eff = dataHistEff.Chi2Test(mcHistEff,
"WW P", residualsEff)
551 print(
"1-CL for NDF = 7 is = ", ROOT.TMath.Prob(chi2Eff, 7))
553 print(
'****************************************************************************************************')
557 qrMCPoints, qrDataBins = np.histogram(mc[:, 0], bins=bins)
558 qrMCPointsStd = np.sqrt(qrMCPoints)
559 qrMCPointsNorm, qrDataBins = np.histogram(mc[:, 0], bins=bins, density=1)
560 qrMCPointsStdNorm = qrMCPointsStd * qrMCPointsNorm / qrMCPoints
562 qrDataMCResiduals = (qrDataPointsNorm - qrMCPointsNorm) \
563 / np.sqrt(qrDataPointsStdNorm * qrDataPointsStdNorm + qrMCPointsStdNorm * qrMCPointsStdNorm)
565 fig3 = plt.figure(3, figsize=(11, 11))
566 ax1 = plt.axes([0.21, 0.37, 0.76, 0.60])
569 p3 = ax1.errorbar(qrBincenters, qrDataPointsNorm, xerr=0.02, yerr=qrDataPointsStdNorm,
570 elinewidth=2, mew=2, ecolor=
'k',
571 fmt=
'o', mfc=
'k', mec=
'#006600', markersize=6, label=
r'${\rm Data}$')
573 ax1.hist(mc[:, 0], bins=bins, density=1, histtype=
'step',
574 edgecolor=
'b', linewidth=4, linestyle=
'dashed', label=
r'${\rm MC}$')
576 p2, = ax1.plot([], label=
r'${\rm MC}$', linewidth=4, linestyle=
'dashed', c=
'b')
578 ax1.set_ylabel(
r'${\rm Fraction\ \ of\ \ Events\ /\ (\, 0.04\, )}$', fontsize=35, labelpad=20)
580 ax1.tick_params(axis=
'y', labelsize=35)
581 ax1.legend([p3, p2], [
r'${\rm Data}$',
r'${\rm MC}$'], prop={
'size': 35}, loc=1, numpoints=1)
584 if qrDataPointsNorm[np.argmax(qrDataPointsNorm)] < 1.4:
586 plt.yticks([0, 0.25, 0.5, 0.75, 1.0, 1.25],
587 [
r'$0.00$',
r'$0.25$',
r'$0.5$',
r'$0.75$',
r'$1.00$',
r'$1.25$'], rotation=0, size=35)
588 ax1.set_xlim(-1.005, 1.005)
589 plt.xticks([-1, -0.75, -0.5, -0.25, 0, 0.25, 0.5, 0.75, 1],
590 [
r'',
r'',
r'',
r'',
r'',
r'',
r'',
r'',
r''], rotation=0, size=35)
592 ax2 = plt.axes([0.21, 0.15, 0.76, 0.2])
593 ax2.errorbar(qrBincenters, qrDataMCResiduals, xerr=0.02, yerr=1, elinewidth=2, mew=2, ecolor=
'k',
594 fmt=
'o', mfc=
'k', mec=
'#006600', markersize=6, label=
r'${\rm Data}$')
595 ax2.set_ylabel(
r'${\rm Normalized}$' +
'\n' +
r'${\rm Residuals}$', fontsize=35, labelpad=20)
596 ax2.set_xlabel(
r'$(q\cdot r)_{\rm ' + methodLabel +
'}$', fontsize=35)
597 plt.xticks([-1, -0.75, -0.5, -0.25, 0, 0.25, 0.5, 0.75, 1],
598 [
r'$-1$',
r'',
r'$-0.5$',
r'',
r'$0$',
r'',
r'$0.5$',
r'',
r'$1$'], rotation=0, size=35)
599 plt.yticks([-5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5],
600 [
r'',
r'$-4$',
r'',
r'$-2$',
r'',
r'$0$',
r'',
r'$2$',
r'',
r'$4$',
r''], rotation=0, size=25)
603 ax2.set_xlim(-1.005, 1.005)
606 plt.axhline(y=3, xmin=-1.005, xmax=1.005, linewidth=2, color=
'#006600', linestyle=
'-')
608 plt.axhline(y=1, xmin=-1.005, xmax=1.005, linewidth=2, color=
'b', linestyle=
'-')
609 plt.axhline(y=-1, xmin=-1.005, xmax=1.005, linewidth=2, color=
'b', linestyle=
'-')
611 plt.axhline(y=-3, xmin=-1.005, xmax=1.005, linewidth=2, color=
'#006600', linestyle=
'-')
613 plt.savefig(PATH +
'/QR_B0_B0bar' + methodLabel +
'WithRes.pdf')
619 def plotWithResiduals(rooFitFrame, rooRealVar, dots, modelCurve, units, nameOfPlot, legend, removeArtifacts):
623 rooFitFrame.GetXaxis().SetTitle(
"")
624 rooFitFrame.GetXaxis().SetLabelSize(0)
626 rooFitFrame.GetYaxis().SetTitleSize(0.072)
627 rooFitFrame.GetYaxis().SetTitleOffset(0.98)
628 rooFitFrame.GetYaxis().SetLabelSize(0.055)
631 yLabelBin =
'{:0.0f}'.format(float(dots.GetErrorXlow(1) + dots.GetErrorXhigh(1)))
633 pointsHist = ROOT.RooHist()
637 xValModel = ROOT.Double(-1.E30)
638 yValModel = ROOT.Double(-1.E30)
639 xValDot = ROOT.Double(-1.E30)
640 yValDot = ROOT.Double(-1.E30)
642 iDotPoint = ROOT.RooRealVar(
"iDotPoint",
"", 0.)
643 iModelPoint = ROOT.RooRealVar(
"iModelPoint",
"", 0.)
644 iDotError = ROOT.RooRealVar(
"iDotError",
"", 0.)
646 iResidual = ROOT.RooFormulaVar(
"yComb",
"",
"(@0 - @1)/@2", ROOT.RooArgList(iDotPoint, iModelPoint, iDotError))
648 while iBin < dots.GetN() - 1:
649 dots.GetPoint(iBin, xValDot, yValDot)
651 iDotPoint.setVal(yValDot)
652 iModelPoint.setVal(modelCurve.interpolate(xValDot))
653 iDotError.setVal(float(dots.GetErrorYlow(iBin) + dots.GetErrorYhigh(iBin)))
655 if np.isnan(iResidual.getVal())
is not True:
657 residualValue = iResidual.getVal()
660 if abs(residualValue) > 3:
667 pointsHist.addBinWithXYError(xValDot, residualValue,
668 dots.GetErrorXlow(iBin),
669 dots.GetErrorXhigh(iBin),
674 pointsHist.SetMarkerStyle(dots.GetMarkerStyle())
675 pointsHist.SetMarkerSize(dots.GetMarkerSize())
677 rooFitFrameRes = rooRealVar.frame()
681 rooFitFrameRes.SetTitle(
"")
682 rooFitFrameRes.GetXaxis().SetTitle(rooRealVar.GetTitle() +
" " + units)
683 rooFitFrameRes.GetXaxis().SetTitleSize(0.2)
684 rooFitFrameRes.GetXaxis().SetTitleOffset(0.9)
686 rooFitFrameRes.GetXaxis().SetTickSize(0.07)
687 rooFitFrameRes.GetYaxis().SetTickSize(0.024)
689 rooFitFrameRes.GetYaxis().SetTitle(
"#splitline{Normalized}{ Residuals}")
690 rooFitFrameRes.GetYaxis().SetTitleSize(0.16)
691 rooFitFrameRes.GetYaxis().SetTitleOffset(0.4)
693 rooFitFrameRes.GetXaxis().SetLabelSize(0.125)
694 rooFitFrameRes.GetYaxis().SetLabelSize(0.120)
696 rooFitFrameRes.SetAxisRange(-5, 5,
"Y")
697 rooFitFrameRes.GetYaxis().SetNdivisions(10)
698 rooFitFrameRes.GetYaxis().ChangeLabel(1, -1, 0.)
699 rooFitFrameRes.GetYaxis().ChangeLabel(3, -1, 0.)
700 rooFitFrameRes.GetYaxis().ChangeLabel(5, -1, 0.)
701 rooFitFrameRes.GetYaxis().ChangeLabel(7, -1, 0.)
702 rooFitFrameRes.GetYaxis().ChangeLabel(9, -1, 0.)
703 rooFitFrameRes.GetYaxis().ChangeLabel(11, -1, 0.)
705 gLine1 = ROOT.TLine(5.2, 3, 5.295, 3)
706 gLine2 = ROOT.TLine(5.2, 1, 5.295, 1)
707 gLine3 = ROOT.TLine(5.2, -1, 5.295, -1)
708 gLine4 = ROOT.TLine(5.2, -3, 5.295, -3)
709 gLine1.SetLineColor(ROOT.kGreen + 3)
710 gLine2.SetLineColor(ROOT.kBlue)
711 gLine3.SetLineColor(ROOT.kBlue)
712 gLine4.SetLineColor(ROOT.kGreen + 3)
714 gLine1.SetLineWidth(2)
715 gLine2.SetLineWidth(2)
716 gLine3.SetLineWidth(2)
717 gLine4.SetLineWidth(2)
721 c1 = ROOT.TCanvas(
"c1",
"c1", 700, 640)
722 c1.SetBottomMargin(0)
725 Pad1 = ROOT.TPad(
"p1",
"p1", 0, 0.277, 1, 1, 0)
726 Pad2 = ROOT.TPad(
"p2",
"p2", 0, 0, 1, 0.276, 0)
730 Pad1.SetLeftMargin(0.15)
731 Pad1.SetBottomMargin(0.02)
732 Pad1.SetTopMargin(0.06)
734 Pad2.SetLeftMargin(0.15)
735 Pad2.SetBottomMargin(0.4)
736 Pad2.SetTopMargin(0.01)
739 rooFitFrameRes.Draw()
744 pointsHist.Draw(
"P SAME")
754 nameOfPlot = nameOfPlot +
"WithResiduals.pdf"
755 c1.SaveAs(nameOfPlot)
759 for treename
in treenames:
763 tdat = ROOT.TChain(treename)
764 tMC = ROOT.TChain(treename)
766 for iFile
in workingFilesData:
769 for iFile
in workingFilesMC:
772 histo_notTaggedEventsMC = ROOT.TH1F(
'notTaggedEventsMC',
773 'Histogram for not tagged events',
776 histo_notTaggedEventsData = ROOT.TH1F(
'notTaggedEventsData',
777 'Histogram for not tagged events',
781 'FBDT_qrCombined>>notTaggedEventsMC',
782 'abs(qrMC) == 0 && isSignal == 1 && FBDT_qrCombined < -1')
785 'FBDT_qrCombined>>notTaggedEventsData',
786 'FBDT_qrCombined < -1')
788 mc_total_notTagged = histo_notTaggedEventsMC.GetEntries()
789 data_total_notTagged = histo_notTaggedEventsData.GetEntries()
791 B0_mbc = ROOT.RooRealVar(
"Mbc",
"#it{m}_{bc}", 0, 5.2, 5.295,
"GeV/#it{c}^{2}")
792 B0_deltaE = ROOT.RooRealVar(
"deltaE",
"#Delta#it{E}", 0, -0.15, 0.15,
"GeV/#it{c}^{2}")
794 B0_FBDT_qrCombined = ROOT.RooRealVar(
"FBDT_qrCombined",
"FBDT_qrCombined", 0, -1.01, 1.01)
795 B0_FANN_qrCombined = ROOT.RooRealVar(
"FANN_qrCombined",
"FANN_qrCombined", 0, -1.01, 1.01)
797 B0_qrMC = ROOT.RooRealVar(
"qrMC",
"qrMC", 0., -100, 100)
799 argSet = ROOT.RooArgSet(B0_mbc, B0_deltaE)
800 argSet.add(B0_FBDT_qrCombined)
801 argSet.add(B0_FANN_qrCombined)
804 cutData =
"abs(FBDT_qrCombined) < 1.01"
805 cutMC =
"abs(qrMC) ==1 && abs(FBDT_qrCombined) < 1.01"
807 data = ROOT.RooDataSet(
814 mc = ROOT.RooDataSet(
821 fitMC = ROOT.RooDataSet(
"fitMC",
"fitMC", ROOT.RooArgSet(B0_mbc, B0_FBDT_qrCombined, B0_FANN_qrCombined, B0_qrMC),
"")
822 fitData = ROOT.RooDataSet(
"fitData",
"fitData", ROOT.RooArgSet(B0_mbc, B0_FBDT_qrCombined, B0_FANN_qrCombined),
"")
827 data_total_Tagged = fitData.numEntries()
828 data_total_notTagged = histo_notTaggedEventsData.GetEntries()
830 maxBo_mbc = tdat.GetMaximum(
"Mbc")
831 print(
"maximum B0_mbc = ", maxBo_mbc)
833 B0_mbc.setRange(
"fitRange", 5.2, maxBo_mbc)
837 mbcMuGauss1 = ROOT.RooRealVar(
"mbcMuGauss1",
"mbcMuGauss1", 5.27943e+00, 5.27, 5.2893,
"GeV")
838 mbcSigmaGauss1 = ROOT.RooRealVar(
"mbcSigmaGauss1",
"mbcSigmaGauss1", 2.62331e-03, 0, 0.1,
"GeV")
839 mbcGauss1 = ROOT.RooGaussModel(
"mbcGauss1",
"mbcGauss1", B0_mbc, mbcMuGauss1, mbcSigmaGauss1)
841 mbcMu2Shift = ROOT.RooRealVar(
"mbcMu2Shift",
"mbcMu2Shift", 7.62270e-03, -0.01, 0.01)
842 mbcSig2Factor = ROOT.RooRealVar(
"mbcSig2Factor",
"mbcSig2Factor", 2.83547e-01, 0, 1)
843 mbcMuGauss2 = ROOT.RooFormulaVar(
"mbcMuGauss2",
"mbcMuGauss2",
"@0-@1", ROOT.RooArgList(mbcMuGauss1, mbcMu2Shift))
844 mbcSigmaGauss2 = ROOT.RooFormulaVar(
"mbcSigmaGauss2",
"mbcSigmaGauss2",
"@0/@1", ROOT.RooArgList(mbcSigmaGauss1, mbcSig2Factor))
845 mbcGauss2 = ROOT.RooGaussModel(
"mbcGauss2",
"mbcGauss2", B0_mbc, mbcMuGauss2, mbcSigmaGauss2)
847 mbcFracGauss2 = ROOT.RooRealVar(
"mbcFracGauss2",
"mbcFracGauss2", 9.93351e-01, 0.0, 1.)
849 mbcSignalModel = ROOT.RooAddModel(
850 "mbcSignalModel",
"mbcSignalModel", ROOT.RooArgList(
851 mbcGauss1, mbcGauss2), ROOT.RooArgList(mbcFracGauss2))
853 mcFit = mbcSignalModel.fitTo(
855 ROOT.RooFit.Minos(ROOT.kFALSE), ROOT.RooFit.Extended(ROOT.kFALSE),
856 ROOT.RooFit.NumCPU(8), ROOT.RooFit.Save())
858 mbcMu2Shift.setConstant()
859 mbcSig2Factor.setConstant()
860 mbcFracGauss2.setConstant()
862 mbcSignalArgset1 = ROOT.RooArgSet(mbcGauss1)
863 mbcSignalArgset2 = ROOT.RooArgSet(mbcGauss2)
865 mbcLogFrameMC = B0_mbc.frame()
866 fitMC.plotOn(mbcLogFrameMC, ROOT.RooFit.LineWidth(4))
867 mbcSignalModel.plotOn(mbcLogFrameMC)
868 mbcSignalModel.plotOn(
870 ROOT.RooFit.Components(mbcSignalArgset1),
871 ROOT.RooFit.LineColor(ROOT.kRed + 2), ROOT.RooFit.LineWidth(4))
872 mbcSignalModel.plotOn(
874 ROOT.RooFit.Components(mbcSignalArgset2),
875 ROOT.RooFit.LineColor(ROOT.kGreen + 3), ROOT.RooFit.LineWidth(4))
877 mbcLogFrameMC.SetTitle(
"")
878 mbcXtitle = B0_mbc.GetTitle() +
" [GeV/#it{c}^{2}]"
879 mbcLogFrameMC.GetXaxis().SetTitle(mbcXtitle)
880 mbcLogFrameMC.GetXaxis().SetTitleSize(0.05)
881 mbcLogFrameMC.GetXaxis().SetLabelSize(0.045)
882 mbcLogFrameMC.GetYaxis().SetTitleSize(0.05)
883 mbcLogFrameMC.GetYaxis().SetTitleOffset(1.5)
884 mbcLogFrameMC.GetYaxis().SetLabelSize(0.045)
886 c1 = ROOT.TCanvas(
"c1",
"c1", 1400, 1100)
888 Pad = ROOT.TPad(
"p1",
"p1", 0, 0, 1, 1, 0, 0, 0)
889 Pad.SetLeftMargin(0.15)
890 Pad.SetBottomMargin(0.15)
896 nPlot = PATH +
"/B0_mbcMC.pdf"
902 mbcSignYield = ROOT.RooRealVar(
"mbcSignYield",
"mbcSignYield", 1.53208e+03, 0.0, 20000)
903 mbcContYield = ROOT.RooRealVar(
"mbcContYield",
"mbcContYield", 1.13456e+03, 0.0, 20000)
905 mbcArgusMax = ROOT.RooRealVar(
"mbcArgusMax",
"mbcArgusMax", maxBo_mbc,
"GeV")
906 mbcArgusC = ROOT.RooRealVar(
"mbcArgusC",
"mbcArgusC", -4.65996e+01, -100, 0)
907 mbcArgusPow = ROOT.RooRealVar(
"mbcArgusPow",
"mbcArgusPow", 5.13442e-01, -2, 2,
"GeV")
908 mbcArgus = ROOT.RooArgusBG(
"mbcArgus",
"mbcArgus", B0_mbc, mbcArgusMax, mbcArgusC, mbcArgusPow)
910 mbcModel = ROOT.RooAddPdf(
911 "mbcModel",
"mbcModel", ROOT.RooArgList(
912 mbcSignalModel, mbcArgus), ROOT.RooArgList(
913 mbcSignYield, mbcContYield))
915 dataFit = mbcModel.fitTo(
917 ROOT.RooFit.Minos(ROOT.kFALSE), ROOT.RooFit.Extended(ROOT.kFALSE),
918 ROOT.RooFit.NumCPU(8), ROOT.RooFit.Save())
920 mbcMuGauss1.setConstant()
921 mbcSigmaGauss1.setConstant()
922 mbcArgusC.setConstant()
923 mbcArgusPow.setConstant()
927 sData = ROOT.RooStats.SPlot(
"sData",
"SPlot using B0_mbc", fitData, mbcModel, ROOT.RooArgList(mbcSignYield, mbcContYield))
932 print(
"Check SWeights:")
934 print(
"Yield of mbcSignYield is ", mbcSignYield.getVal(),
". From sWeights it is ",
935 sData.GetYieldFromSWeight(
"mbcSignYield"))
937 print(
"Yield of mbcContYield is ", mbcContYield.getVal(),
". From sWeights it is ",
938 sData.GetYieldFromSWeight(
"mbcContYield"))
945 dataFitExtended = mbcModel.fitTo(
946 fitData, ROOT.RooFit.Extended(),
947 ROOT.RooFit.Minos(ROOT.kFALSE), ROOT.RooFit.Extended(ROOT.kFALSE),
948 ROOT.RooFit.NumCPU(8), ROOT.RooFit.Save())
952 mbcArgset1 = ROOT.RooArgSet(mbcArgus)
953 mbcArgset2 = ROOT.RooArgSet(mbcSignalModel)
956 mbcFrame = B0_mbc.frame()
957 fitData.plotOn(mbcFrame)
958 mbcModel.plotOn(mbcFrame, ROOT.RooFit.LineWidth(4))
961 ROOT.RooFit.Components(mbcArgset1),
962 ROOT.RooFit.LineColor(
964 ROOT.RooFit.LineStyle(5),
965 ROOT.RooFit.LineWidth(4))
968 ROOT.RooFit.Components(mbcArgset2),
969 ROOT.RooFit.LineColor(
971 ROOT.RooFit.LineStyle(7),
972 ROOT.RooFit.LineWidth(4))
974 mbcFrame.SetTitle(
"")
975 mbcFrame.GetXaxis().SetTitle(mbcXtitle)
976 mbcFrame.SetTitleOffset(0.9,
"X")
977 mbcFrame.GetXaxis().SetTitleSize(0.075)
978 mbcFrame.GetXaxis().SetLabelSize(0.045)
979 mbcFrame.GetYaxis().SetTitleSize(0.06)
980 mbcFrame.GetYaxis().SetTitleOffset(0.95)
981 mbcFrame.GetYaxis().SetLabelSize(0.045)
983 c1 = ROOT.TCanvas(
"c1",
"c1", 10, 15, 700, 500)
984 c1.SetLeftMargin(0.15)
985 c1.SetRightMargin(0.04)
986 c1.SetTopMargin(0.03)
987 c1.SetBottomMargin(0.15)
996 nPlot = PATH +
"/B0_mbc.pdf"
1001 mbcFrameMC = B0_mbc.frame()
1002 mbcFrameMC.SetTitle(
"")
1004 fitMC.plotOn(mbcFrameMC)
1005 mbcSignalModel.plotOn(mbcFrameMC, ROOT.RooFit.LineColor(ROOT.kGreen + 3), ROOT.RooFit.LineStyle(7), ROOT.RooFit.LineWidth(4))
1007 dotsMC = mbcFrameMC.getHist(
"h_fitMC")
1008 modelCurveMC = mbcFrameMC.getCurve(
"mbcSignalModel_Norm[Mbc]")
1010 dotsMC.SetFillColor(ROOT.kWhite)
1011 modelCurveMC.SetFillColor(ROOT.kWhite)
1015 lMC = ROOT.TLegend(0.25, 0.63, 0.5, 0.89)
1016 lMC.AddEntry(dotsMC,
"Belle MC")
1017 lMC.AddEntry(modelCurveMC,
'Signal')
1019 lMC.SetTextSize(0.055)
1021 plotWithResiduals(mbcFrameMC, B0_mbc, dotsMC, modelCurveMC,
"[GeV/#it{c}^{2}]",
"mbcMC", lMC,
True)
1025 dots = mbcFrame.getHist(
"h_fitData")
1026 modelCurve = mbcFrame.getCurve(
"mbcModel_Norm[Mbc]")
1027 signalCurve = mbcFrame.getCurve(
"mbcModel_Norm[Mbc]_Comp[mbcSignalModel]")
1028 continuumCurve = mbcFrame.getCurve(
"mbcModel_Norm[Mbc]_Comp[mbcArgus]")
1030 dots.SetFillColor(ROOT.kWhite)
1031 modelCurve.SetFillColor(ROOT.kWhite)
1032 signalCurve.SetFillColor(ROOT.kWhite)
1033 continuumCurve.SetFillColor(ROOT.kWhite)
1035 l = ROOT.TLegend(0.25, 0.43, 0.5, 0.89)
1036 l.AddEntry(dots,
"Belle data")
1037 l.AddEntry(modelCurve,
'Fit result')
1038 l.AddEntry(signalCurve,
"Signal")
1039 l.AddEntry(continuumCurve,
'Continuum')
1041 l.SetTextSize(0.055)
1043 plotWithResiduals(mbcFrame, B0_mbc, dots, modelCurve,
"[GeV/#it{c}^{2}]",
"Mbc", l,
False)
1045 signalData = ROOT.RooDataSet(
"signalData",
"signalData", fitData, fitData.get(),
"",
"mbcSignYield_sw")
1047 print(
"Variable Names")
1049 qrSigFrame = B0_FBDT_qrCombined.frame()
1050 signalData.plotOn(qrSigFrame, ROOT.RooFit.DataError(ROOT.RooAbsData.SumW2))
1052 qrSigFrame.SetTitle(
"")
1054 qrSigFrame.GetXaxis().SetTitle(mbcXtitle)
1055 qrSigFrame.GetXaxis().SetTitleSize(0.05)
1056 qrSigFrame.GetXaxis().SetLabelSize(0.045)
1057 qrSigFrame.GetYaxis().SetTitleSize(0.05)
1058 qrSigFrame.GetYaxis().SetTitleOffset(1.5)
1059 qrSigFrame.GetYaxis().SetLabelSize(0.045)
1061 c1 = ROOT.TCanvas(
"c1",
"c1", 1400, 1100)
1063 Pad = ROOT.TPad(
"p1",
"p1", 0, 0, 1, 1, 0, 0, 0)
1064 Pad.SetLeftMargin(0.15)
1065 Pad.SetBottomMargin(0.15)
1071 nPlot = PATH +
"/B0_qrSignal.pdf"
1075 continuumData = ROOT.RooDataSet(
"continuumData",
"continuumData", fitData, fitData.get(),
"",
"mbcContYield_sw")
1077 qrContFrame = B0_FBDT_qrCombined.frame()
1078 continuumData.plotOn(qrContFrame, ROOT.RooFit.DataError(ROOT.RooAbsData.SumW2))
1080 qrContFrame.SetTitle(
"")
1082 qrContFrame.GetXaxis().SetTitle(mbcXtitle)
1083 qrContFrame.GetXaxis().SetTitleSize(0.05)
1084 qrContFrame.GetXaxis().SetLabelSize(0.045)
1085 qrContFrame.GetYaxis().SetTitleSize(0.05)
1086 qrContFrame.GetYaxis().SetTitleOffset(1.5)
1087 qrContFrame.GetYaxis().SetLabelSize(0.045)
1089 c1 = ROOT.TCanvas(
"c1",
"c1", 1400, 1100)
1091 Pad = ROOT.TPad(
"p1",
"p1", 0, 0, 1, 1, 0, 0, 0)
1092 Pad.SetLeftMargin(0.15)
1093 Pad.SetBottomMargin(0.15)
1098 nPlot = PATH +
"/B0_qrContinuum.pdf"
1106 sPlotResultsFBDT = np.zeros((fitData.numEntries(), 3))
1107 sPlotResultsFANN = np.zeros((fitData.numEntries(), 3))
1108 for i
in range(fitData.numEntries()):
1109 row = fitData.get(i)
1110 sPlotResultsFBDT[i][0] = row.getRealValue(
'FBDT_qrCombined', 0, ROOT.kTRUE)
1111 sPlotResultsFBDT[i][1] = row.getRealValue(
'mbcSignYield_sw', 0, ROOT.kTRUE)
1112 sPlotResultsFBDT[i][2] = row.getRealValue(
'mbcContYield_sw', 0, ROOT.kTRUE)
1113 sPlotResultsFANN[i][0] = row.getRealValue(
'FANN_qrCombined', 0, ROOT.kTRUE)
1114 sPlotResultsFANN[i][1] = row.getRealValue(
'mbcSignYield_sw', 0, ROOT.kTRUE)
1115 sPlotResultsFANN[i][2] = row.getRealValue(
'mbcContYield_sw', 0, ROOT.kTRUE)
1117 mcSPlotFBDT = np.zeros((fitMC.numEntries(), 2))
1118 mcSPlotFANN = np.zeros((fitMC.numEntries(), 2))
1119 for i
in range(fitMC.numEntries()):
1121 mcSPlotFBDT[i][0] = row.getRealValue(
'FBDT_qrCombined', 0, ROOT.kTRUE)
1122 mcSPlotFBDT[i][1] = row.getRealValue(
'qrMC', 0, ROOT.kTRUE)
1123 mcSPlotFANN[i][0] = row.getRealValue(
'FANN_qrCombined', 0, ROOT.kTRUE)
1124 mcSPlotFANN[i][1] = row.getRealValue(
'qrMC', 0, ROOT.kTRUE)
1126 efficiencyCalculator(mcSPlotFBDT, sPlotResultsFBDT,
"FBDT")
1127 efficiencyCalculator(mcSPlotFANN, sPlotResultsFANN,
"FANN")