C....*...1.........2.........3.........4.........5.........6.........7.*.......8
C     ELAST1       7/14/81
C
C     PURPOSE
C     COMPUTE SUBSTITUTION AND PRICE ELASTICITIES AND THEIR STANDARD
C     ERRORS FROM ESTIMATED COEFFICIENTS OF THE FOURIER COST FUNCTION.
C
C     USAGE
C     CALL FFFCGH(N,KC,IS,KA,IAA,JJA,M,DL,X,CGH,LT,IW)
C     CALL ELAST1(M,LT,CGH,THETA,VAR,SUB,SESUB,PRI,SEPRI,WORK)
C
C     SUBROUTINES CALLED
C     DGMPRD
C
C     ARGUMENTS
C     M     - AS FOR FFFCGH, INPUT.  SET EQUAL TO THE TOTAL NUMBER OF
C             FACTORS.
C     LT    - AS FOR FFFCGH, INPUT.
C     CGH   - AS FOR FFFCGH, INPUT.  NOTE THAT FFFCGH IS CALLED WITH M
C             EQUAL TO THE NUMBER OF FACTOR PRICES NOT THE NUMBER OF
C             FACTOR PRICES LESS ONE.
C     THETA - ESTIMATED COEFFICIENTS OF THE FOURIER COST FUNCTION,
C             INPUT VECTOR OF LENGTH LT.
C             REAL*8
C     VAR   - EXTIMATED VARIANCE-COVARIANCE MATRIX OF THETA, INPUT.
C             MATRIX OF ORDER LT BY LT STORED COLUMNWISE (STORAGE MODE
C             0).
C             REAL* 8
C     SUB   - ESTIMATED ELASTICITIES OF SUBSTITUTION, OUTPUT.  MATRIX OF
C             ORDER M BY M STORED COLUMNWISE (STORAGE MODE 0).
C             REAL*8
C     SESUB - ESTIMATED STANDARD ERRORS OF SUB, OUTPUT.  STORED THE SAME
C             AS SUB.
C             REAL*8
C     PRI   - ESTIMATD PRICE ELASTICITIES, OUTPUT.  MATRIX OF ORDER M
C             BY M STORED COLUMNWISE (STORAGE MODE 0).  ROWS INDEX
C             QUANTITIES AND COLUMNS INDEX PRICES.
C             REAL*8
C     SEPRI - ESTIMATED STANDARD ERRORS OF PRI, OUTPUT.  STORED THE SAME
C             AS SUB.
C             REAL*8
C     WORK  - A WORK VECTOR OF LENGTH 4*LT+1.
C             REAL*8
C
      SUBROUTINE ELAST1(N,LT,CGH,THETA,VAR,SUB,SESUB,PRI,SEPRI,WORK)
      IMPLICIT REAL*8 (A-H,O-Z)
      save
      REAL*8 CGH(1),THETA(1),VAR(LT,LT)
      REAL*8 SUB(N,N),SESUB(N,N),PRI(N,N),SEPRI(N,N),WORK(1)
      INTEGER*4 C0,G0,H0,GI0,GJ0,HIJ0
      C0=0
      G0=C0+1
      H0=G0+N
      GI0=0
      GJ0=GI0+LT
      HIJ0=GJ0+LT
      IE0=HIJ0+LT
      IV0=IE0+1
      IW0=IV0+LT
      LIMIT=IW0+LT
      LCGH=1+N+N*N
      DO 100 I=1,N
      DO 100 J=1,N
      DO 10  K=1,LT
      WORK(GI0+K)=CGH(G0+I+LCGH*(K-1))
      WORK(GJ0+K)=CGH(G0+J+LCGH*(K-1))
10    WORK(HIJ0+K)=CGH(H0+N*(J-1)+I+LCGH*(K-1))
      CALL DGMPRD(WORK(GI0+1),THETA,WORK(IE0+1),1,LT,1)
      GI=WORK(IE0+1)
      CALL DGMPRD(WORK(GJ0+1),THETA,WORK(IE0+1),1,LT,1)
      GJ=WORK(IE0+1)
      CALL DGMPRD(WORK(HIJ0+1),THETA,WORK(IE0+1),1,LT,1)
      HIJ=WORK(IE0+1)
      SUB(I,J)=1.D0+HIJ/(GI*GJ)
      IF(I.EQ.J) SUB(I,J)=SUB(I,J)-1.D0/GI
      DO 20 K=1,LT
      WORK(IV0+K)=WORK(HIJ0+K)/(GI*GJ)
     &           -WORK(GI0+K)*HIJ/(GI*GI*GJ)
     &           -WORK(GJ0+K)*HIJ/(GI*GJ*GJ)
20    IF(I.EQ.J) WORK(IV0+K)=WORK(IV0+K)
     &           +WORK(GI0+K)/(GI*GI)
      CALL DGMPRD(VAR,WORK(IV0+1),WORK(IW0+1),LT,LT,1)
      CALL DGMPRD(WORK(IW0+1),WORK(IV0+1),WORK(IE0+1),1,LT,1)
      SESUB(I,J)=DSQRT(WORK(IE0+1))
      PRI(I,J)=HIJ/GI+GJ
      IF(I.EQ.J) PRI(I,J)=PRI(I,J)-1.D0
      DO 30 K=1,LT
30    WORK(IV0+K)=WORK(HIJ0+K)/GI
     &           -WORK(GI0+K)*HIJ/(GI*GI)
     &           +WORK(GJ0+K)
      CALL DGMPRD(VAR,WORK(IV0+1),WORK(IW0+1),LT,LT,1)
      CALL DGMPRD(WORK(IW0+1),WORK(IV0+1),WORK(IE0+1),1,LT,1)
      SEPRI(I,J)=DSQRT(WORK(IE0+1))
100   CONTINUE
      RETURN
      END