◆ Cfft1f()

Sub Cfft1f	(	N As	Long,
		C() As	Complex,
		Wsave() As	Double,
		Info As	Long,
		Optional Inc As	Long = `1`
	)

One-dimensional complex Fourier transform

Purpose: This routine computes the one-dimensional Fourier transform of a periodic sequence within a complex array. This is referred to as the forward transform or Fourier analysis, transforming the sequence from physical to spectral space.
C(k) = (1/N)ΣC(j)exp(-2πijk/N) (Σ for j = 0 to N-1) (k = 0 to N-1) (i is imaginary unit)

This transform is normalized since a call to Cfft1f followed by a call to Cfft1b (or vice-versa) reproduces the original array subject to algorithmic constraints, roundoff error, etc.

Parameters

[in]	N	The length of the sequence to be transformed. (N >= 1) (The transform is most efficient when N is a product of small primes)
[in,out]	C()	Array C(LC - 1) (LC >= Inc*(N - 1) + 1) [in] The sequence to be transformed. [out] The Fourier forward transformed sequence of data.
[in]	Wsave()	Array Wsave(LWsave - 1) (LWsave >= 2*N + ln(N)/ln(2) + 4) Work data. Its contents must be initialized with a call to Cfft1i before the first call to Cfft1f or Cfft1b for a given transform length N.
[out]	Info	= 0: Successful exit. = -1: The argument N had an illegal value. (N < 1) = -2: The argument C() is invalid. (Array C() is not big enough) = -3: The argument Wsave() is invalid. (Array Wsave() is not big enough) = -5: The argument Inc had an illegal value. (Inc < 1)
[in]	Inc	(Optional) Integer increment between the locations, in array C(), of two consecutive elements within the sequence. (Inc >= 1) (default = 1)

Reference: FFTPACK

Example Program: Compute the Fourier transform and backward transform of 5 random data sequence successively, and compare with the original data sequence.
Sub Ex_Cfft1()

Const N = 5

Dim Wsave() As Double, C(N - 1) As Complex, C0(N - 1) As Complex

Dim LWsave As Long, Info As Long, I As Long

'-- Initialization

LWsave = 2 * N + Log(N) / Log(2) + 4

ReDim Wsave(LWsave - 1)

Call Cfft1i(N, Wsave, Info)

If Info <> 0 Then GoTo Err

'-- Generate test data

For I = 0 To N - 1

C(I) = Cmplx(Rnd(), Rnd())

C0(I) = Cmplx(Creal(C(I)), Cimag(C(I)))

Next

'-- Forward transform

Call Cfft1f(N, C(), Wsave(), Info)

If Info <> 0 Then GoTo Err

'-- Backward transform

Call Cfft1b(N, C(), Wsave(), Info)

If Info <> 0 Then GoTo Err

'-- Print result

For I = 0 To N - 1

Debug.Print "(" & CStr(Creal(C0(I))) & ", " & CStr(Cimag(C0(I))) & ")", "(" & CStr(Creal(C(I))) & ", " & CStr(Cimag(C(I))) & ")", Cabs(Csub(C(I), C0(I)))

Next

Exit Sub

Err:

Debug.Print "Error in Cfft1i/Cfft1f/Cfft1b: Info =", Info

End Sub

Cmplx
Function Cmplx(R As Double, Optional I As Double=0) As Complex
Building complex number

Cabs
Function Cabs(A As Complex) As Double
Absolute value of complex number

Cimag
Function Cimag(A As Complex) As Double
Imaginary part of complex number

Creal
Function Creal(A As Complex) As Double
Real part of complex number

Csub
Function Csub(A As Complex, B As Complex) As Complex
Subtraction of complex number from complex number

Cfft1i
Sub Cfft1i(N As Long, Wsave() As Double, Info As Long)
Initialization of work data for Cfft1f and Cfft1b

Cfft1b
Sub Cfft1b(N As Long, C() As Complex, Wsave() As Double, Info As Long, Optional Inc As Long=1)
One-dimensional complex fast Fourier backward transform

Cfft1f
Sub Cfft1f(N As Long, C() As Complex, Wsave() As Double, Info As Long, Optional Inc As Long=1)
One-dimensional complex Fourier transform

Example Results: (0.961953163146973, 0.871445834636688) (0.961953163146973, 0.871445834636688) 0

(5.62368631362915E-02, 0.949556648731232) (5.62368631362915E-02, 0.949556648731232) 0

(0.364018678665161, 0.524868428707123) (0.364018678665161, 0.524868428707123) 5.55111512312578E-17

(0.767111659049988, 5.35045266151428E-02) (0.767111659049988, 5.35045266151427E-02) 1.11022302462516E-16

(0.592458248138428, 0.468700110912323) (0.592458248138428, 0.468700110912323) 5.55111512312578E-17