Fourier Transform in Image Processing Review Image as

Fourier Transform in Image Processing

Review - Image as a Function • We can think of an image as a function, f, • f: R 2 R – f (x, y) gives the intensity at position (x, y) – Realistically, we expect the image only to be defined over a rectangle, with a finite range: • f: [a, b]x[c, d] [0, 1] • A color image is just three functions pasted together. • We can write this as a “vector-valued” function:

Image as a Function • function, f, • f: R 2 R – f (x, y) gives the intensity at position (x, y)

Linear Transformations on Images

Transformations on Images • Define a new image g in terms of an existing image f – We can transform: 1. either the domain 2. or the range of f • Range transformation: transformation 1. What kinds of operations on the image can this transformation perform? 2. What t transformations are good for what?

Operations that change the domain of image 1. Some operations preserve the range but change the domain of f 2. Here is an example that changes the domain : – What kinds of operations can this perform? 3. Still other operations operate on both the domain and the range of f. 1. 2. In general we are interested in all these transformations Specifically we want them to have some nice (like linear) properties

Linear Shift Invariant Systems (LSIS) Linearity: Shift invariance:

Example of LSIS: ideal lens Linear Shift Invariant Systems (LSIS) Defocused image ( g ) is a processed version of the focused image ( f ) Ideal lens is a LSIS Linearity: Brightness variation Shift invariance: Scene movement (this is not valid for lenses with non-linear distortions)

A Reminder of a Convolution concept

LSIS is Doing a Convolution convolution is linear and shift invariant We negate (mirror) the kernel of convolution h Linear Shift Invariant Systems (LSIS) 1. Continuous convolution defined by integral 2. Discrete convolution from previous lectures defined by finite sums kernel h

Convolution – two more Examples Here we show two more examples of convolution of functions f and g We use two nonidea approximation of Dirac’s Delta function as a signal model From Eric Weinstein’s Math World

Convolution – Another Example 1 -1 1 Here we convolve signals a and b into convolved signal c 1 1 -2 -1 1 2

Convolution Kernel – Impulse Response • We ask question “What h will give us g = f ? ” Dirac Delta Function (Unit Impulse) Sifting property: We convolve f with Dirac Delta: And we get this Therefore g=f when we convolve with delta

Point Spread Function Now we apply information from last slide to optical system scene Optical System image • Ideally, the optical system should be a Dirac delta function. (no distortion) • However, optical systems are never ideal. point source Optical System • Point spread function of Human Eyes point spread function PSF stands for Point Spread Function If the PSF of human eye is different (not a good approximation of delta) then we have eye troubles, related to this PSF type

Point Spread Function PSF stands for Point Spread Function normal vision myopia hyperopia Tell me what your PSF is and I will tell you what is your eye problem. astigmatism Images by Richmond Eye Associates

Properties of Convolution • Commutative • Associative Convolution operator This is unlike matrix multiplication and Kronecker (tensor) products • Cascade system We used in cascading filters in previous lectures. Now we have a theoretical base These are general properties of convolution , not related to signal representation

• So now we know what is convolution and why it is so important. • However signals for convolution can be represented in various ways. • So the question is, “how to represent signals? ” – 1 D – 2 D – 3 D?

How to Represent Signals? • Option 1: Taylor series represents any function using polynomials. • Polynomials are not the best - unstable and not very physically meaningful. (they still have some applications) • Easier to talk about “signals” in terms of its “frequencies” (how fast/often signals change, etc).

Mr. Jean Baptiste Fourier and Fourier Transform idea

Jean Baptiste Joseph Fourier (1768 -1830) • Had crazy idea (1807): • Any periodic function can be rewritten as a weighted sum of Sines and Cosines of different frequencies. • Don’t believe it? – Neither did Lagrange, Laplace, Poisson and other big wigs – Not translated into English until 1878! • But it’s true! – called Fourier Series – Possibly the greatest tool used in Engineering

A Sum of Sinusoids • Our building block: • Add enough of them to get any signal f(x) you want! • How many degrees of freedom? • What does each control? • Which one encodes the coarse vs. fine structure of the signal?

Math Review - Complex numbers • Real numbers: 1 -5. 2 • Complex numbers 4. 2 + 3. 7 i 9. 4447 – 6. 7 i -5. 2 (-5. 2 + 0 i) We often denote in EE i by j

Math Review - Complex numbers • Complex numbers 4. 2 + 3. 7 i 9. 4447 – 6. 7 i -5. 2 (-5. 2 + 0 i) • Amplitude A = | Z | = √(a 2 + b 2) • Phase = Z = tan-1(b/a) • General Form Z = a + bi Re(Z) = a Im(Z) = b Real and imaginary parts

Math Review – Complex Numbers • Polar Coordinate Z = a + bi A • Amplitude A= √(a 2 + b 2) • Phase = tan-1(b/a) a b

Math Review – Complex Numbers and Cosine Waves • Cosine wave has three properties – Frequency – Amplitude – Phase • Complex number has two properties – Amplitude – Wave • Complex numbers to represent cosine waves at varying frequency – Frequency 1: – Frequency 2: – Frequency 3: Z 1 = 5 +2 i Z 2 = -3 + 4 i Z 3 = 1. 3 – 1. 6 i Simple but great idea !!

Fourier Transform Idea • We want to understand the frequency w of our signal. So, let’s reparameterize the signal by w instead of x: Fourier Transform f(x) F(w) • For every w from 0 to infinity, F(w) holds the amplitude A and phase f of the corresponding sine – How can F hold both? F(w) Complex number trick! Inverse Fourier Transform f(x) All info preserved

Math Review Periodic Functions – more illustrations to our formalism

Math Review - Periodic Functions If there is some a, for a function f(x), such that f(x) = f(x + na) then function is periodic with the period a a 0 2 a 3 a

Math Review - Attributes of cosine wave f(x) = cos (x) Amplitude Phase f(x) = 5 cos (x) f(x) = 5 cos (x + 3. 14)

Math Review - Attributes of cosine wave Amplitude Phase f(x) = 5 cos (x) f(x) = 5 cos (x + 3. 14) Frequency f(x) = 5 cos (3 x + 3. 14)

Math Review - Attributes of cosine wave f(x) = cos (x) f(x) = A cos (kx + ) Amplitude, Frequency, Phase This math is base of Fourier Analysis and Fourier Synthesis

Time and Frequency • example : g(t) = sin(2 f t) + (1/3)sin(2 (3 f) t) Our periodic signal is a mixture of two frequencies Before we go to using complex numbers and formulas let us get intuition about time and frequency

Time and Frequency • example : g(t) = sin(2 f t) + (1/3)sin(2 (3 f) t) = +

Frequency Spectra • example : g(t) = sin(2 pf t) + (1/3)sin(2 p(3 f) t) = + In time domain In frequency domain

Frequency Spectra • Usually, frequency is more interesting than the phase 1. 2. We will decompose it to frequencies, we will add one more coefficient based base function in each slide

Frequency Spectra First base function = = second base function + Sum of first and second base functions is already not bad approximation

Frequency Spectra = = + Sum of first , second and third base functions is even better approximation

Frequency Spectra = = +

Frequency Spectra = = +

Frequency Spectra = = +

Frequency Spectra = In time domain In frequency domain

Frequency Spectra Observe the negative values of domain that appear in spectrum

Fourier Transform is Just a change of basis Important concept

Fourier Transform is Just a change of basis M * f(x) = F(w) * . . . M can be a matrix, a decision diagram, a butterfly, a system (physical or simulated) =

IFT: Just a change of basis M-1 * F(w) = f(x) * . . . Inverse Fourier Transform is similarly only a change of basis =

Fourier Transform is an invertible operator Image Fourier Transform FT v 2 will display image or its transform, Matlab examples.

Fourier Transform is an invertible operator Image Fourier Transform Ny ⁄ 2 Ny y F(kx, ky) f(x, y) 0 Nx ⁄ 2 Nx x {F(kx, ky)} = f(x, y)} {f(x, y)} = F(kx, ky) Similarly as for Hadamard-Walsh or Reed -Muller Transforms we can talk about transform pairs.

Continuous Fourier Transform

Continuous Fourier Transform f(x) = F(s) Euler’s Formula Various normalizations of pairs

Some Conventions for Fourier Transforms • Image Domain • Fourier Domain – – – • Forward Transform • f(x, y, z) • g(x) • F Reciprocal space Fourier Space K-space Frequency Space Spectral domain • Reverse Transform, Inverse Transform • F(kx, ky, kz) • G(s) • F

Fourier Analysis Decompose f(x) into a series of cosine waves that when summed reconstruct f(x)

Fourier Analysis in 1 D. Audio signals Amplitude Only 5 10 15 (Hz) 5 10 (Hz) 15

Fourier Analysis in 1 D. Audio signals 5 10 (Hz) Your ear performs fourier analysis. 15

Fourier Analysis in 1 D. Spectrum Analyzer. i. Tunes performs fourier analysis. 1. 2. We know it from our Hi-Fi equipment The same can be done with standard images or radar images or any other images, in any number of dimensions

Fourier Synthesis Summing cosine waves reconstructs the original function

Example: Fourier Synthesis of Boxcar Function Boxcar function Periodic Boxcar Can this function be reproduced with cosine waves?

k=1. One cycle period A 1·cos(2 kx + 1) k=1 1 A ·cos(2 kx + ) k k=1 k

k=2. Two cycles period A 2·cos(2 kx + 2) k=2 2 A ·cos(2 kx + ) k k=1 k

k=3. Three cycles period A 3·cos(2 kx + 3) k=3 3 A ·cos(2 kx + ) k k=1 k

Fourier Synthesis. N Cycles A 3·cos(2 kx + 3) k=3 N A ·cos(2 kx + ) k k=1 k

Example: Fourier Synthesis of a 2 D Function An image is two dimensional data. Intensities as a function of x, y White pixels represent the highest intensities. Greyscale image of iris 128 x 128 pixels

Fourier Synthesis of a 2 D Function F(2, 3) 1. We select coefficients, multiply them by weights and add 2. We can create artificial images

Fourier Filters • Fourier Filters in 1 D data (sound) change the image by changing which frequencies of cosine waves go into the image • Represented by 1 D spectral profile • 2 D Specteral Profile is rotationally symmetrized 1 D profile

Low and High Frequency Image Filters • Low frequency terms – Close to origin in Fourier Space – Changes with great spatial extent (like ice gradient), or particle size • High frequency terms – Closer to edge in Fourier Space – Necessary to represent edges or highresolution features

Frequency-based Filters 1. Low-pass Filter (blurs) – Restricts data to low-frequency components 2. High-pass Filter (sharpens) – Restricts data to high-frequency-components 3. Band-pass Filter – Restrict data to a band of frequencies 4. Band-stop Filter – Suppress a certain band of frequencies We will show many practical examples soon

Cutoff Low-pass Filter Image is blurred Sharp features are lost Ringing artifacts Filter characteristics

Butterworth Low-pass Filter Flat in the pass-band Zero in the stop-band No ringing

Gaussian Low-pass Filter

Butterworth High-pass Filter • Note the loss of solid densities

How the filter looks in 2 D unprocessed bandpass lowpass highpass spectra

Convolution again Our goal is now to approach the most important results of spectral processing

Convolution of some function f(x) with some kernel g(x) Continuous Discrete * =

Convolution in 2 D image kernel result x x = x x x

Microscope Point-Spread. Function is Convolution Point-Spread-Function

Convolution Theorem f g= {FG} f= FG G Convolution in image domain Is equivalent to multiplication in fourier domain Can be derived from above because FG/G = F and f = F-1 (F)

Lowpass Filtering by Convolution f g= • {FG} Camera shake We showed examples of this using standard convolution, but it is done in modern equipment through digitally calculated FFT and IFFT 2 D gaussian is example of kernel 2 D gaussian as a matrix of a kernel

Lens Performs Fourier Transform

Contrast Theory Power spectrum PS = F 2(s) CTF 2(s) Env 2(s) + N 2(s) Incoherant average of transform obs(x) = f(x) psf(x) env(x) + n(x) observed image noise f(x) for true particle envelope function point-spread function

Gibbs Ringing • 5 waves • 25 waves • 125 waves Even with 125 waves we still have some small Gibbs Ringing

Fourier Transform Math Properties

Fourier Transform – more formally Represent the signal as an infinite weighted sum of an infinite number of sinusoids Note: Arbitrary function Single Analytic Expression Spatial Domain (x) Frequency Domain (u) (Frequency Spectrum F(u)) Inverse Fourier Transform (IFT) IFT

Fourier Transform • Also, defined as: Fourier Transform Note: • Inverse Fourier Transform (IFT)

Fourier Transform Pairs (I) Fourier Transform Note that these are derived using angular frequency ( )

Fourier Transform Pairs (I) Fourier Transform Note that these are derived using angular frequency ( )

Fourier Transform Pairs (I) Note that these are derived using angular frequency ( )

Fourier Transform and Convolution

Fourier Transform and Convolution Let Then From definition Here we derive the most important formula of image processing and other applications of Fourier Transform Add and subtract tau Separate domains Convolution in spatial domain Multiplication in frequency domain

MAIN THEOREM OF MODERN IMAGE PROCESSING: Convolution in spatial domain Multiplication in frequency domain

Fourier Transform and Convolution Spatial Domain (x) Frequency Domain (u) So, we can find g(x) by Fourier transform IFT FT FT

Properties of Fourier Transform

Properties of Fourier Transform 1 Spatial Domain (x) Frequency Domain (u) Linearity Scaling Shifting Symmetry Conjugation Convolution Differentiation Note that these are derived using frequency ( )

Properties of Fourier Transform 2 Star will be used for conjugation in complex plane Substitute f instead of g in this formula to get the upper formula For instance, when f is R then F is RE, IO, . etc 1. 2. 3. 4. 5. 6. RO – real odd IE – imaginary even IO – imaginary odd RE – Real even CE – Complex even CO – Complex odd We explained even-odd properties for Walsh, we will explain again for other transforms

Example use of Fourier: Smoothing/Blurring 1. We want a smoothed function of f(x) But we do this in spectral domain • 2. Let us use a Gaussian kernel • 3. Then H(u) attenuates high frequencies in F(u) (Low-pass Filter)!

Our Images are Discrete so discrete Fourier applies

Image as a Discrete Function

Digital Images The scene is – projected on a 2 D plane, – sampled on a regular grid, and each sample is – quantized (rounded to the nearest integer) Image as a matrix

Review Fourier Transform is invertible operator Math Review Periodic functions Amplitude, Phase and Frequency Complex number Amplitude and Phase Fourier Analysis (Forward Transform) Decomposition of periodic signal into cosine waves Fourier Synthesis (Inverse Transform) Summation of cosine waves into multi-frequency waveform Fourier Transforms in 1 D, 2 D, 3 D, ND Image Analysis Image (real-valued) Transform (complex-valued, amplitude plot) Fourier Filters Low-pass High-pass Band-stop Convolution Theorem Deconvolute by Division in Fourier Space All Fourier Filters can be expressed as real-space Convolution Kernels Lens does Fourier transforms Diffraction Microscopy

Where we are in this class? • • • We defined continuous Fourier Transform We defined discrete Fourier Transform DFT is calculated by matrix multiplication It is slow We will need factorizations, Kroneckers, recursions and butterflies to calculate it faster. • Only then FT will become a really useful tool, FFT.

Further Reading • Wikipedia • Mathworld • The Fourier Transform and its Applications. Ronald Bracewell

Filtering with EMAN 2 Low. Pass Filters filtered=image. process(‘filter. lowpass. guass’, {‘sigma’: 0. 10}) filtered=image. process(‘filter. lowpass. butterworth’, {‘low_cutoff_frequency’: 0. 10, ‘high_cutoff_frequency’: 0. 35}) filtered=image. process(‘filter. lowpass. tanh’, {‘cutoff_frequency’: 0. 10, ‘falloff’: 0. 2}) High. Pass Filters filtered=image. process(‘filter. highpass. guass’, {‘sigma’: 0. 10}) filtered=image. process(‘filter. highpass. butterworth’, {‘low_cutoff_frequency’: 0. 10, ‘high_cutoff_frequency’: 0. 35}) filtered=image. process(‘filter. highpass. tanh’, {‘cutoff_frequency’: 0. 10, ‘falloff’: 0. 2}) Band. Pass Filters filtered=image. process(‘filter. bandpass. guass’, {‘center’: 0. 2, ‘sigma’: 0. 10}) filtered=image. process(‘filter. bandpass. butterworth’, {‘low_cutoff_frequency’: 0. 10, ‘high_cutoff_frequency’: 0. 35}) filtered=image. process(‘filter. bandpass. tanh’, {‘cutoff_frequency’: 0. 10, ‘falloff’: 0. 2})

Sources of slides 1. S. Narasimhan 2. Mike Marsh National Center for Macromolecular Imaging Baylor College of Medicine