* C06PXF Example Program Text. * Mark 19 Release. NAG Copyright 1999. * .. Parameters .. INTEGER NIN, NOUT PARAMETER (NIN=5,NOUT=6) INTEGER N1MAX, N2MAX, N3MAX, NMAX, LWORK PARAMETER (N1MAX=16,N2MAX=16,N3MAX=16, + NMAX=N1MAX*N2MAX*N3MAX,LWORK=N1MAX+N2MAX+N3MAX+ + NMAX+45) * .. Local Scalars .. INTEGER IFAIL, N, N1, N2, N3 * .. Local Arrays .. COMPLEX *16 WORK(LWORK), X(NMAX) * .. External Subroutines .. EXTERNAL C06PXF, READX, WRITX * .. Executable Statements .. WRITE (NOUT,*) 'C06PXF Example Program Results' * Skip heading in data file READ (NIN,*) 20 CONTINUE READ (NIN,*,END=40) N1, N2, N3 N = N1*N2*N3 IF (N.GE.1 .AND. N.LE.NMAX) THEN CALL READX(NIN,X,N1,N2,N3) WRITE (NOUT,*) WRITE (NOUT,*) 'Original data values' CALL WRITX(NOUT,X,N1,N2,N3) IFAIL = 0 * * -- Compute transform CALL C06PXF('F',N1,N2,N3,X,WORK,IFAIL) * WRITE (NOUT,*) WRITE (NOUT,*) 'Components of discrete Fourier transform' CALL WRITX(NOUT,X,N1,N2,N3) * * -- Compute inverse transform CALL C06PXF('B',N1,N2,N3,X,WORK,IFAIL) * WRITE (NOUT,*) WRITE (NOUT,*) + 'Original sequence as restored by inverse transform' CALL WRITX(NOUT,X,N1,N2,N3) GO TO 20 ELSE WRITE (NOUT,*) ' ** Invalid value of N1, N2 or N3' END IF 40 CONTINUE STOP END * SUBROUTINE READX(NIN,X,N1,N2,N3) * Read 3-dimensional complex data * .. Scalar Arguments .. INTEGER N1, N2, N3, NIN * .. Array Arguments .. COMPLEX *16 X(N1,N2,N3) * .. Local Scalars .. INTEGER I, J, K * .. Executable Statements .. DO 40 I = 1, N1 DO 20 J = 1, N2 READ (NIN,*) (X(I,J,K),K=1,N3) 20 CONTINUE 40 CONTINUE RETURN END * SUBROUTINE WRITX(NOUT,X,N1,N2,N3) * Print 3-dimensional complex data * .. Scalar Arguments .. INTEGER N1, N2, N3, NOUT * .. Array Arguments .. COMPLEX *16 X(N1,N2,N3) * .. Local Scalars .. INTEGER I, J, K * .. Intrinsic Functions .. INTRINSIC DBLE, DIMAG * .. Executable Statements .. DO 40 I = 1, N1 WRITE (NOUT,*) WRITE (NOUT,99998) 'z(i,j,k) for i =', I DO 20 J = 1, N2 WRITE (NOUT,*) WRITE (NOUT,99999) 'Real ', (DBLE(X(I,J,K)),K=1,N3) WRITE (NOUT,99999) 'Imag ', (DIMAG(X(I,J,K)),K=1,N3) 20 CONTINUE 40 CONTINUE RETURN * 99999 FORMAT (1X,A,7F10.3,/(6X,7F10.3)) 99998 FORMAT (1X,A,I6) END