File: HLSE.FT of Disk: Disks/MyPDP/m8-2-rka1-rkb1

(Source file text) 

---

      SUBROUTINE HLSE(A,X,M,N,L,IER,AUX,IPIV) 
C 
C        PURPOSE
C           TO CALCULATE THE UNNORMALISED COVARIANCE MATRIX OF
C           PARAMETERS AFTER A CALL TO HLS
C 
C        USAGE
C           CALL HLSE (A,X,M,N,L,IER,AUX,IPIV)
C 
C        DESCRIPTION OF PARAMETERS
C           AS BEFORE WITH MATRIX X CONTAINING THE COVARIANCE MATRIX. 
C 
C        REMARKS
C           (1) WITH IER=0, HLSE CAN BE USED INDEPENDENTLY OF HLS 
C 
C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED
C           HHSTEP, HHUK
C 
C        METHOD 
C           CALCULATES THE INVERSE OF (A-TRANSPOSED.A) DIRECTLY 
C           FROM THE UPPER TRIANGULAR MATRIX LEFT AFTER PERFORMING
C           THE HH TRANSFORMATION TO A (SEC 6 OF GOLUB).
C 
C     ..................................................................
C 
      DIMENSION A(2),X(2),AUX(2),IPIV(2)
C 
C 
      K1 = MAX0(N,L)
      IF(IABS(IER).EQ.N) GO TO 5
C 
C     COMPLETE HOUSEHOLDER TRANSFORMATION IF IER LESS THAN N
      IST = IER*(M+1) + 1 
      KS = IER + 1
      DO 4 K=KS,N 
      LV = M - K + 1
C      GENERATION OF VECTOR UK AND PARAMETER BETA 
      CALL HHUK(A(IST),1,LV,SIG,BETA) 
C 
      J = K1 + K
      AUX(J) = -SIG 
      IF (K.EQ.N)  GO TO 4
C     TRANSFORMATION OF MATRIX A
      NK = N - K
      IA = IST + M
      CALL HHSTEP(A(IST),A(IA),1,M,LV,NK,BETA)
C 
      IST = IST + M + 1 
    4 IPIV(K) = K 
C 
C     CALCULATE X  FROM R AND DIAGONAL ELEMENTS OF R-T
    5 DO 40 K=1,N 
      KN=N-K+1
      IA = (KN-1)*(M+1) + 1 
      IX = (KN-1)*(N+1) + 1 
      IEND = (N-1)*M + KN + 1 
      IK = IX 
      IL = IX 
      DO 20 J=K,N 
      JA = IA + M 
      JN=N-J+1
      II=JN+K1
      PIV=1./AUX(II)
      ID=IPIV(JN)-JN
      H=0.
      IF(J.EQ.K) H=PIV
      II=IK 
      IEND = IEND - 1 
      IF(J.EQ.1)  GO TO 15
      DO 10 I=JA,IEND,M 
      II=II+N 
   10 H=H-A(I)*X(II)
      II=IK+ID*N
   15 H=H*PIV 
      X(IL)=H 
      X(IK)=X(II) 
      X(II)=H 
      IK=IK-N 
      IL=IL-1 
   20 IA = IA - M - 1 
C 
C     INTERCHANGE COLUMNS ALREADY FINISHED
      ID=IPIV(KN)-KN
      IF (ID.EQ.0)  GO TO 40
      II=KN+ID
      JA=KN 
      DO 30 J=1,N 
      H=X(II) 
      X(II)=X(JA) 
      X(JA)=H 
      JA=JA+N 
   30 II=II+N 
C     INTERCHANGE REMAINING PART OF THE ROWS
      IF (K.EQ.N)  GO TO 40 
      IA=(KN-1)*( N+1)+1
      KA=IA-KN+1
      KE=IA-1 
      DO 35 J=KA,KE 
      JD= J+ID*N
      H = X(JD) 
      X(JD) = X(J)
   35 X(J) = H
   40 CONTINUE
      RETURN
      END 


---