HW 6

1) Write a simple wrapper function called gauss2d_wrap.pro that accepts as input an array of parameters, P, which has the same format as the parameter array used by gauss2d.pro. Include the following keywords in your function header:

X = X : a 1-D vector of X values

Y = Y : a 1-D vector of Y values

IMAGE = IMAGE : a 2-D array of data values

ERR = ERR : the uncertainties associated with each pixel in the data array

The function should compute a 2-D Gaussian, fit, defined by the input parameters, P, and return the residuals as a 1-D array:

return,reform((IMAGE - fit)/ERR, n_elements(IMAGE))

This function is fairly useless by itself, but necessary for the next part…

2) Write a function called gaussfit2d.pro that fits a 2-D Gaussian to a star image using MPFIT. The input to your function is a 2-D image (or subimage) containing a star. The output should be a 2-D Gaussian with the same dimensions as the input image. Return the best-fit parameters through an optional positional parameter called par.

Download this test image and use restore to read it at the IDL prompt:

IDL> restore,'data.sav'

IDL> help,data

DATA FLOAT = Array[50, 50]

IDL> surface, data, /lego

IDL> fit = gaussfit2d(data, par)

IDL> print,par

109.824 699.876 5.05224 4.98996 28.0130 22.0446 0.00000

IDL> surface, fit, /lego

IDL> surface, data-fit, /lego

HINTS:

a) You’ll need to set up a functargs structure containing fields for X, Y,IMAGE and ERR. Assume the errors are Poisson, or just sqrt(data). The X and Y values should be in pixel units. Inside gaussfit2d, your call to MPFIT should look like this:

par = mpfit(‘gauss2d_wrap’, pguess, functargs=functargs)

b) I basically set this up for you in class. The hardest part is figuring out how to get your code to generate the parameter guesses based on the data.

c) After everything is working properly, you should use the /quiet keyword to MPFIT.