From 030e33897d9939767ab04b50264fed0616d5eb04 Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Wed, 8 May 2024 17:21:38 +0200 Subject: [PATCH 1/7] split dt cuts out of aoe cuts --- src/pygama/pargen/AoE_cal.py | 49 +++++++++++------------------------- 1 file changed, 15 insertions(+), 34 deletions(-) diff --git a/src/pygama/pargen/AoE_cal.py b/src/pygama/pargen/AoE_cal.py index 402edf6ed..fb4e64696 100644 --- a/src/pygama/pargen/AoE_cal.py +++ b/src/pygama/pargen/AoE_cal.py @@ -1536,24 +1536,15 @@ def get_aoe_cut_fit( log.error("A/E cut determination failed") self.low_cut_val = np.nan data[output_cut_param] = False - if self.dt_cut_param is not None and self.dt_cut_hard is True: - self.update_cal_dicts( - { - output_cut_param: { - "expression": f"({aoe_param}>a) & ({self.dt_cut_param})", - "parameters": {"a": self.low_cut_val}, - } - } - ) - else: - self.update_cal_dicts( - { - output_cut_param: { - "expression": f"({aoe_param}>a)", - "parameters": {"a": self.low_cut_val}, - } + + self.update_cal_dicts( + { + output_cut_param: { + "expression": f"({aoe_param}>a)", + "parameters": {"a": self.low_cut_val}, } - ) + } + ) def calculate_survival_fractions_sweep( self, @@ -1777,24 +1768,14 @@ def calibrate( df["AoE_Classifier"] < self.high_cut_val ) - if self.dt_cut_param is not None and self.dt_cut_hard is True: - self.update_cal_dicts( - { - "AoE_High_Side_Cut": { - "expression": f"(a>AoE_Classifier)& ({self.dt_cut_param})", - "parameters": {"a": self.high_cut_val}, - } - } - ) - else: - self.update_cal_dicts( - { - "AoE_High_Side_Cut": { - "expression": "(a>AoE_Classifier)", - "parameters": {"a": self.high_cut_val}, - } + self.update_cal_dicts( + { + "AoE_High_Side_Cut": { + "expression": "(a>AoE_Classifier)", + "parameters": {"a": self.high_cut_val}, } - ) + } + ) self.update_cal_dicts( { From 8d441e4d8e55ace232ab0e4703081bbd1754c106 Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Wed, 8 May 2024 23:13:33 +0200 Subject: [PATCH 2/7] fix name for LQ_Corrected --- src/pygama/pargen/lq_cal.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/pygama/pargen/lq_cal.py b/src/pygama/pargen/lq_cal.py index 4df62276b..bf6f0daea 100644 --- a/src/pygama/pargen/lq_cal.py +++ b/src/pygama/pargen/lq_cal.py @@ -475,7 +475,7 @@ def drift_time_correction( self.update_cal_dicts( { - "LQ_Classifier": { + "LQ_Corrected": { "expression": f"{lq_param} - dt_eff*a - b", "parameters": {"a": self.dt_fit_pars[0], "b": self.dt_fit_pars[1]}, } From 936ab451887f1a544096cf42c1c2b05a21de89c6 Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Thu, 9 May 2024 11:36:34 +0200 Subject: [PATCH 3/7] binning fix --- src/pygama/pargen/energy_cal.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/pygama/pargen/energy_cal.py b/src/pygama/pargen/energy_cal.py index 3959dc781..9d4e3c82f 100644 --- a/src/pygama/pargen/energy_cal.py +++ b/src/pygama/pargen/energy_cal.py @@ -477,8 +477,8 @@ def hpge_cal_energy_peak_tops( (Polynomial(self.pars) - i).roots() for i in (peaks_kev[0] * 0.9, peaks_kev[-1] * 1.1) ) - euc_min = euc_min[0] - euc_max = euc_max[0] + euc_min = np.nanmin(euc_min) + euc_max = np.nanmax(euc_max) if euc_min < 0: euc_min = 0 @@ -831,8 +831,8 @@ def hpge_fit_energy_peaks( (Polynomial(self.pars) - i).roots() for i in (peaks_kev[0] * 0.9, peaks_kev[-1] * 1.1) ) - euc_min = euc_min[0] - euc_max = euc_max[0] + euc_min = np.nanmin(euc_min) + euc_max = np.nanmax(euc_max) if euc_min < 0: euc_min = 0 if euc_max > np.nanmax(e_uncal) * 1.1: From 74bb56c3af956a6f5b475369cb1172e03cee897e Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Thu, 9 May 2024 13:23:52 +0200 Subject: [PATCH 4/7] lq fail fix --- src/pygama/pargen/lq_cal.py | 1 + 1 file changed, 1 insertion(+) diff --git a/src/pygama/pargen/lq_cal.py b/src/pygama/pargen/lq_cal.py index bf6f0daea..dd41b2c99 100644 --- a/src/pygama/pargen/lq_cal.py +++ b/src/pygama/pargen/lq_cal.py @@ -513,6 +513,7 @@ def get_cut_lq_dep(self, df: pd.DataFrame(), lq_param: str, cal_energy_param: st raise (e) log.error("LQ cut determination failed") self.cut_val = np.nan + pars = np.full(2, np.nan) self.update_cal_dicts( { From ffefcb6eb4d23cf0e80c274b05692dd031e51d47 Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Thu, 9 May 2024 14:03:26 +0200 Subject: [PATCH 5/7] finite check on weights --- src/pygama/pargen/energy_cal.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/src/pygama/pargen/energy_cal.py b/src/pygama/pargen/energy_cal.py index 9d4e3c82f..efcd538bd 100644 --- a/src/pygama/pargen/energy_cal.py +++ b/src/pygama/pargen/energy_cal.py @@ -2484,8 +2484,9 @@ def hpge_fit_energy_cal_func( for n in range(len(energy_scale_pars) - 1): d_mu_d_es += energy_scale_pars[n + 1] * mus ** (n + 1) e_weights = np.sqrt(d_mu_d_es * mu_vars) + mask = np.isfinite(e_weights) poly_pars = ( - Polynomial.fit(mus, energies_kev, deg=deg, w=1 / e_weights).convert().coef + Polynomial.fit(mus[mask], energies_kev[mask], deg=deg, w=1 / e_weights[mask]).convert().coef ) if fixed is not None: for idx, val in fixed.items(): @@ -2493,7 +2494,7 @@ def hpge_fit_energy_cal_func( pass else: poly_pars[idx] = val - c = cost.LeastSquares(mus, energies_kev, e_weights, poly_wrapper) + c = cost.LeastSquares(mus[mask], energies_kev[mask], e_weights[mask], poly_wrapper) m = Minuit(c, *poly_pars) if fixed is not None: for idx in list(fixed): From 6d41997c2301bea53ecc1d40ce3f56029d07a262 Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Thu, 9 May 2024 14:03:26 +0200 Subject: [PATCH 6/7] finite check on weights --- src/pygama/pargen/energy_cal.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/src/pygama/pargen/energy_cal.py b/src/pygama/pargen/energy_cal.py index 9d4e3c82f..efcd538bd 100644 --- a/src/pygama/pargen/energy_cal.py +++ b/src/pygama/pargen/energy_cal.py @@ -2484,8 +2484,9 @@ def hpge_fit_energy_cal_func( for n in range(len(energy_scale_pars) - 1): d_mu_d_es += energy_scale_pars[n + 1] * mus ** (n + 1) e_weights = np.sqrt(d_mu_d_es * mu_vars) + mask = np.isfinite(e_weights) poly_pars = ( - Polynomial.fit(mus, energies_kev, deg=deg, w=1 / e_weights).convert().coef + Polynomial.fit(mus[mask], energies_kev[mask], deg=deg, w=1 / e_weights[mask]).convert().coef ) if fixed is not None: for idx, val in fixed.items(): @@ -2493,7 +2494,7 @@ def hpge_fit_energy_cal_func( pass else: poly_pars[idx] = val - c = cost.LeastSquares(mus, energies_kev, e_weights, poly_wrapper) + c = cost.LeastSquares(mus[mask], energies_kev[mask], e_weights[mask], poly_wrapper) m = Minuit(c, *poly_pars) if fixed is not None: for idx in list(fixed): From e1151c00787856c913f6af99e467b6cf60ccea6b Mon Sep 17 00:00:00 2001 From: ggmarshall Date: Thu, 9 May 2024 14:17:22 +0200 Subject: [PATCH 7/7] lq fix and pc --- src/pygama/pargen/energy_cal.py | 10 ++++++++-- src/pygama/pargen/lq_cal.py | 4 +++- 2 files changed, 11 insertions(+), 3 deletions(-) diff --git a/src/pygama/pargen/energy_cal.py b/src/pygama/pargen/energy_cal.py index efcd538bd..bdac157ca 100644 --- a/src/pygama/pargen/energy_cal.py +++ b/src/pygama/pargen/energy_cal.py @@ -2486,7 +2486,11 @@ def hpge_fit_energy_cal_func( e_weights = np.sqrt(d_mu_d_es * mu_vars) mask = np.isfinite(e_weights) poly_pars = ( - Polynomial.fit(mus[mask], energies_kev[mask], deg=deg, w=1 / e_weights[mask]).convert().coef + Polynomial.fit( + mus[mask], energies_kev[mask], deg=deg, w=1 / e_weights[mask] + ) + .convert() + .coef ) if fixed is not None: for idx, val in fixed.items(): @@ -2494,7 +2498,9 @@ def hpge_fit_energy_cal_func( pass else: poly_pars[idx] = val - c = cost.LeastSquares(mus[mask], energies_kev[mask], e_weights[mask], poly_wrapper) + c = cost.LeastSquares( + mus[mask], energies_kev[mask], e_weights[mask], poly_wrapper + ) m = Minuit(c, *poly_pars) if fixed is not None: for idx in list(fixed): diff --git a/src/pygama/pargen/lq_cal.py b/src/pygama/pargen/lq_cal.py index 08d73f25d..669330b2b 100644 --- a/src/pygama/pargen/lq_cal.py +++ b/src/pygama/pargen/lq_cal.py @@ -513,7 +513,9 @@ def get_cut_lq_dep(self, df: pd.DataFrame(), lq_param: str, cal_energy_param: st raise (e) log.error("LQ cut determination failed") self.cut_val = np.nan - self.cut_fit_pars = pars = np.full(2, np.nan) + c = cost.UnbinnedNLL(np.array([0]), gaussian.pdf) + m = Minuit(c, np.full(2, np.nan)) + self.cut_fit_pars = pars = m.values self.update_cal_dicts( {