added figure for fit results

examples/findPeaks.py

@@ -41,19 +41,19 @@ def gauss_plus_bkg(x, Ns=1000, mu=1000, sig=50., Nb=100., s=0., mn=0, mx=1.) :
 parser = argparse.ArgumentParser(description=__doc__)
 parser.add_argument('filename', nargs='?', default='Co60.hst')
 parser.add_argument('-p','--prominence',
-                    help="minimum peak height over baseline for scipy.find_peaks()",
+                    help="minimum peak height over baseline for scipy.find_peaks() (50)",
                     type=float, default=50)
 parser.add_argument('-w','--width',
-                    help="minimum peak width for scipy.find_peaks()",
+                    help="minimum peak width for scipy.find_peaks() (20)",
                     type=float, default=20)
 parser.add_argument('-r','--rel_height',
-                    help=" height at which full width is calculated in scipy.find_peaks()",
+                    help=" height at which full width is calculated in scipy.find_peaks() (0.5)",
                     type=float, default=0.5)
 parser.add_argument('-f','--fit_range_factor',
-                    help = "fit_range_factor: fit range = fit_range_factor * fwhm",
-                    type=float, default=1.5)
+                    help = "fit_range_factor: fit range = fit_range_factor * fwhm (1.75)",
+                    type=float, default=1.75)
 parser.add_argument('-t','--threshold',
-                    help = "threshold, i.e. minimum valid channel nmber",
+                    help = "threshold, i.e. minimum valid channel nmber (51)",
                     type=float, default=51)
 parser.add_argument('-n','--noplot', action='store_true')
 parser.add_argument('-v','--verbose', action='store_true')
@@ -144,9 +144,33 @@ for i, p in enumerate(peaks):
     plot_ranges[i][0] = p-wid
     plot_ranges[i][1] = p+wid+1
 
+# plot results
+sprctrum_fig = plt.figure('Spectrum', figsize=(10, 5))
+xhst = np.linspace(0, len(hst), len(hst), endpoint=False) +0.5
+plt.errorbar(xhst, hst, yerr=np.sqrt(hst),
+                 fmt='.', color='grey', alpha=0.1)
+for i, fit in enumerate(fit_results):
+    xplt = np.linspace(plot_ranges[i][0],plot_ranges[i][1],
+                    10*int((plot_ranges[i][1]-plot_ranges[i][0])))
+    plt.plot(xplt, gauss_plus_bkg(xplt, *fit.parameter_values),
+             linestyle = 'solid', linewidth=3)
+    mu = fit.parameter_values[1]
+    sig = fit.parameter_values[2]
+    mx = gauss_plus_bkg(mu, *fit.parameter_values)
+    h = fit.parameter_values[0]/np.sqrt(2*np.pi)/sig 
+    plt.vlines(mu, mx - h, mx,
+               linewidth=3)
+    plt.vlines(mu, 0, mx - h, 
+               linewidth=2, linestyle='dashed')
+
+plt.show()
+
+
+# kafe2 plots
 if plot:
     fit_results = np.append(fit_results, hdata)
-    plot_ranges[len(peaks+1)][1] = 4000
+    plot_ranges[len(peaks+1)][1] = len(hst)
     kafe2.plot(fits = fit_results, x_range = list(map(tuple, plot_ranges)),
            x_label = "Channel number", y_label = "Counts",
            font_scale = 0.5, save=False)
+