function hi_remove_saturation, im, hdr,silent=silent,saturation_limit=saturation_limit,nsaturated=nsaturated,_extra=ex

;+
; $Id: hi_remove_saturation.pro,v 1.3 2009/06/08 11:03:38 crothers Exp $
;
; Project   : STEREO SECCHI
;
; Name      : hi_remove_sauration
;
; Purpose   : This function masks saturated columns
;
; Explanation: Saturation causes bleeding along the column, since the
;        hi camaras are shutterless this results in inaccurate results
;        in the inversion of data in affected columns.
;
; Use       : IDL> new_image=hi_remove_saturation(im,hdr)
;
; Inputs    : image - level 0.5 image in DN (long)
;             hdr -image header, either fits or SECCHI structure
;
; Outputs   :  new_image - saturated columns are replaced by NaN's
;
;
; Keywords  : saturation_limit [ single pixel saturntion level - default 14000 DN]
;           : nsaturated [ number of pixels saturated before cutoff - default 5]
;
; Common    :
;
; Calls     :
;
; Category    : Prep
;
; Prev. Hist. : None.
;
; Written     : Steve Crothers RAL January 2008
;
; $Log: hi_remove_saturation.pro,v $
; Revision 1.3  2009/06/08 11:03:38  crothers
; updated saturation removal code to allow a column with a small amount of
; saturation to be allowed.
;
; Revision 1.2  2008/02/15 16:18:46  crothers
; *** empty log message ***
;
; Revision 1.1  2008/01/09 16:43:22  crothers
; Initial release
;
;
;-

compile_opt idl2

info="$Id: hi_remove_saturation.pro,v 1.3 2009/06/08 11:03:38 crothers Exp $"

;--Check Header-------------------------------------------------------
IF(DATATYPE(hdr) NE 'STC') THEN hdr=SCC_FITSHDR2STRUCT(hdr)

IF (TAG_NAMES(hdr, /STRUCTURE_NAME) NE 'SECCHI_HDR_STRUCT') THEN $
MESSAGE,'ONLY SECCHI HEADER STRUCTURE ALLOWED'

IF (hdr[0].DETECTOR NE 'HI1' AND hdr[0].DETECTOR NE 'HI2') THEN $
message,'Calibration for HI DETECTOR only'
;----------------------------------------------------------------------


;------------------------------------------------------------------------
;

if n_elements(saturation_limit) eq 0 then saturation_limit=14000
if saturation_limit lt 0 then return,im
if n_elements(nsaturated) eq 0 then nsaturated=5

; compute saturation value, nominal times added images time on board sum factor
n_im=hdr.imgseq+1
sum=hdr.summed
dsatval=saturation_limit*float(n_im)*(2.^(sum-1))^2

ii=where(im gt dsatval, nii)

if nii gt 0 then begin
    mask=im*0
    ans=double(im)
    mask[ii]=1
; determine the columns that contain ANY saturation
    colmask=total(mask,2)
    ii=where(colmask gt nsaturated,ctii)
; set all rows in saturated columns to NAN
    if ctii gt 0 then ans[ii,*]=!values.D_NAN
endif else begin
    ans=im
endelse

RETURN,ans
END
