FFT/Fourier Transforms in IronPython QuickStart Sample

Illustrates how to compute the forward and inverse Fourier transform of a real or complex signal using classes in the Numerics.NET.SignalProcessing namespace in IronPython.

View this sample in: C# Visual Basic F#

```Python
import numerics

from math import *

from Extreme.Mathematics import *
# We'll need real vectors...
from Extreme.Mathematics.LinearAlgebra import *
# and complex vectors...
from Extreme.Mathematics.LinearAlgebra.Complex import *
# The FFT classes reside in the Extreme.Mathematics.SignalProcessing
# namespace.
from Extreme.Mathematics.SignalProcessing import *

# Illustrates the use of the FftProvider and Fft classes for computing 
# the Fourier transform of real and complex signals.

# This QuickStart sample shows how to compute the Fouier
# transform of real and complex signals.

# Some vectors to play with:
r1 = Vector([ 1.0 / (1 + i) for i in xrange(1000)])
c1 = ComplexVector.Create(1000, lambda i: DoubleComplex(sin(0.03 * i), cos(0.07 * i)))

r2 = Vector([ 1.0, 2.0, 3.0, 4.0 ])
c2 = ComplexVector([ 1+2j, 3+4j, 5+6j, 7+8j ])

#
# One-time FFT's
#

# The Vector and ComplexVector classes have static methods to compute FFT's:
c3 = Vector.FourierTransform(r2)
r3 = Vector.InverseFourierTransform(c3)
print "fft(r2) = {0:.3f}".format(c3)
print "ifft(fft(r2)) = {0:.3f}".format(r3)
# The ComplexConjugateSignalVector type represents a complex vector
# that is the Fourier transform of a real signal. 
# It enforces certain symmetry properties:
print "c3[i] == conj(c3[N-i]): {0} == conj({1})".format(c3[1], c3[3])

#
# FFT Providers
#

# FFT's require a fair bit of pre-computation. Using the FftProvider class, # you can get an Fft object that caches these computations.

# Here, we create an FFT implementation for a real signal:
realFft = FftProvider.ManagedProvider.Create1DRealFft(r1.Length)
# For a complex to complex transform:
complexFft = FftProvider.ManagedProvider.Create1DComplexFft(c1.Length)

# You can set the scale factor for the forward transform.
# The default is 1/N.
realFft.ForwardScaleFactor = 1.0 / sqrt(c1.Length)
# and the backward transform, with default 1:
realFft.BackwardScaleFactor = realFft.ForwardScaleFactor

# The ForwardTransform method performs a forward transform:
c4 = realFft.ForwardTransform(r1)
print "First 5 terms of fft(r1):"
for i in range(0,5):
    print "   {0}: {1}".format(i, c4[i])
c4 = complexFft.ForwardTransform(c1)
print "First 5 terms of fft(c1):"
for i in range(0,5):
    print "   {0}: {1}".format(i, c4[i])

# ForwardTransform has many overloads for real to complex and
# complex to complex transforms.

# A one-sided transform returns only the first half of the FFT of
# a real signal. The rest can be deduced from the symmetry properties.
# Here's how to compute a one-sided FFT:
c5 = ComplexVector.Create(r1.Length / 2 + 1)
realFft.ForwardTransform(r1, c5, RealFftFormat.OneSided)

# The BackwardTransform method has a similar set of overloads:
r4 = Vector.Create(r1.Length)
realFft.BackwardTransform(c5, r4, RealFftFormat.OneSided)

# As the last step, we need to dispose the two FFT implementations.
realFft.Dispose()
complexFft.Dispose()

#
# 2D transforms
#

# 2D transforms are handled in a completely analogous way.
m = Matrix.Create(36, 56, lambda i,j: exp(-0.1 * i) * sin(0.01 * (i * i + j * j - i * j)))
mFft = ComplexMatrix.Create(m.RowCount, m.ColumnCount)

with FftProvider.CurrentProvider.Create2DRealFft(m.RowCount, m.ColumnCount) as fft2:
    fft2.ForwardTransform(m, mFft)

    print "First few terms of fft(m):"
    print format(mFft[0:5,0:5], ".3f")

    # and the backward transform:
    fft2.BackwardTransform(mFft, m)

    # Dispose is called automatically.
```