SF_FFT

Description

SciFortran module for Fourier transform

Quick access

Used modules

sf_integrate: SciFortran module for function integration
- simps()
sf_arrays: SciFortran module for array creation and manipulation
- linspace()
sf_constants: contains mathematical, numerical and physical constants
- pi2
- xi
- pi

Subroutines and functions

interface sf_fft_fftpack/ft_direct(ft, t, w)

This function evaluates the direct Fourier transform of a discretized function from time to frequency domain. Takes as input a discretized function \(ft_{i}\) on a set of time points \(t_{i}\) and a set of frequencies \(w_{j}\) and returns a discretized function \(fw_{j}\) = simps \((ft_{i} \cdot e^{-i 2\pi w_{j}t_{i}}, t_{1}, t_{N})\) , where N=size(t)

Parameters:

ft (•) [real/complex, in,required] – Dicretized function to transform (time domain)
t (size(ft)) [real, in] – Discretized time points
w (•) [real, in] – Discretized frequency points

Result:

fw (size(w)) [real/complex] – Discretized Fourier-transform function (frequency domain)

interface sf_fft_fftpack/ft_inverse(fw, t, w)

This function evaluates the inverse Fourier transform of a discretized function from frequency to time domain. Takes as input a discretized function \(fw_{i}\) on a set of time points \(w_{i}\) and a set of frequencies \(t_{j}\) and returns a discretized function \(ft_{j}\) = simps \((fw_{i} \cdot e^{i 2\pi t_{j}w_{i}}, w_{1}, w_{N})\) , where N=size(w)

Parameters:

fw (•) [real/complex, in,required] – Dicretized function to transform (frequency domain)
t (•) [real, in] – Discretized time points
w (size(fw)) [real, in] – Discretized frequency points

Result:

ft (size(t)) [real/complex] – Discretized Fourier-transform function (time domain)

interface sf_fft_fftpack/fft_signal(ft, dt)

This function evaluates the fast Fourier transform of a discretized function ft . It returns fw = dt · tfft ( ft )

Parameters:

ft (•) [real/complex] – Time-domain function
dt [real] – Time step

Result:

fw (size(ft)) [real/complex] – Frequency-domain function

interface sf_fft_fftpack/ifft_signal(fw, dt)

This function evaluates the inverse fast Fourier transform of a discretized function fw . It returns ft = itfft ( fw ) / ft

Parameters:

fw (•) [real/complex] – Frequency-domain function
dt [real] – Time step

Result:

ft (size(fw)) [real/complex] – Time-domain function

interface sf_fft_fftpack/tfft(func_in[, func_out])

This subroutine evaluates the fast Fourier transform of a real or complex array func_in of length N, arranged with the zero of the domain at the center of the array. It performs the sequence

ifftshift of the input
fft (forward transform)
fftshift of the result
multiplication by N

If the optional argument func_out is present it contains the result and func_in is left untouched, otherwise func_in is overwritten. In the real case the imaginary part of the transform is discarded.

Parameters:: func_in (•) [real/complex, in,required] – Array of length N to transform, with the zero of the time domain at the center. Overwritten by the transform if func_out is absent
Options:: func_out (size(func_in)) [real/complex] – If present, receives the transform and func_in is left untouched

interface sf_fft_fftpack/itfft(func_in[, func_out])

This subroutine evaluates the inverse fast Fourier transform of a real or complex array func_in of length N, and is the counterpart of tfft. It performs the sequence

ifft (backward transform)
fftex (alternate the sign of the elements, shifting the output by half a period)
ifftshift of the result
division by N

If the optional argument func_out is present it contains the result and func_in is left untouched, otherwise func_in is overwritten. In the real case the imaginary part of the transform is discarded.

Parameters:: func_in (•) [real/complex, in,required] – Array of length N to transform, with the zero of the frequency domain at the center. Overwritten by the inverse transform if func_out is absent
Options:: func_out (size(func_in)) [real/complex] – If present, receives the inverse transform and func_in is left untouched

interface sf_fft_fftpack/fft(func)

This subroutine evaluates the forward 1-dimensional FFT of an array of lenght N using the FFTPACK 5.1 routines rfft1i + rfft1f for the real case and cfft1i + cfft1f for the complex case. These are called with the parameters

lensav = N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = N
lenr = N
inc = 1

The subroutine modifies the input array to contain

[y(0), y(1), ..., y(N/2), y(-N/2+1), ..., y(-1)] if N is even
[y(0), y(1), ..., y((N-1)/2), y(-(N-1)/2), ..., y(-1)] if N is odd

where

\(y(j) = \sum_{k=0}^{N-1} x(k) \cdot e^{-i 2\pi/N \cdot j \cdot k}\)

for \(j \in [0,N-1]\)

Parameters:: func (•) [real/complex, inout] – Input/output array of size N

interface sf_fft_fftpack/ifft(func)

This subroutine evaluates the backward 1-dimensional FFT of an array of lenght N using the FFTPACK 5.1 routines rfft1i + rfft1b for the real case and cfft1i + cfft1b for the complex case. These are called with the parameters

lensav = N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = N
lenr = N
inc = 1

The subroutine modifies the input array to contain

[y(0), y(1), ..., y(N/2), y(-N/2+1), ..., y(-1)] if N is even
[y(0), y(1), ..., y((N-1)/2), y(-(N-1)/2), ..., y(-1)] if N is odd

where

\(y(j) = \sum_{k=-N/2}^{N/2-1} x(k) \cdot a^{i 2\pi/N \cdot j \cdot k}\)

for \(j \in [0,N-1]\)

Parameters:: func (•) [real/complex, inout] – Input/output array of size N, in the FFTPACK packed form on input, overwritten by the real sequence

interface sf_fft_fftpack/fft2(func)

This subroutine evaluates the forward 2-dimensional FFT of an array of size (L,M) using the FFTPACK 5.1 routines rfft2i + rfft2f for the real case and cfft2i + cfft2f for the complex case. These are called with the parameters

ldim = L
lensav = L + 3*M + int( log(dble(N))/log(2.d0) ) + 2*int( log(dble(M))/log(2.d0) ) + 12 (real)
lensav = 2*(L+M) + int( log(dble(N))/log(2.d0) ) + int( log(dble(M))/log(2.d0) ) + 8 (complex)
lenwrk = M*(L+1) (real), 2*L*M (complex)

The subroutine modifies the input array to contain

\(y(i,j) = \frac{1}{LM} \sum_{l=0}^{L-1} \sum_{m=0}^{M-1} x(l,m) \cdot e^{-i 2\pi (i l/L + j m/M)}\)

for \(i \in [0,L-1]\), \(j \in [0,M-1]\). The real subroutine returns the FFTPACK packed form of the transform of a real array.

Parameters:: func (•, •) [real/complex, inout] – Input/output array of size (L,M), overwritten by the FFTPACK packed form of its transform

interface sf_fft_fftpack/ifft2(func)

This subroutine evaluates the backward 2-dimensional FFT of an array of size (L,M) using the FFTPACK 5.1 routines rfft2i + rfft2b for the real case and cfft2i + cfft2b for the complex case. These are called with the same parameters as fft2.

The subroutine modifies the input array to contain

\(x(i,j) = \sum_{l=0}^{L-1} \sum_{m=0}^{M-1} y(l,m) \cdot e^{i 2\pi (i l/L + j m/M)}\)

for \(i \in [0,L-1]\), \(j \in [0,M-1]\).

Parameters:: func (•, •) [real/complex, inout] – Input/output array of size (L,M), in the FFTPACK packed form on input, overwritten by the real array

interface sf_fft_fftpack/fftn(func, n, lot)

This subroutine evaluates lot forward 1-dimensional FFTs of length N on the sequences stored consecutively in the 1-dimensional array func, using the FFTPACK 5.1 routines rfftmi + rfftmf for the real case and cfftmi + cfftmf for the complex case. These are called with the parameters

lenr = N*lot
lensav = N + int( log(dble(N))/log(2.d0) ) + 4 (real), 2*N + int( log(dble(N))/log(2.d0) ) + 4 (complex)
lenwrk = N*lot (real), 2*N*lot (complex)
jump = N
inc = 1

For the \(l\)-th sequence, \(l \in [0,lot-1]\), the subroutine modifies the input array to contain

\(y(l N + j) = \sum_{k=0}^{N-1} x(l N + k) \cdot e^{-i 2\pi/N \cdot j \cdot k}\)

for \(j \in [0,N-1]\). The program stops if mod(N*lot,size(func))/=0.

Parameters:

func (•) [real/complex, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their transforms (FFTPACK packed form)
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

interface sf_fft_fftpack/ifftn(func, n, lot)

This subroutine evaluates lot backward 1-dimensional FFTs of length N on the sequences stored consecutively in the 1-dimensional array func, using the FFTPACK 5.1 routines rfftmi + rfftmb for the real case and cfftmi + cfftmb for the complex case. These are called with the same parameters as fftn.

For the \(l\)-th sequence, \(l \in [0,lot-1]\), the subroutine modifies the input array to contain

\(x(l N + j) = \sum_{k=0}^{N-1} y(l N + k) \cdot e^{i 2\pi/N \cdot j \cdot k}\)

for \(j \in [0,N-1]\).

Parameters:

func (•) [real/complex, inout] – Input/output array storing lot sequences of length N consecutively, in the FFTPACK packed form on input, overwritten by the real sequences
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

interface sf_fft_fftpack/cosft(func)

This subroutine evaluates the forward 1-dimensional cosine FFT of a real array of lenght N using the FFTPACK 5.1 routines cost1i + cost1f. These are called with the parameters

lensav = 2*N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = N-1
lenr = N
inc = 1

The subroutine modifies the input array in place. The transform is the inverse of icosft, i.e. the normalization is carried by the forward transform.

Parameters:: func (•) [real, inout] – Input/output array of size N, overwritten by its cosine transform

interface sf_fft_fftpack/icosft(func)

This subroutine evaluates the backward 1-dimensional cosine FFT of a real array of lenght N using the FFTPACK 5.1 routines cost1i + cost1b. These are called with the same parameters as cosft.

The subroutine modifies the input array to contain

\(x(i) = y(1) + (-1)^{i-1} y(N) + 2 \sum_{k=2}^{N-1} y(k) \cdot \cos\left(\frac{\pi (k-1)(i-1)}{N-1}\right)\)

for \(i \in [1,N]\).

Parameters:: func (•) [real, inout] – Input/output array of size N, overwritten by its backward cosine transform

interface sf_fft_fftpack/cosftn(func, n, lot)

This subroutine evaluates lot forward cosine FFTs of length N on the sequences stored consecutively in the 1-dimensional array func, using the FFTPACK 5.1 routines costmi + costmf. These are called with the parameters

lenr = N*lot
lensav = 2*N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = lot*(N+1)
jump = N
inc = 1

The program stops if mod(N*lot,size(func))/=0.

Parameters:

func (•) [real, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their cosine transforms
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

interface sf_fft_fftpack/icosftn(func, n, lot)

This subroutine evaluates lot backward cosine FFTs of length N on the sequences stored consecutively in the 1-dimensional array func, using the FFTPACK 5.1 routines costmi + costmb. These are called with the same parameters as cosftn.

Parameters:

func (•) [real, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their backward cosine transforms
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

interface sf_fft_fftpack/sinft(func)

This subroutine evaluates the forward 1-dimensional sine FFT of a real array of lenght N using the FFTPACK 5.1 routines sint1i + sint1f. These are called with the parameters

lensav = N/2 + N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = 2*(N+1)
lenr = N
inc = 1

The subroutine modifies the input array in place. The transform is the inverse of isinft, i.e. the normalization is carried by the forward transform.

Parameters:: func (•) [real, inout] – Input/output array of size N, overwritten by its sine transform

interface sf_fft_fftpack/isinft(func)

This subroutine evaluates the backward 1-dimensional sine FFT of a real array of lenght N using the FFTPACK 5.1 routines sint1i + sint1b. These are called with the same parameters as sinft.

The subroutine modifies the input array to contain

\(x(i) = 2 \sum_{k=1}^{N} y(k) \cdot \sin\left(\frac{\pi k i}{N+1}\right)\)

for \(i \in [1,N]\).

Parameters:: func (•) [real, inout] – Input/output array of size N, overwritten by its backward sine transform

interface sf_fft_fftpack/sinftn(func, n, lot)

This subroutine evaluates lot forward sine FFTs of length N on the sequences stored consecutively in the 1-dimensional array func, using the FFTPACK 5.1 routines sintmi + sintmf. These are called with the parameters

lenr = N*lot
lensav = N/2 + N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = 2*lot*(N+2)
jump = N
inc = 1

The program stops if mod(N*lot,size(func))/=0.

Parameters:

func (•) [real, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their sine transforms
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

interface sf_fft_fftpack/isinftn(func, n, lot)

This subroutine evaluates lot backward sine FFTs of length N on the sequences stored consecutively in the 1-dimensional array func, using the FFTPACK 5.1 routines sintmi + sintmb. These are called with the same parameters as sinftn.

Parameters:

func (•) [real, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their backward sine transforms
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

interface sf_fft_fftpack/fftshift(fin)

This function shifts the zero-frequency component of the output of fft to the center of the array. The input

[y(0), y(1), ..., y(N/2), y(-N/2+1), ..., y(-1)] if N is even
[y(0), y(1), ..., y((N-1)/2), y(-(N-1)/2), ..., y(-1)] if N is odd

is returned as

[y(-N/2), ..., y(-1), y(0), y(1), ..., y(N/2-1)] if N is even
[y(-(N-1)/2), ..., y(-1), y(0), y(1), ..., y((N-1)/2)] if N is odd

Parameters:: fin (•) [real/complex] – Array to shift, of size N, in the aliased order returned by fft
Result:: fout (size(fin)) [real/complex] – Shifted array, of the same size as fin, with the zero-frequency component at the center

interface sf_fft_fftpack/ifftshift(fin)

This function is the inverse of fftshift: it moves the zero-frequency component from the center of the array back to the first position.

Parameters:: fin (•) [real/complex] – Array of size N with the zero-frequency component at the center
Result:: fout (size(fin)) [real/complex] – Array in the aliased order returned by fft, of the same size as fin

interface sf_fft_fftpack/fftex(func)

This subroutine alternates the sign of the elements of an array, \(f(i) \to (-1)^{i-1} f(i)\). It is used by itfft to shift the output of ifft by half a period.

Parameters:: func (•) [real/complex] – Input/output array, overwritten by \((-1)^{i-1} f(i)\)

subroutine sf_fft_fftpack/rfft_1d_forward(func)

This subroutine evaluates the forward 1-dimensional FFT of a real array func of length N, in place, using the FFTPACK 5.1 routines rfft1i + rfft1f. These are called with the parameters

lensav = N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = N
lenr = N
inc = 1

The transform is not normalized. Since the transform of a real sequence is Hermitian, the subroutine returns the FFTPACK packed form of

\(y(j) = \sum_{k=0}^{N-1} x(k) \cdot e^{-i 2\pi/N \cdot j \cdot k}\)

that is, the array contains

[Re y(0), Re y(1), Im y(1), ..., Re y(N/2-1), Im y(N/2-1), Re y(N/2)] if N is even
[Re y(0), Re y(1), Im y(1), ..., Re y((N-1)/2), Im y((N-1)/2)] if N is odd

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•) [real, inout] – Input/output array of size N

subroutine sf_fft_fftpack/cfft_1d_forward(func)

This subroutine evaluates the forward 1-dimensional FFT of a complex array func of length N, in place, using the FFTPACK 5.1 routines cfft1i + cfft1f. These are called with the parameters

lensav = 2*N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = 2*N
lenc = N
inc = 1

The transform is not normalized. The subroutine modifies the input array to contain

\(y(j) = \sum_{k=0}^{N-1} x(k) \cdot e^{-i 2\pi/N \cdot j \cdot k}\)

for \(j \in [0,N-1]\), i.e. the array is returned in the aliased order

[y(0), y(1), ..., y(N/2), y(-N/2+1), ..., y(-1)] if N is even
[y(0), y(1), ..., y((N-1)/2), y(-(N-1)/2), ..., y(-1)] if N is odd

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•) [complex, inout] – Input/output array of size N, overwritten by its transform

subroutine sf_fft_fftpack/rfft_2d_forward(func)

This subroutine evaluates the forward 2-dimensional FFT of a real array func of size (L,M), in place, using the FFTPACK 5.1 routines rfft2i + rfft2f. These are called with the parameters

ldim = L
lensav = L + 3*M + int( log(dble(L))/log(2.d0) ) + 2*int( log(dble(M))/log(2.d0) ) + 12
lenwrk = M*(L+1)

The subroutine computes

\(y(i,j) = \frac{1}{LM} \sum_{l=0}^{L-1} \sum_{m=0}^{M-1} x(l,m) \cdot e^{-i 2\pi (i l/L + j m/M)}\)

for \(i \in [0,L-1]\), \(j \in [0,M-1]\). The transform of a real array is Hermitian, so the array is overwritten with the FFTPACK packed form of \(y\).

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•, •) [real, inout] – Input/output array of size (L,M), overwritten by the FFTPACK packed form of its transform

subroutine sf_fft_fftpack/cfft_2d_forward(func)

This subroutine evaluates the forward 2-dimensional FFT of a complex array func of size (L,M), in place, using the FFTPACK 5.1 routines cfft2i + cfft2f. These are called with the parameters

ldim = L
lensav = 2*(L+M) + int( log(dble(L))/log(2.d0) ) + int( log(dble(M))/log(2.d0) ) + 8
lenwrk = 2*L*M

The subroutine modifies the input array to contain

\(y(i,j) = \frac{1}{LM} \sum_{l=0}^{L-1} \sum_{m=0}^{M-1} x(l,m) \cdot e^{-i 2\pi (i l/L + j m/M)}\)

for \(i \in [0,L-1]\), \(j \in [0,M-1]\).

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•, •) [complex, inout] – Input/output array of size (L,M), overwritten by its transform

subroutine sf_fft_fftpack/rfft_nd_forward(func, n, lot)

This subroutine evaluates lot forward 1-dimensional FFTs of length N on the real sequences stored consecutively in the 1-dimensional array func, in place, using the FFTPACK 5.1 routines rfftmi + rfftmf. These are called with the parameters

lenr = N*lot
lensav = N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = N*lot
jump = N
inc = 1

The transforms are not normalized. For the \(l\)-th sequence, \(l \in [0,lot-1]\), the subroutine computes

\(y(l N + j) = \sum_{k=0}^{N-1} x(l N + k) \cdot e^{-i 2\pi/N \cdot j \cdot k}\)

for \(j \in [0,N-1]\). Since the transform of a real sequence is Hermitian, each sequence is overwritten with the FFTPACK packed form of \(y\), as in rfft_1d_forward.

The program stops if mod(N*lot,size(func))/=0, or if the FFTPACK routines return an error code.

Parameters:

func (•) [real, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their transforms (FFTPACK packed form)
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

subroutine sf_fft_fftpack/cfft_nd_forward(func, n, lot)

This subroutine evaluates lot forward 1-dimensional FFTs of length N on the complex sequences stored consecutively in the 1-dimensional array func, in place, using the FFTPACK 5.1 routines cfftmi + cfftmf. These are called with the parameters

lenc = N*lot
lensav = 2*N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = 2*N*lot
jump = N
inc = 1

The transforms are not normalized. For the \(l\)-th sequence, \(l \in [0,lot-1]\), the subroutine modifies the input array to contain

\(y(l N + j) = \sum_{k=0}^{N-1} x(l N + k) \cdot e^{-i 2\pi/N \cdot j \cdot k}\)

for \(j \in [0,N-1]\). Each sequence is returned in the aliased order described in cfft_1d_forward.

The program stops if mod(N*lot,size(func))/=0, or if the FFTPACK routines return an error code.

Parameters:

func (•) [complex, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their transforms
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

subroutine sf_fft_fftpack/cost_1d_forward(func)

This subroutine evaluates the forward 1-dimensional cosine FFT of a real array func of length N, in place, using the FFTPACK 5.1 routines cost1i + cost1f. These are called with the parameters

lensav = 2*N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = N-1
lenr = N
inc = 1

The transform is the inverse of cost_1d_backward, i.e. the normalization is carried by the forward transform. The subroutine modifies the input array to contain the coefficients \(y(k)\) such that the backward transform cost_1d_backward reconstructs

\(x(i) = y(1) + (-1)^{i-1} y(N) + 2 \sum_{k=2}^{N-1} y(k) \cdot \cos\left(\frac{\pi (k-1)(i-1)}{N-1}\right)\)

for \(i \in [1,N]\).

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•) [real, inout] – Input/output array of size N, overwritten by its cosine transform

subroutine sf_fft_fftpack/sint_1d_forward(func)

This subroutine evaluates the forward 1-dimensional sine FFT of a real array func of length N, in place, using the FFTPACK 5.1 routines sint1i + sint1f. These are called with the parameters

lensav = N/2 + N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = 2*(N+1)
lenr = N
inc = 1

The transform is the inverse of sint_1d_backward, i.e. the normalization is carried by the forward transform. The subroutine modifies the input array to contain the coefficients \(y(k)\) such that the backward transform sint_1d_backward reconstructs

\(x(i) = 2 \sum_{k=1}^{N} y(k) \cdot \sin\left(\frac{\pi k i}{N+1}\right)\)

for \(i \in [1,N]\).

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•) [real, inout] – Input/output array of size N, overwritten by its sine transform

subroutine sf_fft_fftpack/cost_nd_forward(func, n, lot)

This subroutine evaluates lot forward 1-dimensional cosine FFTs of length N on the real sequences stored consecutively in the 1-dimensional array func, in place, using the FFTPACK 5.1 routines costmi + costmf. These are called with the parameters

lenr = N*lot
lensav = 2*N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = lot*(N+1)
jump = N
inc = 1

For the \(l\)-th sequence, \(l \in [0,lot-1]\), the subroutine acts as cost_1d_forward on the elements \(l N + 1, ..., (l+1) N\) of the array, so the transform is the inverse of cost_nd_backward and the normalization is carried by the forward transform.

The program stops if mod(N*lot,size(func))/=0, or if the FFTPACK routines return an error code.

Parameters:

func (•) [real, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their cosine transforms
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

subroutine sf_fft_fftpack/sint_nd_forward(func, n, lot)

This subroutine evaluates lot forward 1-dimensional sine FFTs of length N on the real sequences stored consecutively in the 1-dimensional array func, in place, using the FFTPACK 5.1 routines sintmi + sintmf. These are called with the parameters

lenr = N*lot
lensav = N/2 + N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = 2*lot*(N+2)
jump = N
inc = 1

For the \(l\)-th sequence, \(l \in [0,lot-1]\), the subroutine acts as sint_1d_forward on the elements \(l N + 1, ..., (l+1) N\) of the array, so the transform is the inverse of sint_nd_backward and the normalization is carried by the forward transform.

The program stops if mod(N*lot,size(func))/=0, or if the FFTPACK routines return an error code.

Parameters:

func (•) [real, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their sine transforms
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

subroutine sf_fft_fftpack/rfft_1d_backward(func)

This subroutine evaluates the backward 1-dimensional FFT of a real array func of length N, in place, using the FFTPACK 5.1 routines rfft1i + rfft1b. These are called with the parameters

lensav = N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = N
lenr = N
inc = 1

The transform is not normalized, so applying it after rfft_1d_forward returns the original array multiplied by N. The input array is expected in the FFTPACK packed form produced by rfft_1d_forward, that is

[Re y(0), Re y(1), Im y(1), ..., Re y(N/2-1), Im y(N/2-1), Re y(N/2)] if N is even
[Re y(0), Re y(1), Im y(1), ..., Re y((N-1)/2), Im y((N-1)/2)] if N is odd

and the subroutine modifies the input array to contain the real sequence

\(x(j) = \sum_{k=0}^{N-1} y(k) \cdot e^{i 2\pi/N \cdot j \cdot k}\)

for \(j \in [0,N-1]\), where \(y(k)\) is the Hermitian sequence encoded by the packed input.

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•) [real, inout] – Input/output array of size N, in the FFTPACK packed form on input, overwritten by the real sequence

subroutine sf_fft_fftpack/cfft_1d_backward(func)

This subroutine evaluates the backward 1-dimensional FFT of a complex array func of length N, in place, using the FFTPACK 5.1 routines cfft1i + cfft1b. These are called with the parameters

lensav = 2*N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = 2*N
lenc = N
inc = 1

The transform is not normalized, so applying it after cfft_1d_forward returns the original array multiplied by N. The input array is expected in the aliased order produced by cfft_1d_forward, that is

[y(0), y(1), ..., y(N/2), y(-N/2+1), ..., y(-1)] if N is even
[y(0), y(1), ..., y((N-1)/2), y(-(N-1)/2), ..., y(-1)] if N is odd

and the subroutine modifies the input array to contain

\(x(j) = \sum_{k=0}^{N-1} y(k) \cdot e^{i 2\pi/N \cdot j \cdot k}\)

for \(j \in [0,N-1]\).

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•) [complex, inout] – Input/output array of size N, overwritten by its backward transform

subroutine sf_fft_fftpack/rfft_2d_backward(func)

This subroutine evaluates the backward 2-dimensional FFT of a real array func of size (L,M), in place, using the FFTPACK 5.1 routines rfft2i + rfft2b. These are called with the parameters

ldim = L
lensav = L + 3*M + int( log(dble(L))/log(2.d0) ) + 2*int( log(dble(M))/log(2.d0) ) + 12
lenwrk = M*(L+1)

The input array is expected in the FFTPACK packed form produced by rfft_2d_forward, and the subroutine modifies it to contain the real array

\(x(i,j) = \sum_{l=0}^{L-1} \sum_{m=0}^{M-1} y(l,m) \cdot e^{i 2\pi (i l/L + j m/M)}\)

for \(i \in [0,L-1]\), \(j \in [0,M-1]\), where \(y(l,m)\) is the Hermitian array encoded by the packed input. The transform carries no normalization factor, because the \(1/(LM)\) factor is carried by rfft_2d_forward.

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•, •) [real, inout] – Input/output array of size (L,M), in the FFTPACK packed form on input, overwritten by the real array

subroutine sf_fft_fftpack/cfft_2d_backward(func)

This subroutine evaluates the backward 2-dimensional FFT of a complex array func of size (L,M), in place, using the FFTPACK 5.1 routines cfft2i + cfft2b. These are called with the parameters

ldim = L
lensav = 2*(L+M) + int( log(dble(L))/log(2.d0) ) + int( log(dble(M))/log(2.d0) ) + 8
lenwrk = 2*L*M

The input array is expected in the form produced by cfft_2d_forward, and the subroutine modifies it to contain

\(x(i,j) = \sum_{l=0}^{L-1} \sum_{m=0}^{M-1} y(l,m) \cdot e^{i 2\pi (i l/L + j m/M)}\)

for \(i \in [0,L-1]\), \(j \in [0,M-1]\). The transform carries no normalization factor, because the \(1/(LM)\) factor is carried by cfft_2d_forward.

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•, •) [complex, inout] – Input/output array of size (L,M), overwritten by its backward transform

subroutine sf_fft_fftpack/rfft_nd_backward(func, n, lot)

This subroutine evaluates lot backward 1-dimensional FFTs of length N on the real sequences stored consecutively in the 1-dimensional array func, in place, using the FFTPACK 5.1 routines rfftmi + rfftmb. These are called with the parameters

lenr = N*lot
lensav = N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = N*lot
jump = N
inc = 1

Each sequence is expected in the FFTPACK packed form produced by rfft_nd_forward. The transforms are not normalized, so applying them after rfft_nd_forward returns the original array multiplied by N. For the \(l\)-th sequence, \(l \in [0,lot-1]\), the subroutine modifies the input array to contain the real sequence

\(x(l N + j) = \sum_{k=0}^{N-1} y(l N + k) \cdot e^{i 2\pi/N \cdot j \cdot k}\)

for \(j \in [0,N-1]\), where \(y(l N + k)\) is the Hermitian sequence encoded by the packed input.

The program stops if mod(N*lot,size(func))/=0, or if the FFTPACK routines return an error code.

Parameters:

func (•) [real, inout] – Input/output array storing lot sequences of length N consecutively, in the FFTPACK packed form on input, overwritten by the real sequences
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

subroutine sf_fft_fftpack/cfft_nd_backward(func, n, lot)

This subroutine evaluates lot backward 1-dimensional FFTs of length N on the complex sequences stored consecutively in the 1-dimensional array func, in place, using the FFTPACK 5.1 routines cfftmi + cfftmb. These are called with the parameters

lenc = N*lot
lensav = 2*N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = 2*N*lot
jump = N
inc = 1

Each sequence is expected in the aliased order produced by cfft_nd_forward. The transforms are not normalized, so applying them after cfft_nd_forward returns the original array multiplied by N. For the \(l\)-th sequence, \(l \in [0,lot-1]\), the subroutine modifies the input array to contain

\(x(l N + j) = \sum_{k=0}^{N-1} y(l N + k) \cdot e^{i 2\pi/N \cdot j \cdot k}\)

for \(j \in [0,N-1]\).

The program stops if mod(N*lot,size(func))/=0, or if the FFTPACK routines return an error code.

Parameters:

func (•) [complex, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their backward transforms
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

subroutine sf_fft_fftpack/cost_1d_backward(func)

This subroutine evaluates the backward 1-dimensional cosine FFT of a real array func of length N, in place, using the FFTPACK 5.1 routines cost1i + cost1b. These are called with the parameters

lensav = 2*N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = N-1
lenr = N
inc = 1

The transform is the inverse of cost_1d_forward, i.e. the normalization is carried by the forward transform, so no scaling is applied here. The input array contains the coefficients \(y(k)\) produced by cost_1d_forward, and the subroutine modifies it to contain

\(x(i) = y(1) + (-1)^{i-1} y(N) + 2 \sum_{k=2}^{N-1} y(k) \cdot \cos\left(\frac{\pi (k-1)(i-1)}{N-1}\right)\)

for \(i \in [1,N]\).

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•) [real, inout] – Input/output array of size N, overwritten by its backward cosine transform

subroutine sf_fft_fftpack/sint_1d_backward(func)

This subroutine evaluates the backward 1-dimensional sine FFT of a real array func of length N, in place, using the FFTPACK 5.1 routines sint1i + sint1b. These are called with the parameters

lensav = N/2 + N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = 2*(N+1)
lenr = N
inc = 1

The transform is the inverse of sint_1d_forward, i.e. the normalization is carried by the forward transform, so no scaling is applied here. The input array contains the coefficients \(y(k)\) produced by sint_1d_forward, and the subroutine modifies it to contain

\(x(i) = 2 \sum_{k=1}^{N} y(k) \cdot \sin\left(\frac{\pi k i}{N+1}\right)\)

for \(i \in [1,N]\).

The program stops if the FFTPACK routines return an error code.

Parameters:: func (•) [real, inout] – Input/output array of size N, overwritten by its backward sine transform

subroutine sf_fft_fftpack/cost_nd_backward(func, n, lot)

This subroutine evaluates lot backward 1-dimensional cosine FFTs of length N on the real sequences stored consecutively in the 1-dimensional array func, in place, using the FFTPACK 5.1 routines costmi + costmb. These are called with the parameters

lenr = N*lot
lensav = 2*N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = lot*(N+1)
jump = N
inc = 1

For the \(l\)-th sequence, \(l \in [0,lot-1]\), the subroutine acts as cost_1d_backward on the elements \(l N + 1, ..., (l+1) N\) of the array, so the transform is the inverse of cost_nd_forward and the normalization is carried by the forward transform, so no scaling is applied here.

The program stops if mod(N*lot,size(func))/=0, or if the FFTPACK routines return an error code.

Parameters:

func (•) [real, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their backward cosine transforms
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

subroutine sf_fft_fftpack/sint_nd_backward(func, n, lot)

This subroutine evaluates lot backward 1-dimensional sine FFTs of length N on the real sequences stored consecutively in the 1-dimensional array func, in place, using the FFTPACK 5.1 routines sintmi + sintmb. These are called with the parameters

lenr = N*lot
lensav = N/2 + N + int( log(dble(N))/log(2.d0) ) + 4
lenwrk = 2*lot*(N+2)
jump = N
inc = 1

For the \(l\)-th sequence, \(l \in [0,lot-1]\), the subroutine acts as sint_1d_backward on the elements \(l N + 1, ..., (l+1) N\) of the array, so the transform is the inverse of sint_nd_forward and the normalization is carried by the forward transform, so no scaling is applied here.

The program stops if mod(N*lot,size(func))/=0, or if the FFTPACK routines return an error code.

Parameters:

func (•) [real, inout] – Input/output array storing lot sequences of length N consecutively, overwritten by their backward sine transforms
n [integer, in] – Length of each sequence
lot [integer, in] – Number of sequences

function sf_fft_fftpack/rfft_1d_shift(fin)

This function shifts the zero-frequency component of a real array fin of length N to the center of the array, and returns the result as an array of the same size. It is a circular rotation by floor((N+1)/2) elements, so the input

[y(0), y(1), ..., y(N/2), y(-N/2+1), ..., y(-1)] if N is even
[y(0), y(1), ..., y((N-1)/2), y(-(N-1)/2), ..., y(-1)] if N is odd

is returned as

[y(-N/2), ..., y(-1), y(0), y(1), ..., y(N/2-1)] if N is even
[y(-(N-1)/2), ..., y(-1), y(0), y(1), ..., y((N-1)/2)] if N is odd

where y(N/2) and y(-N/2) are the same element for even N.

Parameters:: fin (•) [real] – Array to shift, of size N, in the aliased order returned by fft
Result:: fout (size(fin)) [real] – Shifted array, of the same size as fin, with the zero-frequency component at the center

function sf_fft_fftpack/cfft_1d_shift(fin)

This function shifts the zero-frequency component of a complex array fin of length N to the center of the array, and returns the result as an array of the same size. It acts as rfft_1d_shift on complex data.

Parameters:: fin (•) [complex] – Array to shift, of size N, in the aliased order returned by fft
Result:: fout (size(fin)) [complex] – Shifted array, of the same size as fin, with the zero-frequency component at the center

function sf_fft_fftpack/rfft_1d_ishift(fin)

This function is the inverse of rfft_1d_shift. It moves the zero-frequency component of a real array fin of length N from the center of the array back to the first position, and returns the result as an array of the same size. It is a circular rotation by N-floor((N+1)/2) elements, so the input

[y(-N/2), ..., y(-1), y(0), y(1), ..., y(N/2-1)] if N is even
[y(-(N-1)/2), ..., y(-1), y(0), y(1), ..., y((N-1)/2)] if N is odd

is returned as

[y(0), y(1), ..., y(N/2-1), y(-N/2), ..., y(-1)] if N is even
[y(0), y(1), ..., y((N-1)/2), y(-(N-1)/2), ..., y(-1)] if N is odd

that is, the aliased order returned by fft.

Parameters:: fin (•) [real] – Array of size N with the zero-frequency component at the center
Result:: fout (size(fin)) [real] – Array in the aliased order returned by fft, of the same size as fin

function sf_fft_fftpack/cfft_1d_ishift(fin)

This function is the inverse of cfft_1d_shift. It moves the zero-frequency component of a complex array fin of length N from the center of the array back to the first position, and returns the result as an array of the same size. It acts as rfft_1d_ishift on complex data.

Parameters:: fin (•) [complex] – Array of size N with the zero-frequency component at the center
Result:: fout (size(fin)) [complex] – Array in the aliased order returned by fft, of the same size as fin

subroutine sf_fft_fftpack/rfft_1d_ex(func)

This subroutine alternates the sign of the elements of a real array func in place,

\(f(i) \to (-1)^{i-1} f(i)\)

for \(i \in [1,N]\), where N=size(func). Applying it twice returns the original array.

Parameters:: func (•) [real] – Input/output array, overwritten by \((-1)^{i-1} f(i)\)

subroutine sf_fft_fftpack/cfft_1d_ex(func)

This subroutine alternates the sign of the elements of a complex array func in place,

\(f(i) \to (-1)^{i-1} f(i)\)

for \(i \in [1,N]\), where N=size(func). Applying it twice returns the original array.

Parameters:: func (•) [complex] – Input/output array, overwritten by \((-1)^{i-1} f(i)\)

function sf_fft_fftpack/fft_tmax(l, dt)

This function returns the maximum time of the time grid used with FFT_signal, for an array of L points with time step dt:

\(t_{max} = L \cdot dt / 2\)

Parameters:

l [integer] – Number of points of the time array
dt [real] – Time step

Result:

fft_tmax [real] – Maximum time \(t_{max}\)

function sf_fft_fftpack/fft_fmax(l, dt)

This function returns the maximum frequency of the frequency grid used with FFT_signal, for a time step dt:

\(\omega_{max} = \pi / dt\)

The argument L is not used, and is kept for symmetry with fft_tmax.

Parameters:

l [integer] – Number of points (not used)
dt [real] – Time step

Result:

fft_fmax [real] – Maximum frequency \(\omega_{max}\)

function sf_fft_fftpack/fft_tarray(l, dt)

This function returns the time array of L equally spaced points in \([-t_{max}, t_{max}]\), with \(t_{max}\) = fft_tmax ( L , dt ) . The endpoint is included, so the actual spacing is \(2 t_{max}/(L-1)\).

Parameters:

l [integer] – Number of points
dt [real] – Time step

Result:

time (l) [real] – Time array of L points in \([-t_{max}, t_{max}]\)

function sf_fft_fftpack/fft_farray(l, dt[, df])

This function returns the frequency array of L equally spaced points in \([-\omega_{max}, \omega_{max})\), with \(\omega_{max}\) = fft_fmax ( L , dt ) . The upper endpoint is excluded, so the spacing is \(2 \omega_{max}/L\). If the optional argument df is present, it returns this spacing.

Parameters:

l [integer] – Number of points
dt [real] – Time step

Options:

df [real] – If present, it returns the frequency spacing

Result:

freq (l) [real] – Frequency array of L points in \([-\omega_{max}, \omega_{max})\)