module fft_mod
  use fft_fftw_mod

  implicit none

  real, private,save :: fract
!  real,allocatable,private    :: centering(:,:)

contains

  subroutine fftw_init(sax,rax)
    complex :: sax(:,:,:),rax(:,:,:)
    integer :: i1,i2

!    allocate( centering(size(sax,1),size(sax,2)))
    !! . . Scaling Factor
    fract=sqrt(real(size(sax,3))/real(size(sax)))
    
    !! . . Set up FFTW plans
    call fftw_init_ser(sax,rax)
  !  centering=fract

  end subroutine fftw_init

  !! ########################################3

  subroutine RAX_Forward(rax,id)
    complex  :: rax(:,:,:)
    integer  :: i1,i2,id
    call RAX_F_fft(rax,id)
    do i2=1,size(rax,2)
       do i1=1,size(rax,1)
          rax(i1,i2,id)=rax(i1,i2,id)*fract !centering(i1,i2)
       end do
    end do

  end subroutine RAX_Forward

  !! ########################################3

  subroutine RAX_Backward(rax,id)
    complex      :: rax(:,:,:)
    integer       :: i1,i2,id
    call RAX_B_fft(rax,id)

    do i2=1,size(rax,2)
       do i1=1,size(rax,1)
          rax(i1,i2,id)=rax(i1,i2,id)*fract !centering(i1,i2)
       end do
    end do

  end subroutine RAX_Backward

  !! ########################################3
  
  subroutine SAX_Forward(sax,id)
    complex  :: sax(:,:,:)
    integer  :: i1,i2,id
    call SAX_F_fft(sax,id)

    do i2=1,size(sax,2)
       do i1=1,size(sax,1)
          sax(i1,i2,id)=sax(i1,i2,id)*fract !centering(i1,i2)
       end do
    end do

  end subroutine SAX_Forward

  !! ########################################3

  subroutine SAX_Backward(sax,id)
    complex      :: sax(:,:,:)
    integer       :: i1,i2,id
    call SAX_B_fft(sax,id)

    do i2=1,size(sax,2)
       do i1=1,size(sax,1)
          sax(i1,i2,id)=sax(i1,i2,id)*fract !centering(i1,i2)
       end do
    end do

  end subroutine SAX_Backward

  !! ########################################3

!  subroutine fftw_deallocate()
!    deallocate( centering )
!  end subroutine fftw_deallocate

end module fft_mod
