program hpftest2
c
c This program calls the MPL library spfft. 2-dim example
c
       integer dimx, dimy
       parameter (dimx=512,dimy=1024)
c  NB: dimx and dimy must be multiples of nxproc and nyproc, respectively
       complex a(0:dimx-1, 0:dimy-1)
       integer ctrl
       integer i, trszx , trszy

c MPL subroutines
CMPF MPL  fft_init, sp_fft, fft_free

c Initialize
       forall (i=0:dimx-1,j=0:dimy-1) a(i,j) = cmplx(i,-j)

C       write(*,*) 'Give transform lengths (x and y):'
c      read(*,*) trszx, trszy
       trszx = 32
       trszy = 32

c      call fft_init (ctrl,mode,numdim,dimx,trszx,dimy,trszy,
c                     dimz,trszz)
       call fft_init (ctrl, 0, 2, dimx, trszx, dimy, trszy, 0, 0) 

       write(*,*) 'Input:'
       write(*,*) 'a(0,1) =', a(0,1)
       write(*,*) 'a(97,17) =', a(97,17)
       write(*,*) 'a(200,250) =', a(200,250)

c      Forward 2-D transform
       write(*,*) 'After FFT:'
       call sp_fft(a,1,ctrl)
       write(*,*) 'a(0,1) =', a(0,1)
       write(*,*) 'a(97,17) =', a(97,17)
       write(*,*) 'a(200,250) =', a(200,250)


c      Inverse 2-D transform
       call sp_fft(a,-1,ctrl)

       write(*,*) 'After inverse FFT:'
       write(*,*) 'a(0,1) =', a(0,1)
       write(*,*) 'a(97,17) =', a(97,17)
       write(*,*) 'a(200,250) =', a(200,250)

c      Free temporary data structures
       call fft_free(ctrl)

       end