Pargen bug fixes (#582)

legend-exp · May 10, 2024 · f4dc729 · f4dc729
1 parent 629e655
commit f4dc729
Show file tree

Hide file tree

Showing 2 changed files with 16 additions and 6 deletions.
diff --git 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:
@@ -2484,16 +2484,23 @@ 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():
                 if val is True or val is None:
                     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):

diff --git a/src/pygama/pargen/lq_cal.py b/src/pygama/pargen/lq_cal.py
@@ -513,6 +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
+            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(
             {