Performance question - PSF photometry

[keflavich]

I'm running some PSF photometry on a relatively small region with not that many stars:

daogroup = DAOGroup(5 * fwhm_pix)
mmm_bkg = MMMBackground()

filtered_errest = stats.mad_std(filtered_data, ignore_nan=True)

daofind_fin = DAOStarFinder(threshold=10 * filtered_errest, fwhm=fwhm_pix, roundhi=1.0, roundlo=-1.0,
                            sharplo=0.30, sharphi=1.40)
phot = IterativelySubtractedPSFPhotometry(finder=daofind_fin, group_maker=daogroup,
                                          bkg_estimator=mmm_bkg,
                                          psf_model=epsf_quadratic_filtered,
                                          fitter=LevMarLSQFitter(),
                                          niters=2, fitshape=(11, 11), aperture_radius=2*fwhm_pix)

>>> len(daofind_fin.find_stars(mask.cutout(filtered_data)))
1655
>>> mask.cutout(filtered_data).shape
(919, 584)

This final step ran for about an hour:

cdata = np.nan_to_num(mask.cutout(filtered_data))
result = phot(cdata)
resid = phot.get_residual_image()

and I now need to run this on the full ~4000x7000 image (NIRCAM)

So I'm wondering if:

• we can add a progress bar to track the progress of either the overall iterations or the subtraction?
• there's any useful advice on how to improve performance. Or, if there are general tips on what significantly harms performance (big DAOGroups are already called out in the docs)

[keflavich]

OK one answer to my question is: avoid including too many sources that are just PSF artifacts; it creates a lot of expensive joint solve groups that are junk anyway. The followup question, then, is how to do that.

[larrybradley]

It's on my to-do list to look into the PSF-fitting code and improve its performance. Some benchmarks in this old issue: https://github.com/astropy/photutils/issues/631.

I'll look into adding progress bars during the PSF fitting and subtraction.

To avoid artifacts, you can run the star finder on the image first, filter any unwanted sources (or adjust threshold/fwhm parameters etc.), and then set the finder=None and input the initial sources init_guesses when calling the finder on your data. With iterative PSF photometry, you must input a finder, which is used after the first iteration (which uses the guesses).

[keflavich]

Thanks. Progressbars would be great.

It's super helpful to know that the non-iterative PSF photometry tools can just take initial source locations. There are plenty of cases where I just need BasicPSFPhotometry, not iterative.

[keflavich]

My temporary hack for a pbar is:

    star_list = filtered_finder(filtered_data)
    group_list = daogroup(star_list)

    pb = ProgressBar(len(group_list))
    lmfitter = LevMarLSQFitter()
    def fitter(*args, **kwargs):
        pb.update()
        return lmfitter(*args, **kwargs)


    phot = BasicPSFPhotometry(finder=None,
                              group_maker=None,
                              bkg_estimator=None,
                              psf_model=epsf_quadratic_filtered,
                              fitter=fitter,
                              fitshape=(11, 11),
                              aperture_radius=2*fwhm_pix)

    result_full = phot(data, init_guesses=star_groups)

EDIT: Nope, this doesn't work, the fitter needs to be a fitter object, not the return from a fitter call

Revised verison:

    pb = tqdm(len(group_list))
    class LevMarLSQFitter_pb(LevMarLSQFitter):
        def __call__(self, *args, **kwargs):
            pb.update()
            return super().__call__(*args, **kwargs)

    phot = BasicPSFPhotometry(finder=None, #filtered_finder,
                              group_maker=None,
                              bkg_estimator=None, #mmm_bkg,
                              #psf_model=psf_modelgrid[0],
                              psf_model=epsf_quadratic_filtered,
                              fitter=LevMarLSQFitter_pb(),
                              fitshape=(11, 11),
                              aperture_radius=2*fwhm_pix)

[keflavich]

And I'll add: part of the motivation for the pbar is that it's hard to tell when my code has frozen vs when it's still working - with these larger JWST images, I've had a lot of cases where the underlying kernel freezes without telling me for several minutes (hours? days?)