;pro Return_Persistence_Mask
;
;PURPOSE:
;
;   To take in a mask of pixels that were saturated in an image and
;   extend the mask radially outward to enclose all of the pixels that
;   have either negative or positive persistence,
;
;INPUTS:
;	FitsFileName - The name of the slopefit that contains the mask
;		       of saturated pixels.  In most cases it should 
;		       be the 0th dither
;
;KEYWORDS:
;
function Return_Persistence_Mask, SatMask, TvFlag = TvFlag, RadMax = RadMax

If N_Elements(TvFlag)  eq 0 then TvFlag = 0 
If N_Elements(RadMask) eq 0 then RadMask = 50
SatRegSize = 2*RadMax+1

Common KeyParams, KeyStr

;The mask to be returned will have 1's in circles around bright stars
PersistenceMask = BytArr(KeyStr.Naxis1, KeyStr.Naxis2)

;Find the regions that were occupied by bright stars as opposed to 
;very high E cosmic rays or messed up pixels
HMin     = 2.		; minimum counts necessary to see a peak for find.pro
FWHM     = 12.		; FWHM cutoff for find.pro
RoundLim = [-2.,2.]	; roundness paramtetrs for find.pro
SharpLim = [0.2,1.]	; sharpness parameters for find.pro
PhotAper = [2.]		; the aperture radius used by aper.pro
RdNoise  = [3.]		; the readnoise used by psf.pro
PhPerADU = [1.]		; the number of photons per adu
PSFRad   = 40.		; the radius of the circular area in which psf will be def
FitRad   = 12.
TReads   = KeyStr.TotalReads

Sky, SatMask, SkyMode, SkySig
Find, SatMask, XStar, YStar, Flux, Sharp, Roundness, HMin, $
	FWHM, RoundLim, SharpLim,/silent 
If N_Elements(XStar) eq 0 then return, PersistenceMask

Aper, SatMask, XStar, YStar, Mag, ErrAp, Sky, SkyErr, 1., PhotAper, -1, $
        [-1., TReads+1], setskyval = SkyMode, /silent

NumStars = N_Elements(XStar)
XSigs = Fltarr(NumStars) & YSigs = Fltarr(NumStars)
print ,'NumStars:' +Strtrim(NumStars,2)
For StarNum = 0, NumStars-1 do begin
  
    XStart = Round(XStar(StarNum)) - 25 > 0 
    XStop  = Round(XStar(StarNum)) + 25 < KeyStr.Naxis1
    YStart = Round(YStar(StarNum)) - 25 > 0 
    YStop  = Round(YStar(StarNum)) + 25 < KeyStr.Naxis1

    Centroid = Find_Centroid(SatMask, XStart, XStop, YStart, YStop, IsStar, $
			     SatRad, TvFlag = TvFlag, Verbose = KeyStr.Verbose)
    XSigs(StarNum) = Centroid(2) & YSigs(StarNum) = Centroid(3)

    ;Check if we'll stay in bounds
    If Centroid(0) + SatRegSize/2 gt KeyStr.Naxis1 or $
       Centroid(0) - SatRegSize/2 lt 0 or $
       Centroid(1) + SatRegSize/2 gt KeyStr.Naxis2 or $
       Centroid(1) - SatRegSize/2 lt 0 then goto, SkipStar

    RelXCoords = Indgen(SatRegSize,SatRegSize) mod SatRegSize
    RelXCoords = RelXCoords - Mean(RelXCoords)
    RelYCoords = reverse(Indgen(SatRegSize,SatRegSize),2) / SatRegSize
    RelYCoords = RelYCoords - Mean(RelYCoords)
    RsqCoords  = RelXCoords^2+RelYCoords^2
    If IsStar then begin
       RelPixInRad = where(sqrt(RsqCoords) lt RadMax)
       AbsPixInds  = RelXCoords(RelPixInRad) + Centroid(0) + $
			(RelYCoords(RelPixInRad) + Centroid(1))* KeyStr.Naxis1
       PersistenceMask(AbsPixInds) = 1

    EndIf	
    SkipStar:  ;The star was too far off

EndFor

return, PersistenceMask

End