diff options
Diffstat (limited to 'utils/plot-fit-results')
| -rwxr-xr-x | utils/plot-fit-results | 16 |
1 files changed, 14 insertions, 2 deletions
diff --git a/utils/plot-fit-results b/utils/plot-fit-results index c90c13a..20c3c68 100755 --- a/utils/plot-fit-results +++ b/utils/plot-fit-results @@ -28,6 +28,7 @@ if __name__ == '__main__': parser = argparse.ArgumentParser("plot fit results") parser.add_argument("filenames", nargs='+', help="input files") + parser.add_argument("-o", "--output", default=None, help="output file") parser.add_argument("--save", action="store_true", default=False, help="save plots") args = parser.parse_args() @@ -110,11 +111,15 @@ if __name__ == '__main__': # 100 bins between 50 MeV and 1 GeV bins = np.arange(50,1000,100) + T = (bins[1:] + bins[:-1])/2 + markers = itertools.cycle(('o', 'v')) fig3, ax3 = plt.subplots(3,1,num=3,sharex=True) fig4, ax4 = plt.subplots(3,1,num=4,sharex=True) + output = pd.DataFrame({'T':T}) + for id in [IDP_E_MINUS, IDP_MU_MINUS]: events = data_true[data_true['mcgn_id'] == id] @@ -128,8 +133,6 @@ if __name__ == '__main__': label = 'Muon' if id == IDP_MU_MINUS else 'Electron' - T = (bins[1:] + bins[:-1])/2 - marker = markers.next() plt.figure(1) @@ -142,6 +145,11 @@ if __name__ == '__main__': ax3[1].errorbar(T,dy['median'],yerr=dy['median_err'],fmt=marker,label=label) ax3[2].errorbar(T,dz['median'],yerr=dz['median_err'],fmt=marker,label=label) + if id == IDP_E_MINUS: + output['e_bias'] = (dT['median']/T).values + else: + output['u_bias'] = (dT['median']/T).values + ax4[0].errorbar(T,dx['iqr_std'],yerr=dx['iqr_std_err'],fmt=marker,label=label) ax4[1].errorbar(T,dy['iqr_std'],yerr=dy['iqr_std_err'],fmt=marker,label=label) ax4[2].errorbar(T,dz['iqr_std'],yerr=dz['iqr_std_err'],fmt=marker,label=label) @@ -152,6 +160,10 @@ if __name__ == '__main__': plt.figure(6) plt.scatter(events['Te'],events['ratio'],marker=marker,label=label) + if args.output: + with h5py.File(args.output,"w") as f: + f.create_dataset('energy_bias',data=output.to_records()) + fig = plt.figure(1) despine(fig,trim=True) plt.xlabel("Kinetic Energy (MeV)") |