- Email: [email protected]
New method for vanishing point detection
CVGIP:
Vol.
IMAGE
54, No.
UNDERSTANDING
2, September,
pp. 289-300,
1991
New Method
for Vanishing
Point Detection
B. BRILLAULT-O'MAHONY Measurement
and Instrumentation
Centre, Department of Electronics and Information, London EC1 V OHB, United Kingdom Received
July 2, 1990; revised
The extraction of 3-dimensionalinformation from an imageis a central problem in computer vision. Location of the vanishing point of a line in the imageprovides its 3-dimensionaldirection. This paperis concernedwith the detectionof the vanishingpoints in an image.A new accumulatorspaceformulation for the whole line accumulationapproachis describedwhich fulfills an important property for robustness:the constancyof the detectability of a vanishing point whatever its location in the image plane. It is comparedwith the Gaussiansphereaccumulator space.Results are presented from interpreting indoor scenesof a nuclear plant. 0 1991 Academic Press, Inc.
INTRODUCTION
The major task in vision research is the extraction of 3D information from one or a number of different 2-D images. When the analysis of a single intensity image is considered, much information can be extracted using models of image formation and geometry. This paper examines an approach that assumes a perspective model, to enable groups of lines parallel in the 3-D scene to be extracted. For the method to be unseable on real tasks, a number of issues have to be addressed: robustness, noise, and uncertainty. The method must work over a wide range of images and produce perceptual groups that have a high certainty of being correct. In an indoor scene environment, directions of some features are more likely than others; e.g., for straight lines, the vertical direction and the horizontal directions parallel to the walls. The search for parallel lines in an image is very useful for the interpretation of the scene. The perspective projection of a set of parallel lines onto an image is a set of lines meeting at a common point, called a vanishing point. The vanishing point coordinates define the direction of the set of lines in the scene relative to the camera coordinate system. A number of methods for the detection of the vanishing points have been developed but for all of them the detectability of a vanishing point depends on its location in the image. The method proposed here is based on the accumulation principle but uses a new mapping that ensures the same detectability
December
City University, Northampton
Square,
26, 1990
of the vanishing point over the space. A further advantage is that the proposed accumulation space implementation is particularly simple. Usually, the vanishing points are detected by projecting lines, or points, each point being the intersection of a pair of straight lines, onto an accumulator space [l-4]. In order to reduce the search to a bounded closed set, Barnard [2] proposed projecting the lines of the image onto a Gaussian sphere centered onto the optic center and to use the Hough transform paradigm to detect the vanishing points. The projection of points onto a Gaussian sphere is equivalent to a resampling of the image plane. The method described here also resamples the image plane but the initial constraint is to keep the uncertainty of the vanishing point location almost constant over the image plane. In this paper the notion of uncertainty is defined, then the equation corresponding to the constraint is solved. The solution also leads to a bounded closed set which is merely a consequence of the method and not an a priori criterion. Once the lines have been accumulated in this new accumulator space, the local maxima are detected as the hypothetical vanishing points. The efficiency of the method is proved on a number of indoor scenes. Finally a comparison is made with methods based on the Gaussian sphere to highlight the main property of the new accumulator space. PREVIOUS WORK
The search for vanishing points consists of finding a small neighborhood in the image plane intersected by a sufficient number of straight lines. The image plane is resampled and bounded to form an accumulator space. The number of resampled straight lines crossing each cell is calculated using the Hough paradigm. Barnard [2] proposed projecting the lines onto a Gaussian sphere centered on the optic center. This projection is a bijection from the image plane to a bounded space. The search for the accumulation points is achieved by sampling the sphere using spherical coordinates. Unfortunately the spherical sampling (elevation, azimuth) is irregular and
289 1049-966019 1 $3 .OO Copyright 0 1991 by Academic Press, Inc. All rights of reproduction in any form reserved.
290
B. BRILLAULT-O’MAHONY
different in the x and y directions. Quan and Mohr [4] use the same method with a dichotomic approach: the sphere is sampled from a coarse resolution to a fine resolution in order to reduce the number of studied cells. The lines are classified by looking for vanishing points. Once a vanishing point is found with its associated lines, the lines are eliminated and the algorithm is performed again. Dickson [5] proposed a triangular sampling of the Gaussian sphere which is very attractive because it is isotropic; however, the computational efficiency has yet to be proved. Magee and Aggarwal [3] accumulate the projection of the intersection points of all pairs of straight lines in the image onto a Gaussian sphere. The accumulation is achieved using the arc distance between two points which leads to an isotropic search. The isotropy of the search for vanishing points using the accumulator spaces previously described can be ensured by an additional cost of complexity and only guarantees isotropy in 0, not necessarily in Y. As a result the probability of detecting a vanishing point depends on its distance from the center of the image. The reason for this dependence with Y is due to the increasing uncertainty of the vanishing point location with r. Kanatani [6] tests hypotheses on parallelism of three or more straight lines by using the projection onto the Gaussian sphere. The uncertainty of the line parameters is taken into account by using a threshold in a concurrency test which depends on the uncertainty of the intersection point coordinates. The author [9] has previously proposed to solve the problem by resampling the accumulator space such that the uncertainty remains approximately constant over the space. The variations of the uncertainty of the intersection point of a pair of lines is statistically estimated, but for Y fixed, the dispersion of this uncertainty is large; also, the solution found is not really isotropic in 13. It will be seen in this paper that a criterion based on distance between intersection points does not seem appropriate to a problem initially defined as the minimization of distances between a point and lines, whatever the representation used. The method described in this paper proposes a solution to the problem of the consistency of the vanishing point detection over the image plane, without increasing the algorithmic complexity. PREPROCESSING
The grey scale image of an indoor scene is digitized using a standard CCD camera. Straight lines are detected in the image using Berthod’s line finding algorithm [71 after the edges are detected using the Canny operator [8]. Each line approximates an edge in the image and is represented in a list by end point Cartesian coordinates.
ACCUMULATOR
SPACE
Notations In the following analysis vectors and matrices are represented in bold. The transposed matrix of M is noted M’. The vector from A to B is denoted AB and AB is its algebraic length. The line supporting the segment AB is denoted (All). The expected value of a real value x is denoted E(x). E(x 1y) is the expected value of x when y is known. The variance of a real value x is noted Vx, and the covariance of the real values x and y is denoted Vxy: Vx = E(x2) - E2(x) Vxy = E(xy) - E(x)E( y). The symbol of differentiation
Uncertainty
is noted 6.
Definition
Under perspective transformation, parallel lines in a 3-D scene are projected onto concurrent lines in the 2-D image. Ideally, in a man-made environment, many lines are parallel, e.g. the edges of a wall and a door frame. In practice this is only approximately true. However, parallelism is a useful concept for the representation of the door frame location. Therefore the first source of uncertainty is the difference between the real scene and the symbolic representation describing this scene, e.g., CAD representations. The second source of uncertainty comes from the image acquisition, digitization and preprocessing. The acquisition of the image is modeled by a pinhole perspective model. This does not take into account the numerous possible distortions of the image caused by the lens. The digitization stage constrains the total error on a pixel location to be at least 0.5 pixels. The Canny edge detector has an accuracy limited by the value of the parameter sigma of the Gaussian filter [S], The line finder may introduce nonnegligible errors in the slope of a segment. To clarify, let us assume that the camera calibration parameters are known (though unnecessary for the vanishing point detection) and the ideal representation of the scene is also known, e.g., in the form of a CAD description. For a segment the error of the end point location is defined as the difference between the observed end point location and the projection of the corresponding ideal end point in the scene using a pinhole model (for which the parameters are known). As a result of the uncertainty, lines assumed parallel in the 3-D scene have their images only approximately concurrent in the image (see Fig. 1).
VANISHING \I P
Lp---- L5 L2 i”\\\ L3
L4
FIG. 1. Image of lines assumed parallel in the scene. They are only approximately concurrent in the image.
Uncertainty
Modelling
The sources of the uncertainty are various and it is difficult and inefficient to try to model them separately. The simplest model commonly used for this type of uncertainty is the Gaussian law, with an expected value equal to zero and a standard deviation equal to cr. Therefore the error on each end point location is assumed to obey the same Gaussian law. Only the component of the error perpendicular to the line is needed as it is responsible for the errors of the line parameters (polar angle and distance from the origin). Let the matrix V = [zY $‘I be the covariance matrix of an end point A. A point M of the uncertainty ellipsoid is defined by AMfV-‘AM
I 1.
In the following u is the unit vector of the line L corresponding to the segment AB and V(u) is the variance of the error in the perpendicular direction to L (Fig. 2), V(u) = (ufV-‘u)*det(V).
(1)
In what follows V is assumed to be constant over the image plane, V(u) = CT& where I is the identity Problem
291
POINT DETECTION
(2)
parallel lines onto the image. The location of P in the image depends on the viewpoint. Searching for P by accumulation is equivalent to counting the real lines in the image passing through a neighborhood of P. A line (L) defined by a line segment in the image is assumed to have its vanishing point in P if D(L, P) < r(P), where D(L, P) is the distance between the line L and the point P and r(P) is the radius of the chosen neighborhood N of P in the perpendicular direction of (L). The constant quality for the detection of P irrespective of its location in the image plane increases the robustness and meaning of the detection. The value of the radius r(P) of the neighborhood N defined above should be proportional to the accuracy with which the line L is knotin in the vicinity of P, that is to say proportional to the uncertainty of the distance D(L, P) from L to P. The neighborhood N has a corresponding neighborhood N’ in the accumulator space. For practical reasons this neighborhood N’ should be constant and have a simple shape. Let (T be the uncertainty of D(L, P). The problem is to jind a transformation Tfrom the image plane to an accumulator space such that the expected value of u’ = T(o) remains constant over the accumulator space. This transformation Twill completely define the new accumulator space (Fig. 3). Statistical
Expression for the Uncertainty
The uncertainty of the end points of the segment AB results in uncertainty u on its distance to a point P. Polar coordinates are used, with pole 0 and axis Ox. Let 0 be the center of the image, r the distance from 0 to P, 6 the polar angle of OP, (Y the polar angle of (AB), and d its distance from 0. Q’ is the intersection of the line (AB) and the circle (C) with center 0 passing through P and CX~ is the polar angle of Q’ (in the ideal case (Y~= 8, see Fig. 4). Let ap be the uncertainty of P along (C). The variance of V~ over all the lines having P for vanishing point is equal to the variance of the error of P in the tangential direction, (r aa!,)*. Let the line (AB) be fixed; the uncertainty cp is equal to
matrix. ai = r*E(&:),
Dejinition
(3)
Let P be the vanishing point to be determined, which is the common point of the ideal projection of the ideal
FIG. 2. Geometric interpretation
of V(u)
FIG. 3. Transformation space.
from the image plane to the accumulator
292
B. BRILLAULT-O’MAHONY
image plane
i/L2
0.002s.
:
.
i
:
:
:
:
0.002
:
i
;
;
’ i ;:
i
0.001s.
j::;::
; :
.!
0.001 0.000s,
FIG. 4.
Line (AB) passing near the point P. 16
d
FIG. 5. Average value of l/L* in function of d. Each dot corresponds to the average value of l/L* over the lines of an image at a distance d from the origin. Twelve images of 6 different scenes have been processed.
where E(Gcu:) is the variance of the noise on (YI. Using the polar equation of the line it follows that r cos(a, - a) = d. This equation is derived to provide the expression 6(arJ, when (r2 - d2)“2 is not too small:
126
64
of
From a scene point of view, the length of a segment of the image depends on its vanishing point location. The distance between a segment and its vanishing point is constrained by
8a1 = 6a - 6d/(r sin(al - ar)) = &Y - 6dl(r* - d2)‘j2.
G’Q’IL
> D,/Df,
(6)
(4) Let L be the length of the segment AB and b the distance between the centroid G of the segment and the projection 0’ of 0 onto (AB), when (r* - d2)1’2 is large enough, the value of cri is given by (see Appendix 1, (A1.4), (AlS), 641 m
l/2 + 2b2/L2 1 - d2/r2 4br
+ L2(1 - d21r2)‘”
VW.
To make the expected value of oi explicit, a statistical model of the a priori repartition of the lines meeting at P within the image is required. In what follows, the image is approximated by a circular disk, in order to simplify the model. The lines are assumed to be distributed at random within the image; i.e., the centroid of the segments AB is uniformly distributed on the image. The length of these segments obeys a probability law of density f(L), which is assumed to be uncorrelated with the centroid location. From an image point of view, these assumptions are not exact because of the segments crossing the boundaries of the image. Therefore, for simplicity, the segments are not constrained to lie entirely within the image, except their centroids, so that the assumptions are consistent. It will be seen that ,?(a*) depends on E(1/L2). A statistical study of l/L* on 12 real images demonstrates the validity of the assumption of the dependence of l/L* from the centroid location (Fig. 5).
where Df is the depth of field and D, the focusing distance. Figure 6 demonstrates that this constraint has a very small effect on the expected value of (l/L2). This is not surprising because the constraint filters the small segments which play a preponderent part in E(1/L2). In fact the distribution of the lengths is mainly due to the line finder algorithm and little to the viewpoint. Let a be E( 1/L2). Using the independence assumptions between l/L2 and (b, d, v) and the fact that E(b 1 d) = 0 for symmetry reasons, the expected value of C; over the set of segments AB at a distance d from the origin is equal to (from (5)) E(a;(d)
= (2ar2 + (l/(1 - d2/r2))/2 + 2aE(b2)d)l(l
- d2/r2))V(u).
(7)
l/L2 I 0.002s. 0.002
1 ,.
0,001s 0.001
.i :*
. : j
. : ! * , : .fi
I
; j
: I.
: 1
: $
i
:i:iii
,,
0.000s. ‘
: 16
64
I12
160
r
FIG. 6. Same graph as in Fig. 5, over the lines filtered by the constraint (6) with D,lDf = 1.3.
VANISHING
POINT
This relation ignores the constraint (6), but still E(b2/ L2/d) gives more weight to large values of b and small values of L corresponding to the lines satisfying (6), which justifies (7). Remark.
minimal
The previous argument only holds when the value allowed for L is small enough.
The centroid G is located at random on a chord of the image disk, the length of which is equal to 2R(l - d2/ R2)‘12; then (Appendix 2): E(b21d))
= R2(1 - d2/R3)/3.
= (2ar2 + (l/(1
The value of the uncertainty of the line location near P depends on the parameters of the corresponding segment in the image. It has been shown that it is possible to express the expected value of this uncertainty by using statistics of the line segment parameters. The expression for the expected values depends only on the distance of the vanishing point from the image center. This result allows almost the same detectability to be guaranteed for any vanishing point, by resampling the image plane proportionally to the expected value of the uncertainty. Transformation Accumulator
Thus E(aj$d)
from the Image Plane to the Space
In
order to simplify the notation, Cr represents in the following. If the transformation T exists it is defined by
- d2/r2))/2
(E((T~))“~
+ 2aR2/3)V(u).
The distance d of the line is smaller than or equal to r, within the uncertainty of r. When d is smaller than r l/(1 - d21r2) = 1 f &(llr2) where c(llr2) is infinitesimal sumption in (2) gives
293
DETECTION
with Ilr2.
Using
the as-
P(r,
P’(x’(r,
19), y’(r,
19)),
x’, y’ being the coordinates of the cell of the accumulator space. Let Q’ the point of the circle (C) located at the distance @,,from P; then T(Q’)
A line having P for its vanishing point satisfies d I r. Equation (8) has been found using (4) under the assumption that (r2 - d2)“2 is not too small. If r > R, d is always smaller than r, and (8) holds and is independent of d. When r I R and d = r, the definition of aP has no sense, as the line (AB) may cross once, twice, or not at all the circle (C) in the vicinity of P. When P is far enough from the origin, a very small proportion of lines are likely to be such that d = r, and (8) is nearly always valid. If r is very small, the radial uncertainty of P has no sense. The uncertainty neighborhood of P is a circle with a constant radius c2, where u2 is defined as the expected value of PQ2. It is shown in appendix 3 that near the origin
f3) 7‘
= T(P(r,
8 + ?F,lr))
= P’(x’(r,
8 + F,lr),
y’(r,
0 + CJr)).
In the accumulator space the uncertainty e/:, corresponding to the uncertainty ZP in the image, is equal to the distance of T(P) to T(Q’). Thus -I2 = (P’(x’(r,
0 + Fp/r), y’(r,
@P
- P’(x’(r,
01, y’(r,
8 + Cp/r))
8>)12
= (x’(r, 8 + ZJr) - x’(r, fl))2 + (y’(r, 13+ @Jr> - y’(r, 13))~.
After linearization
it leads to
-I2 = (6x’/68)2Z~lr2 UP
+ ((6y’/6f3)2Z~lr2.
Thus ((6~‘/68)~
Therefore from (8) and (9) and the previous remarks, uncertainty of P may be approximated by
the
+ (Sy’/619)~) = FL2r2/(ri.
One solution of (11) is y’ = ZLr0lcT
if r is large enough
E(ai)
if r is small
E(u2) = cu;,
where a = E(1/L2)
= (2ar2 + c)oi
and c = l/2 + 2a R2/3.
(10)
(11)
P
x’ = x’(r).
x’ is independent of 8 and may be any bijection from [0, +m] to a bounded interval. For simplicity, x’ has been
B. BRILLAULT-O’MAHONY
Geometrical Interpretation
chosen so that y’
= 8x’lk,
where k is a scale .factor determined resolution. Therefore
by the expected
x’ = (Cf~/aO)krl(2a r2 + c)“~ y’ = (i7-~/cro)rel(2a r2 + c)l12.
(12)
The number of straight lines passing near P(rO, &J is represented by the number of curves crossing the line parallel to y’ axis, y’ = y’(rO, e,), in the range [P’ - iT;j, P’ + Gij] (see Fig. 7). The counting of these lines is performed by accumulating curves dilated by a vertical kernel with a half-width equal to $, . The accumulation of dilated curves prevents the same line from being counted twice. A line passing near P such that d E ]r, r + a[ does not cross the line y’ = y’(rO, &). When P is far from the origin such a configuration is unlikely, but when r is small its probability increases. The problem is solved if the x’ resolution around P’ is less or equal to G. This could be done by choosing the constant k in Eq. (12) such that x’(r) has a resolution equal to (T for small r, but the resolution of x’ has to be consistent with the resolution of y’ (e.g., two curves crossing at x’ = (2n + 1)/2 must cross the same uncertainty neighbhorhood either at x’ = rz or at x’ = n + 1). The two constraints on k are incompatible for large r. The constant k is chosen to ensure the consistency with y’ resolution and the constraint for the small value of r is fulfilled by a pre-accumulation stage in the (d, (Y)parameter space (see later). A transformation T has been found such that any point in the image plane is simply transformed to a single point in the accumulator space and the search for vanishing points is reduced to the search for local maxima in the accumulator space. The sampling is isotropic and homogeneous with respect to the uncertainty criterion. Moreover the accumulator space is bounded.
e=2n
p’ = x’(p)
u’
8= (L) direction
to
FIG. 7. (x’(r),
y’(r,
x
0)) accumulation space.
and
8’ = Blk.
In this sector the neighborhood of uncertainty is a circular arc with arc length FL. Thus, the transformation defined is equivalent to a resampling of p into s(p) followed by a projection onto a cone with vertex 0, axis equal to the optical axis, and angle /3, which is eventually developed onto a plane, so that T=d*p.s
Image plane s P(P, 0) -
Image plane p w-4, 0) -
>
Accumulator Cone d space m(p), 0, p) p’h et)
M is the intersection between a cylinder with radius s(p) and the cone with the same axis and the angle p; M has the Cartesian coordinates (s(p) cos 8, s(p) sin 8, p’ cos p), where p’ is its polar distance. Then
s(p) = p’ sin p
and
8’ = Bs(p)lp’
= 8 sin p.
This gives sin /3 = l/k and s(p) = x’(p)lk. Pre-accumulation
Stage
The lines are first accumulated with respect to d and CK, the resolution of which is determined by E(6d2)‘12 and E(&x~)“~ averaged over all possible b (see Appendix 1, (Al .4), (Al .5)), and for each cell (d, cu)the best estimates of the line parameters are computed. Only one line per cell (d, a) is accumulated in the accumulator space w, Y’). When r is small, d is also small and the expected value of E(6d2) over all b is equal to c, i.e., to (T2. This means that all the lines such that d E [r, r + G[ and (Y’ = 8, are accumulated in the same cell (d, CX). Therefore, for r small, the resolution in x’ is equal to (+ around x’(r) and the lines satisfying d = r cross the uncertainty neighborhood around P’ defined above. The preaccumulation stage avoids a weighting process and thereby increases the significance of the peaks of the accumulator space (x’, y’) described above, as the peak value represents the number N of directions meeting at the same point. For instance, if N = 3 the point may be a corner, if N > 3 the point is likely to be a vanishing point.
liiz&cl A0=n
C = T(L)
e0x,
The accumulator space (x’, y’) described previously is isotropic in 6’ by definition. In the continuous plane, it may be represented as a sector of a disk. In this disk, a point P’ with the polar coordinates (p’, 6’) corresponds to the point P(p, 8) in the image plane and as defined by using (12)
VANISHING
Uncertainty
of the Vanishing
POINT DETECTION
295
Point
It may be important to have a first evaluation of the uncertainty of the vanishing point location, once detected. Let (T; be the variance of the distance Q’P, where Q” is the intersection point of the line (AB) and the line (OP) (see Fig. 4); it is equal to (+z = E(Q’p2)
= r2E(QP2)/d2
= r2c21d2.
(13)
The line such that d = 0 does not provide any information about r and the corresponding uncertainty is infinite; so is the expected value of CT: over all possible d for r large enough. But here the peak in the accumulator space has already been found, which means that the shape of the distribution of Q”P2 around P no longer matters. The information about r is provided by lines passing far from the origin. At least three lines with different directions should cross the neighborhood of P to give it meaning as a possible vanishing point, and therefore at least two lines passing far from the origin, typically at [E(d2)]1’2. Therefore, for large r, using the expression of E(d2) found in Appendix 2,
FIG. 8. Comparison between a “whole line” accumulation proach and an “intersection point” accumulation approach.
ap-
lines classified with it are taken away from the list and the algorithm is performed on the rest of the lines. However, a number of lines in the image may pass through several vanishing points and such a method may ignore important directions in the image. The method used here removes the lines from the list after each candidate vanishing point has been detected, as described above, but only if the number of these lines is high enough to be responsible for several local maxima in the accumulator space. The detected vanishing points correspond to the highest peaks in the accumulator. Once all the candidate vanishing points have been detected, the 1 r2g2 lines are reclassified by a maximum likelihood method in a2 2r2(2ar2 + c)a$R2 cr=--= 2 E(d2) (14) order to take into account the specific uncertainty of each line location with reference to each vanishing point (in where c$~is an estimate of the radial uncertainty on P. contrast with the statistical method described above). If r is small, using the same considerations as for Z,, , c+? This allows the classification of one line with several vancan be replaced by the expression of a2, say c(ri. The ishing points and solution of the ambiguities by appropriexpression (14) for 6, gives a very coarse idea of the ate methods not described here. An additional filter for vanishing point uncertainty and it should be considered discriminating the vanishing point candidates selects the as temporary information. This approximation is refined pairs or triplets of perpendicular directions, which are by using a Kalman filter applied to each line classified hypothesized to be the main directions [9]. The descripwith the point P. tion of such methods is beyond the scope of this paper. For the resolution to be consistent with this uncertainty value, the derivative of r relative to x’ must be RESULTS lower than &,.. The sampling is not regular according to The accumulator space depends on the statistic of l/L2 &,., which does not matter, as c+~does not affect the linewhich depends on the complexity of the image and on the counting in the whole line accumulation approach depreprocessing stages. The size of the image used here is scribed above. This would not have been the case if only 256 and only lines with a length greater than 15 have been the intersection points had been accumulated (see Fig. 8). accumulated (lines with a length below 15 are very nuRemark. It is possible to choose x’(r) in order to have merous and do not provide meaningful information as a sampling approximately regular according to &-, but it their uncertainty is very high and they often correspond slightly increases the complexity. to noise). Thus the pre-accumulation in (d, a) results in a narrow distribution of l/L2 around a = 0.0019, deterVanishing Point Detection mined experimentally over a number of images. The problem posed by local maxima detection in the The method has been applied to various types of inaccumulator space is that a set of close parallel lines, e.g. door scene. Figure 9 shows the initial images. Their rea set of pipes, usually produces several local maxima spective accumulator spaces are displayed in Fig. 10, corresponding to the intersection of this set with various with all the lines accumulated, Figure lOa, and with a lines. number of lines removed, Fig. lob. The final result of the Some of the other methods [3, 41 use the following classification is displayed in Fig. 11. All the main direcstrategy: once a vanishing point has been found, all the tions represented by a sufficient number of segments
B. BRILLAULT-O’MAHONY
296
FIG. 9a. Initial images 1 and 2.
FIG. 9b.
Initial images 1 and 2 and the extracted segments.
have been found. False candidates have been found, corresponding generally to corners. They are filtered out by using the constraint (6) which is usually not satisfied by the corners. It can be noticed that vanishing points have been found corresponding to a small value of Y (18), to an average value of r (135) and to an infinite value of r, with peaks of similar shape along y’ (see Fig. 10).
COMPARISON OTHER METHODS
Comparison with Intersection Accumulation Method
Points
Figure 8 demonstrates that all the lines around Q are counted independently of the value of urr but the count-
VANISHING
POINT
DETECTION
297
FIG. 10. (a) Accumulator spaces with all the lines accumulated. For clarity the image has been scaled from 0 to 255. The axis y = 0 corresponds to the vertical lines. For convenience, the 0 reference has been shifted by a/Z. (b) Image 1 vertical lines parallel to the camera have been removed; image 2 vertical and horizontal lines parallel to the camera have been removed. For clarity the images have been scaled from 0 to 255.
ing of the corresponding intersection points should take or into account. Equation (13) shows that for Y fixed, the disparity of or is very large, and thereby the disparity of the error of the distance between intersection points is very large too. Whatever the representation used, the disparity of the error of the intersection points is very much larger than the disparity of u2, which means that the concurrency test based on the intersection point accumulation approach has a poorer significance than the test used in the whole line accumulation approach. The best estimate of the intersection point of a set of straight lines in the least mean square sense minimizes C(di/af)2. This criterion may be written in terms of intersection point coordinates and the associated covariance matrix (it is complicated because the intersection points are not independent). Because of the disparity of the vari-
ance along the Y direction, an accumulator space based on the intersection point accumulation should necessarily be associated with a weighting process using coefficients equal to the inverse of this variance, so that only significant intersection points would be taken into account. If the dependence of the intersection points is ignored, the value of the peaks obtained would approximately represent the inverse of the corresponding intersection point variance along r. This value is not significant because of the disparity of this value and because it provides no information about the number of lines or directions meeting at this point. Albeit the set of the intersection points and their associated uncertainty matrix theoretically contains sufficient information about the best estimate of the intersection point of a set of straight lines, the corresponding tests
B. BRILLAULT-O’MAHONY
-
--
FIG. 11. Final classification: the top and the bottom left images corresond to the main directions detected; the bottom right images show all the lines classified with a main direction.
have a poor significance. Besides, the information is redundant (5n(n + 1)/2 data in contrast with 5n data in the line accumulation case); i.e., the data are not independent, which is a serious difficulty when reasoning with uncertainty.
cess transforms L into a big circle C, P into P’, and Q into Q’. P is assumed to be near Q, so that the distance D’ from C to P’ can be approximated as D’ = D(P’, Q’), and its uncertainty (T’ approximated by the projection of the uncertainty U. Let h be the distance between 0 and the center of the Gaussian sphere; then
Comparison with Gaussian Sphere Method
If the constraint of the constant detectability is not fulfilled the detection of a vanishing point, i.e., a main direction in the scene, depends on the viewpoint and the orientation of the image. It is important to evaluate the behavior of the more commonly used accumulator space, the Gaussian sphere, with respect to this constraint. (See Fig. 12). Let r be the distance from Q to 0. The projection pro-
d/u
= O’Q’/O’Q
= p/(r2 + h2)l/*,
where p is the radius of the sphere. The constraint U’ = constant implies that (9 + h2) is proportional to (~2. It is shown in Appendix 3 that it is possible to approximate u2 for any ’ by iF2 = (2a r2 + c)cri
which leads, when applied to the gaussian sphere, to O’O=h O”Q=r
h = (cl(2a))1;2. P
FIG. 12. Projection of an image segment onto the Gaussian sphere.
(15)
This relation means that there exists a value for h such that the uncertainty is the same at the origin and at infinity. A large value of h favors the detection of vanishing points located near the center of the image, whereas a small value of h favors the detection of vanishing points located at infinity. Here, the neighborhood used for the counting of the lines is a circular spherical cap with a constant radius a’.
VANISHING
POINT
The Gaussian sphere method appears to be a correct accumulator sphere for the constraint considered here in the context of the whole line accumulation, provided a correct value has been chosen for h and an isotropic sampling is used. This result does not hold when only intersection points are accumulated. But an isotropic sampling of the sphere in the whole line accumulation approach implies an involved process, whereas this problem does not exist for the accumulator space described previously.
299
DETECTION
6d = (b/L + l/2) 6A . u1 - (b/L
- l/2) 6B * ui.
(Al .3) The variance of (SA . ui) and (6B * ur) is V(u) (see (1) in the text), therefore the second moments of &Y and 6d are E(W)
= 2V(u)/L2
(Al .4)
E(6d2)
= (2b2/L2 + 1/2)V(u)
(Al S)
E(6a 6d) = -2b
(Al .6)
V(u)IL2.
Let D be the distance from a fixed point P to the line (Al?); D is equal to
CONCLUSION
This paper has considered the problem of extracting 3D information from a perspective model for the projection of parallel lines onto an image. A new accumulator space for the detection of the vanishing points has been presented. This accumulator space is isotropic and bounded; moreover it ensures a constant detection quality of the vanishing points over the space. This results in a more robust and predictable method that works on a large number of different images of indoor scenes. In addition it is easily implemented. It has been compared with different methods using the Gaussian sphere as an accumulator space. This new method has been used to determine successfully the main directions of indoor scenes. Further work is concerned with the characterization of the signal to noise ratio of the accumulator space, which will be used for scoring each class of lines associated with a vanishing point. APPENDIX 1: VARIANCE OF THE PARAMETERS OF THE LINE (AB) AND THE DISTANCE FROM A POINT P TO (AB)
Let (Ybe the polar angle of the line (AB), d its distance to the origin, u its unit vector, and uI the perpendicular unit vector. Let G be the centroid of the segment All, L its length, 0’ the projection of the origin 0 onto (Ail), and b the algebraic distance between 0’ and G (Fig. 4). The line (AB) is known within the errors of measurement 6A and SB of the end points A and B. These errors results in errors on (Y, d, II, and ul: 6a = (SB * u1 - 6A . u,)lL
(Al.l)
6u = 6au,;
(A1.2)
And 6Ul = -6cY u.
D = AP * ~1. The derivation
of (A1.7) provides 6D, 6D=
-6A*u,
+AP.Sul.
Let Q be the projection of P onto (Ail); using (Al .2) and (A1.3) and then substituting AP . u by AQ and (AQ - L) by BQ, it is found that 6D = @Q/L) E(6D2)
6A . u1 - (AQlL)
6B . uI
Let Y’ be the distance rewritten
from 0’ to Q; (A1.8b) may be
E(6D2) = (2(r’ + b)2/L2 + 1/2)V(u). O”Q’ = r2 - d2 and O’Q = Y’. The vector QQ’ only upon the error on the end points, therefore stitution of Q by Q’ in (Al.ga) only adds terms degree than 1 in 6A * u1 and 6B * ul. Thus the approximation can be made: r’
APPENDIX
(A1.8a) (Al .8b)
= (AQ2 + BQ2)V(u)lL2.
=
(r
2
2: EXPECTED
-
(Al .9) depends the subof higher following
d2)“2.
VALUES OF dZ AND b2
Using the statistical model described in the text, the line (All) may be generated by any centroid G located on the chord (S) of the image disk, (S) being at the distance d from the center of the image, G is assumed to have a uniform density of probability on (S), therefore the density of probability of d is proportional to the length of the string (S): f(d)
Nowd=OA~ul,then6d=6A~ul+OA~6ul=6A~ul - OA * u 6a. Replacing OA . u by (b - L/2) and S(Yby (Al. l), it is found that
(Al .7)
= 2KR(l
- d2/R2)‘12.
The constant K is found by integration (-R
I d I R):
of f(d)
over d
300
B.BRILLAULT-O'MAHONY
K = ll(7fR2).
If r smaller than R then (integration
The expected value of d2 is equal to
by parts)
E(d2) = R2/2(1 - r(1 - r2/R2)3’2/ (R Arcsin (r/R) + (R2 - r2)l12)).
E(d2) = (2/(7rR)) j-R d2(1 - d2/R2)“2 6d = R2i4. When r is small (A2.1) E(d2) = r2/3, The expected value of b2, when d is known, is equal to and substituting E(b*ld)
= (1/(2(R2 - d2P2) \($:;:;:,2
b2 ab (A2 2)
= (R2 - d2)13 APPENDIX
3:
CALCULATION
OF uz
Let u2 be the expected value of 6D2, square of the distance from a given line (AB) and a fixed point P. The variations of 6D are due to the errors of measurement of the end points. Developing (A1.9) by substituting r’ by (r2 - d2)‘12, this gives u2 -_
2r2 I 4br(l - d21r2) I 2b2 i L2 L2 L2 + ; - g
1
E(02) = (IOU r2/9 + 2a R2/3 + 1/2)oi.
(A3.3)
It is more convenient to use an approximation valid for any r. Moreover the comparison with the Gaussian sphere method requires that such an approximation can be made. Using (A3.2) and (A3.3), E(cr2) is approximated for any r by E(a*)
= (2~ r2 + c)cri,
(A3.4)
where c = 2u R2/3 + l/2. ACKNOWLEDGMENTS This research was funded by Electricitt de France and was performed in the Machine Vision Group at the City University, London.
V(u)
If d is fixed and a2 is averaged over all the segments All producing the line (Al?), then by using (A2.2) and the fact that E(b) = 0 for symmetry reasons, it is found that E(u21d) = (2ar 2 - 8ad2/3 + 2aR2/3 + 1/2)V(u), (A3.1) where a = E(1/L2) (for the independence between a and the other parameters, see the paragraph “uncertainty modeling” and Fig. 11). If the line (AB) passes through P then Id ( 5 r. Therefore, to compute the expected value of rJ2 over the set of lines (AB) crossing at the point P, two cases should be considered, whether r is larger or smaller than R: If r is larger than R, E(d21r) = R2/4 (see (A2.1)) and using the assumption V(u) = a;Z (Eq. (2) in the text), the relation (A3.1) may be rewritten E(cr2) = (2~ r* + 1/2)(ri.
in (A3. l),
(A3.2)
REFERENCES 1. D. H. Ballard and C. M. Brown, Computer Vision, Prentice-Hall, Englewood Cliffs, NJ, 1982. 2. S. T. Barnard, Interpreting perspective images, Art(ficial Intelligence 21, 1983,435-462. 3. M. J. Magee and J. K. Aggarwal, Determining vanishing points from perspective images, Comput. Vision Graphics Image Process. 26, 1984, 256-267. 4. L. Quan and R. Mohr, Determining perspective structures using hierarchical Hough transform, Pattern Recognition Lert. 9, 1989, 279-286. 5. W. Dickson, personal communication, Oxford University, June 1989. 6. Kanatani, Reconstruction of consistent shape from inconsistent data, Int. J. Comput. Vision, 3, 1989. 7. M. Berthod, “Un nouvel algorithme d’approximation polygonale,” INRIA internal report, Sophia Antipolis, 1987. 8. J. Canny, A computational approach to edge detection, IEEE Trans. Pattern Anal, Mach. Intelligence 8, 1986, 679-698. 9. B. Brillault, Parallel and perpendicular grouping from a single view, in Proceedinas. Alvey Vision Conference, Reading, UK, 1989.
Vol.
IMAGE
54, No.
UNDERSTANDING
2, September,
pp. 289-300,
1991
New Method
for Vanishing
Point Detection
B. BRILLAULT-O'MAHONY Measurement
and Instrumentation
Centre, Department of Electronics and Information, London EC1 V OHB, United Kingdom Received
July 2, 1990; revised
The extraction of 3-dimensionalinformation from an imageis a central problem in computer vision. Location of the vanishing point of a line in the imageprovides its 3-dimensionaldirection. This paperis concernedwith the detectionof the vanishingpoints in an image.A new accumulatorspaceformulation for the whole line accumulationapproachis describedwhich fulfills an important property for robustness:the constancyof the detectability of a vanishing point whatever its location in the image plane. It is comparedwith the Gaussiansphereaccumulator space.Results are presented from interpreting indoor scenesof a nuclear plant. 0 1991 Academic Press, Inc.
INTRODUCTION
The major task in vision research is the extraction of 3D information from one or a number of different 2-D images. When the analysis of a single intensity image is considered, much information can be extracted using models of image formation and geometry. This paper examines an approach that assumes a perspective model, to enable groups of lines parallel in the 3-D scene to be extracted. For the method to be unseable on real tasks, a number of issues have to be addressed: robustness, noise, and uncertainty. The method must work over a wide range of images and produce perceptual groups that have a high certainty of being correct. In an indoor scene environment, directions of some features are more likely than others; e.g., for straight lines, the vertical direction and the horizontal directions parallel to the walls. The search for parallel lines in an image is very useful for the interpretation of the scene. The perspective projection of a set of parallel lines onto an image is a set of lines meeting at a common point, called a vanishing point. The vanishing point coordinates define the direction of the set of lines in the scene relative to the camera coordinate system. A number of methods for the detection of the vanishing points have been developed but for all of them the detectability of a vanishing point depends on its location in the image. The method proposed here is based on the accumulation principle but uses a new mapping that ensures the same detectability
December
City University, Northampton
Square,
26, 1990
of the vanishing point over the space. A further advantage is that the proposed accumulation space implementation is particularly simple. Usually, the vanishing points are detected by projecting lines, or points, each point being the intersection of a pair of straight lines, onto an accumulator space [l-4]. In order to reduce the search to a bounded closed set, Barnard [2] proposed projecting the lines of the image onto a Gaussian sphere centered onto the optic center and to use the Hough transform paradigm to detect the vanishing points. The projection of points onto a Gaussian sphere is equivalent to a resampling of the image plane. The method described here also resamples the image plane but the initial constraint is to keep the uncertainty of the vanishing point location almost constant over the image plane. In this paper the notion of uncertainty is defined, then the equation corresponding to the constraint is solved. The solution also leads to a bounded closed set which is merely a consequence of the method and not an a priori criterion. Once the lines have been accumulated in this new accumulator space, the local maxima are detected as the hypothetical vanishing points. The efficiency of the method is proved on a number of indoor scenes. Finally a comparison is made with methods based on the Gaussian sphere to highlight the main property of the new accumulator space. PREVIOUS WORK
The search for vanishing points consists of finding a small neighborhood in the image plane intersected by a sufficient number of straight lines. The image plane is resampled and bounded to form an accumulator space. The number of resampled straight lines crossing each cell is calculated using the Hough paradigm. Barnard [2] proposed projecting the lines onto a Gaussian sphere centered on the optic center. This projection is a bijection from the image plane to a bounded space. The search for the accumulation points is achieved by sampling the sphere using spherical coordinates. Unfortunately the spherical sampling (elevation, azimuth) is irregular and
289 1049-966019 1 $3 .OO Copyright 0 1991 by Academic Press, Inc. All rights of reproduction in any form reserved.
290
B. BRILLAULT-O’MAHONY
different in the x and y directions. Quan and Mohr [4] use the same method with a dichotomic approach: the sphere is sampled from a coarse resolution to a fine resolution in order to reduce the number of studied cells. The lines are classified by looking for vanishing points. Once a vanishing point is found with its associated lines, the lines are eliminated and the algorithm is performed again. Dickson [5] proposed a triangular sampling of the Gaussian sphere which is very attractive because it is isotropic; however, the computational efficiency has yet to be proved. Magee and Aggarwal [3] accumulate the projection of the intersection points of all pairs of straight lines in the image onto a Gaussian sphere. The accumulation is achieved using the arc distance between two points which leads to an isotropic search. The isotropy of the search for vanishing points using the accumulator spaces previously described can be ensured by an additional cost of complexity and only guarantees isotropy in 0, not necessarily in Y. As a result the probability of detecting a vanishing point depends on its distance from the center of the image. The reason for this dependence with Y is due to the increasing uncertainty of the vanishing point location with r. Kanatani [6] tests hypotheses on parallelism of three or more straight lines by using the projection onto the Gaussian sphere. The uncertainty of the line parameters is taken into account by using a threshold in a concurrency test which depends on the uncertainty of the intersection point coordinates. The author [9] has previously proposed to solve the problem by resampling the accumulator space such that the uncertainty remains approximately constant over the space. The variations of the uncertainty of the intersection point of a pair of lines is statistically estimated, but for Y fixed, the dispersion of this uncertainty is large; also, the solution found is not really isotropic in 13. It will be seen in this paper that a criterion based on distance between intersection points does not seem appropriate to a problem initially defined as the minimization of distances between a point and lines, whatever the representation used. The method described in this paper proposes a solution to the problem of the consistency of the vanishing point detection over the image plane, without increasing the algorithmic complexity. PREPROCESSING
The grey scale image of an indoor scene is digitized using a standard CCD camera. Straight lines are detected in the image using Berthod’s line finding algorithm [71 after the edges are detected using the Canny operator [8]. Each line approximates an edge in the image and is represented in a list by end point Cartesian coordinates.
ACCUMULATOR
SPACE
Notations In the following analysis vectors and matrices are represented in bold. The transposed matrix of M is noted M’. The vector from A to B is denoted AB and AB is its algebraic length. The line supporting the segment AB is denoted (All). The expected value of a real value x is denoted E(x). E(x 1y) is the expected value of x when y is known. The variance of a real value x is noted Vx, and the covariance of the real values x and y is denoted Vxy: Vx = E(x2) - E2(x) Vxy = E(xy) - E(x)E( y). The symbol of differentiation
Uncertainty
is noted 6.
Definition
Under perspective transformation, parallel lines in a 3-D scene are projected onto concurrent lines in the 2-D image. Ideally, in a man-made environment, many lines are parallel, e.g. the edges of a wall and a door frame. In practice this is only approximately true. However, parallelism is a useful concept for the representation of the door frame location. Therefore the first source of uncertainty is the difference between the real scene and the symbolic representation describing this scene, e.g., CAD representations. The second source of uncertainty comes from the image acquisition, digitization and preprocessing. The acquisition of the image is modeled by a pinhole perspective model. This does not take into account the numerous possible distortions of the image caused by the lens. The digitization stage constrains the total error on a pixel location to be at least 0.5 pixels. The Canny edge detector has an accuracy limited by the value of the parameter sigma of the Gaussian filter [S], The line finder may introduce nonnegligible errors in the slope of a segment. To clarify, let us assume that the camera calibration parameters are known (though unnecessary for the vanishing point detection) and the ideal representation of the scene is also known, e.g., in the form of a CAD description. For a segment the error of the end point location is defined as the difference between the observed end point location and the projection of the corresponding ideal end point in the scene using a pinhole model (for which the parameters are known). As a result of the uncertainty, lines assumed parallel in the 3-D scene have their images only approximately concurrent in the image (see Fig. 1).
VANISHING \I P
Lp---- L5 L2 i”\\\ L3
L4
FIG. 1. Image of lines assumed parallel in the scene. They are only approximately concurrent in the image.
Uncertainty
Modelling
The sources of the uncertainty are various and it is difficult and inefficient to try to model them separately. The simplest model commonly used for this type of uncertainty is the Gaussian law, with an expected value equal to zero and a standard deviation equal to cr. Therefore the error on each end point location is assumed to obey the same Gaussian law. Only the component of the error perpendicular to the line is needed as it is responsible for the errors of the line parameters (polar angle and distance from the origin). Let the matrix V = [zY $‘I be the covariance matrix of an end point A. A point M of the uncertainty ellipsoid is defined by AMfV-‘AM
I 1.
In the following u is the unit vector of the line L corresponding to the segment AB and V(u) is the variance of the error in the perpendicular direction to L (Fig. 2), V(u) = (ufV-‘u)*det(V).
(1)
In what follows V is assumed to be constant over the image plane, V(u) = CT& where I is the identity Problem
291
POINT DETECTION
(2)
parallel lines onto the image. The location of P in the image depends on the viewpoint. Searching for P by accumulation is equivalent to counting the real lines in the image passing through a neighborhood of P. A line (L) defined by a line segment in the image is assumed to have its vanishing point in P if D(L, P) < r(P), where D(L, P) is the distance between the line L and the point P and r(P) is the radius of the chosen neighborhood N of P in the perpendicular direction of (L). The constant quality for the detection of P irrespective of its location in the image plane increases the robustness and meaning of the detection. The value of the radius r(P) of the neighborhood N defined above should be proportional to the accuracy with which the line L is knotin in the vicinity of P, that is to say proportional to the uncertainty of the distance D(L, P) from L to P. The neighborhood N has a corresponding neighborhood N’ in the accumulator space. For practical reasons this neighborhood N’ should be constant and have a simple shape. Let (T be the uncertainty of D(L, P). The problem is to jind a transformation Tfrom the image plane to an accumulator space such that the expected value of u’ = T(o) remains constant over the accumulator space. This transformation Twill completely define the new accumulator space (Fig. 3). Statistical
Expression for the Uncertainty
The uncertainty of the end points of the segment AB results in uncertainty u on its distance to a point P. Polar coordinates are used, with pole 0 and axis Ox. Let 0 be the center of the image, r the distance from 0 to P, 6 the polar angle of OP, (Y the polar angle of (AB), and d its distance from 0. Q’ is the intersection of the line (AB) and the circle (C) with center 0 passing through P and CX~ is the polar angle of Q’ (in the ideal case (Y~= 8, see Fig. 4). Let ap be the uncertainty of P along (C). The variance of V~ over all the lines having P for vanishing point is equal to the variance of the error of P in the tangential direction, (r aa!,)*. Let the line (AB) be fixed; the uncertainty cp is equal to
matrix. ai = r*E(&:),
Dejinition
(3)
Let P be the vanishing point to be determined, which is the common point of the ideal projection of the ideal
FIG. 2. Geometric interpretation
of V(u)
FIG. 3. Transformation space.
from the image plane to the accumulator
292
B. BRILLAULT-O’MAHONY
image plane
i/L2
0.002s.
:
.
i
:
:
:
:
0.002
:
i
;
;
’ i ;:
i
0.001s.
j::;::
; :
.!
0.001 0.000s,
FIG. 4.
Line (AB) passing near the point P. 16
d
FIG. 5. Average value of l/L* in function of d. Each dot corresponds to the average value of l/L* over the lines of an image at a distance d from the origin. Twelve images of 6 different scenes have been processed.
where E(Gcu:) is the variance of the noise on (YI. Using the polar equation of the line it follows that r cos(a, - a) = d. This equation is derived to provide the expression 6(arJ, when (r2 - d2)“2 is not too small:
126
64
of
From a scene point of view, the length of a segment of the image depends on its vanishing point location. The distance between a segment and its vanishing point is constrained by
8a1 = 6a - 6d/(r sin(al - ar)) = &Y - 6dl(r* - d2)‘j2.
G’Q’IL
> D,/Df,
(6)
(4) Let L be the length of the segment AB and b the distance between the centroid G of the segment and the projection 0’ of 0 onto (AB), when (r* - d2)1’2 is large enough, the value of cri is given by (see Appendix 1, (A1.4), (AlS), 641 m
l/2 + 2b2/L2 1 - d2/r2 4br
+ L2(1 - d21r2)‘”
VW.
To make the expected value of oi explicit, a statistical model of the a priori repartition of the lines meeting at P within the image is required. In what follows, the image is approximated by a circular disk, in order to simplify the model. The lines are assumed to be distributed at random within the image; i.e., the centroid of the segments AB is uniformly distributed on the image. The length of these segments obeys a probability law of density f(L), which is assumed to be uncorrelated with the centroid location. From an image point of view, these assumptions are not exact because of the segments crossing the boundaries of the image. Therefore, for simplicity, the segments are not constrained to lie entirely within the image, except their centroids, so that the assumptions are consistent. It will be seen that ,?(a*) depends on E(1/L2). A statistical study of l/L* on 12 real images demonstrates the validity of the assumption of the dependence of l/L* from the centroid location (Fig. 5).
where Df is the depth of field and D, the focusing distance. Figure 6 demonstrates that this constraint has a very small effect on the expected value of (l/L2). This is not surprising because the constraint filters the small segments which play a preponderent part in E(1/L2). In fact the distribution of the lengths is mainly due to the line finder algorithm and little to the viewpoint. Let a be E( 1/L2). Using the independence assumptions between l/L2 and (b, d, v) and the fact that E(b 1 d) = 0 for symmetry reasons, the expected value of C; over the set of segments AB at a distance d from the origin is equal to (from (5)) E(a;(d)
= (2ar2 + (l/(1 - d2/r2))/2 + 2aE(b2)d)l(l
- d2/r2))V(u).
(7)
l/L2 I 0.002s. 0.002
1 ,.
0,001s 0.001
.i :*
. : j
. : ! * , : .fi
I
; j
: I.
: 1
: $
i
:i:iii
,,
0.000s. ‘
: 16
64
I12
160
r
FIG. 6. Same graph as in Fig. 5, over the lines filtered by the constraint (6) with D,lDf = 1.3.
VANISHING
POINT
This relation ignores the constraint (6), but still E(b2/ L2/d) gives more weight to large values of b and small values of L corresponding to the lines satisfying (6), which justifies (7). Remark.
minimal
The previous argument only holds when the value allowed for L is small enough.
The centroid G is located at random on a chord of the image disk, the length of which is equal to 2R(l - d2/ R2)‘12; then (Appendix 2): E(b21d))
= R2(1 - d2/R3)/3.
= (2ar2 + (l/(1
The value of the uncertainty of the line location near P depends on the parameters of the corresponding segment in the image. It has been shown that it is possible to express the expected value of this uncertainty by using statistics of the line segment parameters. The expression for the expected values depends only on the distance of the vanishing point from the image center. This result allows almost the same detectability to be guaranteed for any vanishing point, by resampling the image plane proportionally to the expected value of the uncertainty. Transformation Accumulator
Thus E(aj$d)
from the Image Plane to the Space
In
order to simplify the notation, Cr represents in the following. If the transformation T exists it is defined by
- d2/r2))/2
(E((T~))“~
+ 2aR2/3)V(u).
The distance d of the line is smaller than or equal to r, within the uncertainty of r. When d is smaller than r l/(1 - d21r2) = 1 f &(llr2) where c(llr2) is infinitesimal sumption in (2) gives
293
DETECTION
with Ilr2.
Using
the as-
P(r,
P’(x’(r,
19), y’(r,
19)),
x’, y’ being the coordinates of the cell of the accumulator space. Let Q’ the point of the circle (C) located at the distance @,,from P; then T(Q’)
A line having P for its vanishing point satisfies d I r. Equation (8) has been found using (4) under the assumption that (r2 - d2)“2 is not too small. If r > R, d is always smaller than r, and (8) holds and is independent of d. When r I R and d = r, the definition of aP has no sense, as the line (AB) may cross once, twice, or not at all the circle (C) in the vicinity of P. When P is far enough from the origin, a very small proportion of lines are likely to be such that d = r, and (8) is nearly always valid. If r is very small, the radial uncertainty of P has no sense. The uncertainty neighborhood of P is a circle with a constant radius c2, where u2 is defined as the expected value of PQ2. It is shown in appendix 3 that near the origin
f3) 7‘
= T(P(r,
8 + ?F,lr))
= P’(x’(r,
8 + F,lr),
y’(r,
0 + CJr)).
In the accumulator space the uncertainty e/:, corresponding to the uncertainty ZP in the image, is equal to the distance of T(P) to T(Q’). Thus -I2 = (P’(x’(r,
0 + Fp/r), y’(r,
@P
- P’(x’(r,
01, y’(r,
8 + Cp/r))
8>)12
= (x’(r, 8 + ZJr) - x’(r, fl))2 + (y’(r, 13+ @Jr> - y’(r, 13))~.
After linearization
it leads to
-I2 = (6x’/68)2Z~lr2 UP
+ ((6y’/6f3)2Z~lr2.
Thus ((6~‘/68)~
Therefore from (8) and (9) and the previous remarks, uncertainty of P may be approximated by
the
+ (Sy’/619)~) = FL2r2/(ri.
One solution of (11) is y’ = ZLr0lcT
if r is large enough
E(ai)
if r is small
E(u2) = cu;,
where a = E(1/L2)
= (2ar2 + c)oi
and c = l/2 + 2a R2/3.
(10)
(11)
P
x’ = x’(r).
x’ is independent of 8 and may be any bijection from [0, +m] to a bounded interval. For simplicity, x’ has been
B. BRILLAULT-O’MAHONY
Geometrical Interpretation
chosen so that y’
= 8x’lk,
where k is a scale .factor determined resolution. Therefore
by the expected
x’ = (Cf~/aO)krl(2a r2 + c)“~ y’ = (i7-~/cro)rel(2a r2 + c)l12.
(12)
The number of straight lines passing near P(rO, &J is represented by the number of curves crossing the line parallel to y’ axis, y’ = y’(rO, e,), in the range [P’ - iT;j, P’ + Gij] (see Fig. 7). The counting of these lines is performed by accumulating curves dilated by a vertical kernel with a half-width equal to $, . The accumulation of dilated curves prevents the same line from being counted twice. A line passing near P such that d E ]r, r + a[ does not cross the line y’ = y’(rO, &). When P is far from the origin such a configuration is unlikely, but when r is small its probability increases. The problem is solved if the x’ resolution around P’ is less or equal to G. This could be done by choosing the constant k in Eq. (12) such that x’(r) has a resolution equal to (T for small r, but the resolution of x’ has to be consistent with the resolution of y’ (e.g., two curves crossing at x’ = (2n + 1)/2 must cross the same uncertainty neighbhorhood either at x’ = rz or at x’ = n + 1). The two constraints on k are incompatible for large r. The constant k is chosen to ensure the consistency with y’ resolution and the constraint for the small value of r is fulfilled by a pre-accumulation stage in the (d, (Y)parameter space (see later). A transformation T has been found such that any point in the image plane is simply transformed to a single point in the accumulator space and the search for vanishing points is reduced to the search for local maxima in the accumulator space. The sampling is isotropic and homogeneous with respect to the uncertainty criterion. Moreover the accumulator space is bounded.
e=2n
p’ = x’(p)
u’
8= (L) direction
to
FIG. 7. (x’(r),
y’(r,
x
0)) accumulation space.
and
8’ = Blk.
In this sector the neighborhood of uncertainty is a circular arc with arc length FL. Thus, the transformation defined is equivalent to a resampling of p into s(p) followed by a projection onto a cone with vertex 0, axis equal to the optical axis, and angle /3, which is eventually developed onto a plane, so that T=d*p.s
Image plane s P(P, 0) -
Image plane p w-4, 0) -
>
Accumulator Cone d space m(p), 0, p) p’h et)
M is the intersection between a cylinder with radius s(p) and the cone with the same axis and the angle p; M has the Cartesian coordinates (s(p) cos 8, s(p) sin 8, p’ cos p), where p’ is its polar distance. Then
s(p) = p’ sin p
and
8’ = Bs(p)lp’
= 8 sin p.
This gives sin /3 = l/k and s(p) = x’(p)lk. Pre-accumulation
Stage
The lines are first accumulated with respect to d and CK, the resolution of which is determined by E(6d2)‘12 and E(&x~)“~ averaged over all possible b (see Appendix 1, (Al .4), (Al .5)), and for each cell (d, cu)the best estimates of the line parameters are computed. Only one line per cell (d, a) is accumulated in the accumulator space w, Y’). When r is small, d is also small and the expected value of E(6d2) over all b is equal to c, i.e., to (T2. This means that all the lines such that d E [r, r + G[ and (Y’ = 8, are accumulated in the same cell (d, CX). Therefore, for r small, the resolution in x’ is equal to (+ around x’(r) and the lines satisfying d = r cross the uncertainty neighborhood around P’ defined above. The preaccumulation stage avoids a weighting process and thereby increases the significance of the peaks of the accumulator space (x’, y’) described above, as the peak value represents the number N of directions meeting at the same point. For instance, if N = 3 the point may be a corner, if N > 3 the point is likely to be a vanishing point.
liiz&cl A0=n
C = T(L)
e0x,
The accumulator space (x’, y’) described previously is isotropic in 6’ by definition. In the continuous plane, it may be represented as a sector of a disk. In this disk, a point P’ with the polar coordinates (p’, 6’) corresponds to the point P(p, 8) in the image plane and as defined by using (12)
VANISHING
Uncertainty
of the Vanishing
POINT DETECTION
295
Point
It may be important to have a first evaluation of the uncertainty of the vanishing point location, once detected. Let (T; be the variance of the distance Q’P, where Q” is the intersection point of the line (AB) and the line (OP) (see Fig. 4); it is equal to (+z = E(Q’p2)
= r2E(QP2)/d2
= r2c21d2.
(13)
The line such that d = 0 does not provide any information about r and the corresponding uncertainty is infinite; so is the expected value of CT: over all possible d for r large enough. But here the peak in the accumulator space has already been found, which means that the shape of the distribution of Q”P2 around P no longer matters. The information about r is provided by lines passing far from the origin. At least three lines with different directions should cross the neighborhood of P to give it meaning as a possible vanishing point, and therefore at least two lines passing far from the origin, typically at [E(d2)]1’2. Therefore, for large r, using the expression of E(d2) found in Appendix 2,
FIG. 8. Comparison between a “whole line” accumulation proach and an “intersection point” accumulation approach.
ap-
lines classified with it are taken away from the list and the algorithm is performed on the rest of the lines. However, a number of lines in the image may pass through several vanishing points and such a method may ignore important directions in the image. The method used here removes the lines from the list after each candidate vanishing point has been detected, as described above, but only if the number of these lines is high enough to be responsible for several local maxima in the accumulator space. The detected vanishing points correspond to the highest peaks in the accumulator. Once all the candidate vanishing points have been detected, the 1 r2g2 lines are reclassified by a maximum likelihood method in a2 2r2(2ar2 + c)a$R2 cr=--= 2 E(d2) (14) order to take into account the specific uncertainty of each line location with reference to each vanishing point (in where c$~is an estimate of the radial uncertainty on P. contrast with the statistical method described above). If r is small, using the same considerations as for Z,, , c+? This allows the classification of one line with several vancan be replaced by the expression of a2, say c(ri. The ishing points and solution of the ambiguities by appropriexpression (14) for 6, gives a very coarse idea of the ate methods not described here. An additional filter for vanishing point uncertainty and it should be considered discriminating the vanishing point candidates selects the as temporary information. This approximation is refined pairs or triplets of perpendicular directions, which are by using a Kalman filter applied to each line classified hypothesized to be the main directions [9]. The descripwith the point P. tion of such methods is beyond the scope of this paper. For the resolution to be consistent with this uncertainty value, the derivative of r relative to x’ must be RESULTS lower than &,.. The sampling is not regular according to The accumulator space depends on the statistic of l/L2 &,., which does not matter, as c+~does not affect the linewhich depends on the complexity of the image and on the counting in the whole line accumulation approach depreprocessing stages. The size of the image used here is scribed above. This would not have been the case if only 256 and only lines with a length greater than 15 have been the intersection points had been accumulated (see Fig. 8). accumulated (lines with a length below 15 are very nuRemark. It is possible to choose x’(r) in order to have merous and do not provide meaningful information as a sampling approximately regular according to &-, but it their uncertainty is very high and they often correspond slightly increases the complexity. to noise). Thus the pre-accumulation in (d, a) results in a narrow distribution of l/L2 around a = 0.0019, deterVanishing Point Detection mined experimentally over a number of images. The problem posed by local maxima detection in the The method has been applied to various types of inaccumulator space is that a set of close parallel lines, e.g. door scene. Figure 9 shows the initial images. Their rea set of pipes, usually produces several local maxima spective accumulator spaces are displayed in Fig. 10, corresponding to the intersection of this set with various with all the lines accumulated, Figure lOa, and with a lines. number of lines removed, Fig. lob. The final result of the Some of the other methods [3, 41 use the following classification is displayed in Fig. 11. All the main direcstrategy: once a vanishing point has been found, all the tions represented by a sufficient number of segments
B. BRILLAULT-O’MAHONY
296
FIG. 9a. Initial images 1 and 2.
FIG. 9b.
Initial images 1 and 2 and the extracted segments.
have been found. False candidates have been found, corresponding generally to corners. They are filtered out by using the constraint (6) which is usually not satisfied by the corners. It can be noticed that vanishing points have been found corresponding to a small value of Y (18), to an average value of r (135) and to an infinite value of r, with peaks of similar shape along y’ (see Fig. 10).
COMPARISON OTHER METHODS
Comparison with Intersection Accumulation Method
Points
Figure 8 demonstrates that all the lines around Q are counted independently of the value of urr but the count-
VANISHING
POINT
DETECTION
297
FIG. 10. (a) Accumulator spaces with all the lines accumulated. For clarity the image has been scaled from 0 to 255. The axis y = 0 corresponds to the vertical lines. For convenience, the 0 reference has been shifted by a/Z. (b) Image 1 vertical lines parallel to the camera have been removed; image 2 vertical and horizontal lines parallel to the camera have been removed. For clarity the images have been scaled from 0 to 255.
ing of the corresponding intersection points should take or into account. Equation (13) shows that for Y fixed, the disparity of or is very large, and thereby the disparity of the error of the distance between intersection points is very large too. Whatever the representation used, the disparity of the error of the intersection points is very much larger than the disparity of u2, which means that the concurrency test based on the intersection point accumulation approach has a poorer significance than the test used in the whole line accumulation approach. The best estimate of the intersection point of a set of straight lines in the least mean square sense minimizes C(di/af)2. This criterion may be written in terms of intersection point coordinates and the associated covariance matrix (it is complicated because the intersection points are not independent). Because of the disparity of the vari-
ance along the Y direction, an accumulator space based on the intersection point accumulation should necessarily be associated with a weighting process using coefficients equal to the inverse of this variance, so that only significant intersection points would be taken into account. If the dependence of the intersection points is ignored, the value of the peaks obtained would approximately represent the inverse of the corresponding intersection point variance along r. This value is not significant because of the disparity of this value and because it provides no information about the number of lines or directions meeting at this point. Albeit the set of the intersection points and their associated uncertainty matrix theoretically contains sufficient information about the best estimate of the intersection point of a set of straight lines, the corresponding tests
B. BRILLAULT-O’MAHONY
-
--
FIG. 11. Final classification: the top and the bottom left images corresond to the main directions detected; the bottom right images show all the lines classified with a main direction.
have a poor significance. Besides, the information is redundant (5n(n + 1)/2 data in contrast with 5n data in the line accumulation case); i.e., the data are not independent, which is a serious difficulty when reasoning with uncertainty.
cess transforms L into a big circle C, P into P’, and Q into Q’. P is assumed to be near Q, so that the distance D’ from C to P’ can be approximated as D’ = D(P’, Q’), and its uncertainty (T’ approximated by the projection of the uncertainty U. Let h be the distance between 0 and the center of the Gaussian sphere; then
Comparison with Gaussian Sphere Method
If the constraint of the constant detectability is not fulfilled the detection of a vanishing point, i.e., a main direction in the scene, depends on the viewpoint and the orientation of the image. It is important to evaluate the behavior of the more commonly used accumulator space, the Gaussian sphere, with respect to this constraint. (See Fig. 12). Let r be the distance from Q to 0. The projection pro-
d/u
= O’Q’/O’Q
= p/(r2 + h2)l/*,
where p is the radius of the sphere. The constraint U’ = constant implies that (9 + h2) is proportional to (~2. It is shown in Appendix 3 that it is possible to approximate u2 for any ’ by iF2 = (2a r2 + c)cri
which leads, when applied to the gaussian sphere, to O’O=h O”Q=r
h = (cl(2a))1;2. P
FIG. 12. Projection of an image segment onto the Gaussian sphere.
(15)
This relation means that there exists a value for h such that the uncertainty is the same at the origin and at infinity. A large value of h favors the detection of vanishing points located near the center of the image, whereas a small value of h favors the detection of vanishing points located at infinity. Here, the neighborhood used for the counting of the lines is a circular spherical cap with a constant radius a’.
VANISHING
POINT
The Gaussian sphere method appears to be a correct accumulator sphere for the constraint considered here in the context of the whole line accumulation, provided a correct value has been chosen for h and an isotropic sampling is used. This result does not hold when only intersection points are accumulated. But an isotropic sampling of the sphere in the whole line accumulation approach implies an involved process, whereas this problem does not exist for the accumulator space described previously.
299
DETECTION
6d = (b/L + l/2) 6A . u1 - (b/L
- l/2) 6B * ui.
(Al .3) The variance of (SA . ui) and (6B * ur) is V(u) (see (1) in the text), therefore the second moments of &Y and 6d are E(W)
= 2V(u)/L2
(Al .4)
E(6d2)
= (2b2/L2 + 1/2)V(u)
(Al S)
E(6a 6d) = -2b
(Al .6)
V(u)IL2.
Let D be the distance from a fixed point P to the line (Al?); D is equal to
CONCLUSION
This paper has considered the problem of extracting 3D information from a perspective model for the projection of parallel lines onto an image. A new accumulator space for the detection of the vanishing points has been presented. This accumulator space is isotropic and bounded; moreover it ensures a constant detection quality of the vanishing points over the space. This results in a more robust and predictable method that works on a large number of different images of indoor scenes. In addition it is easily implemented. It has been compared with different methods using the Gaussian sphere as an accumulator space. This new method has been used to determine successfully the main directions of indoor scenes. Further work is concerned with the characterization of the signal to noise ratio of the accumulator space, which will be used for scoring each class of lines associated with a vanishing point. APPENDIX 1: VARIANCE OF THE PARAMETERS OF THE LINE (AB) AND THE DISTANCE FROM A POINT P TO (AB)
Let (Ybe the polar angle of the line (AB), d its distance to the origin, u its unit vector, and uI the perpendicular unit vector. Let G be the centroid of the segment All, L its length, 0’ the projection of the origin 0 onto (Ail), and b the algebraic distance between 0’ and G (Fig. 4). The line (AB) is known within the errors of measurement 6A and SB of the end points A and B. These errors results in errors on (Y, d, II, and ul: 6a = (SB * u1 - 6A . u,)lL
(Al.l)
6u = 6au,;
(A1.2)
And 6Ul = -6cY u.
D = AP * ~1. The derivation
of (A1.7) provides 6D, 6D=
-6A*u,
+AP.Sul.
Let Q be the projection of P onto (Ail); using (Al .2) and (A1.3) and then substituting AP . u by AQ and (AQ - L) by BQ, it is found that 6D = @Q/L) E(6D2)
6A . u1 - (AQlL)
6B . uI
Let Y’ be the distance rewritten
from 0’ to Q; (A1.8b) may be
E(6D2) = (2(r’ + b)2/L2 + 1/2)V(u). O”Q’ = r2 - d2 and O’Q = Y’. The vector QQ’ only upon the error on the end points, therefore stitution of Q by Q’ in (Al.ga) only adds terms degree than 1 in 6A * u1 and 6B * ul. Thus the approximation can be made: r’
APPENDIX
(A1.8a) (Al .8b)
= (AQ2 + BQ2)V(u)lL2.
=
(r
2
2: EXPECTED
-
(Al .9) depends the subof higher following
d2)“2.
VALUES OF dZ AND b2
Using the statistical model described in the text, the line (All) may be generated by any centroid G located on the chord (S) of the image disk, (S) being at the distance d from the center of the image, G is assumed to have a uniform density of probability on (S), therefore the density of probability of d is proportional to the length of the string (S): f(d)
Nowd=OA~ul,then6d=6A~ul+OA~6ul=6A~ul - OA * u 6a. Replacing OA . u by (b - L/2) and S(Yby (Al. l), it is found that
(Al .7)
= 2KR(l
- d2/R2)‘12.
The constant K is found by integration (-R
I d I R):
of f(d)
over d
300
B.BRILLAULT-O'MAHONY
K = ll(7fR2).
If r smaller than R then (integration
The expected value of d2 is equal to
by parts)
E(d2) = R2/2(1 - r(1 - r2/R2)3’2/ (R Arcsin (r/R) + (R2 - r2)l12)).
E(d2) = (2/(7rR)) j-R d2(1 - d2/R2)“2 6d = R2i4. When r is small (A2.1) E(d2) = r2/3, The expected value of b2, when d is known, is equal to and substituting E(b*ld)
= (1/(2(R2 - d2P2) \($:;:;:,2
b2 ab (A2 2)
= (R2 - d2)13 APPENDIX
3:
CALCULATION
OF uz
Let u2 be the expected value of 6D2, square of the distance from a given line (AB) and a fixed point P. The variations of 6D are due to the errors of measurement of the end points. Developing (A1.9) by substituting r’ by (r2 - d2)‘12, this gives u2 -_
2r2 I 4br(l - d21r2) I 2b2 i L2 L2 L2 + ; - g
1
E(02) = (IOU r2/9 + 2a R2/3 + 1/2)oi.
(A3.3)
It is more convenient to use an approximation valid for any r. Moreover the comparison with the Gaussian sphere method requires that such an approximation can be made. Using (A3.2) and (A3.3), E(cr2) is approximated for any r by E(a*)
= (2~ r2 + c)cri,
(A3.4)
where c = 2u R2/3 + l/2. ACKNOWLEDGMENTS This research was funded by Electricitt de France and was performed in the Machine Vision Group at the City University, London.
V(u)
If d is fixed and a2 is averaged over all the segments All producing the line (Al?), then by using (A2.2) and the fact that E(b) = 0 for symmetry reasons, it is found that E(u21d) = (2ar 2 - 8ad2/3 + 2aR2/3 + 1/2)V(u), (A3.1) where a = E(1/L2) (for the independence between a and the other parameters, see the paragraph “uncertainty modeling” and Fig. 11). If the line (AB) passes through P then Id ( 5 r. Therefore, to compute the expected value of rJ2 over the set of lines (AB) crossing at the point P, two cases should be considered, whether r is larger or smaller than R: If r is larger than R, E(d21r) = R2/4 (see (A2.1)) and using the assumption V(u) = a;Z (Eq. (2) in the text), the relation (A3.1) may be rewritten E(cr2) = (2~ r* + 1/2)(ri.
in (A3. l),
(A3.2)
REFERENCES 1. D. H. Ballard and C. M. Brown, Computer Vision, Prentice-Hall, Englewood Cliffs, NJ, 1982. 2. S. T. Barnard, Interpreting perspective images, Art(ficial Intelligence 21, 1983,435-462. 3. M. J. Magee and J. K. Aggarwal, Determining vanishing points from perspective images, Comput. Vision Graphics Image Process. 26, 1984, 256-267. 4. L. Quan and R. Mohr, Determining perspective structures using hierarchical Hough transform, Pattern Recognition Lert. 9, 1989, 279-286. 5. W. Dickson, personal communication, Oxford University, June 1989. 6. Kanatani, Reconstruction of consistent shape from inconsistent data, Int. J. Comput. Vision, 3, 1989. 7. M. Berthod, “Un nouvel algorithme d’approximation polygonale,” INRIA internal report, Sophia Antipolis, 1987. 8. J. Canny, A computational approach to edge detection, IEEE Trans. Pattern Anal, Mach. Intelligence 8, 1986, 679-698. 9. B. Brillault, Parallel and perpendicular grouping from a single view, in Proceedinas. Alvey Vision Conference, Reading, UK, 1989.