*     C06PUF Example Program Text
*     Mark 19 Release. NAG Copyright 1999.
*     .. Parameters ..
      INTEGER          NIN, NOUT
      PARAMETER        (NIN=5,NOUT=6)
      INTEGER          MMAX, NMAX, MNMAX
      PARAMETER        (MMAX=96,NMAX=96,MNMAX=MMAX*NMAX)
*     .. Local Scalars ..
      INTEGER          IFAIL, M, N
*     .. Local Arrays ..
      COMPLEX *16      WORK(MMAX+NMAX+MNMAX+30), X(MNMAX)
*     .. External Subroutines ..
      EXTERNAL         C06PUF, READX, WRITX
*     .. Executable Statements ..
      WRITE (NOUT,*) 'C06PUF Example Program Results'
*     Skip heading in data file
      READ (NIN,*)
   20 CONTINUE
      READ (NIN,*,END=40) M, N
      IF (M*N.GE.1 .AND. M*N.LE.MNMAX) THEN
         CALL READX(NIN,X,M,N)
         WRITE (NOUT,*)
         WRITE (NOUT,*) 'Original data values'
         CALL WRITX(NOUT,X,M,N)
         IFAIL = 1
*
*        -- Compute transform
         CALL C06PUF('F',M,N,X,WORK,IFAIL)
*
         IF (IFAIL.EQ.0) THEN
            WRITE (NOUT,*)
            WRITE (NOUT,*) 'Components of discrete Fourier transform'
            CALL WRITX(NOUT,X,M,N)
*
*           -- Compute inverse transform
            CALL C06PUF('B',M,N,X,WORK,IFAIL)
*
            WRITE (NOUT,*)
            WRITE (NOUT,*)
     +        'Original sequence as restored by inverse transform'
            CALL WRITX(NOUT,X,M,N)
            GO TO 20
         ELSE
            WRITE (NOUT,*)
            WRITE (NOUT,99999) ' ** C06PUF returned with IFAIL = ',
     +        IFAIL
         END IF
      ELSE
         WRITE (NOUT,*) ' ** Invalid value of M or N'
      END IF
   40 CONTINUE
*
99999 FORMAT (1X,A,I5)
      END
*
      SUBROUTINE READX(NIN,X,N1,N2)
*     Read 2-dimensional complex data
*     .. Scalar Arguments ..
      INTEGER          N1, N2, NIN
*     .. Array Arguments ..
      COMPLEX *16      X(N1,N2)
*     .. Local Scalars ..
      INTEGER          I, J
*     .. Executable Statements ..
      DO 20 I = 1, N1
         READ (NIN,*) (X(I,J),J=1,N2)
   20 CONTINUE
      RETURN
      END
*
      SUBROUTINE WRITX(NOUT,X,N1,N2)
*     Print 2-dimensional complex data
*     .. Scalar Arguments ..
      INTEGER          N1, N2, NOUT
*     .. Array Arguments ..
      COMPLEX *16      X(N1,N2)
*     .. Local Scalars ..
      INTEGER          I, J
*     .. Intrinsic Functions ..
      INTRINSIC        AIMAG, DBLE
*     .. Executable Statements ..
      DO 20 I = 1, N1
         WRITE (NOUT,*)
         WRITE (NOUT,99999) 'Real ', (DBLE(X(I,J)),J=1,N2)
         WRITE (NOUT,99999) 'Imag ', (AIMAG(X(I,J)),J=1,N2)
   20 CONTINUE
      RETURN
*
99999 FORMAT (1X,A,7F10.3,/(6X,7F10.3))
      END