!
! Discussion:
!
! In FORTRAN it is the FIRST dimension
! of the complex coefficient array that is "half" the size.
!
!
implicit none
integer ( kind = 4 ), parameter :: Nx = 5
integer ( kind = 4 ), parameter :: Ny = 3
integer ( kind = 4 ), parameter :: Nh = ( Nx / 2 ) + 1
include "fftw3.f"
integer ( kind = 4 ) i,j
real ( kind = 8 ) ux(Nx,Ny)
complex ( kind = 8 ) ukx(Nh,Ny),ukx_dum(Nh,Ny)
integer ( kind = 8 ) plan_forward,plan_backward
integer,parameter :: seed = 99999999
call srand(seed)
write ( *, '(a)' ) ' '
write ( *, '(a)' ) 'TEST04'
write ( *, '(a)' ) ' Demonstrate FFTW3 on a 2D real array'
write ( *, '(a,i8)' ) ' NX = ', nx
write ( *, '(a,i8)' ) ' NY = ', ny
write ( *, '(a)' ) ' '
write ( *, '(a)' ) ' Transform data to FFT coefficients.'
write ( *, '(a)' ) ' Backtransform FFT coefficients to recover'
write ( *, '(a)' ) ' the data.'
write ( *, '(a)' ) ' Compare recovered data to original data.'
!
! Compute the data.
!
do i = 1,Nh
do j = 1,Ny
ukx(i,j) = (1.0d0,0.0d0)*rand() + (0.0d0,1.0d0)*rand()
enddo
enddo
write ( *, '(a)' ) ' '
write ( *, '(a)' ) ' Input Data:'
write ( *, '(a)' ) ' '
do i = 1,Nh
do j = 1,Ny
write(*, *) i,j,ukx(i,j), "ukx"
enddo
enddo
!
! Make a plan for the FFT, and backward transform the data.
!
call dfftw_plan_dft_c2r_2d_ (plan_backward, Nx, Ny, ukx, ux, FFTW_ESTIMATE)
call dfftw_execute_ (plan_backward)
write ( *, '(a)' ) ' '
write ( *, '(a)' ) ' Output FFT Coefficients:'
write ( *, '(a)' ) ' '
do i = 1,Nx
do j = 1,Ny
ux(i,j) = ux(i,j)/(float(Nx)*float(Ny))
write(*, *) i,j,ux(i,j), "ux"
enddo
enddo
!
! Make a plan for the backward FFT, and recover the original data.
!
call dfftw_plan_dft_r2c_2d_ (plan_forward, Nx, Ny, ux, ukx_dum, FFTW_ESTIMATE)
call dfftw_execute_ (plan_forward)
write ( *, '(a)' ) ' '
write ( *, '(a)' ) ' Recovered input data divided by NX * NY:'
write ( *, '(a)' ) ' '
do i = 1,Nh
do j = 1,Ny
write(*, *) i,j,ukx_dum(i,j), "ukx_dum"
enddo
enddo
!
! Discard the information associated with the plans.
!
call dfftw_destroy_plan_ ( plan_forward )
call dfftw_destroy_plan_ ( plan_backward )
! return
end