PRO centroid_thresh_2d, imdata, vardata, bpmdata, sky, thresh, xc, xc_err, yc, yc_err, status, NOSKY = nosky, NOTHRESH = nothresh

; Description: This module calculates the centroid of an image cutout "imdata" with a bad pixel mask
;              "bpmdata" and pixel variances "vardata". The method used is the standard centroid 
;              moment formula with pixel thresholding, meaning that the centroid is only calculated 
;              using good pixels with values greater than the value "sky + thresh". The standard
;              centroid formula is:
;
;              xc = sum_i (x_i * (P_i - S)) / sum_i (P_i - S)
;              yc = sum_i (y_i * (P_i - S)) / sum_i (P_i - S)
;
;              where x_i and y_i are the x and y coordinates of each pixel in the image cutout, P_i
;              are the pixel values in the image cutout, S is the background sky level, and xc
;              and yc are the centroid coordinates. The sum_i represents the sum over all relevant
;              image cutout pixels. Note that the pixel variances are only used when calculating
;              the centroid coordinate uncertainties.
;                The centroid and its uncertainty are returned via the variables "xc" and "xc_err" 
;              for the x coordinate, respectively, and via the variables "yc" and "yc_err" for the 
;              y coordinate, respectively.
;                The option exists for the algorithm to estimate a sky background if not already
;              known by setting the NOSKY keyword. The sky background is estimated as the average 
;              of the four corner pixel values for "imdata", and consequently may not be very reliable. 
;              In this case, the sky estimate used is returned via the "sky" parameter.
;                There is also the option of calculating a non-thresholded centroid by setting the
;              NOTHRESH keyword. In this case, all good pixels in "imdata" are used to calculate the
;              centroid.
;
; Input Parameters:
;
;   imdata - FLOAT/DOUBLE ARRAY - A two-dimensional image array for which the centroid is to be 
;                                 computed.
;   vardata - FLOAT/DOUBLE ARRAY - A two-dimensional array of the same size as "imdata" containing 
;                                  the pixel variances. If this array does not have the same dimensions 
;                                  as "imdata", then the centroid uncertainty is not calculated. 
;   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 - An estimate of the background sky level (ADU). This value is overwritten if 
;                        the keyword NOSKY is set.
;   thresh - FLOAT/DOUBLE - All pixels in "imdata" with values greater than the sky background level 
;                           "sky" by at least the amount "thresh" (ADU) will be included in the 
;                           centroid calculation.
;
; Output Parameters:
;
;   xc - DOUBLE - The x pixel coordinate of the calculated centroid (pix).
;   xc_err - DOUBLE - The uncertainty in the value of "xc" (pix). This parameter is returned with
;                     a value of "-1.0" if the input parameter "vardata" is not of the same
;                     dimensions as "imdata".
;   yc - DOUBLE - The y pixel coordinate of the calculated centroid (pix).
;   yc_err - DOUBLE - The uncertainty in the value of "yc" (pix). This parameter is returned with
;                     a value of "-1.0" if the input parameter "vardata" is not of the same 
;                     dimensions as "imdata".
;   status - INTEGER - If the centroid was calculated successfully, then "status" is returned 
;                      with a value of 1, otherwise it is returned with a value of 0.
;
; Keywords:
;
;   If the keyword NOSKY is set (as "/NOSKY"), then the value of "sky" is ignored and the module 
;   attempts to estimate a local sky background from the average of the four corner pixel values
;   of the image data. Note that this sky estimate will be unreliable if the image data does not 
;   fully encompass the object for which the centroid is to be measured, or if the image data is 
;   contaminated by other nearby objects.
;
;   If the keyword NOTHRESH is set (as "/NOTHRESH"), then the value of "thresh" is ignored and the 
;   module includes all good pixels in the calculation of the image centroid.
;
; Author: Dan Bramich (dmb@ing.iac.es)
;
; History:
;
;   25/07/2008 - Module created (dmb)

;Check that "imdata" is an image

;If "vardata" does not have the same dimensions as "imdata", then do not calculate uncertainties
;on the image centroid coordinates

;If "bpmdata" does not have the same dimensions as "imdata", then assume that all image pixels 
;are good

;If the keyword NOSKY is set

  ;Estimate the sky background as the average of the four corner pixels

;If the keyword NOSKY is not set

  ;Check that "sky" is a number

;If an estimate of the sky background exists

  ;If the keyword NOTHRESH is set

    ;Determine the set of good pixels

  ;If the keyword NOTHRESH is not set

    ;Check that "thresh" is a number

    ;Determine the set of good pixels that have values greater than or equal to the sky 
    ;background plus the input threshold "thresh"

  ;If there is at least one pixel to be used in the centroid calculation 

    ;Prepare the pixel coordinate arrays

    ;Calculate the centroid coordinates

      ;If possible, then calculate the centroid coordinate uncertainties