;+
; NAME:
;	Interpolation under and above tresholds
; PURPOSE:
;	Interpols values under given minimum and/or given 
;	 maximum:
;	that is: replaces them with the linear interpolation
;	of the two nearest values above and/or under the 
;	tresholds, respectively.
; CALLING SEQUENCE:
;	result = TreshInterp( image, min, max [, direction] )
; INPUTS:
;	image: a 2D array
;	min, max: the treshold values
;	direction: 'X' or 'Y'
; RESULT:
;	the interpolated 2D image
; PROCEDURE:
;	1. finds out which pixels ar under/above the tresholds
;	2. use for these ones the function ReplLinInt
;
; MODIFICATION HISTORY:
;	Created in  December, 1991 by A.Csillaghy
;		Inst. of Astronomy, ETH Zurich
;	InterInTresh modified for excluding n=1 cases, 
;		June 93, A.Cs.
;-

FUNCTION InterInTresh, vector, indices, nVect

;
; PURPOSE:
; 	Interpolates the vector values that have
;	the pixel numbers stored in "indices"
;

  n = N_Elements( indices )

  j = 0
  WHILE (j LE n-1) AND (n GT 1) DO BEGIN
    l = (indices( j ) -1) > 0
    WHILE indices(j) EQ indices((j+1)<(n-1))-1 DO j = j+1
    r = (indices(j) + 1) < (nVect-1)
    vector( l:r ) = ReplLinInt( vector( l:r ) )
    j = j+1
  ENDWHILE

  RETURN, vector

END ; interpolate vector


FUNCTION TreshInterp, image, minimum, maximum, direction

  IF N_Elements( direction ) EQ 0 THEN direction = 'X'
  IF direction EQ 'Y' THEN im = Transpose( image ) $
  ELSE im = image
  nx = N_Elements(im(*,0))
  ny = N_Elements(im(0,*))

 
  FOR i = 0, ny-1 DO BEGIN

    lowTresh = Where( im(*,i) LT minimum )
    highTresh = Where( im(*,i) GT maximum )

    im(0,i) = im(0,i) > minimum < maximum
    im(nx-1,i) = im(nx-1,i) > minimum < maximum

    im(*,i) = InterInTresh( im( *,i ), lowTresh, nx )
    im(*,i) = InterInTresh( im(*,i ), highTresh, nx )

  ENDFOR 

  IF direction EQ 'X' THEN RETURN, im $
  ELSE RETURN, Transpose( im )

END ; interpolation inside tresholds