FUNCTION fwhm, imdata, bpmdata, sky, skysig

; Description: This function determines an estimate of the full-width half-maximum (FWHM) and
;              ellipticity of the PSF of an input image "imdata". A reliable estimate of the sky 
;              background "sky" and sky sigma "skysig" are required in order to detect and fit 
;              suspected peaks. A maximum of 100 suspected peaks are used in the FWHM estimation. 
;              The function fits an elliptical Gaussian function to each suspected peak, and 
;              calculates the weighted mean of the FWHMs, where a FWHM measurement is weighted 
;              according to the peak's fitted signal-to-noise.
;
; Input Parameters:
;
;   imdata - FLOAT/DOUBLE ARRAY - The input image data for which the PSF FWHM and ellipticity is 
;                                 to be measured. This image must measure at least 25 pixels along
;                                 the x-axis and at least 25 pixels along the y-axis.
;   bpmdata - BYTE/INTEGER/LONG ARRAY - A bad pixel mask array of the same size as "imdata" which 
;                                       flags bad pixels with a value of "0" and good pixels with
;                                       a value of "1". If this array does not have the same
;                                       dimensions as "imdata", then all image pixels will be
;                                       considered good, and "bpmdata" will be set to an array of
;                                       the same dimensions as "imdata" with all values equal to
;                                       "1".
;   sky - FLOAT/DOUBLE - A reliable estimate of the sky background level in the image (ADU).
;   skysig - FLOAT/DOUBLE - A reliable estimate of the sky background sigma in the image (ADU).
;                           This parameter must be non-negative.
;
; Return value:
;
;   The function returns an IDL structure with the following tags:
;
;   avgfwhm - DOUBLE - The average FWHM of the chosen peaks (pix).
;   avgfwhm_a - DOUBLE - The average FWHM of the major axis of the chosen peaks (pix).
;   avgfwhm_b - DOUBLE - The average FWHM of the minor axis of the chosen peaks (pix).
;   avgell - DOUBLE - The average ellipticity of the chosen peaks.
;   ngood - INTEGER - The number of peaks used to estimate the average FWHM.
;
; Author: Dan Bramich (dan.bramich@hotmail.co.uk)
;
; History:
;
;   06/05/2008 - Module created (dmb)

;Check that "imdata" is an image of the correct number type and dimensions

;If "bpmdata" is not an image of the correct number type or does not have the same dimensions 
;as "imdata", then assume that all image pixels are good

;Check that "sky" and "skysig" are numbers of the correct type

;Calculate the maximum pixel value

;Identify all candidate peaks to be used in the FWHM estimate as good pixels with values greater
;than 5 sigma above the sky background level

;If there is at least one candidate peak

  ;Limit the candidate peaks to those pixels that are a local maximum in a 7x7 pixel square centred
  ;at the peak

;If there is still at least one candidate peak

  ;Limit the candidate peaks to those pixels where a 3x3 pixel square centred at the peak all have
  ;values at least 5 sigma above the sky background level

;If there is still at least one candidate peak

  ;Limit the candidate peaks to those pixels that have no bad pixels in a 21x21 pixel square centred
  ;at the peak

;If there is still at least one candidate peak

  ;Sort the peaks in to ascending order

  ;Repeat until there are no more peak pixels, or until the fitted peak of the current peak pixel
  ;is less than 1% of the fitted peak of the highest peak pixel, or until at least 100 peak pixels
  ;have been processed

    ;Find the centroid of the current peak

    ;If the centroid calculation is successful

      ;Fit the current peak with a 2D elliptical Gaussian with a rotation parameter

      ;If the fit of the 2D elliptical Gaussian was successful

        ;If the fitted peak is high enough, then record the fit parameters

;Calculate the average FWHM and ellipticity, weighted by star signal-to-noise

;Return the fwhm, ellipticity and number of peaks used