;             NAME : mp_dgf_iris
; 
;          PURPOSE : Perform constained/guided double Gaussian fit, see an example in Tian et al. (2011, ApJ, 738, 18)
;
; CALLING SEQUENCE : re = mp_dgf_iris(x, y, e, fit0, fit1c, fit1a, dlambda, good, /double)
; 
;           INPUTS : y - the line profile
;                    
;                    x - wavelength vector
;                    
;                    ee - error vector  
;                    
;                    fit0 - 4-element vecotr containing single Gaussian fit parameters, [background, peak intensity, centroid, Doppler width]), 
;                           used as guess of the parameters for the major (1st) component 
;                    
;                    fit1c - guess of the line centroid of the 2nd component
;                    
;                    fit1a - guess of the peak intensity of the 2nd component 
;                    
;                    dlambda - wavelength pixel size
;                    
;                    good - indices of wavelength pixels at which the intensities are above zero
; 
;         KEYWORDS : double - double precision     
;         
;           OUTPUTS: a structure containing the fitting parameters, tags are descibed in the following:
;                    b: background
;                    i1: peak intensity of the 1st component
;                    p1: centriod position of the 1st component
;                    w1: Doppler width (1/e) of the 1st component
;                    i2: peak intensity of the 2nd component
;                    p2: centriod position of the 2nd component
;                    w2: Doppler width (1/e) of the 2nd component
;                    fit: fiting result
;          
;          EXAMPLES: see an example of usage in dgf_rb_1lp.pro
;
;          HISTORY : Written by Hui Tian at CfA, April 4, 2013
;          
;      RESTRICTIONs: The allowed range of different parameters are fixed in this code. These numbers can be changed for the purpose of general study (see dgf_1lp.pro). 
;          
;; **********************************************************************


function mp_dgf_iris, x, y, e, fit0, fit1c, fit1a, dlambda, good, double = double

; f = x0 + x1*exp((x-x2)/x3)^2 + x4*exp((x-x5)/x6)^2
xlimited=[1,1]
; Background Intensity
x0={limited:xlimited, limits:fit0[0] + 0.5*[-1,1]*fit0[0], value:fit0[0]}  ;0.25

; First Component Peak Intensity
x1={limited:xlimited, limits:fit0[1] + 0.5*[-1,1]*fit0[1], value:fit0[1]}  ;0.25

; First Component Line Center Position
x2={limited:xlimited, limits:fit0[2] + 2*[-dlambda,dlambda], value:fit0[2]} ;2*;0.5*

; First Component exponential Width
x3={limited:xlimited, limits:fit0[3] + 2*[-dlambda,dlambda], value:fit0[3]} ;2*;0.5*

; Second Component Peak Intensity
l =[0.5*fit1a, 1.5*fit1a] ;[0.8*fit1a, 1.2*fit1a] ;[0.25*fit1a, 5*fit1a]; 
x4={limited:xlimited, limits:l, value:fit1a}

; Second Component Line Center Position
x5={limited:xlimited, limits:fit1c + 3*[-dlambda , dlambda], value:fit1c}  ;3*;0.5*

; Second Component exponential Width
x6={limited:xlimited, limits:fit0[3] + 3*[-dlambda,dlambda], value:fit0[3]}  ;3*

parinfo=[x0,x1,x2,x3,x4,x5,x6]
param = parinfo.value

res = mpfitfun('doublegauss', x[good], y[good], e[good], param, parinfo=parinfo, $
	maxiter = 10000, dof = dof, bestnorm = bestnorm, yfit = yfit, double = double,$
	status = status, /quiet, perror=perror)

;PCERROR = PERROR * SQRT(BESTNORM / DOF) 
;help, status

result = {b:res[0], i1:res[1], p1:res[2], w1:res[3], i2:res[4], p2:res[5], w2:res[6], fit:yfit, status:status}

return, result

end