Integral and Discrete Transforms in Image Processing
2nd Generation of Wavelets – Lifting Scheme
David Svoboda
email: svoboda@fi.muni.cz
Centre for Biomedical Image Analysis
Faculty of Informatics, Masaryk University, Brno, CZ
November 6, 2022
David Svoboda (CBIA@FI) Image Transforms autumn 2022 1 / 35
Outline
1 Motivation
2 Analysis of FWT
3 Lifting scheme
4 Integer Wavelet transform
5 Applications
David Svoboda (CBIA@FI) Image Transforms autumn 2022 2 / 35
Motivation
Practical use of DWT or FWT?
Complexity of:
DWT: O(n2
)
FWT: O(cn)
Can we do it faster?
DWT/FWT are computed in floating point arithmetic.
Can we restrict ourselves to integer number?
DWT/FWT are computed out-of-place.
Can we reduce the memory usage?
David Svoboda (CBIA@FI) Image Transforms autumn 2022 4 / 35
Analysis of FWT
Let |f | = N = 8 = 23 = 2J and j0 = 0
f (0) f (1) f (2) f (3) f (4) f (5) f (6) f (7)
∥ ∥ ∥ ∥ ∥ ∥ ∥ ∥
A3(0) A3(1) A3(2) A3(3) A3(4) A3(5) A3(6) A3(7)
⇓
A2(0) A2(1) A2(2) A2(3) D2(0) D2(1) D2(2) D2(3)
⇓
A1(0) A1(1) D1(0) D1(1)
⇓
A0(0) D0(0)
David Svoboda (CBIA@FI) Image Transforms autumn 2022 6 / 35
Analysis of FWT
An example (Haar wavelets)
Let f = Aj be an input signal:
Aj−1 = Aj ∗
√
2
2
,
√
2
2
(↓ 2×)
Dj−1 = Aj ∗
√
2
2
, −
√
2
2
(↓ 2×)
In each sample k, we perform one averaging and one difference:
Aj−1(k) = (Aj (k + 1) + Aj (k))
√
2
2
(↓ 2×)
Dj−1(k) = (Aj (k + 1) − Aj (k))
√
2
2
(↓ 2×)
David Svoboda (CBIA@FI) Image Transforms autumn 2022 7 / 35
Analysis of FWT
Fast Wavelet transform
The filtering is followed by
downsampling
We throw away half of
computed samples!
2x
D
2x
H
L
A
j−1
j
Aj−1
Optimization strategy
Let’s flip the order of filtering and downsampling:
1 we start with downsampling
2 we proceed with filtering
David Svoboda (CBIA@FI) Image Transforms autumn 2022 8 / 35
Analysis of FWT
Downsampling in Z-domain
Let f = [f0, f1, f2, . . . ] be a signal and further let
⇒ feven = [f0, f2, f4, . . . ] and
⇒ fodd = [f1, f3, f5, . . . ] be its even and odd samples, respectively.
The same rule is applied to filter g = [g0, g1, g2, . . . ].
If h = f ∗ g, then in Z-domain:
z0
: h0 = f0g0
z−1
: h1 = f1g0 + f0g1
z−2
: h2 = f2g0 + f1g1 + f0g2
z−3
: h3 = f3g0 + f2g1 + f1g2 + f0g3
z−4
: h4 = f4g0 + f3g1 + f2g2 + f1g3 + f0g4
...
Downsampling H(z) by 2 removes each odd line, which results in:
Heven(z) = Feven(z)Geven(z) + z−1
Fodd (z)Godd (z)
David Svoboda (CBIA@FI) Image Transforms autumn 2022 9 / 35
Analysis of FWT
Downsampling in Z-domain
Downsampling of Aj followed by filtering:
Aj−1(z) = Aj even(z)Leven(z) + z−1
Aj odd (z)Lodd (z)
Dj−1(z) = Aj even(z)Heven(z) + z−1
Aj odd (z)Hodd (z)
In (polyphase) matrix notation:
Aj−1(z)
Dj−1(z)
=
Leven(z) Lodd (z)
Heven(z) Hodd (z)
·
Aj even(z)
z−1Aj odd (z)
P(z) =
Leven(z) Lodd (z)
Heven(z) Hodd (z)
David Svoboda (CBIA@FI) Image Transforms autumn 2022 10 / 35
Analysis of FWT
The use of polyphase matrix
split
even
odd
jA (z)
P(z)
j−1A (z)
j−1D (z)
David Svoboda (CBIA@FI) Image Transforms autumn 2022 11 / 35
Lifting scheme
From Filters to Lifting
Algorithm invented by Wim Sweldens (1996):
1 Input: either ’LoD’ and ’HiD’ filters or φ and ϕ functions
2 Convert both filters ’LoD’ and ’HiD’ to Z-domain
3 Create polyphase matrix (2×2)
4 Factorize matrix into simple (lower and upper diagonal) matrices
Each simple matrix corresponds either to one update or prediction step
5 Convert each matrix from Z-domain to time domain
David Svoboda (CBIA@FI) Image Transforms autumn 2022 13 / 35
Lifting scheme
Example: Haar Filters to Lifting
1 Low (l) and high (h) frequency decomposition filter . . .
l =
√
2
2 ,
√
2
2 h = −
√
2
2 ,
√
2
2
with emphasizing reference positions:
l =
√
2
2
,
√
2
2 h = −
√
2
2
,
√
2
2
2 Z-analysis of both filters:
L(z) =
√
2
2
z0
+
√
2
2
z−1
H(z) = −
√
2
2
z0
+
√
2
2
z−1
David Svoboda (CBIA@FI) Image Transforms autumn 2022 14 / 35
Lifting scheme
Example: Haar Filters to Lifting
2 Z-analysis of both filters (continued):
Using the rule:
Filter(z) = Filtere(z2
) + z−1
Filtero(z2
)
we perform the separation of even and odd part of both filters:
Le(z) =
√
2
2
z0
Lo(z) =
√
2
2
z0
He(z) = −
√
2
2
z0
Ho(z) =
√
2
2
z0
3 Creation of Polyphase matrix:
P(z) =
Le(z) Lo(z)
He(z) Ho(z)
=
√
2
2
√
2
2
−
√
2
2
√
2
2
David Svoboda (CBIA@FI) Image Transforms autumn 2022 15 / 35
Lifting scheme
Example: Haar Filters to Lifting
4 Factorization into lower diagonal matrix ≈ dual step:
P(z) =
√
2
2
√
2
2
−
√
2
2
√
2
2
=
Le(z)
√
2
2
He(z)
√
2
2
·
1 0
Lower(z) 1
We get a set of equations:
√
2
2
= Le(z) + Lower(z)
√
2
2
−
√
2
2
= He(z) + Lower(z)
√
2
2
Solution (divison of polynomials – ambiguous!):
Lower(z) = −1; He(z) = 0; Le(z) =
√
2
P(z) =
√
2
√
2
2
0
√
2
2
·
1 0
−1 1
David Svoboda (CBIA@FI) Image Transforms autumn 2022 16 / 35
Lifting scheme
Example: Haar Filters to Lifting
4 Factorization into upper diagonal matrix ≈ primal step:
P(z) =
√
2
√
2
2
0
√
2
2
=
√
2 Lo(z)
0 Ho(z)
·
1 Upper(z)
0 1
We get a set of equations:
√
2
2
=
√
2 · Upper(z) + Lo(z)
√
2
2
= Ho(z)
Solution:
Upper(z) =
1
2
; Lo(z) = 0
√
2
√
2
2
0
√
2
2
=
√
2 0
0
√
2
2
·
1 1
2
0 1
David Svoboda (CBIA@FI) Image Transforms autumn 2022 17 / 35
Lifting scheme
Example: Haar Filters to Lifting
5 Finally
P(z) =
√
2 0
0
√
2
2
·
1 1
2
0 1
·
1 0
−1 1
√
2 0
0
√
2
2
normalization NA =
√
2, ND =
√
2
2
1 1
2
0 1
primal lifting step Aj ← Aj + 1
2 Dj
1 0
−1 1
dual lifting step Dj ← Dj + (−1)Aj
David Svoboda (CBIA@FI) Image Transforms autumn 2022 18 / 35
Lifting scheme
Adopted idea of FastDWT
The transition from level j to j − 1 is computed efficiently
Basic idea: split → predict → update → normalize
j−1A
j−1D
predictsplit update
even
odd
−
+
jA
N
NA
D
David Svoboda (CBIA@FI) Image Transforms autumn 2022 19 / 35
Lifting scheme
Split step
Also known as lazy wavelet.
The signal Aj is split into odd and even samples:
(Aj−1, Dj−1) ← split(Aj )
j−1A
j−1D
predictsplit update
even
odd
−
+
jA
N
NA
D
David Svoboda (CBIA@FI) Image Transforms autumn 2022 20 / 35
Lifting scheme
Prediction step
j−1A
j−1D
split update
even
odd
−
+
jA
N
NA
D
predict
Dj−1(k) ← Dj−1(k) − predict(Aj−1(k))
Also known as dual lifting step
When using Haar wavelets the neighbouring samples are supposed to
be equal, i.e. the predictor is simple: predictHaar (f (k)) = f (k)
Dj−1(k) ← Dj−1(k) − Aj−1(k)
David Svoboda (CBIA@FI) Image Transforms autumn 2022 21 / 35
Lifting scheme
Update step
j−1A
j−1D
predictsplit update
even
odd
−
+
jA
N
NA
D
Aj−1(k) ← Aj−1(k) + update(Dj−1(k))
Also known as primal lifting step
Update step repairs the wrong estimate of the prediction step.
When using Haar wavelet, we use: updateHaar (f (k)) = 1
2f (k)
Aj−1(k) ← Aj−1(k) + 1
2Dj−1(k)
David Svoboda (CBIA@FI) Image Transforms autumn 2022 22 / 35
Lifting scheme
Normalization step
j−1A
j−1D
predictsplit update
even
odd
−
+
jA
N
NA
D
Aj−1(k) ← Aj−1(k) · NA
Dj−1(k) ← Dj−1(k) · ND
NA · ND = 1
The output signals are normalized to avoid boosting the signal.
David Svoboda (CBIA@FI) Image Transforms autumn 2022 23 / 35
Lifting scheme
Inverse lifting
The forward algorithm can be simply inverted:
Aj−1(k) ← Aj−1(k) − update(Dj−1(k))
Dj−1(k) ← Dj−1(k) + predict(Aj−1(k))
Aj ← merge(Aj−1, Dj−1)
jA
j−1A
j−1D
ND
AN +
−
update predict
even
odd
merge
David Svoboda (CBIA@FI) Image Transforms autumn 2022 24 / 35
Lifting scheme
Inverse lifting (example)
Namely for Haar wavelets we get:
Aj−1(k) ← Aj−1(k) · ND
Dj−1(k) ← Dj−1(k) · NA
Aj−1(k) ← Aj−1(k) − (Dj−1(k)/2)
Dj−1(k) ← Dj−1(k) + Aj−1(k)
Aj = merge(Aj−1, Dj−1)
jA
j−1A
j−1D
ND
AN +
−
update predict
even
odd
merge
David Svoboda (CBIA@FI) Image Transforms autumn 2022 25 / 35
Lifting scheme
From Filters to Lifting (Daubechies-2)
NA
j−1A
ND
j−1D
predictsplit update
even
odd
−
+
jA
−
predict
predict1(f (k)) = (−1.732) f (k)
update1(f (k)) = 0.433 f (k) − 0.067 f (k + 1)
predict2(f (k)) = f (k − 1)
David Svoboda (CBIA@FI) Image Transforms autumn 2022 26 / 35
Lifting scheme
Technical/Implementation notes
Lifting ordering (for N = 8)
f(0) f(1) f(2) f(3) f(4) f(5) f(6) f(7)
A2(0) D2(0) A2(1) D2(1) A2(2) D2(2) A2(3) D2(3)
A1(0) D1(0) A1(1) D1(1)
A0(0) D0(0)
Notice: The computation is performed completedly in-place.
From filters to lifting scheme
conversion ’LoD’,’HiD’ → update(·), predict(·) always exists but is
not unique
conversion is performed in frequency domain via Z-transform
(polyphase matrice factorization)
David Svoboda (CBIA@FI) Image Transforms autumn 2022 27 / 35
Lifting scheme
FWT and DWT comparison (examples)
Price of one decomposition level using DWT (|f | = N)
family of wavelets multiplications additions
Haar 4N 2N
Daubechies-2 8N 6N
Extra memory usage: one memory buffer of size O(N) needed for convolution.
Price of one decomposition level using LS (|f | = N)
family of wavelets multiplications additions
Haar 2N N
Daubechies-2 3N 2N
Extra memory usage: ∅
David Svoboda (CBIA@FI) Image Transforms autumn 2022 28 / 35
Integer Wavelet Transform
Basic idea & Properties
IWT originates from lifting scheme (chain of predictions and updates).
The fixed point arithmetic is guaranteed via ’floor’ function.
The rounding error produced in forward transform is compensated by
mirror ’floor’ in inverse lifting. They cancel out each other.
The lifting is the same as the standard one but
each multiplication is followed but the truncation
no normalization is present
David Svoboda (CBIA@FI) Image Transforms autumn 2022 30 / 35
Integer Wavelet Transform
An example (Haar)
j−1D
j−1A
split
even
odd
−
+
jA predict update
floor floor
& &
(Aj−1, Dj−1) ← split(Aj )
Dj−1(k) ← Dj−1(k) − floor (Aj−1(k))
Aj−1(k) ← Aj−1(k) + floor (Dj−1(k)/2)
David Svoboda (CBIA@FI) Image Transforms autumn 2022 31 / 35
Applications
Fusion of images with different resolution
Image registration
Edge detection
Image compression:
original (979 kB) JPEG (6.21 kB) JPEG2000 (1.83 kB)
David Svoboda (CBIA@FI) Image Transforms autumn 2022 33 / 35
Bibliography
Li H., Manjunath B.S., Mitra S.K. Multisensor Image Fusion Using
the Wavelet Transform, Graphical Models and Image Processing,
Volume 57, Issue 3, May 1995, Pages 235-245, ISSN 1077-3169
Jensen A., La Cour-Harbo A. Ripples in mathematics: the discrete
wavelet transform, Springer, Berlin, 2001, ISBN 3-540-41662-5
Sweldens W. The lifting scheme: A custom-design construction of
biorthogonal wavelets. Applied and Computational Harmonic
Analysis. 1996, Vol 3, Issue 2, pp 186–200
David Svoboda (CBIA@FI) Image Transforms autumn 2022 34 / 35
You should know the answers . . .
Explain three principal phases of lifting scheme.
Compare the time and spatial complexity of FWT and lifting scheme.
Compute one level of wavelet transform (using Haar’s basis) via lifting
scheme for signal f=[3 5 0 -1 4 2].
What does integer wavelet transform mean?
David Svoboda (CBIA@FI) Image Transforms autumn 2022 35 / 35