;+
; NAME:
;	trace_destreak
; PUPROSE:
;	To remove bad pixels in streaks due to cosmic ray hits.
;	The default parameters are optimized for TRACE EUV images.
; CALLING SEQUENCE:
;	result = trace_destreak(image[, sens=sens, thresh=thresh])
; INPUTS:
;	image = 2-dimensional image
; OPTIONAL KEYWORD INPUTS:
;   sens = sensitivity to streaks and spots
;   (sens = 5.0 is default)
;   thresh = threshold for revised pixel map
;   (not used by default)
; PROCEDURE:
;	Streaks at various angles are amplified and located by
;	convolution and thresholding. This information is used to generate
;	a map of good (undamaged) pixels. The pixel map is then revised
;	to include nearest neighbors. Bad pixels are replaced with a 
;	weighted average of nearby good pixels. 
; SIDE-EFFECTS:
;	image is converted to a float array, if it's not already.
; MODIFICATION HISTORY:
;	Originally sxt_destreak by C. Kankelborg
;	1998-Jun-18 CCK Adapted for TRACE SAA/HLZ cleanup
;-

pro convol_enhance, img, ker, goodpix, sens
;Convolve img with kernel ker in all four possible orientations.
;Use this to update the good pixel map.
smoothed=smooth(img,3)
for i=0, 1 do begin  ;used to go to 3, back when z was zero.
	k = rotate(ker,i)/total(abs(ker))
	enhanced = convol(img,k,/edge_wrap)   ;/(smoothed+(5*sens))
	good = float(100.0 GT (sens*enhanced))
	goodpix = goodpix AND good
endfor

end



function trace_destreak, image_in, sens=sens

if keyword_set(sens) then sensitivity = sens $
	else sensitivity = 5.0 ; a reasonable default
	
epsilon = 0.1

image=float(image_in)+epsilon

; (1) Make map of good pixels.

; Find streaks.
a = 0.05
b = 0.0
z = -a/2.0  ;formerly 0.0
good_pixmap = float(image GT 0.0) ;initialize good pixel map

kernel=[[b,b,b,b,b], $
        [z,z,z,z,z], $
        [a,a,a,a,a], $
        [z,z,z,z,z], $
        [b,b,b,b,b]]

convol_enhance, image, kernel, good_pixmap, sensitivity

kernel=[[b,b,b,b,z], $
        [b,z,z,z,a], $
        [z,a,a,a,z], $
        [a,z,z,z,b], $
        [z,b,b,b,b]]

convol_enhance, image, kernel, good_pixmap, sensitivity        

kernel=transpose(kernel)     
          
convol_enhance, image, kernel, good_pixmap, sensitivity          
     
kernel=[[b,b,b,z,z], $
        [b,b,z,a,a], $
        [z,z,a,z,z], $
        [a,a,z,b,b], $
        [z,z,b,b,b]]

convol_enhance, image, kernel, good_pixmap, sensitivity

kernel=transpose(kernel)     
          
convol_enhance, image, kernel, good_pixmap, sensitivity          
     
kernel=[[b,b,b,z,a], $
        [b,z,z,a,z], $
        [b,z,a,z,b], $
        [z,a,z,z,b], $
        [a,z,b,b,b]]
        
convol_enhance, image, kernel, good_pixmap, sensitivity

;print, 'Displaying good pixel map....'
;xtv, good_pixmap

; (2) Revise pixel map, taking nearest neighbors
good_pixmap = float(ck_smooth(temporary(good_pixmap),1) GT 0.9)

; (3) Replace bad pixels
newpix = ck_smooth(good_pixmap*image,4)/ck_smooth(good_pixmap,4) >0 <image
;print,("Displaying newpix....")
;xtv, newpix
return, good_pixmap*image + (NOT good_pixmap)*newpix - epsilon

END