;Mask out circles
;Line mask format: CIRCLE[_IN] ra_cent dec_cent radius [all in
;degrees]
;astrometry keyword for image input use to convert ra,dec mask
;coordinates into pixel using the image WCS
function mask_circle, line, sx, sy, ra_min, dec_min, spix, $
                      astrometry = astrometry, radec_cent = radec_cent, $
                      test = test, image = image
;stop
;get circle parameters
line_tab  = strsplit(line, /extract)
mask_circle = double(line_tab[1 : 3])

;convert in pixels
if (not keyword_set(astrometry)) then begin
    mask_xpix   = long((mask_circle[0] - ra_min)/spix)
    mask_ypix   = long((mask_circle[1] - dec_min)/spix)
    mask_radpix = long(mask_circle[2]/spix)
endif else begin
    circle_cent = deg2radec(radec_cent[0], radec_cent[1], mask_circle[0], $
                            mask_circle[1])
    ad2xy, circle_cent[0], circle_cent[1], astrometry, mask_xpix, mask_ypix
    mask_radpix = long(mask_circle[2]/abs(astrometry.cdelt[0]))
endelse
;select square enclosing the circle (to compute distances on the
;smallest number of pixels)
x = lindgen(sx)
y = lindgen(sy)


squarex = where(x GE mask_xpix-mask_radpix $
                and x LE mask_xpix+mask_radpix, $
                countx)
squarey = where(y GE mask_ypix-mask_radpix $
                and y LE mask_ypix+mask_radpix, $
                county)

if(countx ne 0 and county ne 0) then begin
    dist    = dblarr(countx, county)
    xsquare = lindgen(countx)
    mask_xpix_square = mask_xpix-x[squarex[0]]
    mask_ypix_square = mask_ypix-y[squarey[0]]
    for iy = 0L, county -1L do begin
        dist[*, iy] = ((xsquare - mask_xpix_square)^2 + $
                       (iy -mask_ypix_square)^2)
    endfor
    dist    = sqrt(dist)
    in      = where(dist LT mask_radpix)
    in_cols = in mod countx
    in_rows = in/countx
    
;get pixel coordinate in image frame
    inside_rows = y[squarey[0]]+in_rows
    inside_cols = x[squarex[0]]+in_cols
    
    inside_circle = inside_rows*sx + inside_cols
endif else begin
    inside_circle = -1
endelse
; set_plot, 'x'
; window, 2., xsize = 500, ysize = 500, retain =2
; imtest = intarr(sx, sy)
; imtest[inside_circle] = 1.
; tv, bytscl(congrid(imtest, 500, 500))

; window, 3., xsize = 500, ysize = 500, retain =2

; tv, bytscl(congrid(image, 500, 500))


;stop
; if (keyword_set(test)) then begin
;     imtest = intarr(sx, sy)
;     imtest[inside_circle] = 1.
;     set_plot, 'x'
;     window, 0., xsize = sx, ysize = sy, retain =2
;     tv, bytscl(imtest)
;     stop
; endif
return, inside_circle
end