Integral and Discrete Transforms in Image Processing
Fourier Transform & Spherical Harmonics
David Svoboda
email: svoboda@fi.muni.cz
Centre for Biomedical Image Analysis
Faculty of Informatics, Masaryk University, Brno, CZ
January 9, 2023
David Svoboda (CBIA@FI) Image Transforms autumn 2023 1 / 65
Outline
1 Fourier Transform
Definition
Properties
2 Discrete Fourier Transform
Definition
Properties
Fast Fourier Transform
3 Discrete Fourier transform in 2D
Definition
Properties
4 Spherical Harmonics Transform
5 Hilbert Transform
David Svoboda (CBIA@FI) Image Transforms autumn 2023 2 / 65
Fourier Transform
Motivation
Let
φω(x) = e2πixω
= cos 2πxω + i sin 2πxω
be a orthonormal basis, where x ∈ R defines the index and ω ∈ R defines
the frequency (speed of oscilation). We can perform a projection of any
function f onto the basis φω(x) as follows:
f · φω =
∞
−∞
f (x)φω(x)dx =
∞
−∞
f (x)e−2πixω
dx
David Svoboda (CBIA@FI) Image Transforms autumn 2023 4 / 65
Fourier Transform
Definition
Given 1D integrable function f and a basis (φω, ω ∈ R), let us define:
forward 1D continuous Fourier transform
F(ω) ≡
∞
−∞
f (x)e−2πixω
dx
inverse 1D continuous Fourier transform
f (x) ≡
∞
−∞
F(ω)e2πixω
dω
David Svoboda (CBIA@FI) Image Transforms autumn 2023 5 / 65
Fourier Transform Properties
Oddness and evenness
1 Each function f (x) is sum of its odd and even part:
E(x) =
1
2
[f (x) + f (−x)]
O(x) =
1
2
[f (x) − f (−x)]
f (x) = E(x) + O(x)
2 sin (x) is an odd function
cos (x) is an even function
3 Any FT basis function is composed of sine and cosine waves
Corollary:
FT of even function misses imaginary part (sine waves)
FT of odd function misses real part (cosine waves)
David Svoboda (CBIA@FI) Image Transforms autumn 2023 6 / 65
Fourier Transform Properties
Oddness
−4 −3 −2 −1 0 1 2 3 4
−2
−1.5
−1
−0.5
0
0.5
1
1.5
2
m
f(m)
Original function
−4 −3 −2 −1 0 1 2 3 4
−8
−6
−4
−2
0
2
4
6
8
k
F(k)
Imaginary part
original function (odd) its transform pair
David Svoboda (CBIA@FI) Image Transforms autumn 2023 7 / 65
Fourier Transform Properties
Evenness
−4 −3 −2 −1 0 1 2 3 4
0
0.2
0.4
0.6
0.8
1
1.2
1.4
1.6
1.8
2
m
f(m)
Original function
−4 −3 −2 −1 0 1 2 3 4
−2
0
2
4
6
8
10
12
k
F(k)
Real part
original function (even) its transform pair
David Svoboda (CBIA@FI) Image Transforms autumn 2023 8 / 65
Fourier Transform Propertis
Scaling
Statement:
f (ax) ⊃
1
|a|
F
ω
a
Proof:
F [f (ax)] (ω) =
∞
−∞
f (ax)e−2πixω
dx
=
1
|a|
∞
−∞
f (ax)e−2πi(ax)(ω/a)
d(ax)
=
1
|a|
F
ω
a
□
Notice: Stretch in time domain corresponds to shrinkage in Fourier
domain and vice versa.
David Svoboda (CBIA@FI) Image Transforms autumn 2023 9 / 65
Fourier Transform Properties
Scaling/Reciprocity
David Svoboda (CBIA@FI) Image Transforms autumn 2023 10 / 65
Fourier Transform Properties
Scaling/Reciprocity
David Svoboda (CBIA@FI) Image Transforms autumn 2023 11 / 65
Fourier Transform Properties
Scaling/Reciprocity
David Svoboda (CBIA@FI) Image Transforms autumn 2023 12 / 65
Fourier Transform Properties
Scaling/Reciprocity
David Svoboda (CBIA@FI) Image Transforms autumn 2023 13 / 65
Fourier Transform Properties
Shift
Statement:
f (x − a) ⊃ e−2πiaω
F(ω)
Proof:
F [f (x − a)] (ω) =
∞
−∞
f (x − a)e−2πixω
dx
=
∞
−∞
f (x − a)e−2πi(x−a)ω
e−2πiaω
d(x − a)
= e−2πiaω
F(ω)
□
Notice: Shift affects only phase. The higher the frequency ω is the more
the corresponding cosine wave is affected.
David Svoboda (CBIA@FI) Image Transforms autumn 2023 14 / 65
Fourier Transform Properties
Shift
David Svoboda (CBIA@FI) Image Transforms autumn 2023 15 / 65
Fourier Transform Properties
Shift
David Svoboda (CBIA@FI) Image Transforms autumn 2023 16 / 65
Fourier Transform Properties
Hermitian symmetry of real signal
Statement:
F(−ω) = F(ω)
Proof:
F(−ω) =
∞
−∞
f (x)e−2πix(−ω)
dx
=
∞
−∞
f (x)e2πixω
= F(ω) iff f : R → R
□
David Svoboda (CBIA@FI) Image Transforms autumn 2023 17 / 65
Fourier Transform Properties
Linearity
Statement:
αf (x) + βg(x) ⊃ αF(ω) + βG(ω)
Proof:
F [αf (x) + βg(x)] (ω) =
∞
−∞
[αf (x) + βg(x)] e−2πixω
dx
= α
∞
−∞
f (x)e−2πixω
dx + β
∞
−∞
g(x)e−2πixω
dx
= αF(ω) + βG(ω)
□
David Svoboda (CBIA@FI) Image Transforms autumn 2023 18 / 65
Fourier Transform Properties
Convolution theorem
Statement:
f (x) ∗ g(x) ⊃ F(ω)G(ω)
f (x)g(x) ⊃ F(ω) ∗ G(ω)
Proof:
F [f (x) ∗ g(x)] (ω) =
∞
−∞


∞
−∞
f (x′
)g(x − x′
)dx′

 e−2πixω
dx
=
∞
−∞
f (x′
)


∞
−∞
g(x − x′
)e−2πixω
dx

 dx′
=
∞
−∞
f (x′
)e−2πix′ω
G(ω)dx′
= F(ω)G(ω)
□
David Svoboda (CBIA@FI) Image Transforms autumn 2023 19 / 65
Fourier Transform Properties
Convolution theorem
The meaning:
The convolution in time domain corresponds to point-wise
multiplication in the Fourier domain, and vice versa.
time
domain
domain
Fourier
=
=
f
F G
f*g
FT FT IFT
g
FT(f)FT(g)
faster?
slower?
*
.
David Svoboda (CBIA@FI) Image Transforms autumn 2023 20 / 65
Fourier Transform Properties
Rayleigh’s energy theorem (Parseval’s theorem)
Statement:
∞
−∞
|f (x)|2
dx =
∞
−∞
|F(ω)|2
dω
Proof: ∞
−∞
|f (x)|2
dx =
∞
−∞
f (x)f (x)dx
=
∞
−∞
f (x)f (x)e−2πixω′
dx ω′
= 0
= F(ω′
) ∗ F(−ω′) ω′
= 0
=
∞
−∞
F(ω′
)F(ω − ω′)dω ω′
= 0
=
∞
−∞
F(ω)F(ω)dω □
Notice: Rayleigh’s energy theorem – the integral of the square of a
function is equal to the integral of the square of its transform.
David Svoboda (CBIA@FI) Image Transforms autumn 2023 21 / 65
Discrete Fourier Transform
Motivation
Given a signal f , |f | = N, let
φk(m) =
1
√
N
e
2πimk
N =
1
√
N
cos
2πmk
N
+ i sin
2πmk
N
be an orthonormal basis, where m ∈ {0, . . . , N − 1} defines the index
(position) and k ∈ {0, . . . , N − 1} defines the frequency (speed of
oscilation). We can perform a projection of a signal f onto the basis
φk(m) as follows:
f · φk =
N−1
m=0
f (m)φk(m) =
N−1
m=0
f (m)
1
√
N
e−2πikm
N
David Svoboda (CBIA@FI) Image Transforms autumn 2023 23 / 65
Discrete Fourier Transform
Motivation
An example of 16 sampled basis functions for N = 16:
0 10
−1
0
1
freq=0
0 10
−1
0
1
freq=1
0 10
−1
0
1
freq=2
0 10
−1
0
1
freq=3
0 10
−1
0
1
freq=4
0 10
−1
0
1
freq=5
0 10
−1
0
1
freq=6
0 10
−1
0
1
freq=7
0 10
−1
0
1
freq=8
0 10
−1
0
1
freq=9
0 10
−1
0
1
freq=10
0 10
−1
0
1
freq=11
0 10
−1
0
1
freq=12
0 10
−1
0
1
freq=13
0 10
−1
0
1
freq=14
0 10
−1
0
1
freq=15
0 5 10 15
−1
0
1
freq=0
0 5 10 15
−1
0
1
freq=1
0 5 10 15
−1
0
1
freq=2
0 5 10 15
−1
0
1
freq=3
0 5 10 15
−1
0
1
freq=4
0 5 10 15
−1
0
1
freq=5
0 5 10 15
−1
0
1
freq=6
0 5 10 15
−1
0
1
freq=7
0 5 10 15
−1
0
1
freq=8
0 5 10 15
−1
0
1
freq=9
0 5 10 15
−1
0
1
freq=10
0 5 10 15
−1
0
1
freq=11
0 5 10 15
−1
0
1
freq=12
0 5 10 15
−1
0
1
freq=13
0 5 10 15
−1
0
1
freq=14
0 5 10 15
−1
0
1
freq=15
David Svoboda (CBIA@FI) Image Transforms autumn 2023 24 / 65
Discrete Fourier Transform
Motivation
Properties
periodical:
φk+N (m) =
1
√
N
e
2πim(k+N)
N =
1
√
N
e
2πimk
N = φk (m)
symmetrical:
φN−k (m) =
1
√
N
e
2πim(N−k)
N
=
1
√
N
e− 2πimk
N =
1
√
N
e
2πi(−m)k
N = φk (−m)
orthonormal:
φk (m) · φl (m) =
N−1
m=0
1
√
N
e
2πmki
N
1
√
N
e
2πmli
N =
1
N
N−1
m=0
e
2πm(k−l)i
N = δk,l
David Svoboda (CBIA@FI) Image Transforms autumn 2023 25 / 65
Discrete Fourier Transform
Definition
Given 1D discrete function f of N samples and a basis
(φk(m), k = {0, . . . , N − 1}), let us define:
forward 1D discrete Fourier transform:
F(k) ≡ f · φk =
1
√
N
N−1
m=0
f (m)e−2πimk
N
inverse 1D discrete Fourier transform:
f (m) ≡ F · φm =
1
√
N
N−1
k=0
F(k)e
2πimk
N
f · φk =
1
√
N
N−1
m=0
f (m)e
2πimk
N =
1
√
N
N−1
m=0
f (m)e−2πimk
N
David Svoboda (CBIA@FI) Image Transforms autumn 2023 26 / 65
Discrete Fourier Transform
Matrix notation
If
φk(m) =
1
√
N
e
−2πimk
N =
1
√
N
e−2πi
N
mk
=
1
√
N
ψmk
then
A =
1
√
N







1 1 1 . . . 1
1 ψ1 ψ2 . . . ψN−1
1 ψ2 ψ2·2 . . . ψ(N−1)2
1
...
...
...
...
1 ψN−1 ψ2(N−1) . . . ψ(N−1)(N−1)







and
F = Af ⇒ f = A−1
F = A
T
F
Notice: ψ = e−2πi
N is called the Nth root of unity.
David Svoboda (CBIA@FI) Image Transforms autumn 2023 27 / 65
Discrete Fourier Transform
The obvious question
When we apply FT, we usually say ”let us decompose our signal into the
sine waves . . . ” Why do we use another (so complicated) basis?
Basis function is a ”sine wave”
we avoid complex numbers
more intuitive if basis function is a
simple ”sine wave”
sine waves without phase shift do not
generate the whole space
possible basis function: sin (km − α)
α hidden in the sine function spoils the
linearity; matrix multiplication cannot
be used
Basis function is φk (m)
we have to use complex
numbers
φk (m) functions generate the
whole space (form basis)
this basis is orthonormal
transform is linear, i.e. realized
via matrix multiplication
David Svoboda (CBIA@FI) Image Transforms autumn 2023 28 / 65
Discrete Fourier Transform
The meaning of Fourier coefficients
If you perform inverse FT
f (m) = F · φm =
N−1
k=0
F(k)φk(m) =
1
√
N
N−1
k=0
F(k)ei 2πkm
N
you literally compose the original signal f by combining the individual
basis functions.
Each basis function φk(m) = ei 2πmk
N defines only the frequency.
Fourier coefficient F(k) = |zk| (cos αk + i sin αk) = |zk|eiαk modifies
the corresponding basis function φk(m) by scaling it and shifting it.
F(k)φm(k) = |zk|eiαk
ei 2πmk
N = |zk|eiαk +i 2πmk
N
= |zk| cos
2πk
N
m + αk + i sin
2πk
N
m + αk
David Svoboda (CBIA@FI) Image Transforms autumn 2023 29 / 65
Discrete Fourier Transform
The meaning of Fourier coefficients
F(0) matches the lowest frequency in the signal f and corresponds to
the “mean” of f :
F(0) ≡
1
√
N
N−1
m=0
f (m)e−2πim0
N =
1
√
N
N−1
m=0
f (m)
F(0) is usually called DC term
DC . . . “direct current” (current of zero frequency)
F(1) . . . F(N − 1) are called AC terms
AC . . . “alternating current”
F(N−1
2 ) matches the highest frequency in the signal f
Exercise: Why does the component F(N−1
2 ) correspond to the highest
frequency?
David Svoboda (CBIA@FI) Image Transforms autumn 2023 30 / 65
Discrete Fourier Transform Properties
Some properties can be adopted from continuous Fourier transform
Evenness & Oddness . . . valid
Shift . . . valid
Linearity . . . valid
Rayleigth theorem . . . valid
Hermitian symmetry of real signal . . . valid
Scaling . . . invalid
Convolution theorem . . . modification required
Some new properties are introduced
Periodicity
Stretch
Rearrangement
David Svoboda (CBIA@FI) Image Transforms autumn 2023 31 / 65
Discrete Fourier Transform Properties
Convolution theorem
Let f and g be 1D discrete periodic signals of length N, then:
f ∗ g = IDFT [DFT(f ) · DFT(g)]
√
N
Proof:
f (m) ∗ g(m) =
N−1
k=0
f (k)g(m − k)
=
N−1
k=0
1
√
N
N−1
n=0
F(n)e2πikn/N 1
√
N
N−1
l=0
G(l)e2πi(m−k)l/N
=
1
N
N−1
n=0
F(n)
N−1
l=0
G(l)
N−1
k=0
e2πikn/N
e2πi(m−k)l/N
=
1
N
N−1
n=0
F(n)
N−1
l=0
G(l)e2πiml/N
N−1
k=0
e2πik(n−l)/N
=
1
N
N−1
n=0
F(n)
N−1
l=0
G(l)e2πiml/N
δ(n − l)N =
N−1
n=0
[F(n) · G(n)] e2πimn/N
□
David Svoboda (CBIA@FI) Image Transforms autumn 2023 32 / 65
Discrete Fourier Transform Properties
Stretch
If f (m) is a 1D function of length N, p ∈ N, and stretchp{f } = {g}, where
g(n) =
f (n/p) n = 0, p, 2p, . . . , (N − 1)p
0 otherwise
then
G(k) =
1
√
p



F(k) k = 0, . . . , N − 1
F(k − N) k = N, . . . , 2N − 1
...
...
F(k − (p − 1)N) k = (p − 1)N, . . . , pN − 1
Notice: Stretch by a factor p in the time domain results in p-fold
repetition of F(k) in the frequency domain.
David Svoboda (CBIA@FI) Image Transforms autumn 2023 33 / 65
Discrete Fourier Transform Properties
Stretch
An example of stretch for p = 2
0 1 2 3 4 5 6 7
0
2
4
6
8
10
m
f(m)
time domain
1 2 3 4 5 6 7 8
0
2
4
6
8
10
u
|F(u)|
Fourier (frequency) domain
0 5 10 15
0
2
4
6
8
10
m
g(m)
time domain
0 5 10 15
0
2
4
6
8
10
u
|G(u)|
Fourier (frequency) domain
David Svoboda (CBIA@FI) Image Transforms autumn 2023 34 / 65
Discrete Fourier Transform Properties
Rearrangement
Let f be a 1D discrete signal, where |f | = N, and let a ∈ N be number
that satisfies the condition gcd(a, N) = 1. The rearrangement of signal f
determines the rearangement of Fourier coefficients in the following way:
F(k) = F(f (m))
⇓
F(⟨bk⟩N) = F(f (⟨am⟩N))
where
⟨m⟩N = m mod N
b = ⟨aEuler(N)−1⟩N
David Svoboda (CBIA@FI) Image Transforms autumn 2023 35 / 65
Discrete Fourier Transform Properties
Why scaling does not work
David Svoboda (CBIA@FI) Image Transforms autumn 2023 36 / 65
Discrete Fourier Transform Properties
Periodicity
Statement:
F(k + N) = F(k)
Proof:
F(k + N) =
1
√
N
N−1
m=0
f (m)e−
2πim(k+N)
N
=
1
√
N
N−1
m=0
f (m)e−2πimk
N e−2πimN
N
=
1
√
N
N−1
m=0
f (m)e−2πimk
N
= F(k)
□
David Svoboda (CBIA@FI) Image Transforms autumn 2023 37 / 65
Discrete Fourier Transform Properties
Periodicity and Symmetry ⇒ Domain centering
⇔
David Svoboda (CBIA@FI) Image Transforms autumn 2023 38 / 65
Fast Fourier Transform
Idea: N-point signal (N = 2m, m ∈ N) is decomposed into two N/2-point
signals:
one with all odd samples
one with all even samples
Example:
1 input signal: [ f0 f1 f2 f3 f4 f5 f6 f7 ] ⊃ [ ? ? ? ? ? ? ? ? ]
2 2× simpler DFT: [ f0 f2 f4 f6 ] ⊃ [ A B C D ]
2× simpler DFT: [ f1 f3 f5 f7 ] ⊃ [ P Q R S ]
3 stretch: [ f0 0 f2 0 f4 0 f6 0 ] ⊃ 1√
2
[ A B C D A B C D ]
4 stretch: [ f1 0 f3 0 f5 0 f7 0 ] ⊃ 1√
2
[ P Q R S P Q R S ]
shift: [ 0 f1 0 f3 0 f5 0 f7 ] ⊃ 1√
2
[ P ψQ ψ2
R ψ3
S ψ4
P ψ5
Q ψ6
R ψ7
S ]
5 linearity: [ f0 f1 f2 f3 f4 f5 f6 f7 ] = [ f0 0 f2 0 f4 0 f6 0 ] + [ 0 f1 0 f3 0 f5 0 f7 ]
⇓
[ f0 f1 f2 f3 f4 f5 f6 f7 ] ⊃ 1√
2
[ (A + P) (B + ψQ) (C + ψ2
R) . . . ]
David Svoboda (CBIA@FI) Image Transforms autumn 2023 39 / 65
Fast Fourier Transform
Derivation:
F(k) =
1
√
N
N−1
m=0
f (m)e−2πimk
N
=
1
√
N
N/2−1
m=0
f (2m)e−
2πi(2m)k
N +
1
√
N
N/2−1
m=0
f (2m + 1)e−
2πi(2m+1)k
N
=
1
√
N
N/2−1
m=0
f (2m)e
−2πimk
N/2 + e−2πik
N
1
√
N
N/2−1
m=0
f (2m + 1)e
−2πimk
N/2
= Fe
(k) + ψk
Fo
(k)
Notice: ψ = e−2πi
N
David Svoboda (CBIA@FI) Image Transforms autumn 2023 40 / 65
Fast Fourier Transform
Idea: While it is possible, repeat the division.
F (k) → Fe
(k), Fo
(k)
→ Fee
(k), Feo
(k), Foe
(k), Foo
(k)
→ Feee
(k), Feeo
(k), Feoe
(k), Feoo
(k), Foee
(k), Foeo
(k),
Fooe
(k), Fooo
(k)
→ . . .
After log2(n) divisions we have Feeeoe...eoeooee(k) which is just one point
long signal in Fourier domain.
You should know that: {X} ⊃ {X}
Exercise: What is the complexity of FFT?
David Svoboda (CBIA@FI) Image Transforms autumn 2023 41 / 65
Fast Fourier Transform
One more illustration
150 1 2 3 4 5 6 7 10 11 12 13 1498
0 2 4 14121086 1 3 5 7 9 11 13 15
0 4 8 12 2 6 10 14 1 5 9 13 3 7 1511
0 8 4 12 2 10 6 14 1 9 5 13 11 7 153
0 8 4 12 2 10 6 14 1 9 5 13 3 11 15716 signals (1 point each)
2 signals (8 points each)
8 signals (2 points each)
4 signals (4 points each)
1 signal (16 points)
David Svoboda (CBIA@FI) Image Transforms autumn 2023 42 / 65
Discrete Fourier Transform in 2D
Definition
Given 2D discrete function f of (M, N) samples and two bases
(φk, k = {0, . . . , M − 1}) and (φl , l = {0, . . . , N − 1}), let us define:
forward 2D discrete Fourier transform:
F(k, l) ≡
1
√
MN
M−1
m=0
N−1
n=0
f (m, n)e−2πi(mk
M
+nl
N )
inverse 2D discrete Fourier transform:
f (m, n) ≡
1
√
MN
M−1
k=0
N−1
l=0
F(k, l)e2πi(mk
M
+nl
N )
David Svoboda (CBIA@FI) Image Transforms autumn 2023 44 / 65
Discrete Fourier Transform in 2D
Properties
Some properties are adopted from lower-dimensional case
Shift
Periodicity
Convolution theorem
Stretch
Some new properties are introducted
Separability
Rotation
David Svoboda (CBIA@FI) Image Transforms autumn 2023 45 / 65
Discrete Fourier Transform in 2D
Properties
Separability
The evaluation of 2D-(D)FT can be decomposed into two simpler tasks:
F(k, l) =
1
√
MN
M−1
m=0
N−1
n=0
f (m, n)e−2πi(mk
M
+nl
N )
=
1
√
N
N−1
n=0
1
√
M
M−1
m=0
f (m, n)e−2πimk
M e−2πinl
N
=
1
√
N
N−1
n=0
F(k, n)e−2πinl
N
David Svoboda (CBIA@FI) Image Transforms autumn 2023 46 / 65
Discrete Fourier Transform in 2D
Stretch
f(x,y)=x+y
5 10 15
5
10
15
Fourier Transform (amplitude)
5 10 15
5
10
15
f(x,y) stretched by factor 2
5 10 15
5
10
15
Fourier Transform (amplitude)
5 10 15
5
10
15
David Svoboda (CBIA@FI) Image Transforms autumn 2023 47 / 65
Discrete Fourier Transform in 2D
Rotation
Input image "I" |FT(I)|
rotated I |FT(rotated I)|
David Svoboda (CBIA@FI) Image Transforms autumn 2023 48 / 65
Discrete Fourier Transform in 2D
Rotation
Let us introduce the polar coordinates:
x = r cos ϕ
y = r sin ϕ
ωx = R cos ψ
ωy = R sin ψ
Then
f (x, y) → f (r, ϕ)
F(ωx , ωy ) → F(R, ψ)
David Svoboda (CBIA@FI) Image Transforms autumn 2023 49 / 65
Discrete Fourier Transform in 2D
Rotation
It is now clear to see that:
f (r, ϕ + ϕ0) ⊃ F(R, ψ + ϕ0)
Conclusion: Rotating f (x, y) by an angle ϕ0 rotates F(ωx , ωy ) by the
same angle, and vice versa.
David Svoboda (CBIA@FI) Image Transforms autumn 2023 50 / 65
Spherical Harmonics Transform
Motivation
How does the time/input domain influences the transform selection
Discrete-time unit step sequence:
u1(n) =
1 if n ≥ 0
0 otherwise
suitable transform: DFT in 1D
Discrete-time unit step function:
u2(m, n) = u1(m) · u1(n)
suitable transform: DFT in 2D
Discrete-time regularly sampled unit sphere:
suitable transform: SHT (Spherical Harmonics Transform)
Notice: Keep in mind, that Fourier transform requires the input domain to
be regularly sampled!
David Svoboda (CBIA@FI) Image Transforms autumn 2023 52 / 65
Spherical Harmonics Transform
Motivation
Let us define a new basis Y m
ℓ (θ, φ) on a unit sphere as follows:
Y m
ℓ (θ, φ) = kℓ,mPm
ℓ (cos θ)eimφ
where
ℓ and m are respectively the degree and order of the harmonic
kℓ,m is the normalization coefficient:
kℓ,m =
2ℓ + 1
4π
(ℓ − m)!
(ℓ + m)!
θ and φ are respectively the azimuth and the elevation
Pm
ℓ (x) is the associated Legendre polynomial:
Pm
ℓ (x) =
(−1)m
2ℓℓ!
(1 − x2
)m/2 dℓ+m
dxℓ+m
(x2
− 1)ℓ
Notice: Associated Legendre polynomials are mutualy orthogonal
polynomials for fixed m or ℓ.
David Svoboda (CBIA@FI) Image Transforms autumn 2023 53 / 65
Spherical Harmonics Transform
Motivation
An example of some items (their magnitudes only) from the Y m
ℓ basis
David Svoboda (CBIA@FI) Image Transforms autumn 2023 54 / 65
Spherical Harmonics Transform
Definition
Using this Y m
ℓ basis, any spherical scalar function f (θ, φ) can be expressed
as follows:
f (θ, φ) =
∞
ℓ=0
ℓ
m=−ℓ
H(ℓ, m)Y m
ℓ (θ, φ),
where H(ℓ, m) is a harmonic coefficient, given by:
H(ℓ, m) =
π
0
2π
0
f (θ, φ)Y m
ℓ (θ, φ) sin θdφdθ
Notice: Compare the terms Harmonic coefficient and Fourier coefficient.
David Svoboda (CBIA@FI) Image Transforms autumn 2023 55 / 65
Spherical Harmonics Transform
Application
3D cell shape modeling
ℓ = 21
H(0, 0) =


8.892
8.282
9.269


H(1, −1) =


0.609 − 1.354i
1.228 + 1.264i
−1.797 + 0.119i


H(1, 0) =


1.251
2.045
3.078


H(1, 1) =


−0.609 − 1.354i
−1.228 + 1.264i
1.797 + 0.119i


courtesy: images prepared by D. Wiesner
Notice: For more details how to get regular sampling over a sphere, see
http://lishenlab.com/spharm/SPHARM-MAT.pdf
David Svoboda (CBIA@FI) Image Transforms autumn 2023 56 / 65
Hilbert Transform
Definition
Let f be a real-valued signal. The Hilbert transform of f is defined as
a convolution:
H(f ) = f ∗ h,
where h(t) = 1
πt is called a Hilbert kernel.
Alternatively, we can define Hilbert transform via frequency domain:
F(H(f ))(ω) = −i sign(ω)F(f )(ω),
where
sign(x) =



1 x > 0
0 x = 0
−1 x < 0
David Svoboda (CBIA@FI) Image Transforms autumn 2023 58 / 65
Hilbert Transform
Properties
Hilbert transform is further used to form so called Analytic Signal:
fA = f + i H(f )
0.0 0.1 0.2 0.3 0.4 0.5
1
0
1
Original signal
0.0 0.1 0.2 0.3 0.4 0.5
1
0
1
Components of Analytic signal
Real part
Imaginary part
Notice: Hilbert transform acts as a phase shift.
David Svoboda (CBIA@FI) Image Transforms autumn 2023 59 / 65
Hilbert Transform
Properties
Analytic signal has the following properties:
real(fA) = f
imag(fA) = H(f )
fA is a complex signal, i.e. can be interpreted in polar form:
fA(x) = A(x)eiφ(x)
where
A(x) = f (x)2 + H2(f )(x) . . . amplitude/envelope
φ(x) = arctan H(f )(x)
f (x) . . . phase
David Svoboda (CBIA@FI) Image Transforms autumn 2023 60 / 65
Hilbert Transform
Application
Positive versus Negative frequencies
Fourier coefficent F(k) corresponds to:
positive frequency . . . bk(m) = e
2πimk
N if 0 ≤ k ≤ N/2
negative frequency . . . b−k(m) = e−2πimk
N otherwise
Real-valued signal
Fourier transform produces:
double-sided spectrum with
redundant items (Fourier
coefficients)
pairs of identical items with
positive and negative
frequencies
Analytic signal
Fourier transform produces:
single-sided spectrum (no
redundancies)
spectrum contains only
positive frequencies
David Svoboda (CBIA@FI) Image Transforms autumn 2023 61 / 65
Hilbert Transform
Application
Double-sided versus single-sided spectrum
David Svoboda (CBIA@FI) Image Transforms autumn 2023 62 / 65
Hilbert Transform
Application
Signal decomposition
0 100 200 300 400 500 600
0.5
1.0
1.5
sin(x)
0 100 200 300 400 500 600
1
0
1
Chirpsignal
0 100 200 300 400 500 600
1
0
1
Mixedsignal
0 100 200 300 400 500 600
1
0
1
Envelope
0 100 200 300 400 500 600
2
0
2
Phase
David Svoboda (CBIA@FI) Image Transforms autumn 2023 63 / 65
Bibliography
Bracewell, R. N., Fourier transform and its applications / 2nd ed.
New York: McGraw-Hill, 474 pages, ISBN 0070070156
Gonzalez, R. C., Woods, R. E., Digital image processing / 2nd ed.,
Upper Saddle River: Prentice Hall, 2002, pages 793, ISBN
0201180758
Veit, J., Integr´aln´ı transformace, Praha, 1983
Smith, Steven W. Digital signal processing: A practical guide for
engineers and scientists; Amsterdam: Newnes, 2003, 650 pages;
on-line version: http://www.dspguide.com/
Talwalkar, S.A. and Marple S. L., ”Time-frequency scaling property of
Discrete Fourier Transform (DFT),” 2010 IEEE International
Conference on Acoustics, Speech and Signal Processing, 2010, pp.
3658-3661
David Svoboda (CBIA@FI) Image Transforms autumn 2023 64 / 65
You should know the answers . . .
Express the discrete Fourier transform as a matrix multiplication.
Derive this matrix.
How many Fourier basis functions do we need if we transform the
signal of length N into the frequency domain?
Formulate the forward discrete Fourier transform. Explain all the
variables and constants.
What does DC and AC terms mean?
Explain the meaning of one particular Fourier coefficient in inverse
Fourier transform.
What is the product of the projection of an even function into a sin
wave?
Why are the wide functions in time domain transformed into their
narrow counterparts in frequency domain and vice versa?
Derive FFT for a signal of length 3m, m ∈ N.
David Svoboda (CBIA@FI) Image Transforms autumn 2023 65 / 65