diff --git a/pes_to_spec/model.py b/pes_to_spec/model.py
index 32430745e3ad060335abe81771c96b34c2da8dac..e198b82acc550d38e2c0ee4351aaec45133934fd 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 a199e8b663e5a0d9c813c9c967afbdd5729da4c6..3911d0f630ed8f4239d528be1ff1bc288ea7ceae 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 52c6caea532d442d801645206df116635d714086..d2a128265c895c6d9ac93e49ad067f8a0b9dc8bc 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 75ba648acf981744ba1d4e7033d07e6cde185cb5..3f6ecf64f0b7199a1b69bc22b0293cb649e016a7 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')