Merge pull request #206 from desihub/fix-3.1.1-bugs

fix miscellaneous bugs present in 3.1.1

doc/changes.rst

@@ -5,7 +5,14 @@ Change Log
 3.2.0 (not released yet)
 ------------------------
 
-* 
+*
+
+3.1.2 (2025-01-07)
+------------------
+
+* Miscellaneous bug fixes revealed in Loa/Y3 processing [`PR #206`_].
+
+.. _`PR #206`: https://github.com/desihub/fastspecfit/pull/206
 
 3.1.1 (2025-01-05)
 ------------------

py/fastspecfit/emlines.py

@@ -641,14 +641,13 @@ def optimize(self, linemodel,
             fit_info = least_squares(objective, initial_guesses, jac=jac, args=(),
                                      max_nfev=5000, xtol=1e-10, ftol=1e-5, #x_scale='jac' gtol=1e-10,
                                      tr_solver='lsmr', tr_options={'maxiter': 1000, 'regularize': True},
-                                     method='trf', bounds=bounds,) # verbose=2)
+                                     method='trf', bounds=bounds)#, verbose=2)
             free_params = fit_info.x
 
             if not fit_info.success:
                 errmsg = 'least_squares optimizer failed' + \
                     (f' for {self.uniqueid}' if self.uniqueid is not None else '')
                 log.critical(errmsg)
-                raise RuntimeError(errmsg)
             elif fit_info.status == 0:
                 log.warning('optimizer failed to converge')
 
@@ -1283,8 +1282,8 @@ def test_broad_model(EMFit,
     # get the pixels of the broad Balmer lines
     IBalmer = EMFit.isBalmerBroad_noHelium_Strong
 
-    balmer_linesigmas = broad_values[line_params[IBalmer, ParamType.SIGMA ] ]
-    balmer_linevshifts = broad_values[line_params[IBalmer, ParamType.VSHIFT] ]
+    balmer_linesigmas = broad_values[line_params[IBalmer, ParamType.SIGMA]]
+    balmer_linevshifts = broad_values[line_params[IBalmer, ParamType.VSHIFT]]
 
     balmerpix = LineMasker.linepix_and_contpix(
         emlinewave, emlineivar, EMFit.line_table[IBalmer],
@@ -1293,6 +1292,15 @@ def test_broad_model(EMFit,
     balmerlines =  [EMFit.line_map[ln] for ln in balmerpix['linepix']]
     balmerpixels = [px for px in balmerpix['linepix'].values()]
 
+    # balmerlines and balmerpixels can be an empty set when a camera is fully
+    # masked; if so, politely quit here! Example: loa/main/backup/21126/2305843031363822582
+    if len(balmerlines) == 0:
+        log.debug(f'Dropping broad-line model: no good data.')
+        adopt_broad = False
+        delta_linechi2_balmer = 0
+        delta_linendof_balmer = np.int32(0)
+        return adopt_broad, delta_linechi2_balmer, delta_linendof_balmer
+
     # Determine how many lines (free parameters) are in wavelengths in and
     # around the Balmer lines, with and without broad lines.
     balmer_nfree_broad = 0
@@ -1364,7 +1372,6 @@ def test_broad_model(EMFit,
         log.info(f'  sigma_broad={Habroad:.1f} km/s, sigma_narrow={Hanarrow:.1f} km/s')
         if _broadsnr:
             log.info(f'  {_broadsnr} > {minsnr_balmer_broad:.0f}')
-
     else:
         adopt_broad = False
         if dchi2test == False:

py/fastspecfit/linemasker.py

@@ -47,16 +47,23 @@ def _get_linepix(zlinewave, sigma):
             I = (wave > (zlinewave - nsigma * sigma)) * (wave < (zlinewave + nsigma * sigma)) * (ivar > 0.)
             return np.where(I)[0]
 
-        def _get_contpix(zlinewaves, sigmas, nsigma_factor=1., linemask=None, zlyawave=[]):
+        def _get_contpix(zlinewaves, sigmas, nsigma_factor=1., linemask=None,
+                         use_ivar=True, zlyawave=[]):
             # never use continuum pixels blueward of Lyman-alpha
             minwave = np.min(zlinewaves - nsigma_factor * nsigma * sigmas)
             if len(zlyawave) > 0 and minwave < zlyawave:
                 minwave = zlyawave
             maxwave = np.max(zlinewaves + nsigma_factor * nsigma * sigmas)
             if linemask is None:
-                J = (wave > minwave) * (wave < maxwave) * (ivar > 0.)
+                if use_ivar:
+                    J = (wave > minwave) * (wave < maxwave) * (ivar > 0.)
+                else:
+                    J = (wave > minwave) * (wave < maxwave)
             else:
-                J = (wave > minwave) * (wave < maxwave) * (ivar > 0.) * ~linemask
+                if use_ivar:
+                    J = (wave > minwave) * (wave < maxwave) * (ivar > 0.) * ~linemask
+                else:
+                    J = (wave > minwave) * (wave < maxwave) * ~linemask
             return np.where(J)[0]
 
         if linevshifts is None:
@@ -133,7 +140,10 @@ def _get_contpix(zlinewaves, sigmas, nsigma_factor=1., linemask=None, zlyawave=[
                     # line.
                     J = _get_contpix(zlinewave, sigma, nsigma_factor=FACTOR_DEFAULT,
                                      linemask=None, zlyawave=zlyawave)
-                    J = np.delete(J, np.isin(J, pix['linepix'][linename]))
+                    J2 = np.delete(J, np.isin(J, pix['linepix'][linename]))
+                    # extreme check for, e.g., Lya on load/sv1/bright/6682/39632989667198594
+                    if len(J2) > 0:
+                        J = J2
 
                 if len(J) > 0:
                     pix['contpix'][linename] = J
@@ -171,6 +181,14 @@ def _get_contpix(zlinewaves, sigmas, nsigma_factor=1., linemask=None, zlyawave=[
                     if len(J) >= mincontpix:
                         break
 
+                # desperate measures; drop the linemask condition and then drop
+                # the ivar condition
+                if len(J) == 0:
+                    J = _get_contpix(zlinewaves_patch, sigmas_patch, nsigma_factor=factor,
+                                     linemask=None, zlyawave=zlyawave)
+                if len(J) == 0:
+                    J = _get_contpix(zlinewaves_patch, sigmas_patch, nsigma_factor=factor,
+                                     linemask=None, zlyawave=zlyawave, use_ivar=False)
                 if len(J) == 0:
                     errmsg = f'Unable to measure the continuum pixels for patch {patchid}'
                     log.critical(errmsg)