From c1523ca023c0a304c6795e961c020540a8bf617e Mon Sep 17 00:00:00 2001
From: Danilo Ferreira de Lima <danilo.enoque.ferreira.de.lima@xfel.de>
Date: Mon, 28 Aug 2023 18:51:33 +0200
Subject: [PATCH] Updated plots.

---
 pes_to_spec/model.py                         |  1 +
 pes_to_spec/test/offline_analysis.py         |  8 +-
 pes_to_spec/test/plot_channel_sensitivity.py | 15 ++--
 pes_to_spec/test/prepare_plots.py            | 92 +++++++++++++++++---
 4 files changed, 95 insertions(+), 21 deletions(-)

diff --git a/pes_to_spec/model.py b/pes_to_spec/model.py
index 3243074..e198b82 100644
--- a/pes_to_spec/model.py
+++ b/pes_to_spec/model.py
@@ -989,6 +989,7 @@ class Model(TransformerMixin, BaseEstimator):
                     unc=unc.reshape((B, P, -1)),
                     pca=pca_unc,
                     total_unc=total_unc.reshape((B, P, -1)),
+                    nopca_unc=unc.reshape((B, P, -1)),
                     expected_pca=high_pca.reshape((B, P, -1)),
                     expected_pca_unc=high_pca_unc.reshape((B, P, -1)),
                     )
diff --git a/pes_to_spec/test/offline_analysis.py b/pes_to_spec/test/offline_analysis.py
index a199e8b..3911d0f 100755
--- a/pes_to_spec/test/offline_analysis.py
+++ b/pes_to_spec/test/offline_analysis.py
@@ -364,12 +364,18 @@ def main():
 
         # rmse
         rmse = np.sqrt(np.mean((spec_smooth[:, np.newaxis, :] - spec_pred["expected"])**2, axis=(-1, -2)))
+        nopca_unc = np.sqrt(np.mean(spec_pred["nopca_unc"]**2, axis=(-1, -2)))
+        total_unc = np.sqrt(np.mean(spec_pred["total_unc"]**2, axis=(-1, -2)))
+        median_unc = np.median(spec_pred["total_unc"], axis=(-1, -2))
 
         q = dict(chi2_prepca=chi2_prepca,
                  ndof=spec_smooth_pca.shape[-1]*np.ones_like(chi2_prepca),
                  xgm_flux_t=xgm_flux_t[:,0],
                  rmse=rmse,
-                 root_mean_squared_pca_unc=np.sqrt((spec_pred["expected_pca_unc"][:, 0, :]**2).sum(axis=-1))
+                 nopca_unc=nopca_unc,
+                 total_unc=total_unc,
+                 median_unc=median_unc,
+                 root_mean_squared_pca_unc=np.sqrt((spec_pred["expected_pca_unc"][:, 0, :]**2).mean(axis=-1))
                  )
         q.update({f'res_prepca_{k}': res_prepca[:, k]
                   for k in range(res_prepca.shape[1])
diff --git a/pes_to_spec/test/plot_channel_sensitivity.py b/pes_to_spec/test/plot_channel_sensitivity.py
index 52c6cae..d2a1282 100755
--- a/pes_to_spec/test/plot_channel_sensitivity.py
+++ b/pes_to_spec/test/plot_channel_sensitivity.py
@@ -44,22 +44,23 @@ def plot_unc(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.number_channels, df.unc, c='tab:blue', lw=3)
+    ax.plot(df.number_channels, 2*df.unc, c='tab:blue', lw=3, alpha=0.5)
     ax.spines['top'].set_visible(False)
     ax.set(
            xlabel="Number of channels",
-           ylabel="Uncertainty [a.u.]",
+           ylabel="Average uncertainty [a.u.]",
            )
-    maxY = 1.1*df.unc.max()
-    ax.set_ylim((None, maxY))
+    #maxY = 1.1*df.unc.max()
+    #ax.set_ylim((None, maxY))
     rax = ax.twinx()
-    rax.plot(df.number_channels, df.rmse, c='tab:red', lw=3)
+    rax.plot(df.number_channels, df.rmse, c='tab:red', lw=3, alpha=0.5)
     rax.spines['right'].set_color('tab:red')
+    rax.spines['top'].set_visible(False)
     rax.tick_params(axis='y', colors='tab:red')
     rax.set_ylabel("Root-mean-squared error [a.u.]")
     rax.yaxis.label.set_color("tab:red")
-    maxY = 1.1*df.rmse.max()
-    rax.set_ylim((None, maxY))
+    #maxY = 1.1*df.rmse.max()
+    #rax.set_ylim((None, maxY))
     plt.tight_layout()
     fig.savefig(filename)
     plt.close(fig)
diff --git a/pes_to_spec/test/prepare_plots.py b/pes_to_spec/test/prepare_plots.py
index 75ba648..3f6ecf6 100755
--- a/pes_to_spec/test/prepare_plots.py
+++ b/pes_to_spec/test/prepare_plots.py
@@ -71,19 +71,6 @@ def plot_rmse(df: pd.DataFrame, filename: str):
     fig.savefig(filename)
     plt.close(fig)
 
-def plot_rmse_intensity(df: pd.DataFrame, filename: str):
-    fig = plt.figure(figsize=(12, 8))
-    gs = GridSpec(1, 1)
-    ax = fig.add_subplot(gs[0, 0])
-    ax.scatter(df.rmse, df.xgm_flux_t, c='r', s=30)
-    ax = plt.gca()
-    ax.set(title=f"",
-           xlabel=r"Root-mean-squared error",
-           ylabel="Beam intensity [uJ]",
-           )
-    fig.savefig(filename)
-    plt.close(fig)
-
 def plot_residue(df: pd.DataFrame, filename: str):
     cols = [k for k in df.columns if "res_prepca" in k]
     df_res = df.loc[:, cols]
@@ -143,6 +130,83 @@ def plot_chi2_intensity(df: pd.DataFrame, filename: str):
     fig.savefig(filename)
     plt.close(fig)
 
+def plot_rmse_intensity(df: pd.DataFrame, filename: str):
+    fig = plt.figure(figsize=(12, 8))
+    gs = GridSpec(1, 1)
+    ax = fig.add_subplot(gs[0, 0])
+    sns.kdeplot(x=df.rmse, y=df.xgm_flux_t*1e-3,
+                fill=True,
+                ax=ax)
+    sns.scatterplot(x=df.rmse, y=df.xgm_flux_t*1e-3,
+                    s=5,
+                    alpha=0.4,
+                    c="tab:red",
+                    #size=df.root_mean_squared_pca_unc,
+                    #sizes=(20, 200),
+                    ax=ax)
+    ax = plt.gca()
+    ax.set(title=f"",
+           xlabel=r"Root-mean-squared error [a.u.]",
+           ylabel="Beam intensity [mJ]",
+           ylim=(0, df.xgm_flux_t.mean()*1e-3 + 3*df.xgm_flux_t.std()*1e-3)
+           )
+    ax.spines['top'].set_visible(False)
+    ax.spines['right'].set_visible(False)
+    plt.tight_layout()
+    fig.savefig(filename)
+    plt.close(fig)
+
+def plot_unc_intensity(df: pd.DataFrame, filename: str):
+    fig = plt.figure(figsize=(12, 8))
+    gs = GridSpec(1, 1)
+    ax = fig.add_subplot(gs[0, 0])
+    sns.kdeplot(x=df.total_unc, y=df.xgm_flux_t*1e-3,
+                fill=True,
+                ax=ax)
+    sns.scatterplot(x=df.total_unc, y=df.xgm_flux_t*1e-3,
+                    s=5,
+                    alpha=0.4,
+                    c="tab:red",
+                    #size=df.root_mean_squared_pca_unc,
+                    #sizes=(20, 200),
+                    ax=ax)
+    ax = plt.gca()
+    ax.set(title=f"",
+           xlabel=r"Uncertainty [a.u.]",
+           ylabel="Beam intensity [mJ]",
+           ylim=(0, df.xgm_flux_t.mean()*1e-3 + 3*df.xgm_flux_t.std()*1e-3)
+           )
+    ax.spines['top'].set_visible(False)
+    ax.spines['right'].set_visible(False)
+    plt.tight_layout()
+    fig.savefig(filename)
+    plt.close(fig)
+
+def plot_unc_rmse(df: pd.DataFrame, filename: str):
+    fig = plt.figure(figsize=(12, 8))
+    gs = GridSpec(1, 1)
+    ax = fig.add_subplot(gs[0, 0])
+    sns.kdeplot(x=2*df.total_unc, y=df.rmse,
+                fill=True,
+                ax=ax)
+    sns.scatterplot(x=2*df.total_unc, y=df.rmse,
+                    s=5,
+                    alpha=0.4,
+                    c="tab:red",
+                    #size=df.root_mean_squared_pca_unc,
+                    #sizes=(20, 200),
+                    ax=ax)
+    ax = plt.gca()
+    ax.set(title=f"",
+           xlabel=r"Root-mean-squared unc. [a.u.]",
+           ylabel="Root-mean-squared error [a.u.]",
+           )
+    ax.spines['top'].set_visible(False)
+    ax.spines['right'].set_visible(False)
+    plt.tight_layout()
+    fig.savefig(filename)
+    plt.close(fig)
+
 def pca_variance_plot(df: pd.DataFrame, filename: str, max_comp_frac: float=0.99):
     """
     Plot variance contribution.
@@ -359,6 +423,8 @@ if __name__ == '__main__':
 
     plot_chi2(pd.read_csv(f'{indir}/quality.csv'), f'chi2_prepca.pdf')
     plot_chi2_intensity(pd.read_csv(f'{indir}/quality.csv'), f'intensity_vs_chi2_prepca.pdf')
+    plot_unc_intensity(pd.read_csv(f'{indir}/quality.csv'), f'intensity_vs_unc.pdf')
+    plot_unc_rmse(pd.read_csv(f'{indir}/quality.csv'), f'rmse_vs_unc.pdf')
     plot_rmse(pd.read_csv(f'{indir}/quality.csv'), f'rmse.pdf')
     plot_rmse_intensity(pd.read_csv(f'{indir}/quality.csv'), f'intensity_vs_rmse.pdf')
 
-- 
GitLab