PRO line_centroid, x, y, cent, wrange=wrange, min=min, backg=backg, $
                   silent=silent, err=err

;+
; NAME
;
;     LINE_CENTROID
;
; EXPLANATION
;
;     Finds the centroid of an emission or absorption line by minmizing the
;     function SUM( (wavel_i-ref) * flux_i ), for the centroid wavelength,
;     ref. Note that the accuracy of the derived wavelength is improved by
;     having the input wavelength vector (str.wavel) on a fine scale.
;
; INPUTS
;
;     X     Wavelength vector
;
;     Y     Flux vector.
;
; OUTPUTS
;
;     CENT   The centroid.
;
; OPTIONAL INPUTS
;
;     WRANGE A wavelength range over which centroiding is applied. By default
;            it is set to the minimum and maximum values of str.wavel.
;
; KEYWORD
;
;     MIN    This needs to be set if you are finding the centroid for an
;            absorption line.
;
; USE
;
;     IDL> line_centroid, wvl, flx, [1030,1035]
;
; HISTORY
;
;     Ver. 1, 14-Mar-04, Peter Young
;     Ver. 2, 6-Apr-04, Peter Young
;        changed to give x, y as vector inputs, rather than structure
;     Ver. 3, 15-Apr-04, Peter Young
;        changed method so that the centroid is a correctly minimised
;        quantity
;        corrected bug with wrange
;     Ver. 4, 20-May-04, Peter Young
;        Added error check for X.
;     Ver. 5, 1-Jul-04, Peter Young
;        Removed error array, and added /silent keyword.
;-

IF n_params() LT 2 THEN BEGIN
  print,'Use:  IDL> line_centroid, xx, yy, cent, wrange=wrange, min=min'
  return
ENDIF

IF datatype(x) EQ 'STR' THEN BEGIN
  print,'%LINE_CENTROID: XX must be an array of wavelengths.'
  return
ENDIF


IF n_elements(backg) EQ 0 THEN backg=0.
IF n_elements(wrange) EQ 0 THEN wrange=[min(x),max(x)]

ind=where( (x GE wrange[0]) AND (x LE wrange[1]))

xx=x[ind]
yy=y[ind]-backg

yy=double(yy)

cent=total(yy^2*xx)/total(yy^2)

IF NOT keyword_set(silent) THEN BEGIN
  print,format='(f15.3)',cent
  plot,xx,yy,psym=10,/xsty
  oplot,cent*[1,1],[-1d10,1d10]
ENDIF

END