fft - fortran complex to real fftw issue -

i working on project need implement fourier transform , inverse transform. testing program modified online example; print or write commands debugging purposes:

  program testit   include 'fftw3.f'   double complex out!, in   real in   parameter (n=100)   dimension in(n), out(n)   integer*8 p,p2   integer i,j   real x   real fact   write(*,*)"stuff in data"   open(unit=12, file="input.txt", action="write", status="replace")   open(unit=20, file="dftoutput.txt", action="write", status="replace")   x=0   in = 0   i=1,n/2     in(i)=1   enddo   i=1,n   write(*,"(f10.2,1x,f10.2)")in(i)   write(12,*)real(in(i))   enddo   write(*,*)"create plans"   call dfftw_plan_dft_r2c_1d(p ,n,in,out,fftw_estimate)   call dfftw_plan_dft_c2r_1d(p2,n,in,out,fftw_estimate)   write(*,*)"do it"   call dfftw_execute_dft_r2c(p,in,out)   i=1,n     write(*,"(f12.4,1x,f12.4)")out(i)     write(20,*)abs(out(i))   enddo   write(*,*)"undo it"   call dfftw_execute_dft_c2r(p2,in,out)   fact=1.0/n   i=1,n     write(*,)in(i)     write(*,)out(i)   enddo   write(*,*)"clean up"   call dfftw_destroy_plan(p,in,out)   call dfftw_destroy_plan(p2,in,out)   end program 

the real complex transformation works fine. inverse transformation gives wrong values , somehow modifies both in , out variables. not know problem , not able figure out answers online. appreciated.

thanks in advance!

chad w. freer

edit: wondering if there similar function fftshift() , ifftshift() matlab in fftw package.

there precision issue : on computers, real refers single precision 32bit floating point numbers. real*8 can used specify double precision floating point numbers, consistent double complex , dfftw_....

for fftw real complex transform, if size of input n real*8, size of output n/2+1 double complex. since input real, coefficients of negative frequencies (higher n/2+1) conjugate of positive frequencies, , computation avoided.

according documentation of fftw regarding planner flags ,

fftw_preserve_input specifies out-of-place transform must not change input array. ordinarily default, except c2r , hc2r (i.e. complex-to-real) transforms fftw_destroy_input default...

if wish save coefficients in fourier space out, either copy array, or try use fftw_preserve_input. first solution seems best way go if have enough memory.

order of arguments fftw functions origin,destination. backward transform performed out in :

call dfftw_plan_dft_c2r_1d(p2,n,out,in,fftw_estimate) ... call dfftw_execute_dft_c2r(p2,out,in) 

here code based on yours, performs forward , backward fft. compiled gfortran main.f90 -o main -lfftw3 -lm -wall :

program testit include 'fftw3.f' double complex out!, in real*8 in parameter (n=100) dimension in(n), out((n/2+1)) integer*8 p,p2 integer real x real fact write(*,*)"stuff in data" open(unit=12, file="input.txt", action="write", status="replace") open(unit=20, file="dftoutput.txt", action="write", status="replace") x=0 in = 0 i=1,n/2   in(i)=1 enddo i=1,n   write(*,"(f10.2,1x/)")in(i)   write(12,*)real(in(i)) enddo write(*,*)"create plans" call dfftw_plan_dft_r2c_1d(p ,n,in,out,fftw_estimate) call dfftw_plan_dft_c2r_1d(p2,n,out,in,fftw_estimate) write(*,*)"do it" call dfftw_execute_dft_r2c(p,in,out) i=1,n/2+1   write(*,"(f12.4,1x,f12.4/)")out(i)   write(20,*)abs(out(i)) enddo write(*,*)"undo it" call dfftw_execute_dft_c2r(p2,out,in) fact=1.0/n write(*,*)"input back, except scaling" i=1,n   write(*,"(f10.2,1x)")in(i) enddo write(*,*)"output may modified c2r...except if fftw_preserve_input set" i=1,n/2+1   write(*,"(f10.2,1x,f10.2/)")out(i) enddo write(*,*)"clean up" call dfftw_destroy_plan(p,in,out) call dfftw_destroy_plan(p2,out,in) end program 

after forward fft in->out , backward fft out->in , in identical initial value, except scale factor n.