# Fast Fourier Transformation (FFT)

## Filtering images

There are two main approaches to filtering an image. The first is the convolution of an image and kernel in the spatial domain. The second is the multiplication of an image's Fourier transform with a filter in the frequency domain. The Fast Fourier transformation (FFT) algorithm, which is an example of the second approach, is used to obtain a frequency-filtered version of an image. Processing images by filtering in the frequency domain is a three-step process:

1. Perform a forward fast Fourier transform to convert a spatial image to its complex fourier transform image. The complex image shows the magnitude frequency components. It does not display phase information.
2. Enhance selected frequency components of the image and attenuate other frequency components of the image by multiplying by a low-pass, high-pass, band-pass, or band-stop filter. In MIPAV, you can construct frequency filters using one of three methods: finite impulse response filters constructed with Hamming windows, Gaussian filters, and Butterworth filters. However, for the Gaussian filters, only low-pass and high-pass filters are available.
3. Perform an inverse fast fourier transform to reconvert the image from the frequency domain to the spatial domain.

You can also use Algorithms > Filters (frequency) to obtain a frequency-filtered version of an image in a single step. A contrast of the FFT and Filter (frequency) commands on the Algorithms menu is helpful in determining which one to select.

Use . . .
To . . .
Algorithms > Filters (frequency)
Obtain a frequency-filtered version of an image in a single step. All three steps in the process are performed internally and only the resultant filtered image is displayed.
Algorithms > FFT
Produce pictures of the fourier transform of the image and the filtered fourier transform of the image.

### Background

Fast Fourier Transformation (FFT) algorithm:background;

Since all dimensions of an n-dimensional dataset must be powers of 2, datasets with arbitrary dimensions are zero padded to powers of 2 before applying the Fast Fourier Transform and stripped down to the original dimensions after applying the Inverse Fast Fourier Transform.

If real source data is used, then the output obtained after the final Inverse Fast Fourier Transform step is purely real, no matter what filter construction method is used. Any non-zero imaginary part present after the Inverse Fast Fourier Transform is due to Randolph error. Hence, the imaginary part is discarded and only the real part is imported into the final image. Figure 1 shows the application of the FFT algorithm to a 2D image.

In MIPAV, frequency filters may be constructed using one of three methods: finite impulse response filters constructed with Hamming windows, Gaussian filters, and Butterworth filters. You select the method when you enter the parameters for the algorithm.

#### Finite Impulse Response (FIR) filters with Hamming windows

If you use FIR filters with Hamming windows, an ideal filter kernel is constructed in the spatial domain. A Hamming window kernel is also constructed in the spatial domain. Then, the ideal filter kernel and the Hamming window kernel are multiplied together in the spatial domain. The resultant kernel is a circle or odd diameter sphere with its center at (kernel diameter - 1)/2 in every dimension. The kernel is zero padded up to the same dimensions as the original image data.

The two fourier transforms (image and filter) are multiplied, and the inverse fourier transform is obtained. The result is a filtered version of the original image shifted by (kernel diameter - 1)/2 toward the end of each dimension. The data is shifted back by (kernel diameter - 1)/2 to the start of each dimension before the image is stripped to the original dimensions.

In every dimension the length of the resulting fourier image must be greater than or equal to the length of non-zero image data kernel diameter - 1 so that aliasing does not occur. In aliasing, higher frequency components incorrectly shift to lower frequencies. Thus, if the image is not cropped, each frequency space dimension must be greater than or equal to the minimum power of two that equals or exceeds the original dimension size kernel diameter 1. If cropping is selected, the fourier image is only increased to the minimum power of 2 that equals or exceeds the original dimension size rather than the minimum power 2 that equals or exceeds the original dimension size kernel diameter - 1. With cropping, a band of pixels at the boundaries of the image is set equal to 0. A power of 2 increase in dimension length is avoided at the cost of losing edge pixels.

The following formula defines the discrete-space fourier transform pair:

$X(\omega_1, \omega_2) = \sum_{n_{1} =-\infty}^\infty \sum_{n_2=-\infty}^\infty x(n_1,n_2)e^{-j\omega_1 n_1} e^{-j\omega_2n_2}$

$x(n_1n_2) = \frac {1}{(2\pi)^2}\int_{\omega_1}^\pi=-\pi\int_{w_2}^\pi=-\pi^{X(\omega_1,\omega_2) e^{-j\omega_1 n_1}e^{-j\omega_2 n_2}}\partial\omega_1\partial\omega_2$

While centering is invoked in FFT to provide a clearer display, centering is not used in frequency filter where no fourier images are produced.

#### Gaussian and Butterworth filters

In the frequency filtering step, both the real and imaginary parts of the fourier transform are scaled by the same formula, which depends only on the distance from the center of image.

The FIR filter must create two fourier transforms: one for the original image and one for the kernel. The Butterworth and Gaussian filters only need to create one fourier transform for the image since frequency scaling is done with a formula. Therefore, the FIR filter uses about twice as much memory as the Gaussian or Butterworth filters.

Centering is needed for the Gaussian and Butterworth filters since the scaling formulas are based on frequency magnitude increasing from the center of the image.

##### Butterworth 2D and 3D filters
• Butterworth filters:2D;
• Butterworth filters:3D;

The following are the formulas for the Butterworth 2D and 3D filters, respectively.'

$distsq = \frac {(x-xcenter)^2}{(xcenter)^2} + \frac {(y-ycenter)^2}{(ycenter)^2}$

$distsq = \frac {(x-xcenter)^2}{(xcenter)^2} + \frac {(y-ycenter)^2}{(ycenter)^2}+ \frac {(z-zcenter)^2}{(zcenter)^2}$

The real data is the product of the real data and the coefficient. The imaginary data is the product of the imaginary data and the coefficient.

##### Butterworth low-pass, high-pass, band-pass, and band-stop filters
• Butterworth filters:high-pass;
• Butterworth filters:low-pass;
• Butterworth filters:band-stop;
• Butterworth filters:band-stop;

The following are the equations for the Butterworth low-pass, high-pass, band-pass, and band-stop filters.

###### Butterworth low-pass filter

$coefficient = \frac {1}{1 +\left ( \frac {distsq}{(f_1)} \right )^{2butterworthorder}}$

###### High-pass filter

$coefficient = \frac {1}{1 +\left ( \frac {(f_1)}{distsq} \right )^{2butterworthorder}}$

###### Band-pass filter

$width = f_2-f_1centersq = f_1f_2$

$num = (distsq-centersq)^{2butterworthorder}$

$coeff = \frac {num} {(num + (distsq-centersq)^{2butterworthorder})}$

###### Band-stop filter

$width = f_2=f_1,centersq= f_1f_2$

$num = ({\sqrt distsq}{width})^{2butterworthorder}$

$coeff= \frac {num} {\left ( num + ({\sqrt distsq}{width})^{2butterworthorder} \right )}$

##### Gaussian low- and high-pass filters

The following equations are for the Gaussian low- and high-pass filters.

###### Gaussian low-pass (2D)

$coefficientforgaussianhighpass= 1-gaussianlowpass$

###### Gaussian low-pass (3D)

$coefficient = {e^{-\frac {(x-xcenter)^2}{(2f^2xcenter^2)}}}{e^{-\frac {(y-ycenter)^2} {(2f^2ycenter^2)}}} {e^{-\frac {(z-zcenter)^2} {(2f^2zcenter^2)}}}$

 Note: The routine bessj1 for determining a bessel function of the first kind and first order is taken from "Numerical Recipes in C."
• i.routine bessj1;

#### Image types

• Fast Fourier Transformation (FFT) algorithm:image types;
• image types:Fast Fourier Transformation (FFT) algorithm;

You cannot apply a forward FFT to complex data. However, you can apply the frequency filter or perform an inverse FFT on conjugate, symmetric, and complex data. You cannot apply this algorithm to Alpha Red Green Blue (ARGB) datasets, but you can apply it to other gray-scale 2D and 3D images.

#### Special notes

By default, the resultant image is a float-type image. Intermediate fourier transform images and filtered fourier transform images obtained in FFT are complex images.