cffti1 Subroutine

subroutine cffti1(n, Wa, Ifac)

Uses

Arguments

Type IntentOptional Attributes Name
integer :: n
real(kind=rk) :: Wa
integer :: Ifac

Contents

Variables

Source Code

cffti1

Variables

Type Visibility Attributes Name Initial
real(kind=rk), public :: arg
real(kind=rk), public :: argh
real(kind=rk), public :: argld
real(kind=rk), public :: fi
integer, public :: i
integer, public :: i1
integer, public :: ib
integer, public :: ido
integer, public :: idot
integer, public :: ii
integer, public :: ip
integer, public :: ipm
integer, public :: j
integer, public :: k1
integer, public :: l1
integer, public :: l2
integer, public :: ld
integer, public :: nf
integer, public :: nl
integer, public :: nq
integer, public :: nr
integer, public :: ntry
integer, public, parameter, dimension(4) :: ntryh = [3, 4, 2, 5]
real(kind=rk), public, parameter :: tpi = 2.0_rk*acos(-1.0_rk)

Source Code

      subroutine cffti1(n,Wa,Ifac)
      use fftpack_kind
      implicit none
      real(rk) :: arg , argh , argld , fi , Wa
      integer :: i , i1 , ib , ido , idot , Ifac , ii , ip , ipm , j , k1, &
                 l1 , l2 , ld , n , nf , nl , nq , nr , ntry
      dimension Wa(*) , Ifac(*)  
      integer,dimension(4),parameter :: ntryh = [3 , 4 , 2 , 5]
      real(rk),parameter :: tpi = 2.0_rk * acos(-1.0_rk) ! 2 * pi
      nl = n
      nf = 0
      j = 0
 100  j = j + 1
      if ( j<=4 ) then
         ntry = ntryh(j)
      else
         ntry = ntry + 2
      endif
 200  nq = nl/ntry
      nr = nl - ntry*nq
      if ( nr/=0 ) goto 100
      nf = nf + 1
      Ifac(nf+2) = ntry
      nl = nq
      if ( ntry==2 ) then
         if ( nf/=1 ) then
            do i = 2 , nf
               ib = nf - i + 2
               Ifac(ib+2) = Ifac(ib+1)
            enddo
            Ifac(3) = 2
         endif
      endif
      if ( nl/=1 ) goto 200
      Ifac(1) = n
      Ifac(2) = nf
      argh = tpi/real(n, rk)
      i = 2
      l1 = 1
      do k1 = 1 , nf
         ip = Ifac(k1+2)
         ld = 0
         l2 = l1*ip
         ido = n/l2
         idot = ido + ido + 2
         ipm = ip - 1
         do j = 1 , ipm
            i1 = i
            Wa(i-1) = 1.0_rk
            Wa(i) = 0.0_rk
            ld = ld + l1
            fi = 0.0_rk
            argld = real(ld, rk)*argh
            do ii = 4 , idot , 2
               i = i + 2
               fi = fi + 1.0_rk
               arg = fi*argld
               Wa(i-1) = cos(arg)
               Wa(i) = sin(arg)
            enddo
            if ( ip>5 ) then
               Wa(i1-1) = Wa(i-1)
               Wa(i1) = Wa(i)
            endif
         enddo
         l1 = l2
      enddo
      end subroutine cffti1