; Copyright (c) 1998-2012 A.P. Hitchcock  All rights reserved
;+
;NAME:
;	CLIP
;
;LAST CHANGED: ----------------------------------- 04-Aug-12
;
;PURPOSE:
;	This function selects a subset of intensities of an image.
;
;CATEGORY:
;	STAND ALONE: image processing
;
;CALLING SEQUENCE:
;	Result = CLIP(tmp)
;
;CALLED FROM AXIS:
;	->Images->Clip

;INPUTS:
;	TMP - an AXIS 2d structure (image)
;
;KEYWORDS:
;	SAVE_HIST  - if set, write out histogram
;
;OUTPUTS:
;	The result is an AXIS 2d structure with modified intensities.
;
;COMMON BLOCKS:
;	@AXIS_COM	standard set of common blocks
;
;SIDE EFFECTS: none
;
;PROCEDURE:
;	A histogram of the intensities of an image are displayed as a spectrum on
; MainImg (if in AXIS) or on the current or a generated window.
; Two cursor clicks then identify lower and upper bounds of interest.
; Lower intensities are set to the lower bound.
; Higher intensities are set to the upper bound.
; This is an alternate to DEGLITCH. DEGLITCH works nicely for localized,
; single pixel 'glitches' but CLLP is better for extended artefacts

;
;EXAMPLE:
; Create a sample nxn image, display, clip, and redisplay
;	n = 100 & t = findgen(n) & d = dist(n)
;	a = {t:'2d', x:t, y:t, xl:'x', yl:'y', d: d, dl:'test'}
;	splot2d, a   ; display input
;	b = CLIP(a)	 ; clip
;	splot2d, b
;
;MODIFICATION HISTORY:
; (31-dec-99 aph) isolated from AXIS_C as a stand-alone
;				  AXIS standard documentation
; (09-jan-00 aph) fixed up error from histogram with big images (out of virtual memory)
; (19-dec-00 cgz) substituted ax_white_color_index with ax_plot_axes_color_index
;			substituted ax_black_color_index with ax_plot_bkgd_color_index
; (08-may-01 aph) tried to fix problem with flakey clipping - NOT FIXED !!
;		output of histogram somehow seems to get mis-identified as a structure
;		suspect variable confusion in histogram.pro, but this is not an open file
;		REPLACING FLOAT with a LONG seems to have fixed this - values in array is the issue
; (30-sep-01 aph) removed 'LONG' - this does not seem to do anything useful
; (13-mar-02 aph) add save_hist keyword
; (10-apr-02 aph) replace IDL clip ruotine with where command
; (02-may-02 aph) adopt RSI suggestion to make clip work reliably
; (10-apr-04 aph) add annotation to header
; (20-oct-04 aph) simplify tag
; (15-Jul-12 aph) add numerical clipping approach; replace NaN values with minimum
; (04-Aig-12 aph) save values of clip limits in common (aparam/xlock(1), ylock(1) )
;-

FUNCTION CLIP, tmp, save_hist=save_hist, numeric=numeric
@axis_com
on_error,2

dt = tmp.d
n_bins = 400
;if widget_info(/active) EQ 1 then begin
;		n_bins = get_num(Prompt='Number of bins',Val=n_bins, group = axis_id)
;	endif else 	n_bins = get_num(Prompt='Number of bins',Val=n_bins)
xt = findgen(n_bins) & xmn = min(dt, max=xmx, /nan)
xt = xmn + xt*(xmx-xmn)/n_bins
;
; from RSI - 2-may-02
; dtz = histogram(double(dt), binsize = double((xmx-xmn)/n_bins), /nan)
;
; print, "xmx, xmn, n_bins ", xmx, xmn, n_bins
; print, 'binsize = ', double((xmx-xmn)/n_bins)
dtz = histogram(double(dt), binsize = double((xmx-xmn)/n_bins), /nan)
t = size(dtz)
dtzn = make_array(t(1), /long, value=xmn)
for i = 0, t(1)-1 do dtzn(i) = dtz(i)
if keyword_set(save_hist) then begin
  	tmph = {t:'1d', x:xt, d:dtzn,dn:dtzn, xl:'value', yl:'number', dl:'hist '+tmp.dl}
	t = spc_save(tmph)
endif

if keyword_set(numeric) then begin
	clo = min(dt, /nan) & chi = max(dt, /nan)   ; start with min / max of data
	clo = xlock(1)  & chi = ylock(1)		;  restore last values (will be 0, at start)
	clo = get_num(prompt='lower limit', val = clo, group=axis_ID)
	chi = get_num(prompt='upper limit', val = chi, group=axis_ID)
endif else begin
	if widget_info(/active) EQ 1 then $
		plot,xt,dtzn,color=ax_plot_axes_color_index, background=ax_plot_bkgd_color_index $
		else plot,xt,dtzn
	if widget_info(/active) EQ 1 then $
	    WIDGET_CONTROL, Uprompt, BAD_ID = BadID, SET_VALUE='Clip: select lower limit'
	cursor, x1, y1, /data, /down   ; select lower limit
	if widget_info(/active) EQ 1 then $
	   WIDGET_CONTROL, Uprompt, BAD_ID = BadID, SET_VALUE='Clip: select upper limit'
	cursor, x2, y2, /data, /down   ; select upper limit
	zzz = fltarr(2)   &  zzz(0) = x1(0)   &  zzz(1) = x2(0)
	clo = min(zzz)   & chi = max(zzz)
endelse
; print, 'CLIP: lower limit = ', clo, '  upper limit = ', chi

; ----- save values to xlock(1) and ylock(1)
xlock(1) = clo  & ylock(1) = chi

; --------- pre 4-jul-99 - prompt for average
;ddd = moment(tmp.d)  &  avg = ddd(0)  ; find average of image
;ddd = get_num(prompt='Replace by', val=avg, group=axis_ID)
;avg = ddd(0)
; --------- post 4-jul-99 - replace out-of-bounds with limits
cnt = long(0) & index = lonarr(1) & nval = long(0) ; declare long for big images
index = where(dt LT clo(0),cnt)
if cnt GT 0 THEN dt(index) = clo
nval = cnt
index = where(dt GT chi(0),cnt)
if cnt GT 0 THEN dt(index) = chi
nval = nval + cnt
cnt = size(dt)

; ----- scan for NaN values and replace with minimum
index_nan = where(~finite(dt), count)
if count GT 0 then dt(index_nan)=min(dt, /Nan)
;

print, 'CLIP kept', STRCOMPRESS(string(cnt(4) - nval)), $
       ' values between', STRCOMPRESS(string(clo(0))),' and ', STRCOMPRESS(string(chi(0)))
print, ' and replaced', STRCOMPRESS(string(nval)),' points with limits'
; check if energy defined; if not, add it
names = tag_names(tmp)  & test_E=0
for  i = 0, n_tags(tmp)-1 do if strupcase(names(i)) EQ 'E' then test_E=1
if test_E EQ 1 then energy = tmp.E else energy = 0.
tmp.dl = tmp.dl + 'Cl '+ strcompress(string(clo, format='(G8.2)') + '-' + string(chi, format='(G8.2)'), /remove_all)
tmpN = {t:'2d', x:tmp.x, y:tmp.y, d:dt, xl:tmp.xl, yl:tmp.yl, dl:tmp.dl, E:energy}
return, tmpN
end