SPHERICAL IMAGE DENOISING AND ITS APPLICATION TO
OMNIDIRECTIONAL IMAGING
St
´
ephanie Bigot
1
, Djemaa Kachi
2
, Sylvain Durand
1
and El Mustapha Mouaddib
2
1
L.A.M.F.A., UMR CNRS 6140 - U.P.J.V. 33 rue Saint Leu 80000 Amiens Cedex 1, France
2
C.R.E.A., EA 3299 - U.P.J.V. 7 rue du moulin neuf 80039 Amiens Cedex 1, France
Keywords:
Omnidirectional imaging, Spherical filters, Spherical harmonics, Wiener filter, Tikhonov regularization.
Abstract:
This paper addresses the problem of spherical image processing. Thanks to projective geometry, the omnidi-
rectional image can be presented as a function on sphere S
2
. The target application includes omnidirectional
image smoothing. We describe a new method of smoothing for spherical images. For that purpose, we in-
troduce a suitable Wiener filter and we use the Tikhonov method to these images. In order to compare their
performances, we present the most used classical spherical kernels.
We present several examples for filtering real and synthetical spherical images.
a
.
a
This project was financed by the D.G.A.
1 INTRODUCTION
Images defined on a spherical surface rather than in
the image plane arise in various domains like en-
vironment mapping, and medical imaging (Chung,
2006). . . . Omnidirectional images appear with distor-
sion (because of the used mirror) in the image plane
and are naturally defined on the sphere as shown in the
model of Geyer (Geyer and Daniilidis, 2001). Most
approaches analysing omnidirectional images applied
classical operators. However the analysis or the fil-
tering of these images by traditional methods pro-
duces erroneous results. In (Daniilidis et al., 2002),
Daniilidis proposed of the very first approaches of
adaptation of the operations of filtering to the omni-
directional images. The proposal of the authors is to
project the image plane to a non-deformed space (unit
sphere) (Geyer and Daniilidis, 2000).
The spherical images are obtained by stere-
ographic projection of omnidirectional image, as
shown in Fig.1.Daniilidis then define kernels in the
space of the spherical coordinates and propose to
carry out the convolution on the surface of the sphere.
In (Strauss and Comby, 2005), the authors present
new morphological operators that use the projective
property of omnidirectional sensor. These operators
Figure 1: Omnidirectional image of Office and the corre-
sponding spherical image.
make use of a structural element with a varying shape.
Another work (Ieng et al., 2003), advocates the prin-
ciple of kernel with varying shape to match indices
from two omnidirectional images. Except a few at-
tempts to adapt the method of processing omnidirec-
tional images, no normalization method has been es-
tablished for these images.
A simple technique for smoothing planar images is
the convolution of the image with a Gaussian kernel
as
g(x;t) =
1
4πkt
exp(
x
2
4kt
)
Rather then just using this kernel for noise-removal, a
linear scale-space theory can be built upon the appli-
cation of this convolution kernel. This idea has been
101
Bigot S., Kachi D., Durand S. and Mustapha Mouaddib E. (2007).
SPHERICAL IMAGE DENOISING AND ITS APPLICATION TO OMNIDIRECTIONAL IMAGING.
In Proceedings of the Second International Conference on Computer Vision Theory and Applications - IFP/IA, pages 101-108
Copyright
c
SciTePress
introduced by Witkin (Witkin, 1984) and Lindeberg
(Lindeberg, 1994). The linear scale space of an im-
age is defined as the set of smoothed images created
by convolving the image with Gaussian functions of
different scales.
The Gaussian function is the Green function of the
linear diffusion equation, i.e. the solution of equation
(1) with the δ-impulse as initial condition.
u(x;t) =
1
k
t
u(x;t) (1)
In the case of omnidirectional images, the major dif-
ficulties of these approaches have to deal with explic-
itly include shrinkage, dependency of the smoothing
results on the mesh connectivity. In order to be able
to deal with well defined globally parameterized do-
main we restrict the class of omnidirectional images
to those which can be described as functions on the
sphere, as defined in the model of Geyer.
Diffusion smoothing of surfaces then corresponds
to convolution of the surface with the spherical
Gaussian.
We propose in this paper to adapt the classical Wiener
filter and the Tikhonov regularization to the spheri-
cal images. We compare the denoising performances
of these methods with the classical spherical kernels.
The remainder of this paper is structured as follows.
In the next section, we recall mathematical framework
related to spherical harmonics, which will be used for
the solution of the spherical diffusion equation, and
convolution on the sphere. In Section 3, we define
the most important classical spherical kernels, and we
describe the spherical form of Wiener and Tikhonov
filtering. Section 4 shows results of the spherical fil-
tering process applied to synthetical and real indoor
omnidirectional images. Finally we give some con-
clusions and future work.
2 SPHERICAL THEORY
Filtering is based on convolution operators. We
present in this section the various definitions of this
operation and the mathematical tools used for its im-
plementation.
2.1 Spherical Fourier Transform
We parametrize the unit sphere, embedded in R
3
, by
using the spherical coordinates η S
2
: η(θ,ϕ) =
(cos(ϕ)sin(θ), sin(ϕ)sin(θ),cos(θ)) with ϕ [0,2π[,
angle of longitude and θ [0,π], angle of colatitude
(latitude + π/2).
The effect of planar diffusion smoothing can be well
understood in the frequency domain as a low-pass fil-
ter. Since we are going to carry out the corresponding
analysis on the sphere we need a spherical analog of
the Fourier transform. Such a tool exists in the expan-
sion of a function into a series of spherical harmonic
functions.
Notice by Y
lm
the spherical harmonic of l(N) degree
and order m as follows
Y
lm
(θ,ϕ) =
s
2l + 1
4π
(l m)!
(l + m)!
P
m
l
(cos(θ))e
imϕ
for m 0
where the P
m
l
(x) are the polynomials of Legendre as-
sociated with l degree and order m. We can notice that
the spherical harmonics of l degree form a subspace
of L
2
(S
2
) of dimension 2l +1 which is invariant under
rotations of the sphere. Since the spherical harmonics
form an orthonormal basis for L
2
(S
2
), we have
b
f
lm
=
b
f (l,m) =
h
f ,Y
lm
i
where the scalar product on the sphere is defined as
h
f ,h
i
=
2π
0
π
0
f (θ, ϕ)h(θ,ϕ)sin(θ)dθdϕ
The set of the coefficients
b
f
lm
is called spherical
Fourier transform or spectrum of f. For the imple-
mentation, we will use the sampling theorem (Healy
et al., 1998).
Theorem: Let f L
2
(S
2
) be a bandlimited function
of bandwith B, then:
b
f
lm
=
2π
2B
2B1
j=0
2B1
k=0
a
(B)
j
f (θ
j
,ϕ
k
)Y
lm
(θ
j
,ϕ
k
)
for
|
m
|
l < B. The sampling grid is the equiangular
or lat-lon grid with θ
j
=
π(2 j+1)
4B
et ϕ
k
=
πk
B
.
2.2 The Convolution on the Sphere
Two definitions were proposed to carry out a prod-
uct of Convolution on the sphere : that introduced
by Driscoll-Healy (Driscoll and Healy, 1994) and that
used, by Daniilidis (Daniilidis et al., 2002) and Wan-
delt (Wandelt and G
´
orski, 2001).
We begin by introducing some notations.
We represent the sphere S
2
as the quotient
SO(3)/SO(2) where S O(3) is the group of rotations
which acts on the sphere (Vilenkin, 1969). The rota-
tion of a function f L
2
(S
2
) by an element g SO(3)
is then defined with the operator Λ
g
such as
Λ
g
f (η) = f (g
1
η)
VISAPP 2007 - International Conference on Computer Vision Theory and Applications
102
Let f L
2
(S
2
) and g L
1
(S
2
)
The convolution product between f et g is the func-
tion on SO(3) defined by
( f
e
g)(R) =
S
2
f (R
1
η)g(η)dη
We have
f
e
g L
2
(SO(3))
with L
2
(SO(3)) L
2
(SO(3),dR) where dR is the
Haar measure on SO(3).
An important problem arises; whereas the functions f
and g are defined on the sphere, the product of con-
volution is defined on the group of rotations. Conse-
quently, it is clear that the product of convolution is
not associative. Moreover
f
e
δ
0
(R) = f (R
1
n
0
)
that implies a loss of symmetry in this definition of
convolution.
Now let us see another definition of the convolution,
introduced by Driscoll and Healy. Let f , h L
2
(S
2
),
we have
( f h)(η) =
SO(3)
f (Rn
0
)h(R
1
η)dR
where η S
2
, n
0
= (0,0, 1) represent the North Pole
of the unit sphere.
The effectiveness of this second definition becomes
obvious thanks to the theorem of convolution, shown
by Driscoll and Healy. For f , h L
2
(S
2
), we have that
\
( f h)
lm
= 2π
r
4π
2l + 1
b
f
lm
b
h
l0
(2)
We can observe that only the coefficients of spheri-
cal transform of the filter
b
h
lm
with m = 0 are used in
expression (2). These coefficients correspond to the
zonal harmonics Y
l0
and thus represent the invariant
part by rotation of the filter h. So the filtering opera-
tor is invariant under the actions of rotations.
3 SPHERICAL FILTERS
3.1 Spherical Kernels
The question which we put is to know to define
Gaussian on the sphere. B
¨
ulow (Bulow, 2004) pro-
posed to determine the Green function, as a solu-
tion of the spherical diffusion equation. The spherical
equation of diffusion is :
S
2
u(θ,φ,t) =
1
k
t
u(θ,φ,t)
where
S
2
=
1
sin(θ)
∂θ
sin(θ)
∂θ
+
1
sin
2
(θ)
2
2
ϕ
Let us recall that the spherical harmonics are the
eigenfunctions of the spherical Laplace operator :
S
2
Y
l,m
= l(l + 1)Y
l,m
Consequently, it is checked easily that the functions
u
l,m
(θ,φ,t) = Y
l,m
(θ,φ)exp(l(l + 1)kt)
are solutions of the spherical equation of diffusion.
To obtain the function of Green G, we force like initial
condition :
u(θ,φ, 0) = G(θ,φ, 0) := δ
S
2
(θ,φ)
where δ is defined by f (n
0
) =
S
2
f (η)δ
S
2
(η)dη.
By using the decomposition of spherical Dirac in the
basis of the spherical harmonics, we obtain
G(.,0) = δ
S
2
=
lN
r
2l + 1
4π
Y
l,0
=
lN
r
2l + 1
4π
u
l,0
(.,0)
From where finally for the Green function
G(θ,φ,t) =
lN
2l + 1
4π
P
0
l
(cos(θ))e
l(l+1)kt
we carry out then the spherical Fourier transform to
find
\
G(.,t)(l, m) =
(
q
2l+1
4π
e
l(l+1)kt
if m = 0
0 if not
Let us define another functions candidates to build
Gaussian on the sphere.
3.2 Spherical Gaussian
If we start with the definition of the Gaussian on a
plane, we can obtain the spherical Gaussian using the
inverse stereographic projection
G
s
(θ,φ,t) = e
(tan(θ/2)/t)
2
3.3 Poisson Kernel
This Kernel is given by the function
P(θ,φ,t) =
nN
(2n + 1)t
n
P
0
n
(cos(θ))
=
1 t
2
(1 2t cos(θ) + t
2
)
3/2
for t [0, 1]
This equation is the solution of the equation
S
2
P(θ,ϕ,t) = t
2
t
2
(hP(θ,ϕ,t))
SPHERICAL IMAGE DENOISING AND ITS APPLICATION TO OMNIDIRECTIONAL IMAGING
103
Let us calculate the Fourier transform of this function.
\
P(.,t)(l, m) =
(
p
(2l + 1)4π t
l
if m = 0
0 else
For more details, we refer to (Bulow, 2004).
3.4 Spherical Form of the Wiener Filter
We present in this section a Wiener filter adapted to
spherical functions.
We consider our original image f L
2
(S
2
), corrupted
by an additive white noise n. We denotes x = f + n,
the data, i.e. the corrupted image.
We seek to obtain the best possible estimate g of f
starting from our data x.
For that purpose we proceed to obtain the maximum
of the SNR defined as
SNR = 10log
10
E(
k
f
k
2
)
E(
k
f g
k
2
)
by minimizing the average quadratic error
e = E(
k
f g
k
2
)
The image f and the estimator g are related to L
2
(S
2
),
whose spherical harmonics form a basis. As a conse-
quence, we can use the theorem of Riesz-Fischer to
obtain
e = E(
lN
|
m
|
l
b
f (l,m)
b
g(l, m)
2
)
We propose to search g as the result of a filter (of
impulse response h) applied to the data x, which is
equivalent to write g as
g =
1
2π
x h
We apply the convolution theorem of Driscoll-Healy,
we obtain :
e = E
lN
|
m
|
l
b
f (l,m)
r
4π
2l + 1
b
x(l, m)
b
h(l, 0)
2
(3)
We must thus find the filter h which minimizes e, de-
fined in (3). We wish that this filter to be invariant
by rotation, i.e. of form
h(η) =
kN
b
h(k, 0)Y
k,0
(η)
That amounts determining the
b
h(k, 0) which mini-
mizes e, and so solve the equation
e
b
h(k, 0)
= 0
By considering that noise is a white noise with
standard deviation σ, we have E|
b
n(k, m)|
2
= σ
2
.
We assume
E|
b
f (k,m)|
2
=
c
k
2
where c is a constant to be chosen, we are led with
b
h(k, 0) =
r
2k + 1
4π
1
1 +
σ
2
k
2
c
3.5 Regularization Process of Tikhonov
In 1977, Arsenin and Tikhonov (Tikhonov and Ars-
enin, 1977) have proposed a method for determining
the best possible estimator.
The idea is to find
min
f H
1
(S
2
)
S
2
|
f
|
2
+ λ
k
f x
k
2
where as above,
k
.
k
represents the standard L
2
(S
2
)
norm. We recall that for a smooth function f with
compact support, we have
|
f
|
2
= f f
Therefore :
S
2
|
f
|
2
+ λ
k
f x
k
2
=
S
2
f f + λ
k
f x
k
2
However the spherical harmonics, which form an or-
thonormal basis of L
2
(S
2
), are the eigenfunctions of
Laplace-Beltrami operator. Consequently, we obtain
S
2
|
f
|
2
+ λ
k
f x
k
2
=
lN
|
m
|
l
b
f
2
(l,m)l(l + 1) + λ
k
f x
k
2
Then by using the Riesz-Fischer theorem, we have
J =
S
2
|
f
|
2
+ λ
k
f x
k
2
=
lN
|
m
|
l
b
f
2
(l, m)l(l +1)+λ
lN
|
m
|
l
b
f (l,m)
b
x(l, m)
2
Tikhonov methods aims to find the function f where
the minimum is achieved. However by using the ex-
pansion of f in the basis of the spherical harmonics,
f (η) =
lN
|
m
|
l
b
f (l,m)Y
l,m
(η) η S
2
We have to estimate the function J according to
b
f (k, p).
J
b
f (k, p)
= 2k(k + 1)
b
f (k, p) + 2λ(
b
f (k, p)
b
x(k, p))
VISAPP 2007 - International Conference on Computer Vision Theory and Applications
104
Then
J
b
f (k, p)
= 0
b
f (k, p) =
λ
λ + k(k + 1)
b
x(k, p)
Let g be the estimator which minimizes the function
J. We infer from the computations above that
b
g(k, p) =
λ
λ + k(k + 1)
b
x(k, p)
As for the Wiener filter, we suppose that the estimator
is the result of a filtering of the data, i.e.,
g =
1
2π
x h
with a h that is the impulse response of the filter of
restoration. We thus find
b
h(k, 0) =
r
2k + 1
4π
1
1 +
k(k + 1)
λ
The filter h of restoration has a form similar to the
Wiener filter, that is a classical result in the framework
of the planar images denoising.
4 EXPERIMENTAL RESULTS
We have compared our denoising process on synthet-
ical, real indoor and outdoor omnidirectional images.
The used mirror is parabolic and verifies the prop-
erty of single effective viewpoint. So we can ap-
(a) Office image (b) Office image noised
Figure 2: Noised Office image with additive noise
N(σ, µ) = N(40,0).
ply the Geyer model and process the filtering on
spherical images. For the implementation, omnidi-
rectional images are projected on S
2
and redefined on
the lat-lon grid. The spherical transform of images
and kernels are computed by using the yawtb tool-
box (http://rhea.tele.ucl.ac.be/yawtb/). The noise is
an additive white noise with varying standard devi-
ation σ and is added directly to the spherical images
Table 1: Optimal SNR for Office image N(σ,µ) = N(40,0).
Kernel Scale t SNR
Green 2.6650 10
4
19.7852
Gaussian spherical 0.0164 19.7842
Wiener 2903000 19.1174
Tikhonov 1900 19.1013
Poisson 0.9850 18.6118
as seen in Fig.2. The scale is an important parameter
in the process, and the result of filtering depends on
this value. Fig. 3 shows an example of filtering results
(Green filter) with different scale value. For each fil-
ter, we search (numerically) the optimal scale, which
provides the highest SNR. This SNR is calculated on
the whole sphere, while the projected omnidirectional
image is defined on the hemisphere. The optimal SNR
is approximately the same whatever the filter and the
type of images (indoor or outdoor) as proved in Fig.4,
Fig.5, and Table 1, even if the Green function and the
gaussian kernel provide often ligthly higher SNR. We
have obtained similar results for other types of noise
(like exponential noise. . . ).
However if we carry out a zoom in some images,
see Fig.6 and Fig.7 we can notice (visually) that the
Wiener filter and the method of Tikhonov give further
precision on edges. This precision is independent of
position on the sphere, (close or far from the pole of
the sphere). The SNR calculated at different regions
on the sphere are similar, see Fig. 6. We can also
show this conclusion for various images. These re-
sults have a considerable benefit if we want to detect
edges, or extract important primitive from images.
5 CONCLUSION AND FUTURE
WORK
We presented a panorama on the tools and the
methods to treat the spherical images. We adapted the
Wiener filter and Tikhonov method in that case and
compared their performances to classical spherical
gaussian filters. The SNR obtained for various types
of images prove the effectiveness of these filters for
the denoising application. Moreover the treatments
can be carried out in real time. The filters are
invariant by rotation, the results are thus independent
of the position on the sphere. We also noticed that
the Wiener filter and Tikhonov regularization are
relatively better, if we want to analyse precisely some
parts of images for edge detection for example.
For the continuation of our work, we propose to
SPHERICAL IMAGE DENOISING AND ITS APPLICATION TO OMNIDIRECTIONAL IMAGING
105