photutils: Fitting/Modeling/Sigma Clipping Crashes with NANs, whereas before it just gave a warning

I’m working on PSF Fitting using Photutils. My script previously worked (astropy == 5.0.4 ) worked. I upgraded astropy, and now it fails. Although I’m not exactly sure why. Error message below. Nothing has changed.

Now (astropy == 5.2.dev75+ga3c27114b):

---------------------------------------------------------------------------
NonFiniteValueError                       Traceback (most recent call last)
Input In [26], in <cell line: 22>()
     15 tic = time.perf_counter()
     17 phot = IterativelySubtractedPSFPhotometry(finder=daofind, group_maker=daogroup,
     18                                             bkg_estimator=mmm_bkg, psf_model=psf_model,
     19                                             fitter=fitter,
     20                                             niters=10, fitshape=(7, 7), aperture_radius=ap_radius, 
     21                                             extra_output_cols=('sharpness', 'roundness2'))
---> 22 psf_phot_results = phot(data)
     23 toc = time.perf_counter()
     25 print('Time needed to perform photometry:', '%.2f' % ((toc - tic) / 3600), 'hours')

File ~/miniconda3/envs/hackphot/lib/python3.8/site-packages/photutils/psf/photometry.py:228, in BasicPSFPhotometry.__call__(self, image, init_guesses)
    223 def __call__(self, image, init_guesses=None):
    224     """
    225     Perform PSF photometry. See `do_photometry` for more details
    226     including the `__call__` signature.
    227     """
--> 228     return self.do_photometry(image, init_guesses)

File ~/miniconda3/envs/hackphot/lib/python3.8/site-packages/photutils/psf/photometry.py:758, in IterativelySubtractedPSFPhotometry.do_photometry(self, image, init_guesses)
    755     if self.aperture_radius is None:
    756         self.set_aperture_radius()
--> 758     output_table = self._do_photometry()
    760 output_table.meta = {'version': _get_version_info()}
    762 return QTable(output_table)

File ~/miniconda3/envs/hackphot/lib/python3.8/site-packages/photutils/psf/photometry.py:817, in IterativelySubtractedPSFPhotometry._do_photometry(self, n_start)
    810         init_guess_tab.add_column(
    811             Column(name=param_tab_name,
    812                    data=(getattr(self.psf_model,
    813                                  param_name) *
    814                          np.ones(len(sources)))))
    816 star_groups = self.group_maker(init_guess_tab)
--> 817 table, self._residual_image = super().nstar(
    818     self._residual_image, star_groups)
    820 star_groups = star_groups.group_by('group_id')
    821 table = hstack([star_groups, table])

File ~/miniconda3/envs/hackphot/lib/python3.8/site-packages/photutils/psf/photometry.py:425, in BasicPSFPhotometry.nstar(self, image, star_groups)
    419 for row in star_groups.groups[n]:
    420     usepixel[overlap_slices(large_array_shape=image.shape,
    421                             small_array_shape=self.fitshape,
    422                             position=(row['y_0'], row['x_0']),
    423                             mode='trim')[0]] = True
--> 425 fit_model = self.fitter(group_psf, x[usepixel], y[usepixel],
    426                         image[usepixel])
    427 param_table = self._model_params2table(fit_model,
    428                                        star_groups.groups[n])
    429 result_tab = vstack([result_tab, param_table])

File ~/miniconda3/envs/hackphot/lib/python3.8/site-packages/astropy/modeling/fitting.py:269, in fitter_unit_support.<locals>.wrapper(self, model, x, y, z, **kwargs)
    264         raise NotImplementedError("This model does not support being "
    265                                   "fit to data with units.")
    267 else:
--> 269     return func(self, model, x, y, z=z, **kwargs)

File ~/miniconda3/envs/hackphot/lib/python3.8/site-packages/astropy/modeling/fitting.py:1222, in _NonLinearLSQFitter.__call__(self, model, x, y, z, weights, maxiter, acc, epsilon, estimate_jacobian)
   1219 model_copy.sync_constraints = False
   1220 farg = (model_copy, weights, ) + _convert_input(x, y, z)
-> 1222 init_values, fitparams, cov_x = self._run_fitter(model_copy, farg,
   1223                                                  maxiter, acc, epsilon, estimate_jacobian)
   1225 self._compute_param_cov(model_copy, y, init_values, cov_x, fitparams, farg)
   1227 model.sync_constraints = True

File ~/miniconda3/envs/hackphot/lib/python3.8/site-packages/astropy/modeling/fitting.py:1283, in LevMarLSQFitter._run_fitter(self, model, farg, maxiter, acc, epsilon, estimate_jacobian)
   1281     dfunc = self._wrap_deriv
   1282 init_values, _, _ = model_to_fit_params(model)
-> 1283 fitparams, cov_x, dinfo, mess, ierr = optimize.leastsq(
   1284     self.objective_function, init_values, args=farg, Dfun=dfunc,
   1285     col_deriv=model.col_fit_deriv, maxfev=maxiter, epsfcn=epsilon,
   1286     xtol=acc, full_output=True)
   1287 fitter_to_model_params(model, fitparams)
   1288 self.fit_info.update(dinfo)

File ~/miniconda3/envs/hackphot/lib/python3.8/site-packages/scipy/optimize/_minpack_py.py:410, in leastsq(func, x0, args, Dfun, full_output, col_deriv, ftol, xtol, gtol, maxfev, epsfcn, factor, diag)
    408 if not isinstance(args, tuple):
    409     args = (args,)
--> 410 shape, dtype = _check_func('leastsq', 'func', func, x0, args, n)
    411 m = shape[0]
    413 if n > m:

File ~/miniconda3/envs/hackphot/lib/python3.8/site-packages/scipy/optimize/_minpack_py.py:24, in _check_func(checker, argname, thefunc, x0, args, numinputs, output_shape)
     22 def _check_func(checker, argname, thefunc, x0, args, numinputs,
     23                 output_shape=None):
---> 24     res = atleast_1d(thefunc(*((x0[:numinputs],) + args)))
     25     if (output_shape is not None) and (shape(res) != output_shape):
     26         if (output_shape[0] != 1):

File ~/miniconda3/envs/hackphot/lib/python3.8/site-packages/astropy/modeling/fitting.py:1077, in _NonLinearLSQFitter.objective_function(self, fps, *args)
   1074     value = np.ravel(weights * (model(*args[2: -1]) - meas))
   1076 if not np.all(np.isfinite(value)):
-> 1077     raise NonFiniteValueError("Objective function has encountered a non-finite value, "
   1078                               "this will cause the fit to fail!")
   1080 return value

NonFiniteValueError: Objective function has encountered a non-finite value, this will cause the fit to fail!`