C
C-----------------------------------------------------------------------
C  MAIN PROGRAM:  TEST PROGRAM FOR LPTRN
C  AUTHORS:       A. H. GRAY, JR.  AND J. D. MARKEL
C       GRAY -    UNIV. OF CALIF., SANTA BARBARA, CA 93106
C       BOTH -    SIGNAL TECHNOLOGY INC., 15 W. DE LA GUERRA AVE.
C                 SANTA BARBARA, CA 93101
C-----------------------------------------------------------------------
C
      DIMENSION A(21), C(21), R(21), RC(21), ALAR(21)
      DIMENSION AREA(21)
      COMMON IOUTD
      IOUTD = I1MACH(2)
      DO 10 J=4,21
        A(J) = 0.
  10  CONTINUE
      A(1) = 1.
      A(2) = -.45
      A(3) = .81
      M = 3
      ALPHA = 10.
      DO 30 I=1,6
        CALL LPTRN(I, M, A, C, R, RC, ALAR, AREA, ALPHA)
        MP = M + 1
        WRITE (IOUTD,9999) I, ALPHA
9999    FORMAT (/I10, F11.6/)
        DO 20 J=1,MP
          WRITE (IOUTD,9998) A(J), C(J), R(J), RC(J), ALAR(J), AREA(J)
  20    CONTINUE
  30  CONTINUE
9998  FORMAT (5X, 6F11.6)
      STOP
      END
C
C-----------------------------------------------------------------------
C SUBROUTINE:  LPTRN
C THIS SUBROUTINE CARRIES OUT THE TRANSFORMATIONS BETWEEN
C VARIOUS PARAMETER SETS USED IN LINEAR PREDICTION
C-----------------------------------------------------------------------
C
      SUBROUTINE LPTRN(I, M, A, C, R, RC, ALAR, AREA, ALPHA)
C
C  INPUTS:  I  -  VARIABLE IDENTIFYING WHICH ARE INPUTS
C           M  -  FILTER ORDER (M<51, SEE NOTE BELOW*)
C
C     I=1     INPUT=A(.) & ALPHA
C     I=2     INPUT=C(.)
C     I=3     INPUT=R(.)
C     I=4     INPUT=RC(.) & R(1)
C     I=5     INPUT=ALAR(.) & R(1)
C     I=6     INPUT=AREA(.) & R(1)
C
C        A(.)    =  FILTER COEF.
C        ALPHA   =  GAIN
C        C(.)    =  CEPSTRAL COEF.
C        R(.)    =  AUTOCORRELATION COEF.
C        R(1)    =  FIRST AUTO. COEF. (ENERGY)
C        RC(.)   =  REFLECTION COEF.
C        ALAR(.) =  LOG AREA RATIOS
C        AREA(.) =  AREA FUNCTIONS
C
C    * PROGRAM LIMITED TO M=50 BY DIMENSION SA(50),
C      A TEMPORARY STORAGE FOR FILTER COEFFICIENTS
C
      DIMENSION A(1), C(1), R(1), RC(1), ALAR(1), AREA(1)
      DIMENSION SA(50)
      COMMON IOUTD
C
C TEST FOR M OUT OF RANGE
C
      MTEST = (M-1)*(50-M)
      IF (MTEST) 340, 340, 10
  10  MP = M + 1
      IF (I-2) 50, 20, 50
C
C ..GENERATES A(.) ,ALPHA, FROM C(.)
C
  20  ALPHA = EXP(C(1))
      A(1) = 1.
      DO 40 K=1,M
        KP = K + 1
        SUM = 0.
        DO 30 J=1,K
          JB = K - J + 2
          SUM = SUM + A(J)*C(JB)*FLOAT(JB-1)
  30    CONTINUE
        A(KP) = -SUM/FLOAT(K)
  40  CONTINUE
  50  GO TO (110, 110, 160, 60, 200, 220), I
C
C ..GENERATES SA(.),R(.), & ALPHA FROM RC(.) & R(1)
C
  60  DO 70 J=1,M
        SA(J) = RC(J)
  70  CONTINUE
      R(2) = -RC(1)*R(1)
      ALPHA = R(1)*(1.-RC(1)*RC(1))
      DO 100 J=2,M
        MH = J/2
        Q = RC(J)
        ALPHA = ALPHA*(1.-Q*Q)
        DO 80 K=1,MH
          KB = J - K
          AT = SA(K) + Q*SA(KB)
          SA(KB) = SA(KB) + Q*SA(K)
          SA(K) = AT
  80    CONTINUE
        SUM = 0.
        DO 90 L=1,J
          LB = J + 1 - L
          SUM = SUM + SA(L)*R(LB)
  90    CONTINUE
        R(J+1) = -SUM
 100  CONTINUE
      IF (I-4) 240, 240, 260
 110  DO 120 J=1,M
        SA(J) = A(J+1)
 120  CONTINUE
C
C ..GENERATES RC(.),R(1), FROM A(.) & ALPHA
C
      DO 130 J=1,M
        RC(J) = SA(J)
 130  CONTINUE
      ALT = ALPHA
      DO 150 J=2,M
        JB = M + 2 - J
        MH = JB/2
        RCT = RC(JB)
        DR = 1./(1. - RCT*RCT)
        DO 140 K=1,MH
          KB = JB - K
          Q = (RC(K)-RCT*RC(KB))*DR
          RC(KB) = (RC(KB)-RCT*RC(K))*DR
          RC(K) = Q
 140    CONTINUE
        ALT = ALT*DR
 150  CONTINUE
      R(1) = ALT/(1.-RC(1)*RC(1))
      GO TO 60
C
C ..GENERATES RC(.),SA(.),& ALPHA FROM R(.)
C
 160  RC(1) = -R(2)/R(1)
      SA(1) = RC(1)
      ALPHA = R(1)*(1.-RC(1)*RC(1))
      DO 190 J=2,M
        MH = J/2
        JM = J - 1
        Q = R(J+1)
        DO 170 L=1,JM
          LB = J + 1 - L
          Q = Q + SA(L)*R(LB)
 170    CONTINUE
        Q = -Q/ALPHA
        RC(J) = Q
        DO 180 K=1,MH
          KB = J - K
          AT = SA(K) + Q*SA(KB)
          SA(KB) = SA(KB) + Q*SA(K)
          SA(K) = AT
 180    CONTINUE
        SA(J) = Q
        ALPHA = ALPHA*(1.-Q*Q)
 190  CONTINUE
      GO TO 240
C
C ..GENERATES RC(.) & AREA(.) FROM ALAR(.)
C
 200  AREA(MP) = 1.
      DO 210 J=1,M
        JB = M + 1 - J
        AR = EXP(ALAR(JB))
        RC(JB) = (1.-AR)/(1.+AR)
        AREA(JB) = AREA(JB+1)/AR
 210  CONTINUE
      GO TO 60
C
C ..GENERATES ALRAR(.) & RC(.) FROM AREA(.)
C
 220  DO 230 J=1,M
        AR = AREA(J+1)/AREA(J)
        ALAR(J) = ALOG(AR)
        RC(J) = (1.-AR)/(1.+AR)
 230  CONTINUE
      GO TO 60
C
C ..GENERATES AREA(.) & ALAR(.) FROM RC(.)
C
 240  AREA(MP) = 1.
      DO 250 J=1,M
        JB = M + 1 - J
        AR = (1.-RC(JB))/(1.+RC(JB))
        ALAR(JB) = ALOG(AR)
        AREA(JB) = AREA(JB+1)/AR
 250  CONTINUE
      IF (I-2) 280, 310, 260
 260  DO 270 J=2,MP
        A(J) = SA(J-1)
 270  CONTINUE
      A(1) = 1.
C
C ..GENERATE C(.) FROM A(.) & ALPHA
C
 280  C(1) = ALOG(ALPHA)
      C(2) = -A(2)
      DO 300 L=2,M
        LP = L + 1
        SUM = FLOAT(L)*A(LP)
        DO 290 J=2,L
          JB = L - J + 2
          SUM = SUM + A(J)*C(JB)*FLOAT(JB-1)
 290    CONTINUE
        C(LP) = -SUM/FLOAT(L)
 300  CONTINUE
 310  DO 330 J=1,M
        IF (ABS(RC(J))-1.) 330, 320, 320
 320    WRITE (IOUTD,9999)
 330  CONTINUE
9999  FORMAT (19H FILTER IS UNSTABLE)
      RETURN
 340  WRITE (IOUTD,9998) M
9998  FORMAT (3H M=, I6, 13H OUT OF RANGE)
      RETURN
      END