      SUBROUTINE BINMRQMI(COVAR,ALPHA,ALAMDA,REPEAT)
C----------------------------------------------------------------------
C This is adapted from MRQMIN in Numerical Recipes, p. 527, for BINFIT.
C     J.T. Armstrong      7 Feb 1991
C----------------------------------------------------------------------

      INCLUDE 'BINFIT.INC'

      REAL*8    ALPHA(MXPARM,MXPARM), BETA(MXPARM), COVAR(MXPARM,MXPARM)
      REAL*8    DA(MXPARM), DAA(MXPARM,1)
      REAL*8    ALAMDA
      INTEGER*4 KK, IHIT
      LOGICAL*1 REPEAT

      EQUIVALENCE ( DA(1), DAA(1,1) )

C Initialization
D     WRITE(OUTC,*) ' Entering BINMRQMI'
      IF (ALAMDA.LT.0.) THEN
          KK=NVPARTOT+1
          DO J=1,MXPARM
              IHIT=0
              DO K=1,NVPARTOT
                  IF(LISTA(K).EQ.J) IHIT=IHIT+1
              END DO
              IF(IHIT.EQ.0) THEN
                  LISTA(KK)=J
                  KK=KK+1
              ELSE IF (IHIT.GT.1) THEN
                  WRITE(OUTC,'(A,I3)') 
     1                ' Improper permutation in LISTA; IHIT = ',IHIT
                  RETURN
              END IF
          END DO
          IF (KK.NE.(MXPARM+1)) THEN
              WRITE(OUTC,'(A,I3,A,I3)') 
     1            ' Improper permutation in LISTA; KK = ',KK,
     1            '; MXPARM+1 = ',MXPARM+1
              RETURN
          END IF
          ALAMDA=0.001
D         WRITE(OUTC,*) ' Callin BINMRQCO (initialization)'
          CALL BINMRQCO(ALPHA,BETA)
          TCHISQ0=TCHISQ
          DO J=1,NFILT
              FCHISQ0(J)=FCHISQ(J)
          END DO
          DO J=1,MXPARM
              PARMSV1D(J)=PARM1D(J)
          END DO  
      ENDIF
C Alter the linearized fitting matrix by augmenting diagonal elements
      REPEAT=.FALSE.
      DO J=1,NVPARTOT
          DO K=1,NVPARTOT
              COVAR(J,K)=ALPHA(J,K)
          END DO
          COVAR(J,J)=ALPHA(J,J)*(1.+ALAMDA)
          DAA(J,1)=BETA(J)       ! Because CGAUSSJ wants a 2-dim array here
      END DO  
D     WRITE(OUTC,*) ' Calling CGAUSSJ'
      CALL CGAUSSJ(COVAR,NVPARTOT,MXPARM,DAA,1,1)    ! Matrix solution
      DO J=1,NVPARTOT
          DELTAP1D(LISTA(J))=DAA(J,1)
C         WRITE(OUTC,'(2(A,I1,A,G11.4))') 
C    1        ' PARMSV1D(LISTA(',J,')) = ',PARMSV1D(LISTA(J)),
C    1        '; DELTAP1D(LISTA(',J,')) = ',DELTAP1D(LISTA(J))
      END DO
      IF(ALAMDA.EQ.0.) THEN           ! Convergence
          CALL BINCVSRT(COVAR,MXPARM,MXPARM,LISTA,NVPARTOT)
          RETURN
      ENDIF
      DO J=1,NVPARTOT
          PARM1D(LISTA(J))=PARMSV1D(LISTA(J))+DELTAP1D(LISTA(J))
      END DO  
D     WRITE(OUTC,*) ' Calling BINMRQCO'
      CALL BINMRQCO(COVAR,DA)
      IF(TCHISQ.LE.TCHISQ0)THEN
          ALAMDA=0.1*ALAMDA
          TCHISQ0=TCHISQ
          DO J=1,NFILT
              FCHISQ0(J)=FCHISQ(J)
          END DO
          DO J=1,NVPARTOT
              DO K=1,NVPARTOT
                  ALPHA(J,K)=COVAR(J,K)
              END DO  
              BETA(J)=DA(J)
              PARMSV1D(LISTA(J))=PARM1D(LISTA(J))
          END DO  
      ELSE
          ALAMDA=10.*ALAMDA
D         WRITE(OUTC,'(A,2(F12.9,A))') ' BINMRQMI: TCHISQ = ',TCHISQ,
D    1       '; TCHISQ0 = ',TCHISQ0,'; try smaller step.'
          TCHISQ=TCHISQ0
          REPEAT=.TRUE.
          DO J=1,NFILT
              FCHISQ(J)=FCHISQ0(J)
          END DO
          DO J=1,MXPARM
              PARM1D(J)=PARMSV1D(J)
          END DO
      ENDIF

      RETURN
      END

