Fitting a single PSF model to multiple stars

[entaylor]

I have not been able to figure out how to fit a single PSF model to a group of stars across an image. The only indication i have of the mode of failure is the message:

Final Chi^2/DoF: 66.63856259729535, L: 1.0. Converged: fail. Could not find step to improve Chi^2
Unable to update uncertainty due to non finite covariance matrix

As i describe below, i have a structure that works for a handful of stars, but fails when the number of stars is too big. I haven't been able to figure out what 'too big' means in any precise way. Some groups of 10 or even 15 stars i can fit without too much of an issue, but some groups of 5 will fail. I have verified that i can fit each individual star; it's when they are grouped that i run into issues.

Below is a block of code that will (i hope) show the issue. The code block uses example data that you can grab from this dropbox link.

image = fits.getdata('example_image.fits').astype(np.float32)
variance = fits.getdata('example_variance.fits').astype(np.float32)

# create crude image mask
okay = (np.isfinite(image) & (variance > 0))
mask = ~okay

# create the astrophot target image object
target_image = astrophot.image.Target_Image(
    data=image, variance=variance, mask=mask, pixelscale=1 )

# load an empirical PSF image as basis for target model
psf_image = fits.getdata( 'example_psf.fits' ).astype(np.float32)
# note need for conversion to float32 to manage fits big v little-endian? 
halfsize=psf_target.shape[0]//2

# create the astrophot psf image object
psf_target = astrophot.image.PSF_Image( data=psf_image, pixelscale=1 )
psf_target.normalize()

# define the basic (Moffat) psf model
psf_model_pars = {
    "n":{"value":2., "limits":[0.5,5.]}, 
    "Rd":{"value:":2., "limits":[0.5,10]}
    }
psf_model = astrophot.models.AstroPhot_Model(
    model_type = "moffat psf model",
    parameters = psf_model_pars, 
    target = psf_target
)
psf_model.initialize()

# load the ra/dec and x/y positions of a set of bright-ish stars
stars = Table.read( 'example_stars.fits' )

# now we loop over each star to create a point source model
all_models, all_psfs = [], []
for ii, star in enumerate( stars[0:10], start=1 ):

    # define fitting window based on centroid positions    
    xcen, ycen = star['x_centroid_pixel'], star['y_centroid_pixel']
    # use pixel centroids for astrophot windows
    xlo = max(int(xcen)-halfsize, 0)
    xhi = min(int(xcen)+halfsize+1, imaging.image.shape[0])
    ylo = max(int(ycen)-halfsize, 0)
    yhi = min(int(ycen)+halfsize+1, imaging.image.shape[1])
    window = [[xlo, xhi], [ylo, yhi]]

    # create the star model, including target window
    star_model = astrophot.models.AstroPhot_Model(
        name=f'star{ii:03d}', 
        model_type='point model', 
        parameters={"center":[xcen, ycen], "flux":1},
        psf=psf_model,
        psf_mode='full',
        target=target, 
        pixelscale=1,
        window=window,                    
        )

    all_models.append( star_model )

model = astrophot.models.AstroPhot_Model(
    name='psf star models',
    model_type='group model',
    psf_mode='full',
    models=all_models,
    target=target
)

model.initialize()
print( model.parameters )

result = astrophot.fit.LM( model, verbose=1 ).fit()
result.update_uncertainty()

print( model.parameters )

If i fit each star independently and average the results, the Moffat parameters should be something like n=1.32 and Re=1.89. When i try to stars[0:10], as above, the fit fails with the non finite covariance matrix error. When i try to fit stars[0:5], the fit similarly fails.

When i try to fit stars[0:4], then that works with approximately the expected result, and when i try to fit stars[4:10], then that also works. So it seems that i can get simultaneous fits to groups of stars, but some stars do not play well with some others stars?

This is working in pycharm on a Mac, with astrophot 0.15.0

I'd be very grateful for any advice or direction on the right way to structure this problem!

And thank you, thank you, thank you for this awesome tool!! I am having a lot of fun wrapping my mind around its capabilities, and formulating grand plans!! : D

[ConnorStoneAstro]

Hi @entaylor, I've been playing around with this this morning and I think I've identified the issue, and it's not with your setup which seems perfectly reasonable to me.

There are in fact two problems and both are quite solvable. One has to do with the fact that these stars take up a small part of the much larger image, so when AstroPhot is trying to update the chi^2 it sees a very small change despite large adjustments to the parameters. It isn't something you should need to worry about and so I'll make a small patch to fix it. Similarly, there was one spot where AstroPhot wasn't properly accounting for the mask that was also causing problems. I'll have a patch put together today to fix this, so be on the lookout for AstroPhot v0.15.3 and everything should start working as expected :)

[ConnorStoneAstro]

HI @entaylor ok I think everything should be ready to go. Please let me know if you run into any other problems!

Also, one minor point about the script you sent. When setting the values for a parameter you have an extra colon in the Rd parameter value:

psf_model_pars = {
    "n":{"value":2., "limits":[0.5,5.]}, 
    "Rd":{"value:":2., "limits":[0.5,10]}
    }

should be

psf_model_pars = {
    "n":{"value":2., "limits":[0.5,5.]}, 
    "Rd":{"value":2., "limits":[0.5,10]}
    }

Closed by #175

[entaylor]

This is just phenomenal to have such a quick and effective response. Thank you, Connor!

I can confirm that i am now able to run through the larger set of stars for this test image. Woohoo!

In the way of things, though, just one more question ... I am getting an error when trying to update the uncertainties for the larger fits: "Unable to update uncertainty due to non finite covariance matrix"

Is there any chance the inability to update the uncertainties is related to a similar issue?

PS. Thank you for spotting the extra colon!! Never in a million years would have spotted it. Even with your message, i had to cut and paste into a text editor and do a character by character comparison to find it!!

[ConnorStoneAstro]

You are absolutely right that there was a similar issue in the uncertainty update. I've fixed that now, so you should be able to get full uncertainties for your fits :) Glad you like the project, feel free to let me know if you spot any other issues!

[entaylor]

Works!! Thank you again! You are amazing!!!