I'm trying to implement a Fast Fourier Transform (Radix-2) in MS's Excel VBA. The code I'm using pulls data from a range in the worksheet, does the calculations, then dumps the results in the adjacent columns. What I'm having trouble with is 1) know what to do with the resulting X[k] arrays, and 2) matching these results with the results from Excel's built in FFT (they do not currently match). The code is shown below. Thanks in advance for your help.
Sub Enforce_DecimationInTime()
On Error GoTo ERROR_HANDLING
Dim SubName As String
SubName = "Enforce_DecimationInTime()"
Dim WS As Worksheet
Dim n As Long, v As Long, LR As Long, x As Long
Set WS = Worksheets("FFT")
LR = WS.Range("A" & Rows.Count).End(xlUp).Row
n = LR - 1
Do Until 2 ^ x <= n And 2 ^ (x + 1) > n 'locates largest power of 2 from size of input array
x = x + 1
Loop
n = n - (n - 2 ^ x) 'calculates n using the largest power of 2
If n + 1 <> WS.Range("A" & Rows.Count).End(xlUp).Row Then
WS.Range("A" & 2 ^ x + 2 & ":A" & LR).Delete xlUp 'deletes extra input data
End If
v = WorksheetFunction.Log(n, 2) 'calculates number of decimations necessary
Application.ScreenUpdating = False
For x = 1 To v
Call Called_Core.DecimationInTime(WS, n, 2 ^ x, x) 'calls decimation in time subroutine
Next x
Application.ScreenUpdating = True
Exit Sub
ERROR_HANDLING:
MsgBox "Error encountered in " & SubName & ": exiting subroutine." _
& vbNewLine _
& vbNewLine & "Error description: " & Err.Description _
& vbNewLine & "Error number: " & Err.Number, vbCritical, Title:="Error!"
End
End Sub
The above subroutine calls the below subroutine through a For/Next loop to the count of "v".
Sub DecimationInTime(WS As Worksheet, n As Long, Factor As Integer, x As Long)
On Error GoTo ERROR_HANDLING
Dim SubName As String
SubName = "DecimationInTime()"
Dim f_1() As Single, f_2() As Single
Dim i As Long, m As Long, k As Long
Dim TFactor_N1 As String, TFactor_N2 As String, X_k() As String
Dim G_1() As Variant, G_2() As Variant
ReDim f_1(0 To n / Factor - 1) As Single
ReDim f_2(0 To n / Factor - 1) As Single
ReDim G_1(0 To n / 1 - 1) As Variant
ReDim G_2(0 To n / 1 - 1) As Variant
ReDim X_k(0 To n - 1) As String
TFactor_N1 = WorksheetFunction.Complex(0, -2 * WorksheetFunction.Pi / (n / 1)) 'twiddle factor for N
TFactor_N2 = WorksheetFunction.Complex(0, -2 * WorksheetFunction.Pi / (n / 2)) 'twiddle factor for N/2
For i = 0 To n / Factor - 1
f_1(i) = WS.Range("A" & 2 * i + 2).Value 'assign input data
f_2(i) = WS.Range("A" & 2 * i + 3).Value 'assign input data
Next i
WS.Cells(1, 1 + x).Value = "X[" & x & "]" 'labels X[k] column with k number
For k = 0 To n / 2 - 1
For m = 0 To n / Factor - 1
G_1(m) = WorksheetFunction.ImProduct(WorksheetFunction.ImPower(TFactor_N2, k * m), WorksheetFunction.Complex(f_1(m), 0)) 'defines G_1[m]
G_2(m) = WorksheetFunction.ImProduct(WorksheetFunction.ImPower(TFactor_N2, k * m), WorksheetFunction.Complex(f_2(m), 0)) 'defines G_2[m]
Next m
X_k(k) = WorksheetFunction.ImSum(WorksheetFunction.ImSum(G_1), WorksheetFunction.ImProduct(WorksheetFunction.ImSum(G_2), WorksheetFunction.ImPower(TFactor_N1, k))) 'defines X[k] for k
If k <= n / 2 Then X_k(k + n / 2) = WorksheetFunction.ImSum(WorksheetFunction.ImSum(G_1), WorksheetFunction.ImProduct(WorksheetFunction.ImSum(G_2), WorksheetFunction.ImPower(TFactor_N1, k), WorksheetFunction.Complex(-1, 0))) 'defines X[k] for k + n/2
WS.Cells(k + 2, 1 + x).Value = X_k(k)
WS.Cells(k + 2 + n / 2, 1 + x).Value = X_k(k + n / 2)
Next k
Exit Sub
ERROR_HANDLING:
MsgBox "Error encountered in " & SubName & ": exiting subroutine." _
& vbNewLine _
& vbNewLine & "Error description: " & Err.Description _
& vbNewLine & "Error number: " & Err.Number, vbCritical, Title:="Error!"
End
End Sub
The fast Fourier transform has become a major DSP tool since being popularized by Cooley and Tuckey in 1965. In the FFT formula, the DFT equation X(k) = ∑x(n)WNnk is decomposed into a number of short transforms and then recombined.
Octave uses the FFTW libraries to perform FFT computations. When Octave starts up and initializes the FFTW libraries, they read a system wide file (on a Unix system, it is typically /etc/fftw/wisdom ) that contains information useful to speed up FFT computations. This information is called the wisdom.
Y = fft( X ) computes the discrete Fourier transform (DFT) of X using a fast Fourier transform (FFT) algorithm. If X is a vector, then fft(X) returns the Fourier transform of the vector. If X is a matrix, then fft(X) treats the columns of X as vectors and returns the Fourier transform of each column.
I went back through the process and determined my problem was that I had assigned the wrong values to the twiddle factors, TFactor_N1 and TFactor_N2. After fixing this problem and adjusting which values are displayed, I was able to get the same results as Excel's built in FFT. The fixed code is show below.
Sub Enforce_DecimationInTime()
On Error GoTo ERROR_HANDLING
Dim SubName As String
SubName = "Enforce_DecimationInTime()"
Dim WS As Worksheet
Dim n As Long, v As Long, LR As Long, x As Long
Dim TFactor_N1 As String, TFactor_N2 As String
Set WS = Worksheets("FFT")
LR = WS.Range("A" & Rows.Count).End(xlUp).Row
n = LR - 1
Do Until 2 ^ x <= n And 2 ^ (x + 1) > n 'locates largest power of 2 from size of input array
x = x + 1
Loop
n = n - (n - 2 ^ x) 'calculates n using the largest power of 2
If n + 1 <> WS.Range("A" & Rows.Count).End(xlUp).Row Then
WS.Range("A" & 2 ^ x + 2 & ":A" & LR).Delete xlUp 'deletes extra input data
End If
v = WorksheetFunction.Log(n, 2) 'calculates number of decimations necessary
TFactor_N1 = WorksheetFunction.ImExp(WorksheetFunction.Complex(0, -2 * WorksheetFunction.Pi / (n / 1))) 'twiddle factor for N
TFactor_N2 = WorksheetFunction.ImExp(WorksheetFunction.Complex(0, -2 * WorksheetFunction.Pi / (n / 2))) 'twiddle factor for N/2
Application.ScreenUpdating = False
For x = 1 To v
Call Called_Core.DecimationInTime(WS, n, 2 ^ x, x, TFactor_N1, TFactor_N2) 'calls decimation in time subroutine
Next x
Application.ScreenUpdating = True
Exit Sub
ERROR_HANDLING:
MsgBox "Error encountered in " & SubName & ": exiting subroutine." _
& vbNewLine _
& vbNewLine & "Error description: " & Err.Description _
& vbNewLine & "Error number: " & Err.Number, vbCritical, Title:="Error!"
End
End Sub
Sub DecimationInTime(WS As Worksheet, n As Long, Factor As Integer, x As Long, TFactor_N1 As String, TFactor_N2 As String)
On Error GoTo ERROR_HANDLING
Dim SubName As String
SubName = "DecimationInTime()"
Dim f_1() As String, f_2() As String
Dim i As Long, m As Long, k As Long
Dim X_k() As String
Dim G_1() As Variant, G_2() As Variant
ReDim f_1(0 To n / Factor - 1) As String
ReDim f_2(0 To n / Factor - 1) As String
ReDim G_1(0 To n / 1 - 1) As Variant
ReDim G_2(0 To n / 1 - 1) As Variant
ReDim X_k(0 To n - 1) As String
For i = 0 To n / Factor - 1
f_1(i) = WS.Cells(2 * i + 2, 1).Value 'assign input data
f_2(i) = WS.Cells(2 * i + 3, 1).Value 'assign input data
Next i
For k = 0 To n / 2 - 1
For m = 0 To n / Factor - 1 'defines G_1[m] and G_2[m]
G_1(m) = WorksheetFunction.ImProduct(WorksheetFunction.ImPower(TFactor_N2, k * m), f_1(m))
G_2(m) = WorksheetFunction.ImProduct(WorksheetFunction.ImPower(TFactor_N2, k * m), f_2(m))
Next m 'defines X[k] for k and k + n/2
X_k(k) = WorksheetFunction.ImSum(WorksheetFunction.ImSum(G_1), WorksheetFunction.ImProduct(WorksheetFunction.ImSum(G_2), WorksheetFunction.ImPower(TFactor_N1, k)))
If k <= n / 2 Then X_k(k + n / 2) = WorksheetFunction.ImSub(WorksheetFunction.ImSum(G_1), WorksheetFunction.ImProduct(WorksheetFunction.ImSum(G_2), WorksheetFunction.ImPower(TFactor_N1, k)))
If x = 1 Then
WS.Cells(k + 2, 1 + x).Value = X_k(k)
WS.Cells(k + 2 + n / 2, 1 + x).Value = X_k(k + n / 2)
End If
Next k
Exit Sub
ERROR_HANDLING:
MsgBox "Error encountered in " & SubName & ": exiting subroutine." _
& vbNewLine _
& vbNewLine & "Error description: " & Err.Description _
& vbNewLine & "Error number: " & Err.Number, vbCritical, Title:="Error!"
End
End Sub
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With