FUNCTION xrt_unspike, image_in, sens=sens, cleanjpg=cleanjpg

; =========================================================================
;+
; PROJECT:
;       Solar-B / XRT 
;
; NAME:
;       XRT_UNSPIKE
;
; CATEGORY:
;       Data calibration
;
; PURPOSE:
;       To remove bad pixels due to cosmic ray hits.
;
; CALLING SEQUENCE:
;       result = XRT_UNSPIKE(image_in [,sens=sens] [,/cleanjpg])
;
; INPUTS:
;
;       IMAGE_IN  - [Mandatory] (2-dim number array, [Nx,Ny])
;                   The input image array.
;
; KEYWORDS:
;
;       SENS      - [Optional] (integer scalar)
;                   Sensitivity to spots. (sens = 1.5 is default; 
;                   for UV, 0.4 is better)
;       /CLEANJPG - [Optional] (Boolean)
;                   If set, call <xrt_cleanjpg.pro>. This cleans up
;                   the residue JPEG artifacts from the particle hits.
;                   
; OUTPUTS:
;
;       return    - (2-dim number array, [Nx,Ny])
;                   The output image array, with hits removed.
;
; PROCEDURE:
;       Spikes, i.e., bright pixels, are amplified and located by
;       convolution and thresholding. This information is used to
;       generate a map of bad pixels. The bad pixel map is then
;       revised to include nearest neighbors.
;
; EXAMPLES:
;
;       Basic usage:
;       IDL> image_out = xrt_unspike(image_in, /cleanjpg)
;
;
; COMMON BLOCKS:
;       None.
;
; NOTES:
;
;       1) Image is converted to a float array, if it's not already.
;
; MODIFICATION HISTORY:
;
progver = 'v2006.Jun.12' ; --- (JWC) Written. (Based heavily on
;                              <trace_unspike.pro>.)
progver = 'v2007.Feb.16' ; --- (MW) Cleaned up.
;
;-
; =========================================================================

;===== Preparation

  if keyword_set(sens) then sensitivity = sens $
                       else sensitivity = 1.5  ; a reasonable default

  epsilon = 0.1

  image=float(image_in)+epsilon


;===== (1) Make map of good pixels.  

  kernel=[[-0.125,-0.125,-0.125], $
          [-0.125, 1.0,  -0.125], $
          [-0.125,-0.125,-0.125]]

  enhanced = convol(image,kernel,/edge_truncate)


;===== The sensitivity below must be verified.

  good_pixmap = float(100.0 GT (sensitivity*enhanced))   

;===== (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,3)/ck_smooth(good_pixmap,3) >0 <image
  result= good_pixmap*image + (NOT good_pixmap)*newpix - epsilon

;===== Clean artifacts, if requested.
  if keyword_set(cleanjpg) then result=xrt_cleanjpg(temporary(result))

;===== Finish.
  return, result

END ;======================================================================