Added further resolution analysis.

pes_to_spec/test/offline_analysis.py

@@ -246,23 +246,8 @@ def main():
     for k in pes_raw.keys():
         pes_raw[k] = pes_raw[k][train_idx]
 
-    # fit model with no smearing to estimate maximum resolution
-    print("Fitting model only to detect resolution")
-    start = time_ns()
-    model = Model(channels=channels, model_type=args.model_type)
-    model.uniformize(xgm_flux[train_idx])
-    model.fit(pes_raw,
-              spec_raw_int[train_idx],
-              spec_raw_pe[train_idx],
-              pulse_energy=xgm_flux[train_idx],
-              )
-    t += [time_ns() - start]
-    t_names += ["Resolution estimate"]
-    resolution = model.resolution
-    print(f"Resolution detected: {resolution} eV")
-
     # do it again with smearing, but now with the knowledge of the resolution
-    model = Model(channels=channels, model_type=args.model_type, high_res_fwhm=resolution)
+    model = Model(channels=channels, model_type=args.model_type)
 
     model.debug_peak_finding(pes_raw, os.path.join(args.directory, "test_peak_finding.pdf"))
     if len(args.model) == 0:
@@ -352,8 +337,8 @@ def main():
     showSpec = False
     if len(args.model) == 0:
         showSpec = True
-        # always smooth by 0.2, even if no smoothing was requested, for visualization
-        smoother = HighResolutionSmoother(0.2)
+        # always smoothen by the resolution, even if no smoothing was requested, for visualization
+        smoother = HighResolutionSmoother(model.resolution)
         spec_smooth = smoother.fit_transform(spec_raw_int_t, energy=spec_raw_pe_t)
         # use whatever smoothing was used in the model for the calculations
         spec_data = model.preprocess_high_res(spec_raw_int_t)
@@ -413,6 +398,13 @@ def main():
         df = pd.DataFrame(q)
         df.to_csv(os.path.join(args.directory, "quality.csv"))
 
+    df = pd.DataFrame(dict(auto_corr_hr=model.auto_corr_hr,
+                           auto_corr_virt=model.auto_corr_virt,
+                           energy=model.wiener_energy
+                          )
+                     )
+    df.to_csv(os.path.join(args.directory, "resolution.csv"))
+
     first, last = model.get_low_resolution_range()
     # plot
     high_int_idx = np.argsort(xgm_flux_t[:,0])

pes_to_spec/test/prepare_plots.py

@@ -50,6 +50,23 @@ def plot_final(df: pd.DataFrame, smooth: bool, filename: str):
     fig.savefig(filename)
     plt.close(fig)
 
+def plot_autocov(df: pd.DataFrame, filename: str):
+    fig = plt.figure(figsize=(12, 8))
+    gs = GridSpec(1, 1)
+    ax = fig.add_subplot(gs[0, 0])
+    ax.plot(df.energy, df.auto_corr_hr, c='b', lw=3, label="Grating spectrometer")
+    ax.plot(df.energy, df.auto_corr_virt, c='r', lw=3, label="Virtual spectrometer")
+    ax.legend(frameon=False, borderaxespad=0, loc='upper left')
+    ax.spines['top'].set_visible(False)
+    ax.spines['right'].set_visible(False)
+    ax.set(
+           xlabel="Photon energy [eV]",
+           ylabel="Auto-covariance [a.u.]",
+           ylim=(0, 1.0))
+    plt.tight_layout()
+    fig.savefig(filename)
+    plt.close(fig)
+
 def plot_chi2(df: pd.DataFrame, filename: str):
     fig = plt.figure(figsize=(12, 8))
     gs = GridSpec(1, 1)
@@ -529,3 +546,5 @@ if __name__ == '__main__':
     pca_variance_plot(pd.read_csv(f'{indir}/pca_spec.csv'), f'pca_spec.pdf', max_comp_frac=0.99)
     pca_variance_plot(pd.read_csv(f'{indir}/pca_pes.csv'), f'pca_pes.pdf', max_comp_frac=0.95)
 
+    plot_autocov(pd.read_csv(f'{indir}/autocov.csv'), f'autocov.pdf')
+