diff options
Diffstat (limited to 'utils/plot-energy')
| -rwxr-xr-x | utils/plot-energy | 76 |
1 files changed, 23 insertions, 53 deletions
diff --git a/utils/plot-energy b/utils/plot-energy index 5c33969..f082b8c 100755 --- a/utils/plot-energy +++ b/utils/plot-energy @@ -56,7 +56,9 @@ def plot_hist2(df, muons=False): plt.hist(np.log10(df_id.ke.values/1000), bins=np.linspace(0,4.5,100), histtype='step') plt.xlabel("log10(Energy (GeV))") else: - plt.hist(df_id.ke.values, bins=np.linspace(20,10e3,100), histtype='step') + bins = np.logspace(np.log10(20),np.log10(10e3),21) + plt.hist(df_id.ke.values, bins=bins, histtype='step') + plt.gca().set_xscale("log") plt.xlabel("Energy (MeV)") plt.title('$' + ''.join([particle_id[int(''.join(x))] for x in grouper(str(id),2)]) + '$') @@ -103,6 +105,7 @@ if __name__ == '__main__': import h5py from sddm.plot_energy import * from sddm.plot import despine + from sddm import setup_matplotlib parser = argparse.ArgumentParser("plot fit results") parser.add_argument("filenames", nargs='+', help="input files") @@ -110,46 +113,9 @@ if __name__ == '__main__': parser.add_argument("--save", action='store_true', default=False, help="save corner plots for backgrounds") args = parser.parse_args() - if args.save: - # default \textwidth for a fullpage article in Latex is 16.50764 cm. - # You can figure this out by compiling the following TeX document: - # - # \documentclass{article} - # \usepackage{fullpage} - # \usepackage{layouts} - # \begin{document} - # textwidth in cm: \printinunitsof{cm}\printlen{\textwidth} - # \end{document} - - width = 16.50764 - width /= 2.54 # cm -> inches - # According to this page: - # http://www-personal.umich.edu/~jpboyd/eng403_chap2_tuftegospel.pdf, - # Tufte suggests an aspect ratio of 1.5 - 1.6. - height = width/1.5 - FIGSIZE = (width,height) - - import matplotlib.pyplot as plt - - font = {'family':'serif', 'serif': ['computer modern roman']} - plt.rc('font',**font) - - plt.rc('text', usetex=True) - else: - # on retina screens, the default plots are way too small - # by using Qt5 and setting QT_AUTO_SCREEN_SCALE_FACTOR=1 - # Qt5 will scale everything using the dpi in ~/.Xresources - import matplotlib - matplotlib.use("Qt5Agg") - - import matplotlib.pyplot as plt - - # Default figure size. Currently set to my monitor width and height so that - # things are properly formatted - FIGSIZE = (13.78,7.48) + setup_matplotlib(args.save) - # Make the defalt font bigger - plt.rc('font', size=22) + import matplotlib.pyplot as plt ev = pd.concat([pd.read_hdf(filename, "ev") for filename in args.filenames],ignore_index=True) fits = pd.concat([pd.read_hdf(filename, "fits") for filename in args.filenames],ignore_index=True) @@ -339,11 +305,11 @@ if __name__ == '__main__': plot_corner_plot(noise,"Noise") plot_corner_plot(signal,"Signal") - fig = plt.figure(figsize=FIGSIZE) + fig = plt.figure() plot_hist2(flashers) despine(fig,trim=True) plt.suptitle("Flashers") - fig = plt.figure(figsize=FIGSIZE) + fig = plt.figure() plot_hist2(muon,muons=True) despine(fig,trim=True) plt.suptitle("Muons") @@ -460,15 +426,19 @@ if __name__ == '__main__': muon_best_fit = muons.sort_values('fmin').groupby(['run','gtid']).nth(0) muons = muons[muons.id == 22].sort_values('fmin').groupby(['run','gtid'],as_index=False).nth(0).reset_index(level=0,drop=True) - # require r < 6 meters + # require (r/r_psup)^3 < 0.9 prompt = prompt[prompt.r_psup < 0.9] atm = atm[atm.r_psup < 0.9] print("number of events after radius cut = %i" % len(prompt)) - # Note: Need to design and apply a psi based cut here + # require psi < 6 + prompt = prompt[prompt.psi < 6] + atm = atm[atm.psi < 6] + + print("number of events after psi cut = %i" % len(prompt)) - fig = plt.figure(figsize=FIGSIZE) + fig = plt.figure() plot_hist2(prompt) despine(fig,trim=True) if args.save: @@ -476,7 +446,7 @@ if __name__ == '__main__': plt.savefig("prompt.eps") else: plt.suptitle("Without Neutron Follower") - fig = plt.figure(figsize=FIGSIZE) + fig = plt.figure() plot_hist2(atm) despine(fig,trim=True) if args.save: @@ -484,7 +454,7 @@ if __name__ == '__main__': plt.savefig("atm.eps") else: plt.suptitle("With Neutron Follower") - fig = plt.figure(figsize=FIGSIZE) + fig = plt.figure() plot_hist2(michel_best_fit) despine(fig,trim=True) if args.save: @@ -492,7 +462,7 @@ if __name__ == '__main__': plt.savefig("michel_electrons.eps") else: plt.suptitle("Michel Electrons") - fig = plt.figure(figsize=FIGSIZE) + fig = plt.figure() plot_hist2(muon_best_fit,muons=True) despine(fig,trim=True) if len(muon_best_fit): @@ -504,7 +474,7 @@ if __name__ == '__main__': plt.suptitle("External Muons") # Plot the energy and angular distribution for external muons - fig = plt.figure(figsize=FIGSIZE) + fig = plt.figure() plt.subplot(2,1,1) plt.hist(muons.ke.values, bins=np.logspace(3,7,100), histtype='step') plt.xlabel("Energy (MeV)") @@ -533,7 +503,7 @@ if __name__ == '__main__': stopping_muons['T_dx'] = dx_to_energy(stopping_muons.dx) stopping_muons['dT'] = stopping_muons['energy1'] - stopping_muons['T_dx'] - fig = plt.figure(figsize=FIGSIZE) + fig = plt.figure() plt.hist((stopping_muons['energy1']-stopping_muons['T_dx'])*100/stopping_muons['T_dx'], bins=np.linspace(-100,100,200), histtype='step') despine(fig,trim=True) plt.xlabel("Fractional energy difference (\%)") @@ -554,7 +524,7 @@ if __name__ == '__main__': dT = stopping_muons.groupby(pd_bins)['dT'].agg(['mean','sem','std',std_err,median,median_err,iqr_std,iqr_std_err]) - fig = plt.figure(figsize=FIGSIZE) + fig = plt.figure() plt.errorbar(T,dT['median']*100/T,yerr=dT['median_err']*100/T) despine(fig,trim=True) plt.xlabel("Kinetic Energy (MeV)") @@ -566,7 +536,7 @@ if __name__ == '__main__': else: plt.title("Stopping Muon Energy Bias") - fig = plt.figure(figsize=FIGSIZE) + fig = plt.figure() plt.errorbar(T,dT['iqr_std']*100/T,yerr=dT['iqr_std_err']*100/T) despine(fig,trim=True) plt.xlabel("Kinetic Energy (MeV)") @@ -578,7 +548,7 @@ if __name__ == '__main__': else: plt.title("Stopping Muon Energy Resolution") - fig = plt.figure(figsize=FIGSIZE) + fig = plt.figure() bins=np.linspace(0,100,100) plt.hist(michel.ke.values, bins=bins, histtype='step', label="Dark Matter Fitter") if michel.size: |