!+ Subroutine to perform a least squares fit for an ESFT.
!
! Method:
!	A matrix equation is set up. This is solved using an SVD
!	decomposition.
!
! Current Owner of Code: J. M. Edwards
!
! History:
!	Version		Date			Comment
!	1.0		12-04-95		First Version under RCS
!						(J. M. Edwards)
!
! Description of Code:
!   FORTRAN 77  with extensions listed in documentation.
!
!- ---------------------------------------------------------------------
      SUBROUTINE LEAST_SQUARES(IERR, N_TERM, N_DATA
     &   , TRANS, WEIGHT, EXP_KU
     &   , A
     &   )
!
!
      IMPLICIT NONE
!
!
!     INCLUDE HEADER FILES.
      INCLUDE 'incl/prec_machine.ucf'
      INCLUDE 'incl/precision_pp.tcf'
      INCLUDE 'incl/dimensions_spec.ucf'
      INCLUDE 'incl/dimensions_pp.ucf'
      INCLUDE 'incl/error.pcf'
!
!     DUMMY ARGUMENTS.
      INTEGER 	!, INTENT(OUT)
     &     IERR                             
!             ERROR FLAG
      INTEGER 	!, INTENT(IN)
     &     N_TERM                           
!             NUMBER OF TERMS
     &   , N_DATA                           
!             NUMBER OF DATA POINTS
      REAL*8 	!, INTENT(IN)
     &     TRANS(NPD_AMOUNT+1)          
!             TRANSMISSIVITIES
     &   , WEIGHT(NPD_AMOUNT+1)         
!             WEIGHTINGS
     &   , EXP_KU(NPD_ESFT_TERM, NPD_AMOUNT+1)
!             PRE-CALCULATED EXPONENTIAL TERMS
      REAL*8 	!, INTENT(OUT)
     &     A(NPD_ESFT_TERM)                   
!             ESFT WEIGHTS
!
!     LOCAL ARGUMENTS.
      INTEGER
     &     I                                
!             LOOP VARIABLE
     &   , J                                
!             LOOP VARIABLE
     &   , K                                
!             LOOP VARIABLE
      REAL*8
     &     WMAX                             
!             MAXIMUM TERM OF DIAGONAL MATRIX
     &   , THRESHOLD                        
!             THRESHOLD FOR ZEROING
     &   , COEFF_MATRIX(NPD_ESFT_TERM, NPD_ESFT_TERM)
!             COEFFICIENT MATRIX
     &   , B(NPD_ESFT_TERM)                   
!             RHS OF EQUATION
     &   , W(NPD_ESFT_TERM)                   
!             DIAGONAL MATRIX
     &   , V(NPD_ESFT_TERM, NPD_ESFT_TERM)      
!             V-MATRIX OF SVD DECOMPOSITION
     &   , WRK(NPD_ESFT_TERM)                 
!             WORKING ARRAY
!
!     SUBROUTINES CALLED:
      EXTERNAL
     &     SVD_DECOMPOSE, BACK_SUBSTITUTE
!
!
!
!     SET THE COEFFICIENTS OF THE MATRICES COEFF_MATRIX AND B.
      DO I=1, N_TERM
         DO J=1, N_TERM
            COEFF_MATRIX(I, J)=0.0D+00
            DO K=1, N_DATA
               COEFF_MATRIX(I, J)=COEFF_MATRIX(I, J)
     &            +WEIGHT(K)*EXP_KU(J, K)*EXP_KU(I, K)
            ENDDO
         ENDDO
         B(I)=0.0D+00
         DO K=1, N_DATA
            B(I)=B(I)+WEIGHT(K)*TRANS(K)*EXP_KU(I, K)
         ENDDO
      ENDDO
!
!     MAKE THE SINGULAR VALUE DECOMPOSITION OF THE MATRIX.
      CALL SVD_DECOMPOSE(IERR
     &   , COEFF_MATRIX, N_TERM, N_TERM
     &   , NPD_ESFT_TERM, NPD_ESFT_TERM
     &   , W, V, WRK
     &   )
      IF (IERR.NE.I_NORMAL) RETURN
!
!     ZERO VERY SMALL TERMS OF W TO AVOID ILL-CONDITIONING.
      WMAX=0.0D+00
      DO J=1, N_TERM
         IF (W(J).GT.WMAX) THEN
            WMAX=W(J)
         ENDIF
      ENDDO
      THRESHOLD=TOL_CONDITION*WMAX
      DO J=1, N_TERM
         IF (W(J).LT.THRESHOLD) THEN
            W(J)=0.0D+00
         ENDIF
      ENDDO
!
!     SOLVE BY BACKSUBSTITUTION.
      CALL BACK_SUBSTITUTE(COEFF_MATRIX, W, V
     &   , N_TERM, N_TERM
     &   , NPD_ESFT_TERM, NPD_ESFT_TERM
     &   , B, A, WRK
     &   )
!
!
      RETURN
      END
