Function:     LSQFIT

Purpose:      LSQFIT is a routine for making a least-squares fit of a
              function to a set of data points. The method used is
              described in: Marquardt, J.Soc.Ind.Appl.Math. 11, 431 (1963).
              This  method is a mixture of the steepest descent method and
              the Taylor method.

Category:     MATH

File:         lsqfit.c

Author:       K.G. Begeman

Use:          INTEGER LSQFIT( XDAT ,    Input      REAL ARRAY
                              XDIM ,    Input      INTEGER
                              YDAT ,    Input      REAL ARRAY
                              WDAT ,    Input      REAL ARRAY
                              NDAT ,    Input      INTEGER
                              FPAR ,   In/Output   REAL ARRAY
                              EPAR ,    Output     REAL ARRAY
                              MPAR ,    Input      INTEGER ARRAY
                              NPAR ,    Input      INTEGER
                              TOL  ,    Input      REAL
                              ITS  ,    Input      INTEGER
                              LAB  ,    Input      REAL
                              FOPT )    Input      INTEGER

              LSQFIT   Returns number of iterations needed to achieve
                       convergence according to TOL. When this
                       number is negative, the fitting was not
                       continued because a fatal error occurred:
                       -1 Too many free parameters, maximum is 32.
                       -2 No free parameters.
                       -3 Not enough degrees of freedom.
                       -4 Maximum number of iterations too small to
                          obtain a solution which satisfies TOL.
                       -5 Diagonal of matrix contains elements which
                          are zero.
                       -6 Determinant of the coefficient matrix is zero.
                       -7 Square root of negative number.
              XDAT     Contains coordinates of data points.
                       XDAT is two-dimensional: XDAT(XDIM,NDAT)
              XDIM     Dimension of fit.
              YDAT     Contains data points.
              WDAT     Contains weigths for data points.
              NDAT     Number of data points.
              FPAR     On input contains initial estimates of the
                       parameters for non-linear fits, on output the
                       fitted parameters.
              EPAR     Contains estimates of errors in fitted
                       parameters.
              MPAR     Logical mask telling which parameters are
                       free (MPAR(J)=non-zero) and which parameters
                       are fixed (MPAR(J)=0).
              NPAR     Number of parameters (free+fixed).
              TOL      Relative tolerance. LSQFIT stops when
                       successive iterations fail to produce a
                       decrement in reduced chi-squared less than
                       TOL. If TOL is less than the minimum tolerance
                       possible, TOL will be set to this value. This
                       means that maximum accuracy can be obtained by
                       setting TOL=0.0.
              ITS      Maximum number of iterations.
              LAB      Mixing parameter, LAB determines the initial
                       weight of steepest descent method relative to
                       the Taylor method. LAB should be a small
                       value (i.e. 0.01). LAB can only be zero when
                       the partial derivatives are independent of
                       the parameters. In fact in this case LAB
                       should be exactly equal to zero.
              FOPT     A value which is passed unmodified to FUNC
                       and DERV (see below).

Notes:        The following routines have to be defined by the user:
              
              REAL FUNC( XDAT ,   Input    REAL ARRAY
                         FPAR ,   Input    REAL ARRAY
                         NPAR ,   Input    INTEGER
                         FOPT )   Input    INTEGER

              FUNC    Returns the function value of the function to
                      be fitted.
              XDAT    Coordinate(s) of data point.
              FPAR    Parameter list.
              NPAR    Number of parameters.
              FOPT    A user defined option.
              
              CALL DERV( XDAT ,   Input    REAL ARRAY
                         FPAR ,   Input    REAL ARRAY
                         DPAR ,   Output   REAL ARRAY
                         NPAR ,   Input    INTEGER
                         FOPT )   Input    INTEGER

              XDAT    Coordinate(s) of data point.
              FPAR    Parameter list.
              EPAR    Partial derivatives to the parameters of the
                      function to be fitted.
              NPAR    Number of parameters.
              FOPT    A user defined option.

              In FORTRAN applications you need to specify the
              following f2cvv syntax for the C to Fortran
              interface somewhere in your source code:
              
              C@ real function func( real, real, integer, integer )
              C@ subroutine derv( real, real, real, integer, integer )
              
Example:      Fitting y(x) = a + b * x 

              REAL FUNCTION FUNC( XDAT, FPAR, NPAR, FOPT )
              ...
              FUNC = FPAR(1) + FPAR(2) * XDAT
              RETURN
              END
              
              SUBROUTINE DERV( XDAT, FPAR, DPAR, NPAR, FOPT )
              ...
              DPAR(1) = 1.0
              DPAR(2) = XDAT
              RETURN
              END

Updates:      May  7, 1990: KGB, Document created.
              May 14, 1990: MXV, Document refereed.