G-Class of Medical Image Information

Fourier Transform Ⅱ

September 22, 2009

Shigehiko Katsuragawa, Ph.D.

 Wave propagating in Real Domain

Signal in an image:
Optical density on radiograph
Brightness on display monitor

Spatial Frequency is defined as

the number of waves in a unit
length
Signal intensity

Distance (unit of spatial freq.):

mm (cycles/mm)
Distance cm (cycles/cm)
Expression of Function in Spatial Freq. Domain

Real Domain Spatial Freq. Domain

Periodic function Fourier Series

Non-periodic function Fourier Transform

f (x)

Signal
f (x)
Signal

x (mm) x (mm)

Medical images are 2D non-periodic functions.

 Complex Fourier Series

The complex Fourier series of f(x) is defined as follows;

∞
f ( x) = ∑
n = −∞
cn e inω0 x (2.19)

The complex Fourier coefficient, cn, is obtained as follows;

1 L/2
cn = ∫ f ( x)e −inω0 x dx
L −L / 2 (2.20)
2π
where ω0 = , L is a period of periodic function, f ( x).
L
f(x) L (cm)

0 x (cm)
 Fourier Transform

∞ 1 ∞ 1 ∞
∫ ∫ ∫
−iωξ iωx
f ( x) = [ f (ξ )e dξ ]e dω = F (ω )e iωx dω
− ∞ 2π −∞ 2π −∞

∞
F (ω ) = ∫ −∞
f ( x)e −iωx dx Fourier Transform (2.24)

1 ∞
f ( x) =
2π ∫
−∞
F (ω )e iωx dω Inverse Fourier Transform (2.25)

ω = 2πu , dω = 2πdu, ω : angular freq., u : spatial freq.

∞
F (u ) = ∫ f ( x)e −i 2πux dx
−∞
∞
f ( x) = ∫ F (u )e i 2πux du
−∞
 Fourier Transform

f(x) F(u)

FT

Inverse FT

x (cm) u (cycles/cm)

Function in real domain Spectrum in spatial frequency domain

Function, f(x), includes the distribution of spatial frequency components

(spectrum) as shown by F(u).
 Fourier Transform of Image

FT

Inverse FT

Image in real domain Spectrum in spatial frequency domain

 Parseval’s Theorem and Convolution Integral

∞ 1 ∞
∫ f ( x) dx = ∫ F (ω ) dω ; Parseval' s Theorem
2 2
−∞ 2π −∞

F (ω )
2
; Power Spectrum

The convolution integral of f(x) and g(x) is defined as follows,

∞
f ( x) * g ( x) = ∫−∞
f (τ ) g ( x − τ )dτ

If F (ω ) = ℑ[ f ( x)] and G (ω ) = ℑ[ g ( x)] ,

ℑ[ f ( x) * g ( x)] = F (ω )G (ω ) convolution integral theorem

FT of convolution of f(x) and g(x) is equal to the multiplication of each FT.

 Fourier Transform of Rectangle Pulse 1

f(x)

1/2d ⎧1 / 2d ( x ≤ d)
f ( x) = ⎨
⎩0 ( x > d)
-d d x
∞
From (2.24), F (ω ) = ∫
−∞
f ( x)e −iωx dx
1 −iωx
d 1 1 −iωx d
= ∫
− d 2d
e dx = −
2 d iω
[e ]− d

1 e i ωd − e − i ωd
=
ωd 2i
sin( ωd ) e iθ − e −iθ
= (sin(θ ) = )
ωd 2i
= sin c(ωd )
 Fourier Transform of Rectangle Pulse 2

f(x) F(ω)
1
1/2d
FT

-d d ω
x
−2π/d −π/d
π/d 2π/d

⎧1 / 2d ( x ≤ d) sin( ωd )
f ( x) = ⎨ F (ω ) = = sin c(ωd )
⎩0 ( x > d) ωd

Rectangle Pulse Sinc Function

Relation between Width of Rectangle Pulse and Sinc Function

f(x) F(ω)
1

FT

ω
x
wide width narrow width
(low freq. component)
F(ω)
f(x)
1

FT

ω
x
narrow width wide width
(high freq. component)
 Fourier Transform of Sinc Function
F(x) f(ω)
1
1/2d
FT

x -d d ω
−2π/d −π/d
π/d 2π/d
sin( xd ) ⎧1 / 2d (ω ≤ d)
F ( x) = = sin c( xd ) f (ω ) = ⎨
xd ⎩0 (ω > d)

Sinc Function Rectangle Pulse

From symmetric property of Fourier transform,

if F (ω ) = ℑ[ f ( x)],
ℑ[ F ( x)] = f (−ω ) = f (ω ) ( f (ω ) : even function)
sin( xd ) ⎧1 / 2d (ω ≤ d)
Therefore, ℑ[ ] = f (−ω ) = f (ω ) = ⎨
xd ⎩ 0 (ω > d)
 Fourier Transform of Delta δ-Function 1

δ(x) δ(x-a) Definition of δ-function

∞
∫
−∞
δ ( x)dx = 1
∞

0 a
x ∫
−∞
f ( x)δ ( x)dx = f (0)
∞
∫
−∞
f ( x)δ ( x − a )dx = f (a )

Fourier Transform of δ-function

∞
∫
−iωx − iω 0
F (ω ) = δ ( x )e dx = e = e =1
0
−∞
 Fourier Transform of δ-Function 2

F(ω )
δ(x)
FT
1

x 0 ω
0

The δ-function includes all frequency component.

The spectrum including all frequency component
with constant value is called “white noise”.

δ(x-x0)
∞
ℑ[δ ( x − x0 )] = ∫
−∞
δ ( x − x0 )e −iωx dx = e −iωx0
0 x0 x
 Fourier Transform of δ-Function Sequence 1
f(x)

δ-function sequence
(or comb function):

Periodic function of δ-function

with period, x0
0 x
x0

f ( x) = L + δ ( x + 2 x0 ) + δ ( x + x0 ) + δ ( x) + δ ( x − x0 ) + +δ ( x − 2 x0 )L
∞
= ∑ δ ( x − nx )
n = −∞
0
 Fourier Transform of δ-Function Sequence 2

f(x) can be expressed by a Fourier series, because f(x) is periodic.

∞ ∞
f ( x) = ∑
n = −∞
cn e inω0 x
= ∑
n = −∞
cn e i 2πnu0 x (ω0 = 2πu0 )

where, u0 = 1 / x0
1 x0 / 2
cn =
x0 ∫− x0 / 2
f ( x)e −i 2πnu0 x dx

1 x0 / 2 1 0
=
x0 ∫ − x0 / 2
δ ( x)e −i 2πnu0 x dx =
x0
e

1
= = u0
x0

Therefore,
∞ ∞
f ( x) = ∑u e
n = −∞
0
i 2πnu0 x
= u0 ∑e
n = −∞
i 2πnu0 x
 Fourier Transform of δ-Function Sequence 3
∞
F (u ) = ℑ[u0 ∑e
n = −∞
i 2πnu0 x
]

∞
= u0 ∑ ℑ[e
n = −∞
i 2πnu0 x
]

ℑ[δ ( x − u0 )] = e −iωu0 = e −i 2πu0u

ℑ[δ ( x + u0 )] = e i 2πu0u
ℑ[δ ( x + nu0 )] = e i 2πnu0u

From symmetric property of FT, if F (u ) = ℑ[ f ( x)], ℑ[ F ( x)] = f (−u )

Therefore,
ℑ[e i 2πnu0 x ] = δ (−u + nu0 ) = δ (u − nu0 ) (even function: δ (U ) = δ (−U ))

Finally, ∞
F (u ) = u0 ∑ δ (u − nu )
n = −∞
0
 Fourier Transform of δ-Function Sequence 4

∞ ∞
ℑ[ ∑ δ ( x − nx )] = u ∑ δ (u − nu ),
n = −∞
0 0
n = −∞
0 where u0 = 1 / x0

f(x) F(u)

FT

x u
0 0
x0 u0 =1/x0
∞ ∞
f ( x) = ∑ δ ( x − nx )
n = −∞
0 F (u ) = u0 ∑ δ (u − nu )
n = −∞
0

δ-function sequence δ-function sequence

 System for Image Formation

X-ray X-ray X-ray

Object Object Object

Screen
Imaging Plate Flat Panel Detector
Film
Digitization Digitization
Film Processing
system system system
Image Image Image

Screen/Film CR DR
 Input-Output System for Image Formation

input output
System
transmitted x-ray image
f(x) g(x)

The system includes screen/film or imaging plate or flat panel detector.

g ( x ) = L[ f ( x )]
where L indicate a conversion rule from f(x) to g(x).

If we can know the conversion rule, L, we can extrapolate input

from output accurately.

It is very difficult to determine L of a general system.

However, L of a linear system can be determined easily.
 Linear System

input output
System
f(x) g(x)

Additivity

if g1 ( x) = L[ f1 ( x)], and g 2 ( x) = L[ f 2 ( x)],

L[a1 f1 ( x) + a2 f 2 ( x)] = a1 g1 ( x) + a2 g 2 ( x)

Steadiness

L[ f ( x − x1 )] = g ( x − x1 )

If the system has additivity and steadiness,

the system is called a linear system.
 Is the screen/film system a linear system?

Film Characteristic Curve

The screen/film system is not linear,

optical density
because film has not additivity.

However, if the optical density can be

converted to exposure by using a film
characteristic curve, we can handle the
screen/film system as a linear system.

exposure
 Linear System Response 1

f(x) g(x)
Linear System

By using δ-function, input f(x) can be expressed as follows,

∞ ∞
f ( x) = ∫ f (τ )δ ( x − τ )dτ because ∫ f ( x)δ ( x)dx = f (0)
−∞ −∞

δ(x) h(x)

Consider the impulse response,

h( x) = L[δ ( x)]
x x
Linear System

Impulse Response
 Linear System Response 2

f(x) g(x)
Linear System

g ( x) = L[ f ( x)]
∞
= L[ ∫ f (τ )δ ( x − τ )dτ ]
−∞
∞
= ∫ f (τ ) L[δ ( x − τ )]dτ
−∞
∞
= ∫ f (τ )h( x − τ )dτ
−∞

= f ( x ) * h( x )
 Linear System Response 3

In summary, output of a linear system can be obtained from convolution

integral of input and impulse response as follows,

g ( x ) = f ( x ) * h( x )
if F (ω ) = ℑ[ f ( x)], G (ω ) = ℑ[ g ( x)], and H (ω ) = ℑ[h( x)],
from convolution integral theorem,
G (ω ) = F (ω ) H (ω )

h(x)
f (x ) g ( x) = f ( x) * h( x ) :in real domain
Linear System
F (ω ) G (ω ) = F (ω ) H (ω ) :in frequency domain
H (ω )
 Linear System Response in Real Domain

h(x)
f (x ) g ( x) = f ( x ) * h( x)
Linear System

input：f(x) output：g(x)
f(τ)δ(x-τ) f(τ)h(x-τ)

x x

The output of a linear system is superposition of many impulse responses.

System Function (Frequency Response Function) 1

δ(x) h(x)

x x
Linear System h(x): impulse response

The Fourier transform of impulse response, h(x), is called

the system function, H(ω).

The impulse response of screen/film is obtained by a slit

image (line spread function) with very narrow width.

|H(ω)| is called MTF or the response function to evaluate

the frequency response of contrast on medical images .
 Impulse Response of Screen/Film

slit slit image

δ(x) h(x)

line spread function

x x
Linear System

ℑ[δ ( x)]
H (u ) MTF (modulation
transfer function)
u
u
 Modulation Transfer Function (MTF) 1

δ (x) h(x)
Linear System
H (u ) = ℑ[ h( x )]
ℑ[δ ( x)]
1.0
1.0
0.8
u | H (u ) |

MTF
0.4

0.2
0.0
0 1 2 3 4 u(cycles/mm)
 Modulation Transfer Function (MTF) 2

A B

System System

1.0 0.8 1.0 0.4

1 mm
1 mm
1 mm u=1 cycles/mm 1 mm u=2 cycles/mm

C 1.0
A
0.8
| H (u ) |
System

MTF
B
0.4
1.0 0.2
C
1 mm 0.2
1 mm u=4 cycles/mm 0.0
0 1 2 3 4 u(cycles/mm)
 Two Dimensional (2D) Fourier Transform

∞ ∞
1D FT F (u ) = ∫ f ( x)e −i 2πux dx, f ( x) = ∫ F (u )ei 2πux du
−∞ −∞

2D FT can be obtained from the expansion of 1D FT.

F (u , v) = ℑ[ f ( x, y )]
∞ ∞
F (u , v) = ∫ ∫ f ( x, y )e −i 2π ( ux+ vy ) dxdy
−∞ −∞
∞ ∞
f ( x, y ) = ∫ ∫ F (u , v)e i 2π (ux+ vy ) dudv
−∞ −∞
 Two Dimensional (2D) Discrete Fourier Transform (DFT)

2D DFT of a digital image with M x N matrix size

M −1 N −1 − i 2π (
ux vy
+ )
F (u , v) = ∑ ∑ f ( x, y )e M N x
x =0 y =0 0 1 2 3 M-1
0
M −1 N −1 ux vy
i 2π ( + )
f ( x, y ) = ∑ ∑ F (u, v)e M N 1
2
u =0 v =0 y
3
x, u = 0, 1, 2, LL , M − 1
y, v = 0, 1, 2, LL , N − 1

Q e i 2πk = 1 (k = ±1, ± 2, ± 3, L) N-1

f ( x + k x M , y + k y N ) = f ( x, y ) Digital Image
k x , k y = ±1, ± 2, ± 3, L

FFT (fast Fourier transform) is usually used for calculation.

 Rearrangement Frequency Domain after FT

origin of an image
u

4 3
v

2 1
u
3 4

1 2

FFT is calculated for image with an origin at upper left.

Therefore, we must rearrange the frequency domain to move
an origin from upper left to center.
Thank you for your attention!!