Identifier     flux_conserv

 Purpose        Perform linear interpolation conserving flux 


 Synopsis       status = flux_conserv(xarr, yarr, sarr, refarr, yarr_out,
                         sarr_out)


 Arguments      Name      I/O Type       Description
                -----------------------------------------------------
                xarr       I  Fltarr     x-scale (e.g. wavelength array) that
                                         is to be resampled
                yarr       I  Fltarr     y-scale (e.g. flux array) that
                                         is to be resampled
                sarr       I  Fltarr     sigma array (stdev) for yarr
                refarr     I  Fltarr     reference scale to which x-scale
                                         is to be resampled
                sarr_out   O  Fltarr     the propagated errors of yarr_out
                yarr_out   O  Fltarr     the rebinned yarr (interpolated
                                         flux array)


 Returns        Error status : 1 for fatal, 0 otherwise

 Description    xarr has to be strictly monotonically increasing .
                only process refarr which are within xarr limits.

                A middle point array is generated from the given
                ref scale i.e
                  xmid(i) = (refarr(i+1) - refarr(i) )/ 2

                The flux at ref(i) is then the total integral of
                xmid(i) <= xarr <= xmid(i+1)  divided by 
                (xmid(i) - xmid(i-1)). The integral is obtained 
                using the formula (trapezoidal rule)

                sum of( (yarr(k) + yarr(k-1)) * 
                        (xarr(k) - xarr(k-1) / 2 )

                for all k such that xmid(i) <= xarr(k) <= xmid(i+1)
                wavelength range.

                For the first element ref(0), use all the points in 
                xarr within the range ref(0) to xmid(0)
                For mid(0) < xarr(0) set mid(0) = xarr(0)

                The reference scale can be locally coarser or finer than the
                data. Both may be unequally sampled.

                The error for each refmid point is computed as follows.
                TBD

 Comment        
               
 Example        1) flux_conserv,lsanwav, lsanflx, lsanstdv, ref, yout, sout
                   new grid is given by lsanwav2 (the wavelength scale
                   of a second spectrum). By default linear interpolation is
                   performed.

 Dependencies   CALLS: KWD_NOTSET, POLYINT, SAP_ERROR
                CALLED FROM: SAP_REBIN

 Category       ISAP

 Filename       flux_conserv.pro

 Author         Iffat R. Khan (irk@ipac.caltech.edu)

 Version        0.3

 History        0.1  25-10-95  definition of requirements --> mazz@ipac
                0.1  02-11-95  Design and code            --> irk@ipac
                0.2  22-11-95  revised code               --> irk@ipac
                0.3  28-11-95  revised to correct range   --> irk@ipac
                0.4  02-05-95  revised to terminate if    --> irk@ipac
                               ref array is too short (<= 1 point)    
                1.0  24-05-96  change usage to FUNCTION from PRO   -->irk@ipac
                               to handle error status flag


;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
Copyright (C) 1995, California Institute of Technology.  ;
U.S. Government Sponsorship under NASA Contract NAS7-918 ;
is acknowledged.                                         ;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;