in: M. van Ginkel, P.W. Verbeek, and L.J. van Vliet, Improved Orientation Selectivity for Orientation Estimation,
in: M. Frydrych, J. Parkkinen, A. Visa (eds.), SCIA’97, Proc. 10th Scandinavian Conference on Image Analysis (Lappeenranta, Finland, June 9-11), 1997, 533-537.
Improved Orientation Selectivity for Orientation Estimation
M. van Ginkel, P.W. Verbeek, L.J. van Vliet Faculty of Applied Physics,
Delft University of Technology, Lorentzweg 1, 2628 CJ Delft, The Netherlands, e-mail:{michael,piet,lucas}@ph.tn.tudelft.nl
Abstract
Filtering of an image with rotated versions of an orien-tation selective filter yields a set of images which can be stacked to form an orientation space. Orientation space provides a means of analyzing overlapping and touching anisotropic textures. A set of rotated kth
or-der directional or-derivatives yields a discrete orientation space, which allows interpolation. Next we apply a de-convolution scheme that results in improved orientation selectivity. This scheme allows decomposition of noisy multi-orientation patterns.
1
Introduction
Images are composed of various types of structures, such as lines, edges and textures. These in turn can be characterized by various properties, the most impor-tant being intensity, scale and orientation. In this paper we focus on improved orientation selectivity for multi-orientation estimation.
An image can often locally be modeled as a weighted sum of translation invariant patterns or (paintbrush) strokes. Each stroke has a one-dimensional intensity profile and a typical orientation: the profile orientation across the stroke, see figure 1.
x
φ
Figure 1: An oriented pattern
We define directions as angles in the interval (0, 2π), thus making a distinction between an arrow pointing to the left or right. When we refer to orientation (0, π) we make no such distinction. For many applications a reliable estimator for a single locally dominant orienta-tion suffices, Kass & Witkin [4], Haglund [3] and Van Vliet & Verbeek [9]. There are a number of applications
for which an estimate of the local orientation will not be sufficient. Textures can often be characterized by a number of overlapping patterns with a different orien-tation. At boundaries between two single orientation regions, an estimator for a locally dominant orientation will give the wrong answer. In such cases a more ad-vanced analysis is required.
For another important property, scale, an extensive framework has been created over the last decade, ini-tiated by Witkin [11] and Koenderink [6]. In the scale space paradigm, scale is explicitly dealt with by embed-ding an image in a new image with one extra dimension, the scale dimension. This image is generated from the original by iteratively applying an isotropic (e.g. Gaus-sian) filter.
In order to deal with anisotropy an image or scale space can be replaced by a stack of directionally filtered versions of the original. The resulting orientation space or orientation+scale space is periodic along the orien-tation axis with a period of π. Orienorien-tation space and scale space differ in two aspects; the original image is part of scale space, but not of orientation space; scale space has a hierarchical structure, whereas orientation space does not. As shown in Van Vliet and Verbeek [10], such an orientation space can also be used to perform a segmentation of overlapping objects. The main ob-jective of this paper is to establish how to generate the orientation space representation of an image.
2
Constructing
an
orientation
space
Orientation space is created by applying rotated ver-sions of some orientation selective filter to the image. Initially the only constraint on the choice of filter is that it allows proper sampling of orientation space. The filters we use are all two dimensional and operate in the (x, y) plane. Yet, most of the ensuing (Fourier) analy-sis takes place in the orientation dimension and will be one-dimensional.
-1.5 -1 -0.5 0.5 1 1.5 -0.2 0.2 0.4 0.6 0.8 1 φ ω ωφ −π− − 2 π2 0 0 0 0 252 210 120 45 10 1 (a) (b) (c)
Figure 2: a.) Fourier coefficientsea(10)(ω
φ). The numbers show relative amplitude. b.) Fourier coefficients eb(10)(ωφ). c.)
The kernels a(10)(φ), dashed, and b(10)(φ), solid line
2.1
Directional derivatives
We will consider the possibility of using a Gaussian di-rectional derivative of order k as orientation selective filter. The Gaussian is only used for regularization of the derivative operators. Therefore the σ, which speci-fies the width of the Gaussian, is not of interest in what follows. The Gaussian derivative filter of order k in the φ direction is given by:
G(k)(x, y, φ) = (cos φ ∂ ∂x + sin φ ∂ ∂y) kexp( −x 2+ y2 2σ2 ) (1) Next consider an image I(x, y) that consists of a sum of strokes fi, with corresponding orientations θi. Thus
I(x, y) can locally be described by: I(x, y) =X
i
fi(x cos θi+ y sin θi) (2)
The response of the kth order directional Gaussian
derivative to the stroke fi when φ = θi, is denoted by
fik,max. The response of the directional derivative to the image I(x, y) can then be written as:
I(k)(x, y, φ) = G(k)(x, y, φ)∗(x,y)I(x, y) =X i cosk(φ− θi)f k,max i (x, y) (3)
Where ∗(x,y) denotes convolution in the (x, y) plane.
This can be rewritten as: I(k)(x, y, φ) =X
i
coskφ∗(φ)fk,max
i (x, y)δ(φ− θi)
(4) Convolution in the φ-dimension is denoted by∗(φ). The
last equation has a simple interpretation. In orientation space each stroke gives an impulse in the φ-dimension convolved with the following kernel:
a(k)(φ) = coskφ (5)
A property of this kernel is that the orientation selectiv-ity increases with increasing order k. This means that a specific selectivity can be chosen for a particular ap-plication.
A less desirable property is that not only the orienta-tion selectivity changes with the order k, but the radial frequency sensitivity as well. Therefore we will follow Knutsson [5] and decompose the filter into an angular and a radial part. The Fourier transform of a directional Gaussian derivative using polar coordinates ω and θ is given by: e G(k)(ω, θ, φ) = (jω)kexp(−1 2σ 2ω2) cosk(φ − θ) (6) Although separable, both the angular response, cosk(φ
− θ), and the radial response depend on k. Since we wish to vary only the orientation selectivity of the filter, we will keep the radial response fixed at k = 2 by multi-plying eG(k)by ω−k+2. This is a rather arbitrary choice,
but will suffice for the purpose of this paper. The ac-tual choice of radial function depends on the applica-tion. The Fourier transform of the resulting filter A is as follows: e A(k)(ω, θ, φ) = ω−k+2Ge(k)(ω, θ, φ) = jkω2exp(−1 2σ 2ω2) cosk(φ − θ) (7)
2.2
Sampling the orientation space
Given a certain filter, Freeman & Adelson [2] and Per-ona [8] address the problem of creating a set of filters that allows interpolation and minimizes the error be-tween this set and the filter response obtained by rota-tion of the original filter. They also present the filter constraints required to allow interpolation.
The following analysis is essentially identical to the derivations in the references given above and shows that our filter allows interpolation. The Fourier coeffients of 534
(a) A(9)(x, y) (b) eA(9)(ωx, ωy) (c) B(9)(x, y) (d) eB(9)(ωx, ωy)
Figure 3: Spatial and frequency domain versions of the filters A and B
cos(φ) are zero, except for ωφ=±1. Repeated
convolu-tion (k− 1 times) of the Fourier series representation of cos φ by itself, yields the Fourier coefficientsea(k)(ω
φ) of
a(k)(φ), which are then given by Pascal’s triangle:
e a(k)(ωφ) = 1 2 k k 1 2(k + ωφ) k odd, ωφ odd ≤ k k even, ωφ even≤ k 0 elsewhere (8)
The coefficients are depicted in figure 2a for order 10. The filter is bandlimited and the coefficients that cor-respond to frequencies higher than kωφ are zero. The
Nyquist sampling theorem states that we need 2k + 1 samples on the interval (0, 2π) to allow reconstruction, or equivalently 2k + 1 filters. The number of filters needed can be reduced to k + 1, by noting the following symmetry:
a(k)(φ + π) = (−1)ka(k)(φ) (9)
2.3
Improving the angular resolution
The kernel a(k)(φ) can be well approximated by a
Gaus-sian shape with a variance that decreases as k1. There-fore the peak width decreases as √1
k. Since the kernels
are periodic the notion of peak width only makes sense on the interval (−π 2, π 2). The 1/ √ k behaviour is dis-appointing, because one would wish that doubling the number of filters would halve the peak width. We will now show how this can be achieved.
Equation 4 shows that each stroke gives an impulse convolved with the kernel a(k)(φ). To improve the orien-tation selectivity we can simply deconvolve the resulting signal, because the blurring kernel is exactly known. As equation 8 shows, some Fourier series coefficients are zero, meaning that only a partial deconvolution can be performed. For order k the amplitude of the Fourier coefficients can be flattened (deconvolved) to:
eb(k)(ωφ) = 1 k+1 k odd, ωφ odd ≤ k k even, ωφ even≤ k 0 elsewhere (10)
In figure 2b the coefficients eb(k)(ω
φ) are depicted for
k = 10. The resulting response in orientation space is given by:
b(k)(φ) = sin((k + 1)φ)
sin φ (11)
In figure 2c the kernels a(k)(φ) and b(k)(φ) are depicted for k = 10. It is clear that the peak width is indeed re-duced, but at the cost of adding some side lobes. These lobes will cause problems if the amplitude of the differ-ent strokes differs too much. The first zerocrossing of b(k)(φ) lies at π/(k +1). This shows that the peak width
indeed decreases (nearly) linearly as 1/k.
The response along the φ axis in orientation space to a stroke is equivalent to the angular response of the filter that is used. So instead of performing the deconvolution, it is also possible to use a filter that uses equation 11 as its angular response. The Fourier transform of the filter created by replacing the angular part of equation 7 by b(k)(φ) is given by: e B(k)(ω, θ, φ) = jkω2exp(−1 2σ 2ω2)sin((k + 1)(φ− θ)) sin(φ− θ) (12)
In figure 3 the A and B filters are shown in the the frequency domain as well as in the spatial domain.
2.4
Measuring the rms orientation
in-tensity
The filters A and B do not measure the rms orienta-tion intensity, but respond to either locally symmetric or antisymmetric structures. For texture analysis we do not wish to discriminate between the two. The rms intensity is measured by the following procedure; First the filter result is squared, yielding an estimate of the local orientation energy, but with an image containing high frequency signals superimposed on it. By adding
φ
y
φ
x
x
y
C
A
D
B
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
slice CD
slice AB
orientation
(a)
φ
x
y
A
D
C
B
image
original
space
k=19
Figure 4: a) Schematic representation of orientation space. An x, φ slice as indicated is used in (d-f). b) A noise free image containing two superimposed patterns. The intersection of the slice with the (x, y) plane is indicated by line AB. The intersection with a second slice is indicated by line CD. c) Noisy version of (b). d) The slice indicated in (a) for the rms orientation intensity measured by filter A in image (b). e) Same as (d), but using filter B. f) Same as (e), but applied to image (c). (g-i) Same as (d-f) but using the other slice as indicated by line CD. Note the structures next to (f) indicating where in the slice each pattern gives a response.
the square of its quadrature counterpart [5] the high fre-quency signals will cancel. Creating a quadrature coun-terpart for filter B results in a filter that is discontinuous for even k and non-differentiable for odd k. Therefore it is not possible to use a quadrature filter approach. Instead the image will be smoothed by a Gaussian to remove the high frequency signals. The square root is taken of the smoothed image to yield the final estimate of the rms orientation intensity. Any interpolation must take place before the squaring operation, otherwise the necessary conditions to allow interpolation will no longer hold. All these operations take place in the (x, y) planes of orientation space.
3
Experiments
In figure 4 we show the results of applying the filters pro-posed in this paper to an image consisting of two over-lapping patterns. Each of the patterns has been created by calculating the Euclidean distance to some point in the image and computing the sine of the result. The amplitude of each individual pattern is one. The second test image was created by adding Gaussian noise with a variance of one to the first test image (SNR=6dB). We have computed the orientation space for both im-ages using both filter A and filter B for k = 19. We show two slices of the orientation space for the follow-536
ing combinations: noise free image and filter A, noise free image and filter B, noisy image and filter B. The slices show the rms orientation intensity, using σ = 5 for the smoothing.
The resulting slices show that the orientation selectiv-ity is indeed improved. It is also clear that the scheme is quite robust with respect to noise. Furthermore, close examination of the slices shows that there is indeed sub-pixel information present along the orientation axis, and it is therefore possible to interpolate in (the original) orientation space to get improved accuracy.
4
Conclusions
We have shown that the peak width of the response in orientation space can be made to depend linearly on the number of filters k. This occurs at the cost of intro-ducing some side lobes. Simple experiments with this scheme show that it works well for noisy images.
The final challenge is the interpretation of orientation space. Two touching orientation fields will indeed give separate responses in orientation space, but it will still be difficult to use this information to accurately determine the boundary. Interesting related work on junction classification has been done by Andersson [1] and Michaelis [7] and might provide a starting point for further analysis. We are currently investigating several possible applications of orientation space.
Acknowledgements
This work was partially supported by the Rolling Grants program 94RG12 of the Netherlands Organization for Fundamental Research of Matter (FOM) and by the Royal Dutch Academy of Sciences (KNAW).
References
[1] M. Andersson, Controllable Multidimensional Fil-ters and Models in Low Level Computer Vision, PhD thesis, Link¨oping University, Sweden, 1992. [2] W.T. Freeman and E.H. Adelson, The Design and
Use of Steerable filters, IEEE transactions on Pat-tern Analysis and Machine Intelligence, vol. 13, no. 9, September 1991, pp. 891-906.
[3] L. Haglund, Adaptive Multidimensional Filtering, PhD thesis, Link¨oping University, Sweden, 1992. [4] M. Kass and A. Witkin, Analyzing Oriented
Pat-terns, Computer Vision, Graphics and Image Pro-cessing, vol. 37, 1987, pp. 362-385
[5] H. Knutsson, Filtering and Reconstruction in Im-age Processing, PhD thesis, Link¨oping University, Sweden, 1982.
[6] J.J. Koenderink, The Structure of Images, Biologi-cal Cybernetics, vol. 50, 1984, pp. 363-370.
[7] M. Michaelis and G. Sommer, Junction Classifica-tion by multiple orientaClassifica-tion detecClassifica-tion, in: Jan-Olof Eklundh (ed.) ECCV ’94, Third European Con-ference on Computer Vision, Stockholm, Sweden, May 1994, pp. 101-108.
[8] P. Perona, Deformable Kernels for Early Vision, IEEE transactions on Pattern Analysis and Ma-chine Intelligence, vol. 17, no. 5, May 1995, pp. 488-499.
[9] L.J. van Vliet and P.W. Verbeek, Estimators for Orientation and Anisotropy in Digitized Images, in: J. van Katwijk, J.J. Gerbrands, M.R. van Steen, J.F.M. Tonino (eds.), ASCI’95, Proc. First Annual Conference of the Advanced School for Computing and Imaging (Heijen, NL, May 16-18), ASCI, Delft, 1995, pp. 442-450.
[10] L.J. van Vliet and P.W. Verbeek, Segmentation of overlapping objects, in: L.J. van Vliet, I.T. Young (eds.), Abstracts of the ASCI Imaging Workshop 1995, Venray, The Netherlands, October 1995, pp. 5-6.
[11] A. Witkin, Scale space filtering, Proc. Int. Joint Conf. on Artif. Intell., Karlsruhe, Germany, 1983, pp. 1019-1023.
!["Ordine di Malta : punti fermi per un rilancio", Roma 2008 : [recenzja]](data:image/gif;base64,R0lGODlhAQABAIAAAP///wAAACH5BAEAAAAALAAAAAABAAEAAAICRAEAOw==)