- Email: [email protected]
An SIMD algorithm for range image segmentation
Pattern Recoonition, Vol. 28, No. 2, pp. 255 267, 1995 Elsevier Science Ltd Copyright © 1995 Pattern Recognition Society Printed in Great Britain. All rights reserved 0031 3203/95 $9.50 +.00
Pergamon
0031-3203(94)00093-X
AN S I M D A L G O R I T H M F O R R A N G E IMAGE SEGMENTATION PRABIR KR. BISWAS, t S. S. BISWAS and B. N. C H A T T E R J I Department of Electronics and Electrical Communication Engineering, Indian Institute of Technology, Kharagpur - 721 302, India (Received 10 March 1993; in revised form 26 July 1994; received for publication 5 August 1994)
Abstract--Understanding 3D scenes from range images need the segmentation of 3D surfaces into approximate planar surface patches from which curved surfaces can be constructed quickly for high level vision purpose [Biswas et al., Qualitative description of three-dimensional scenes, Intl. J. Pattern Recognition and Artificial Intelligence 6(4), 651-672 (1992)]. In this paper we have presented a new parallel algorithm for 3D surface segmentation wherein the problem of surface segmentation is modelled as a quantization problem. The surface normals are quantized to some predefined directions, or stated otherwise, the surface regions are approximated by planar surface patches which are parallel to some predefined planes. The novelty of the algorithm lies in the fact that, though surface segmentation is achieved by quantization of surface normals, the algorithm does not compute the surface normals explicitly. Rather the quantization is achieved using simple operations of shift, subtract and threshold. This technique suitably avoids the computations of the differential properties of the surfaces or the surface fitting expressions which are used in most of the other existing techniques. Hence this approach is computationally attractive. This algorithm is easily implementable on an SIMD array computer. Another advantage of this technique is that it is robust to noise present in the image. The algorithm has been explained with a number of examples. Experimental results with synthetic as well as real range images are cited in this paper to highlight distinctive features of the algorithm. Segmentation
Range image
Parallel algorithms
I. INTRODUCTION Processing and recognition of three-dimensional objects from range images has gained importance in recent years. Various approaches towards this problem have been reported in literature. A range image is a twodimensional array of numbers, where each array element denotes the distance between the image plane and the object surface at that point. Since an image array contains a huge a m o u n t of sensed values, it is often difficult for a computer to interpret a scene directly from raw image data. Thus a preliminary processing of the raw sensor data is needed which partitions (segments) an image array into groups of pixels with particular properties that can be conveniently used by the later stages of representation, description and recognition of the scene. Thus every segment is a collection of pixels which have some uniform property quite distinct from the property of pixels belonging to other segments. It is often difficult to give a strict mathematical definition to this property because the requirement of segmentation is often guided by the high-level description/recognition stage. Most of the attempt# t) for range image segmentation either fit the surface points to a surface function or calculate a local surface property such as normal, curvature etc. These tAuthor to whom all correspondence should be addressed.
SIMD array architecture
algorithms are often computation intensive and not amenable to easy parallelization. In this paper we present a new parallel segmentation algorithm which segments range images into approximate planar surface patches. Segmenting range images into approximate planar regions provides a basis from which curved surfaces can be constructed quickly for high-level vision purposes. ~2~ A n u m b e r of approaches for range image processing have been reported. Among those using surface based properties some extract only planar faces. Milgrim and Bjorklund t3) have suggested a method in which planar surface is extracted using a 5 × 5 window and surfaces are fitted using least square method. Henderson 14) organized the points into a three-dimensional binary tree from which he obtained a three-dimensional spatial proximity graph. Then a sequential region growing algorithm is used to create convex planar faces. A method known as 3-point seed method was suggested in reference (5) where region growing started from three neighboring points. Planar, cylinderical, conical and spherical primitives are handled in references (6) and (7). Spheres and ellipsoids were dealt with by Sethi and Jayaramamurthy Im Besl and Jain tg) developed an algorithm based on variable order surface fitting. In addition to the above there are segmentation techniques using local properties of a pixel. Besl and Jain t 1o) proposed a method where each pixel is labelled 255
256
P.K. BISWAS et al. Plane Set (NPS) values indicative of the surface normal at every range pixel. We have modified an extended the concept of digital neighborhood plane (DNP) and neighborhood plane set (NPS), introduced in reference (13) for segmentating surfaces in range images. In this paper we have presented a new parallel algorithm for segmenting surfaces in range images. The problem of surface segmentation in range image is modelled as a quantization problem, where the normal to the surface at a point is quantized
by the sign of mean and Gaussian curvatures at that point. Brady et al. ill) discussed a method of segmentation where features such as curvatures, asymptotes, boundary contours, certain surface patches etc. are calculated using principal curvatures, principal directions and a surface primal sketch. Taylor et al. (12) used a split and merge algorithm with a homogeneity criteria based on a 3-parameter surface description technique. Mukherjee et al. (13) have proposed a technique where a range image is segmented based on Neighborhood
1
2
-i
- - t"~,
II
li
71
,,,
; ,'-l--J "
= 0 ,180 °
~
O = 90 °
:
0
=
'I,
!
o
0
--
m
--..i,"
0 ° ~ ~6 ~< 3 6 0 °
90 °
O= 0 °
,'
,
,/''~
I
, ',Xi
,/I
,IL--+I,J 7 7,"
',
I, ~ ,, , / t - - - - - , :-i ,.
I
I
-)--)
_1,"
I
:
/+5°, 2 2 5 °
¢, :
135°315 °
@ = 90 °
:
90 °
0
90 °
0
:
7
J
,
+,-t-'l l
N / N/ I
I
1 "i---I-~
'U_ "_ __ d ; : @ =
"_~_ _ _ u,"
2700
~
:
0°
£~5 °
6
=
&5 °
10
: 0 =
: ,+-y':-r, ,> I
,
~__~_ _-,
"_>'C__?',"
180 ° ~5°
12
11
4
I
~5 °
,'~,~
J
'~
=
I
9
e
r0
i
Z
6
(-- -I..: - ,7/( i
li
~ L,,-----!,-I~
I,
II~ -0
L_
E:
II
$
,--..(k: -~- .-.~,-f,
,,_J < z
,,
90 , 2 7 0
4
cl
III
kC.__
0
cr~ L£J Z <
+I '- ---"I
i
0
Y
x,'F-- - :--'i
,,,
-+-r-r,/l '
i
Id_ _ k'_ _I,"
,_I
3
,,4"
I
~,\l
13
,4"~-
#~
i
i
, / - . ~ "u"+
,~, = 135°
~ : 31 5 ° '
~. : 2 2 5 °
•
0 = $/.,.73 =
0 = 5/4.73 °
0 :
EXTENDED
I
IZ
O : 5&.73 °
(b)
-- --7l
=
45 ° 5&.73 °
PLANES
Fig. 1. Orientations of the normals of digital neighborhood planes.
An SIMD algorithm for range image segmentation to the normal of a given plane based on a proximity measure. The angular difference between the surface normal at a given point and the normal to the given plane is taken as a measure of proximity between them. The new algorithm has got the advantage of being simple in design, efficient in execution and robust to noise. The complexity of the algorithm is 0(1) on a mesh connected computer. A preliminary version of this work has been reported in reference (14). The following two notations are used extensively in this paper.
Notation 1. A two-dimensional digital space Z 2 is defined as Z 2 = { (Xl, x 2 ) l X l , x 2
are integers}.
Notation 2. A three-dimensional digital space Z 3 is defined as Z 3 = {(xl, x2, x3)Pxl, x2, x3 are integers}. 2. DIGITAL N E I G H B O R H O O D P L A N E
A point p E Z 3 can have 26 different neighbors in its 3 × 3 x 3 neighborhood which are however classified into three categories: 6, 18, 26 neighbors, according to their relative coordinates. A point q is a 6 (18, 26) neighbor ofp ifp and q differ in at most one (two, three) coordinates by unity. For a 3 x 3 x 3 digital cube at p, Mukherjee et al. (13) have used nine planes containing p. These planes have different orientations and their normals through p are the 18-connected straight lines in the cube. These planes are termed as the digital neighborhood planes (DNP) of the point p. The planes are illustrated in Fig. l(a). The neighboring condition of p determines in which of these planes p will lie. The set of these planes is called the neighborhood plane set (NPS) of point p. The NPS value is computed for all the points in three-dimensional space and all the connected points having same NPS value form a segment.
and the extended four DNPs is that, the normals to the original nine DNPs are 18 connected straight lines, whereas the normal to the extended DNPs are 26 connected straight lines in the cube. 4. M O T I V A T I O N BEHIND T H E A L G O R I T H M
In this work the problem of range image segmentation is modelled as a quantization problem, where the surface normals are quantized to some predefined directions; or stated otherwise, the surface regions are approximated by planer surface patches which are parallel to some predefined planes. Range images provide explicit information about the visible surface of a scene. We can associate a discrete function r with a range image such that r(i, j) will give the depth value at (i, j) where (i, j ) e Z 2 and (i, j, r(i, j))eZ 3. Our interest is to extract the local orientation at a point. This orientation will be reflected in the relation between (i, j, r(i, j)) and (i', j', r(i', j') where (i', j) is a neighboring point to (i, j). Though the segmentation problem is now converted to a problem of quantization of surface normals, we now establish with the help of following examples that the quantization can be achieved using simple operations of shift, subtract and thresholding.
Example 1. Consider a curve in two dimensions represented by a function v = f(u) (refer Fig. 2a) where f(u) is a single valued function of u. It represents an analogous case of a range image in-two dimensions. The problem here is to extract the portion of the curve f(u) lying along the direction of a vector t. Let R be the set of points on f(u) lying along t. Consider any other
J iI
I
-T,~
/'
3. EXTENDED DIGITAL N E I G H B O R H O O D PLANES
Range images represent actual three-dimensional information about an object's visible surface. The normal at a point on the visible surface can have any orientation in the range (0 < 4~< 2n, 0 < 0 < n/2), where the angles q~ and 0 are the azimuthal and elevation angles, respectively, in spherical coordinate system. A hemispherical surface contains surface elements of all these orientations. Experiments with the hemispherical surface in a synthetic range image indicated that the nine DNP's introduced in reference (13) are not sufficient to segment all possible orientations in the visible range. Hence four more digital neighborhood planes are introduced as shown in Fig. l(b). These additional four planes along with the nine planes defined in reference (13) are now sufficient to represent all possible directions of normals on any visible surface in a range image. The difference between the original nine DNPs
257
:=
K
R
IW U
(a)
R
~
R1 . - - ~ ~ IL~ M
(b)
Fig. 2. Extraction of the portion of a curve v = F(u) lying along the direction of vector t.
P.K. BISWAS et al.
258
vector tl which is deviated from t by an angle V~,, and let R 1 be the set of points on f(u) lying along the direction of t~. Let the curve f(u) be shifted by the vector t and let the shifted curve be represented by ft(u) and let R' ~ R1 as shown in Fig. 2(b). Then the increment (positive or negative) measured along a direction perpendicular to t is given by IAI = 0 ,
Vp~R
< Itl*tan(V~b),
Vp~R'
= Itl*tan(V~),
V p s R ~ - R',
(1)
where Itl represents the magnitude of the vector t, and 121represents the absolute value of the increment 2 and p is a point on the curve f(u). For fixed Itl, the maximum variation in the increment 2 with V~ is shown in Fig. 3. The absolute value of the increment is very small for small values of V~b, and increases monotonically for larger values of V~,. If certain deviation from the direction of t is permitted then by applying a suitable threshold on V~, the portion of the curve closely aligned with the direction of t can be extracted. Let this threshold be V~bt. Applying the threshold V~, on V~, is equivalent to applying a threshold 6 = Itl* tan (V~,t) on I;~1-Thus the portion of the curve for which the increment (resulting from shift) measured in a direction perpendicular to t is IAI < 6, is taken to be aligned with the direction of t. The angular deviation V~kt is also termed as angular tolerance. If the normal to t makes an angle ~k with the v-axis, then the increment along v-axis is given by Vv = 2/cos (~), when Itl and V~b are small. Thus applying a threshold 6 on the increment measured along a direction perpendicular to t is equivalent to applying a threshold 6/cos(~b) on the increment Vv, where Vv = f t ( t t ) - f ( u ) . Note that this technique can not be applied if @ = 90 °, i.e., if shift vector is parallel to v-axis.
A quantitative relation between the threshold 6 and the magnitude of shift vector Itl can be deduced as follows. = Itl*tan(V~,,)~-- = tan (V~,,).
(2)
Itl
[] Example 1 illustrates the techniques of segmenting smooth curves in two dimensions into piecewise linear segments using shift, subtract and threshold operations. However, shift only along the positive direction of t alongwith the subsequent operations of subtraction and thresholding is not sufficient to extract all the desired points on the curve v = f(u) as shown in the following example.
Example 2. Consideracurvev = f(u) where in addition to deviation from t, certain portion of the curve f(u) deviates from the direction of - t also. Such a situation is depicted in Fig. 4(a). If f(u) is shifted by the vector t (Fig. 4b) a portion of the shifted curve ft(u) overlaps with f(u), where the deviation of f,(u) from - t is V~b and the deviation of f(u) from - t is less than V~'t, V~b > VSt. The set of points on f(u) in such region is represented by R~. Thus increment at a point p in such overlapped region is given by
Vv = If,(u)-f(u)l > [tl*tan(V~kt), Vp~R1.
(3)
Thus if R is the set of points lying along the direction oft, shift t extracts the set of points (R - R1). Similarly
v,~t • '~flu)
o
(a)
-IS
f ITt
.ta, iv~,)
~",-~i--4 ~ . ~ % ~- -. - .
-----
....
-.
~-ft(u)
,<
''% ,
[
!
(b)
-v~'t
v~
-fib
nil
--
-$
V~'
I v
~
R ....
,l
~R;~
-10
(c) •-15
Fig. 3. Variation of 2 with V~k for fixed Itl.
Fig. 4. Shift along t and - t to extract portion of a curve v = f ( u ) lying a l o n g the direction of t.
An SIMD algorithm for range image segmentation shift - t extracts the set of points ( R - R'I) as shown in Fig. 4(c). Taking union of these two sets (R- R,)u(R-
R'~)= R - ( R ~ c~R'~)
259
Z
(4)
r~... I
-x
If It[ and V~k, are so selected that R1 ~ R ' , = ~b, then (R - R O w ( R - g',) = R.
; tX" - .
(5) []
Thus Example 2 shows that union of the set of points extracted by shifts t and - t can extract the set of points R on f(u) lying along the direction of t. 5. F R A M E W O R K
Segmentation of surfaces in range images into piecewise planar regions, where a planar segment is parallel to one of the 13 D N P ' s is essentially a quantization of the surface normals such that the direction of the normal at a point on the surface is quantized to the direction of the normal of either of the 13 DNP's. Since all possible surface orientations in a range image are also present in a hemispherical surface, in the following discussions references will be made to hemispherical surface only unless otherwise mentioned. It is also assumed that the surfaces are piecewise smooth and the curvature at any point on the surface (except at the edge points between two smooth surface patches) is very small, i.e. the radius of curvature at such a point is sufficiently large compared to the a m o u n t of shift. Quantization of a surface point p to a D N P ~ will mean that the direction of the normal to the surface at point p is quantized to the direction of the normal N= to the D N P ~. Each of the 13 D N P ' s can be represented by its normal vector. The orientations of the normal vectors (in terms of q5 and 0) for the 13 D N P ' s are shown in Fig. 1. Let a D N P ot (1 < ct < 13) be represented by its unit normal vector N~ whose orientation is given by (~b~,0~). A point p on the surface is quantized to the D N P ct if the angular difference between the normal to the surface at point p (denoted by Np) and N~ is less than some value ft. The angular difference between the two normal vectors Np and N= is denoted by (Np, N~>. The value of fl actually depends upon the amount of tolerance given in the definition of DNP's. Thus a conical space is defined around the vector N~ such that all the surface points with normals lying in the cone with the centre of curvature of the hemisphere as apex, N= as asix and fl as semiapex angle, are quantized to the D N P ct (refer Fig. 5). This cone is denoted as CONE(~, fl). For designing a quantization strategy which quantizes the surface orientation at a point to that of one of the 13 DNP's, the following two points are to be considered. (1) It must be able to quantize all possible surface orientations in the visible range, i.e. each of the points on the hemispherical surface must be quantized to one of the 13 DNP's.
t~%='7., •~::-:.~.:.7' :?." c>
]~ [~ \.- ~ I, \
I I
Fig. 5. Intersection of a hemisphere with CONE (ct,33.75°) is a circle C. (2) The quantized direction of the normal to the surface at a given point on it should be a close approximation of the actual normal direction. The discrepency between the actual direction and the quantized direction is determined by the value of the semiapex angle ft. Considering these two points we have chosen the value of fl as 33.75 °. Choosing an semiapex angle (fl) of 33.75 ° provides a tolerance in angular deviation (Vfft) of 33.75 °. It is to be noted that with a tolerance of 33.75 ° the D N P ' s can overlap. Such overlapping of D N P ' s is necessary because if no overlapping is permitted then some of the orientations can not be quantized to any of the 13 DNP's. The physical interpretation of overlapping of D N P ' s is explained later in this section. In the following lemma we show that the chosen directions of the 13 D N P ' s alongwith the tolerance in angular deviation (V~t = 33.75 °) associated with each of them is sufficient to quantize all possible orientations in the visible range. Lemma 1. The chosen directions of the normals to the 13 DNP's alongwith the tolerance in angular deviation Vfft = 33.75 ° are sufficient to quantize all possible surface orientations in a range image. Proof. The direction of a vector in three-dimensional space can be mapped to a point in q~- 0 space, when q~- 0 is expressed as a two-dimensional coordinate system. Since the normal to a visible surface can have any orientation in the range 0 < ~ _< 2H, 0 < 0 _< I-I/2, the directions of all such surface normals are mapped to points lying in a rectangular area given by 0 < ~b < 360 °, 0 < 0 < 90 °. Consider a D N P ~twhose normal N, is mapped to the point (qS,,0,) in q~- 0 space. The intersection of the hemisphere and C O N E (cq 33.75 °) is a circle C as shown in Fig. 5. Take a point p on the circle C. Let/5 be the mapping (in ~b - 0 space) of the direction of the surface normal at point p. It is clear from the geometric construction that as p traces the perimeter of the circle C, the locus of 15 (in ~ b - 0 space) forms a closed contour. The mapping of all other normal directions which are quantized to D N P ct lies in this closed contour. All
260
P.K. BISWAS et al.
--
9O
g O) 0 o
so
loo
1so
c~
2oo
zso
300
35o
(DEGREE)
Fig. 6. Loci of the directions of surface normals at points lying on C for different DNP's.
these loci corresponding to the 13 DNP's are plotted (with V~'t = 33.75 °) in Fig. 6. The set of points bounded by the locus of/~ corresponding to a particular DNP is marked by the serial number of the DNP. Note that the area covered by DNP 3 is a rectangle in the range 0 < ~b< 360 °, 0 _<0 < 33.75 °. It is found from Fig. 6 that all the points in ~b- 0 space (0 < q~_<360 °, 0 < 0 < 90 °) is covered by all these loci. This shows that the set of 13 DNP's along with their angular orientations and tolerance is sufficient to quantize all possible surface orientations in the visible range. Q.E.D. Lemma 1 states a sufficient condition, but it is not necessary because some other choice of the orientations of the DNP's and the value of V~bt could have been used to quantize all possible orientations. But the DNP orientations and the tolerance in angular deviation used in this case is a reasonably good choice. It is to be noted that due to overlapping of DNP's certain points are quantized to more than one DNP's. The neighborhood plane set of a point p e Z 3 can now be defined as follows. Definition 1. The neighborhood plane set (NPS) of a point p on a surface in range image is defined as [p] = {~tl(N,, Np) _< V~t },
(6)
where ct, 1 _<~ < 13 denotes the serial number of the DNP's, N= is the unit normal vector to the DNP ~, Np is the normal to the surface at point p and (N=,Np) is the angular difference between the two vectors N, and Np. Thus [p] actually lists the serial numbers of those neighborhood planes whose normals have angular deviations from Np (the normal to the surface at point p) less than the tolerance V~'t. If II[p] H= 1 (where [I'll represents the cardinality of a finite set), then [p] represents the DNP to which the point p is quantized. However if I[[P] [I > 1, confusion arises regarding to which of the DNP's contained in [p] the point p is actually quantized. In such cases the point p is quantized to a plane whose orientation is given by the vector sum of the unit normal vectors to all the DNP's contained in [p]. Such situations occur in the regions where more than one DNP's overlap. Thus though there are only 13 neighborhood planes (DNP's) but due to the overlapping of neighbor-
hood planes the number of quantized orientations which can actually be obtained is more than 13.
6. THE ALGORITHM
With this definition and the physical interpretation of neighborhood plane set (NPS) of a point p on a surface in a range image, a parallel algorithm for mesh connected computer (MCC) is now presented, which computes the NPS values of all the points in a range image. In a mesh connected computer the processors are arranged in a two-dimensional array of size n × n, where each processor is connected to its eight neighbors (except for edge processors which have fewer connections). The input image is an n x n array of depth values, stored one pixel per processor. Each processor computes the NPS value of its own pixel. We have shown in Example 1, the points along a vector t on a curve v = f(u), can be extracted by shifting the curve by vector t. Since two nonparallel intersecting vectors hi and h 2 define a plane in three dimensions, the same concept can be extended in three dimensions to segment piecewise smooth surfaces in range images into piecewise planar regions. Let hi and h Ebe two nonparallel intersecting vectors. Let R~ and R 2 be the sets of points (R~, R E~_Z 3, Z 3 is three-dimensional digital space) lying along the directions of h~ and hE, respectively. Then R~ c~ R 2 is the set of points lying on the plane containing both h I and h 2. If the tolerance in angular deviation given for computing DNP's in nil (Vet = 0 °) then shift along any two nonparallel vectors hi and h 2 would have been sufficient to ensure that all the extracted surface points lie on the plane containing h i and h 2. The following example shows that when V~t > 0 ° shifts in only two directions extract some extra points the normais at which deviate from the orientation of a DNP ct by an angle more than the tolerance in angular deviation V0t. Example 3. Consider the illustration in Fig. 7(a). The normal to a DNP ct is represented by N=. Surface points are extracted using shifts along two vectors h t and h2, both lying in DNP ct. With a shift along h x and V~t as tolerance in angular deviation, the set of surface points R 1lying within two planes A and B are extracted,
An SIMD algorithm for range image segmentation
26
,z
z
.--~ ....
///'P4. ++ I ((-" v~'t i -..~
----.~,
x
(i)
(b)
Fig. 7. (a) Portion of a hemispherical surface extracted by applying shift along two intersecting h I and both lying in DNP ~. (b) Expanded view of the extracted surface patch.
0
/
,T-T.R--~ I,~A/I
/ I
,
,4 . . . . I
,
I
13 ,/\,
I
I
13 71"~'-.l---f "71 ~
,
1/i4-t-~
1
L#
III
2
5
I
I
0/I '/
0
"
t,'.____
~"
I
0
0
//~i~
....
I
I
6
71
I '
i
- ,
3
0 __
T
/'~
'
I \ / i
2
6
o
,
r-~t--,--f
'.l/l~l,
~
(----'l~"-X
'
,
2/c,,,/'
,,~
¢--,~fdfl-( i IklV! I
.4
74
,,/"~
,,"
7
8
9
0
0
F - ; ~ , S ~ ,
! /I
,I . . . .
/ i1
211
0~-4,
1/=.j~"
~'~,~
"1~
1
,zf
-~
~
-
-
, \ i _Ji-1 t
,
z
i
;.. 10
11
!
i
2
Jr
_
<.',," .... 12
13
Fig. 8. Straight lines on different digital neighborhood planes passing through the point p.
h2
P.K. BlSWAS et al.
262
where both the planes A and B pass through the centre of curvature of the hemisphere, perpendicular to the plane containing h a and N,, and are inclined at an angle of Vq6 with N,. Similarly, shift along h2 extracts the set of surface points R 2 lying within two planes C and D, where both C and D pass through the centre of curvature of the hemisphere, perpendicular to the plane containing h 2 and N,, and make angle V~kt with N,. Then R 1 n R 2 gives the set of points lying on the portion of the surface enclosed by the planes A, B, C and D and is shown in Fig. 7(b) in expanded form. Consider a point p at a corner of the surface element. Clearly ( N ~ , N p ) > V~t. Thus shifts along two nonparallel vectors only include some points the angular deviations at which are more than the tolerance V~bt. [] It follows from Example 3, if R i is the set of points extracted with a shift along a vector h i (h i lying in the plane of ct), then the point p is quantized to D N P ~ if
direction i, 0 < i < k (k = 3 for D N P ' s 1-9 and k = 2 for D N P ' s 10-13) on D N P ~t, then the set of points lying along the direction of t,~ is given by I(rh7 (x) + zio) - r(x)[n~w I(r_ g(x) - ze) - r(x) ln', (9) where 6, represents the threshold for D N P ct (recall that 6, = 6/cos 0n). The set of points quantized to D N P ~t is given by (~
(l(rh;(x) + z e) - r ( x ) l ~"
O
w ](r _ hT(x) - z e) - r(x)]~),
(10)
where k is a constant and its value is different for different neighborhood planes (k = 3 for DNP's 1 9 and k = 2 for D N P ' s 10-13). Ediuation 10 gives the set of points quantized to D N P ct as a binary image B _~ Z 2 such that Vx~B, b(x) = 1 [-where b(x) represents the binary value of x] iff (x,r(x))EZ 3 is quantized to D N P ~t.
p~Ri Yh i lying on ct.
(7)
Since there are infinite n u m b e r of straight lines, all passing through the point p and lying in a plane parallel to the D N P ~, the above method ideally requires infinite number of shifts. Fortunately, in digital domain the number of shifts in each of the D N P ' s are finite. In a 3 x 3 x 3 digital neighborhood of a point p, there are four straight lines in each of the D N P ' s 1 9 (these are represented by directions 0, I, 2 and 3) and three straight lines in each of the D N P ' s 10-13 (represented by directions 0, 1 and 2), all passing through the point p. These straight lines in all the D N P ' s are shown in Fig. 8. Thus a D N P ~ can also be defined by a set of straight lines {PI 0 < i < k}; where k is a constant (k = 3 for 1 < c t < 9 and k = 2 for 10 < ct < 13). As we have already mentioned, a range image is modelled as a function r(x), where x e Z 2 and (x, r(x)) Z 3. The function r(x) represents the range value at point x e Z 2. Thus the information stored in a range image is essentially a three-dimensional surface map. A shift t in three-dimension is denoted as [h, z], where the component of the shift on the image plane is h and on the depth axis is z. Note that z can assume both positive and negative values. Such a shift in threedimension can be realized by applying a shift by a vector h to the image r(x) along the image plane and addition of the scalar quantity z. An image r(x) shifted by a vector t is represented by rt(x), and the image r(x) thresholded by 6 is represented by r(x) ~. Thus a picture r(x) shifted by a vector t can be represented by rt(x ) = (rh(x) + Z). Similarly a picture r(x) shifted by a vector - t can be represented as r_t(x ) =(r_h(X)--z ). Thus the set of surface points lying along the direction of t is given by
[(rh(x)+z)--r(x)16ul(r_h(X)--z)--r(x)l '~,
(8)
where [.[ represents the absolute value. If t~, = [hzo,z~o] represents the shift vector along the
7. S M O O T H I N G
In the above discussion it was assumed that the image is noise free. However in practice, due to the presence of noise in the range data, discreteness in the directions of shift vectors and the amount of threshold, discreteness in the orientations of D N P ' s and the error occuring from the tolerance given in the definition of the DNP's, the points belonging to a D N P are obtained in the form of clusters of white pixels. In between clusters of such white pixels, n u m b e r of black pixels may be found. Similarly some isolated or small patches of white pixels may also occur which do not actually belong to that DNP. Hence the binary image needs to be smoothed to remove the undesired white or black patches (if the size of such patches are not too large) to produce segments of definite boundary and uniform interior. The size of the patches to be smoothed depends upon the noise level and also on the size of the objects present in the image. If this size is too large then some of the valid planar segments may also be smoothed out. A patch size of 5 × 5 gives good result with the images used in the present work. Smoothing is a two stage process. In the first stage all the black pixels inside a cluster of white pixels are removed. The second stage removes all the undesired white pixels. This is done using morphological opening and closing operations. Equation (10) gives an algorithmic method to extract the set of points in Z 3 which can be quantized to a D N P ~. Ifa point p e Z 3 is quantized to more than one neighborhood planes, then the set of all such neighborhood planes determines the neighborhood plane set (NPS) of point p denoted by [p]. As pointed out in example 1, that this technique can not be used to extract line segments which are parallel to v-axis, by similar reason since the normals to the DNP's 1, 2, 4 and 5 are perpendicular to (i.e., these D N P ' s are parallel to) the depth axis (Z-axis), points quantized to
An SIMD algorithm for range image segmentation these planes can not be determined by this method. An algorithmic description of the technique is given below. A loorithm S E G M E N T Input: Range image I _ Z 2
Output: Array of NPS values N / / A t the end of the algorithm entry N(x) contains the NPS value of the point (x, r(X))EZ 3, VX~I. Initially N(x) contains N U L L for all x . / / / / A r r a y s rl, r2, r3, r4, r l and f2 are used to store intermediate results.// for e = 1 to 13 do / / E x c e p t e = 1, 2, 4 and 5. denotes the serial number of D N P . / / (1) Vx e I do in parallel r4(x ) , - 1 (2) for i = l to k do / / i denotes the direction of the shift vector on DNP ~.// (a) V x e l do in parallel r x ( x ) ~ r (x)+ z i, ' hl (b) Vx e l do in'parallel if Ir1(x) - r(x) l < 6, then ra (x) = 1 else fl (x) = 0 (c) V x e l do in parallel r ~ ( x ) ~ - r_ ~;(x) - x i,
(d) ¥ x e l do in parallel iflrz(x) - r(x)l < 6, then f2(x) = 1 else f2(x) = 0 (e) V x e l do in parallel r3(x)*--rx(x) v f2(x) (f) ¥ x e l do in parallel r4(x)~-r4(x) ^ r3(x) end for (3) V x e l do in parallel if r4(x) = 1 then N(x) ~- N(x) u {~} end for End S E G M E N T Theorem 1. Algorithm S E G M E N T takes constant time on a mesh connected computer. Proof. In algorithm SEGMENT, Step 1 takes constant time. Each of the Steps 2(a) through 2(f) takes constant time and each of them is repeated for (k + 1) times where the value ofk is bounded (2 or 3 depending upon the D N P being computed). Thus the time taken by Step 2 is also bounded and is constant. The time taken by Step 3 is also constant. Thus each of the Steps I - 3 takes constant time and since each of them is repeated nine times (DNP's 1, 2, 4 and 5 being excluded), the algorithm S E G M E N T also takes constant time. Q.E.D. Computation of NPS value alone does not guarantee proper segmentation, since two adjacent parallel planes separated by a jump boundary will also have the same NPS value, and thus will form a single segment. Hence some restriction must be imposed on the algorithm S E G M E N T to prevent such segment merging. For
263
every pixel 'x' both its NPS value [x] and actual range value r(x) are considered for this purpose. For two neighboring pixels 'x' and 'y' if [x] = [y] and It(x) r(y) l > 6, where 6 is some threshold (in this case J is taken as 5), then 'x' and 'y' are boundary pixels of two adjacent parallel segments separated by a jump boundary. The NPS values of both 'x' and 'y' are set to NULL. This shrinks the boundaries of such segments by one pixel on all sides, but avoids the problem of segment merging. S. RESULTS AND DISCUSSIONS
Experimental results are shown for a synthetic range image of a sphere (Fig. 9a) and three real range images referred to as 'Poly' (Fig. 10a), 'Quad' (Fig. 11a)and 'Ring' (Fig. 12a). The images of Poly and Quad are of size 256 x 256 and are obtained from PRIMULA range data base of Osaka University, Japan. The other image Ring is of size 128 x 128 and is obtained from Environmental Research Institute of Michigan (ERIM), USA. The image Sphere contains only quantization noise, whereas the SNR value in P01Y and Quad is around 20dB and in Ring the SNR value is around 30dB. The results are shown with a shift of 3. The value of threshold was calculated accordingly. For the synthetic range image of sphere, the sets of points quantized to different digital neighborhood planes are shown in Fig. 9(b)-(j). In each of these figures the region of points belonging to a D N P is indicated by the closed boundary inside an outer circle, where the outer circle represents the boundary of the hemisphere. The positions of these regions with respect to the centre of the hemisphere give an indication of the range of the surface normals which are quantized. to the corresponding DNP. For example, Fig. 9(b) indicates that all the surface points where the directions of the surface normals are closed to the positive direction of the Z-axis are quantized to the D N P 3. Closed boundaries in Fig. 9(k) represent the regions of points with the same NPS value. Each such region is marked with an integer. The NPS values for these marked regions are listed in Table 1. Though region 1 is background, the points in it assume the NPS value of 3. Points belonging to region 2 are the points at which the directions of the surface normals are nearly perpendicular to the Z-axis. These points do not assume any NPS value because the DNP's 1, 2, 4 and 5 whose normals are perpendicular to Z-axis cannot be realized using this technique. For the image Poly the segmentation result without smoothing is shown in Fig. 10(b). Fig. 10(c) and (d) illustrate the segmented output when smoothing is carried out over a neighborhood of 3 x 3 and 5 x 5, respectively. It can be seen that without the smoothing operation the segment boundaries are not clearly defined (Fig. 10b). The segmentation results with range images Quad and Ring are illustrated in Figs 1 l(b) and 12(b) where the smoothing is carried out over a neighborhood of 5x5.
264
P.K. BISWAS et al.
(a)
(b) DNP 3
(c) DNP 6
(d) DNP 7
(e) DNP 8
(f) DNP 9
(g) DNP 10
(h) DNP 11
(i) DNP 12
(j) DNP 13
Fig. 9. (a) Synthetic ramge image of sphere. (b)-(j) Sets of surface points on synthetic sphere quantized to (k) different DNP's.
An SIMD algorithm for range image segmentation
265
Table 1. Regions of uniform NPS value in a hemispherical surface
Fig. 9. (Continued.)
(a) rt _
~'-' ~ t
-
-!
Region no.
NPS value
Region no.
NPS value
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
Background Nil
19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35
{3} {3,9} {3,8} {9,12} {8,11} {3,7,9,12} {3,7,8, 11} {3,7,9} {3,7} {3,7,8} {7,9,12} {7,8,12} {12} ill} {12,7} {7,11} {7}
{6}
{10} {13} {6,10} {6,13} {3,6} {3,6,9,10} {3,6,8,13} {6,9, 10} {6,8,13} {9,10} {8,13} {9} {8} {3,6,9} {3,6,8}
(b)
¢
(c)
CO)
Fig. 10. (a) Range image poly. (b) Segmentation output without smoothing. (c)Segmentation output with smoothing over 3 x 3 neighborhood. (d) Segmentation output with smoothing over 5 x 5 neighborhood.
266
P . K . BISWAS et al.
.bC
~--i
I;
I
I
(
(a)
(b)
Fig. 11. (a) Range image quad. (b) Segmentation output with smoothing over 5 × 5 neighborhood.
• o . . . . . . . . .
. . . . . . . . . . . . . . . . . .
.
. . . . .
. . o
~¢~,~,,,v......st ~tv*,*0ta...,a~,.,_
f.(t" "t
;to I
t
i,"
Z'I ii
. . . . . . .
: P
oN
~,
a
1
I
p
.", IO
•
s,II
JI
@
!
..~,.w,.,._.,.~.,..,..,.:
(a)
. . . . . . . . .
(b)
Fig. 12. (a) Range image Ring. (b) Segmentation output with smoothing over 5 x 5 neighborhood.
9. CONCLUSION
In this paper the problem of segmentation of surfaces in range images is modelled as a quantization problem, where the orientations of the surface normals are quantized to the orientations of normals to some predefined planes (DNP's). Though the algorithm quantizes the surface normals it does not compute the normals at all the points on a surface. Rather the segmentation is achieved by simple shift, subtract and threshold operations. Hence this approach for surface segmentation in range images is computationally more efficient than the existing range image segmentation algorithms(3,6-10.12) which normally attempt either to fit the surface points to a surface function or to calculate a local property (such as normal, curvature etc.). The algorithm described in this paper is parallel and its time complexity is shown to be of 0(1) on a mesh connected computer. The algorithm presented in this paper produces excellent segmentation results with polyhedral ob-
jects (Fig. 10d) and reasonably good results with curved objects (Fig. l i b and 12b). The segmentation produced by our algorithm is comparable to the ones in Mukherjee e t al. (13) Morphological operations of opening and closing are conveniently incorporated in our technique to smooth the regions quantized to different DNP's, thus making our technique robust to noise. Smoothing takes care of the noise present in the range data and the noise arising out of the discreteness in shift vectors, discreteness in the threshold and the tolerance given in the DNP's. Excellent segmentation results show the power of our technique. Acknowledgement--The authors are grateful to the anonymous reviewer for his/her valuable comments on an earlier
version of this paper. REFERENCES
1. P.J. Besl and R.'C. Jain, Three-dimensional object recognition, A C M Comp. Surv. 17(1), 75-145 (1985).
An SIMD algorithm for range image segmentation
2. P.K. Biswas, J. Mukherjee, B.N. Chatterji and P.P. Chakraborty, Qualitative description of three-dimensional scenes, Int. J. Pattern Recognition and Artificial Intelligence 6(4), 65t-672 (1992). 3. D. L. Milgrim and C. M. Bjorklund, Range image processing planar surface extraction, Proc. 5th Int. Conf. Pattern Recognition, pp. 912-919, Miami, Florida (1-4 December 1980). 4. T. C. Henderson, Efficient 3-D object representation for industrial vision systems, IEEE Trans. Pattern Anal. Mach. InteU. 5(6), 609 617 (1983). 5. T. C. Henderson and B. Bhanu, Three-point seed method for the extraction of planar faces from range data, Proc. Workshop Industrial Applications of Machine Vision Research, pp. 181 186, Triangle Park, North Corolina (May 1982). 6. J. H. Han, R. A. Volz and T. H. Mudge, Range image segmentation and surface parameter extraction for 3-D object recognition of industrial parts, Proc. on 1987 IEEE Int. Conf. on Robotics and Automation, 1,380-386 (March 1987). 7. M. Herbert and J. Ponce, A new Method for segmenting 3-D scenes into primitives, Proc. Int. Conf. Pattern Recognition, pp. 836-838, Munich, West Germany, (21-22 October 1982).
267
8. I. K. Sethi and S. N. Jayaramamurthy, Surface classification using characteristic contours, Proc. 7th Int. Conf. Pattern Recognition, pp. 438-440, Montreal, P. Q., Canada (30 July-2 August 1984). 9. P. J. Besl and R. C. Jain, Segmentation through variableorder surface fitting, IEEE Trans. Pattern Anal. Much. Intell. 10(2), 167-192 (1988). 10. P. J. Besl and R. C. Jain, Invariant Surface Characteristics for 3-D object recognition in range images, Computer Vision Graphics and Image Processing, 33, 33-80 (1986). 11. M. Brady, J. Ponce, A. Yuille and H. Asada, Describing surfaces, Comput. Vision Graphics Image Process. 32, 1-28 (1985). 12. Russel, W. Taylor, M. Savini and A. P. Reeves, Fast segmentation of range imagery into planar regions, Computer Vision Graphics and Image Process. 45, 42-60 (1989). 13. Jayanta Mukherjee, P. P. Das and B. N. Chatterji, Segmentation of Range Images, Pattern Recognition 25(10), 1141-1156 (1992). 14. P. K. Biswas, S. S. Biswas, A. K. Roy and B. N. Chatterji, Segmentation of range images using morphological operators, Proc. IEEE Region 10 Intl. Conf. Energy, Computer, Communication and Control Systems, pp. 338343, New Delhi (28 30 August 1991).
About the Author--PRABIR KR. BISWAS completed his B.Tech(Hons), M.Tech and Ph.D. in electronics and electrical communication eng. from Indian Institute of Technology, Kharagpur, in the years 1985, 1988 and 1991, respectively. From 1985 to 1987 he was with Bharat Electronics Ltd., Ghaziabad, India. Currently he is an Assistant Professor in the Dept. of Electronics and Electrical Communication Eng., Indian Institute of Technology, Kharagpur. His areas of interest include pattern recognition, computer vision, parallel processing, distributed control and computer network.
About the Author--S. S. BISWAS completed his B.Tech (hons), M. Tech and Ph.D. in electronics and electrical communication engineering from Indian Institute of Technology, Kharagpur in the years, 1985, 1986 and 1993, respectively. Currently he is a Software Engineer in Tata Consultancy Services Ltd., India. His areas of interest include pattern recognition, computer vision and systems engineering.
About the Author--B. N. CHATTERJI completed his B.Tech(hons) and Ph.D. in electronics and electrical communication eng. from Indian Institute of Technology in the years 1965 and 1970. He was with Telerad Pvt. Ltd. Bombay in 1965 and Central Electronics Engineering Research Institute, Pilani in 1966. Since 1967 he has worked as a faculty member in the Dept. of Electronics and Electrical Communication Eng., Indian Institute of Technology, Kharagpur where he has been Professor since 1980 and was Head of the Department during 1987-1991. He is a fellow of several professional societies like the Institution of Engineers (India), the Institution of Electronics and Telecommunication Engineers (India) etc. Dr Chatterji has published extensively in journals and is the author of two text books. He is the recipient of Hari Om Ashram Prerit Dr Vikram Sarabhai Award for the year 1983 for his contributions in the area of Electronics and Telecommunications. His areas of interest are pattern recognition, image processing, computer vision, parallel processing, system identification and digital signal processing.
Pergamon
0031-3203(94)00093-X
AN S I M D A L G O R I T H M F O R R A N G E IMAGE SEGMENTATION PRABIR KR. BISWAS, t S. S. BISWAS and B. N. C H A T T E R J I Department of Electronics and Electrical Communication Engineering, Indian Institute of Technology, Kharagpur - 721 302, India (Received 10 March 1993; in revised form 26 July 1994; received for publication 5 August 1994)
Abstract--Understanding 3D scenes from range images need the segmentation of 3D surfaces into approximate planar surface patches from which curved surfaces can be constructed quickly for high level vision purpose [Biswas et al., Qualitative description of three-dimensional scenes, Intl. J. Pattern Recognition and Artificial Intelligence 6(4), 651-672 (1992)]. In this paper we have presented a new parallel algorithm for 3D surface segmentation wherein the problem of surface segmentation is modelled as a quantization problem. The surface normals are quantized to some predefined directions, or stated otherwise, the surface regions are approximated by planar surface patches which are parallel to some predefined planes. The novelty of the algorithm lies in the fact that, though surface segmentation is achieved by quantization of surface normals, the algorithm does not compute the surface normals explicitly. Rather the quantization is achieved using simple operations of shift, subtract and threshold. This technique suitably avoids the computations of the differential properties of the surfaces or the surface fitting expressions which are used in most of the other existing techniques. Hence this approach is computationally attractive. This algorithm is easily implementable on an SIMD array computer. Another advantage of this technique is that it is robust to noise present in the image. The algorithm has been explained with a number of examples. Experimental results with synthetic as well as real range images are cited in this paper to highlight distinctive features of the algorithm. Segmentation
Range image
Parallel algorithms
I. INTRODUCTION Processing and recognition of three-dimensional objects from range images has gained importance in recent years. Various approaches towards this problem have been reported in literature. A range image is a twodimensional array of numbers, where each array element denotes the distance between the image plane and the object surface at that point. Since an image array contains a huge a m o u n t of sensed values, it is often difficult for a computer to interpret a scene directly from raw image data. Thus a preliminary processing of the raw sensor data is needed which partitions (segments) an image array into groups of pixels with particular properties that can be conveniently used by the later stages of representation, description and recognition of the scene. Thus every segment is a collection of pixels which have some uniform property quite distinct from the property of pixels belonging to other segments. It is often difficult to give a strict mathematical definition to this property because the requirement of segmentation is often guided by the high-level description/recognition stage. Most of the attempt# t) for range image segmentation either fit the surface points to a surface function or calculate a local surface property such as normal, curvature etc. These tAuthor to whom all correspondence should be addressed.
SIMD array architecture
algorithms are often computation intensive and not amenable to easy parallelization. In this paper we present a new parallel segmentation algorithm which segments range images into approximate planar surface patches. Segmenting range images into approximate planar regions provides a basis from which curved surfaces can be constructed quickly for high-level vision purposes. ~2~ A n u m b e r of approaches for range image processing have been reported. Among those using surface based properties some extract only planar faces. Milgrim and Bjorklund t3) have suggested a method in which planar surface is extracted using a 5 × 5 window and surfaces are fitted using least square method. Henderson 14) organized the points into a three-dimensional binary tree from which he obtained a three-dimensional spatial proximity graph. Then a sequential region growing algorithm is used to create convex planar faces. A method known as 3-point seed method was suggested in reference (5) where region growing started from three neighboring points. Planar, cylinderical, conical and spherical primitives are handled in references (6) and (7). Spheres and ellipsoids were dealt with by Sethi and Jayaramamurthy Im Besl and Jain tg) developed an algorithm based on variable order surface fitting. In addition to the above there are segmentation techniques using local properties of a pixel. Besl and Jain t 1o) proposed a method where each pixel is labelled 255
256
P.K. BISWAS et al. Plane Set (NPS) values indicative of the surface normal at every range pixel. We have modified an extended the concept of digital neighborhood plane (DNP) and neighborhood plane set (NPS), introduced in reference (13) for segmentating surfaces in range images. In this paper we have presented a new parallel algorithm for segmenting surfaces in range images. The problem of surface segmentation in range image is modelled as a quantization problem, where the normal to the surface at a point is quantized
by the sign of mean and Gaussian curvatures at that point. Brady et al. ill) discussed a method of segmentation where features such as curvatures, asymptotes, boundary contours, certain surface patches etc. are calculated using principal curvatures, principal directions and a surface primal sketch. Taylor et al. (12) used a split and merge algorithm with a homogeneity criteria based on a 3-parameter surface description technique. Mukherjee et al. (13) have proposed a technique where a range image is segmented based on Neighborhood
1
2
-i
- - t"~,
II
li
71
,,,
; ,'-l--J "
= 0 ,180 °
~
O = 90 °
:
0
=
'I,
!
o
0
--
m
--..i,"
0 ° ~ ~6 ~< 3 6 0 °
90 °
O= 0 °
,'
,
,/''~
I
, ',Xi
,/I
,IL--+I,J 7 7,"
',
I, ~ ,, , / t - - - - - , :-i ,.
I
I
-)--)
_1,"
I
:
/+5°, 2 2 5 °
¢, :
135°315 °
@ = 90 °
:
90 °
0
90 °
0
:
7
J
,
+,-t-'l l
N / N/ I
I
1 "i---I-~
'U_ "_ __ d ; : @ =
"_~_ _ _ u,"
2700
~
:
0°
£~5 °
6
=
&5 °
10
: 0 =
: ,+-y':-r, ,> I
,
~__~_ _-,
"_>'C__?',"
180 ° ~5°
12
11
4
I
~5 °
,'~,~
J
'~
=
I
9
e
r0
i
Z
6
(-- -I..: - ,7/( i
li
~ L,,-----!,-I~
I,
II~ -0
L_
E:
II
$
,--..(k: -~- .-.~,-f,
,,_J < z
,,
90 , 2 7 0
4
cl
III
kC.__
0
cr~ L£J Z <
+I '- ---"I
i
0
Y
x,'F-- - :--'i
,,,
-+-r-r,/l '
i
Id_ _ k'_ _I,"
,_I
3
,,4"
I
~,\l
13
,4"~-
#~
i
i
, / - . ~ "u"+
,~, = 135°
~ : 31 5 ° '
~. : 2 2 5 °
•
0 = $/.,.73 =
0 = 5/4.73 °
0 :
EXTENDED
I
IZ
O : 5&.73 °
(b)
-- --7l
=
45 ° 5&.73 °
PLANES
Fig. 1. Orientations of the normals of digital neighborhood planes.
An SIMD algorithm for range image segmentation to the normal of a given plane based on a proximity measure. The angular difference between the surface normal at a given point and the normal to the given plane is taken as a measure of proximity between them. The new algorithm has got the advantage of being simple in design, efficient in execution and robust to noise. The complexity of the algorithm is 0(1) on a mesh connected computer. A preliminary version of this work has been reported in reference (14). The following two notations are used extensively in this paper.
Notation 1. A two-dimensional digital space Z 2 is defined as Z 2 = { (Xl, x 2 ) l X l , x 2
are integers}.
Notation 2. A three-dimensional digital space Z 3 is defined as Z 3 = {(xl, x2, x3)Pxl, x2, x3 are integers}. 2. DIGITAL N E I G H B O R H O O D P L A N E
A point p E Z 3 can have 26 different neighbors in its 3 × 3 x 3 neighborhood which are however classified into three categories: 6, 18, 26 neighbors, according to their relative coordinates. A point q is a 6 (18, 26) neighbor ofp ifp and q differ in at most one (two, three) coordinates by unity. For a 3 x 3 x 3 digital cube at p, Mukherjee et al. (13) have used nine planes containing p. These planes have different orientations and their normals through p are the 18-connected straight lines in the cube. These planes are termed as the digital neighborhood planes (DNP) of the point p. The planes are illustrated in Fig. l(a). The neighboring condition of p determines in which of these planes p will lie. The set of these planes is called the neighborhood plane set (NPS) of point p. The NPS value is computed for all the points in three-dimensional space and all the connected points having same NPS value form a segment.
and the extended four DNPs is that, the normals to the original nine DNPs are 18 connected straight lines, whereas the normal to the extended DNPs are 26 connected straight lines in the cube. 4. M O T I V A T I O N BEHIND T H E A L G O R I T H M
In this work the problem of range image segmentation is modelled as a quantization problem, where the surface normals are quantized to some predefined directions; or stated otherwise, the surface regions are approximated by planer surface patches which are parallel to some predefined planes. Range images provide explicit information about the visible surface of a scene. We can associate a discrete function r with a range image such that r(i, j) will give the depth value at (i, j) where (i, j ) e Z 2 and (i, j, r(i, j))eZ 3. Our interest is to extract the local orientation at a point. This orientation will be reflected in the relation between (i, j, r(i, j)) and (i', j', r(i', j') where (i', j) is a neighboring point to (i, j). Though the segmentation problem is now converted to a problem of quantization of surface normals, we now establish with the help of following examples that the quantization can be achieved using simple operations of shift, subtract and thresholding.
Example 1. Consider a curve in two dimensions represented by a function v = f(u) (refer Fig. 2a) where f(u) is a single valued function of u. It represents an analogous case of a range image in-two dimensions. The problem here is to extract the portion of the curve f(u) lying along the direction of a vector t. Let R be the set of points on f(u) lying along t. Consider any other
J iI
I
-T,~
/'
3. EXTENDED DIGITAL N E I G H B O R H O O D PLANES
Range images represent actual three-dimensional information about an object's visible surface. The normal at a point on the visible surface can have any orientation in the range (0 < 4~< 2n, 0 < 0 < n/2), where the angles q~ and 0 are the azimuthal and elevation angles, respectively, in spherical coordinate system. A hemispherical surface contains surface elements of all these orientations. Experiments with the hemispherical surface in a synthetic range image indicated that the nine DNP's introduced in reference (13) are not sufficient to segment all possible orientations in the visible range. Hence four more digital neighborhood planes are introduced as shown in Fig. l(b). These additional four planes along with the nine planes defined in reference (13) are now sufficient to represent all possible directions of normals on any visible surface in a range image. The difference between the original nine DNPs
257
:=
K
R
IW U
(a)
R
~
R1 . - - ~ ~ IL~ M
(b)
Fig. 2. Extraction of the portion of a curve v = F(u) lying along the direction of vector t.
P.K. BISWAS et al.
258
vector tl which is deviated from t by an angle V~,, and let R 1 be the set of points on f(u) lying along the direction of t~. Let the curve f(u) be shifted by the vector t and let the shifted curve be represented by ft(u) and let R' ~ R1 as shown in Fig. 2(b). Then the increment (positive or negative) measured along a direction perpendicular to t is given by IAI = 0 ,
Vp~R
< Itl*tan(V~b),
Vp~R'
= Itl*tan(V~),
V p s R ~ - R',
(1)
where Itl represents the magnitude of the vector t, and 121represents the absolute value of the increment 2 and p is a point on the curve f(u). For fixed Itl, the maximum variation in the increment 2 with V~ is shown in Fig. 3. The absolute value of the increment is very small for small values of V~b, and increases monotonically for larger values of V~,. If certain deviation from the direction of t is permitted then by applying a suitable threshold on V~, the portion of the curve closely aligned with the direction of t can be extracted. Let this threshold be V~bt. Applying the threshold V~, on V~, is equivalent to applying a threshold 6 = Itl* tan (V~,t) on I;~1-Thus the portion of the curve for which the increment (resulting from shift) measured in a direction perpendicular to t is IAI < 6, is taken to be aligned with the direction of t. The angular deviation V~kt is also termed as angular tolerance. If the normal to t makes an angle ~k with the v-axis, then the increment along v-axis is given by Vv = 2/cos (~), when Itl and V~b are small. Thus applying a threshold 6 on the increment measured along a direction perpendicular to t is equivalent to applying a threshold 6/cos(~b) on the increment Vv, where Vv = f t ( t t ) - f ( u ) . Note that this technique can not be applied if @ = 90 °, i.e., if shift vector is parallel to v-axis.
A quantitative relation between the threshold 6 and the magnitude of shift vector Itl can be deduced as follows. = Itl*tan(V~,,)~-- = tan (V~,,).
(2)
Itl
[] Example 1 illustrates the techniques of segmenting smooth curves in two dimensions into piecewise linear segments using shift, subtract and threshold operations. However, shift only along the positive direction of t alongwith the subsequent operations of subtraction and thresholding is not sufficient to extract all the desired points on the curve v = f(u) as shown in the following example.
Example 2. Consideracurvev = f(u) where in addition to deviation from t, certain portion of the curve f(u) deviates from the direction of - t also. Such a situation is depicted in Fig. 4(a). If f(u) is shifted by the vector t (Fig. 4b) a portion of the shifted curve ft(u) overlaps with f(u), where the deviation of f,(u) from - t is V~b and the deviation of f(u) from - t is less than V~'t, V~b > VSt. The set of points on f(u) in such region is represented by R~. Thus increment at a point p in such overlapped region is given by
Vv = If,(u)-f(u)l > [tl*tan(V~kt), Vp~R1.
(3)
Thus if R is the set of points lying along the direction oft, shift t extracts the set of points (R - R1). Similarly
v,~t • '~flu)
o
(a)
-IS
f ITt
.ta, iv~,)
~",-~i--4 ~ . ~ % ~- -. - .
-----
....
-.
~-ft(u)
,<
''% ,
[
!
(b)
-v~'t
v~
-fib
nil
--
-$
V~'
I v
~
R ....
,l
~R;~
-10
(c) •-15
Fig. 3. Variation of 2 with V~k for fixed Itl.
Fig. 4. Shift along t and - t to extract portion of a curve v = f ( u ) lying a l o n g the direction of t.
An SIMD algorithm for range image segmentation shift - t extracts the set of points ( R - R'I) as shown in Fig. 4(c). Taking union of these two sets (R- R,)u(R-
R'~)= R - ( R ~ c~R'~)
259
Z
(4)
r~... I
-x
If It[ and V~k, are so selected that R1 ~ R ' , = ~b, then (R - R O w ( R - g',) = R.
; tX" - .
(5) []
Thus Example 2 shows that union of the set of points extracted by shifts t and - t can extract the set of points R on f(u) lying along the direction of t. 5. F R A M E W O R K
Segmentation of surfaces in range images into piecewise planar regions, where a planar segment is parallel to one of the 13 D N P ' s is essentially a quantization of the surface normals such that the direction of the normal at a point on the surface is quantized to the direction of the normal of either of the 13 DNP's. Since all possible surface orientations in a range image are also present in a hemispherical surface, in the following discussions references will be made to hemispherical surface only unless otherwise mentioned. It is also assumed that the surfaces are piecewise smooth and the curvature at any point on the surface (except at the edge points between two smooth surface patches) is very small, i.e. the radius of curvature at such a point is sufficiently large compared to the a m o u n t of shift. Quantization of a surface point p to a D N P ~ will mean that the direction of the normal to the surface at point p is quantized to the direction of the normal N= to the D N P ~. Each of the 13 D N P ' s can be represented by its normal vector. The orientations of the normal vectors (in terms of q5 and 0) for the 13 D N P ' s are shown in Fig. 1. Let a D N P ot (1 < ct < 13) be represented by its unit normal vector N~ whose orientation is given by (~b~,0~). A point p on the surface is quantized to the D N P ct if the angular difference between the normal to the surface at point p (denoted by Np) and N~ is less than some value ft. The angular difference between the two normal vectors Np and N= is denoted by (Np, N~>. The value of fl actually depends upon the amount of tolerance given in the definition of DNP's. Thus a conical space is defined around the vector N~ such that all the surface points with normals lying in the cone with the centre of curvature of the hemisphere as apex, N= as asix and fl as semiapex angle, are quantized to the D N P ct (refer Fig. 5). This cone is denoted as CONE(~, fl). For designing a quantization strategy which quantizes the surface orientation at a point to that of one of the 13 DNP's, the following two points are to be considered. (1) It must be able to quantize all possible surface orientations in the visible range, i.e. each of the points on the hemispherical surface must be quantized to one of the 13 DNP's.
t~%='7., •~::-:.~.:.7' :?." c>
]~ [~ \.- ~ I, \
I I
Fig. 5. Intersection of a hemisphere with CONE (ct,33.75°) is a circle C. (2) The quantized direction of the normal to the surface at a given point on it should be a close approximation of the actual normal direction. The discrepency between the actual direction and the quantized direction is determined by the value of the semiapex angle ft. Considering these two points we have chosen the value of fl as 33.75 °. Choosing an semiapex angle (fl) of 33.75 ° provides a tolerance in angular deviation (Vfft) of 33.75 °. It is to be noted that with a tolerance of 33.75 ° the D N P ' s can overlap. Such overlapping of D N P ' s is necessary because if no overlapping is permitted then some of the orientations can not be quantized to any of the 13 DNP's. The physical interpretation of overlapping of D N P ' s is explained later in this section. In the following lemma we show that the chosen directions of the 13 D N P ' s alongwith the tolerance in angular deviation (V~t = 33.75 °) associated with each of them is sufficient to quantize all possible orientations in the visible range. Lemma 1. The chosen directions of the normals to the 13 DNP's alongwith the tolerance in angular deviation Vfft = 33.75 ° are sufficient to quantize all possible surface orientations in a range image. Proof. The direction of a vector in three-dimensional space can be mapped to a point in q~- 0 space, when q~- 0 is expressed as a two-dimensional coordinate system. Since the normal to a visible surface can have any orientation in the range 0 < ~ _< 2H, 0 < 0 _< I-I/2, the directions of all such surface normals are mapped to points lying in a rectangular area given by 0 < ~b < 360 °, 0 < 0 < 90 °. Consider a D N P ~twhose normal N, is mapped to the point (qS,,0,) in q~- 0 space. The intersection of the hemisphere and C O N E (cq 33.75 °) is a circle C as shown in Fig. 5. Take a point p on the circle C. Let/5 be the mapping (in ~b - 0 space) of the direction of the surface normal at point p. It is clear from the geometric construction that as p traces the perimeter of the circle C, the locus of 15 (in ~ b - 0 space) forms a closed contour. The mapping of all other normal directions which are quantized to D N P ct lies in this closed contour. All
260
P.K. BISWAS et al.
--
9O
g O) 0 o
so
loo
1so
c~
2oo
zso
300
35o
(DEGREE)
Fig. 6. Loci of the directions of surface normals at points lying on C for different DNP's.
these loci corresponding to the 13 DNP's are plotted (with V~'t = 33.75 °) in Fig. 6. The set of points bounded by the locus of/~ corresponding to a particular DNP is marked by the serial number of the DNP. Note that the area covered by DNP 3 is a rectangle in the range 0 < ~b< 360 °, 0 _<0 < 33.75 °. It is found from Fig. 6 that all the points in ~b- 0 space (0 < q~_<360 °, 0 < 0 < 90 °) is covered by all these loci. This shows that the set of 13 DNP's along with their angular orientations and tolerance is sufficient to quantize all possible surface orientations in the visible range. Q.E.D. Lemma 1 states a sufficient condition, but it is not necessary because some other choice of the orientations of the DNP's and the value of V~bt could have been used to quantize all possible orientations. But the DNP orientations and the tolerance in angular deviation used in this case is a reasonably good choice. It is to be noted that due to overlapping of DNP's certain points are quantized to more than one DNP's. The neighborhood plane set of a point p e Z 3 can now be defined as follows. Definition 1. The neighborhood plane set (NPS) of a point p on a surface in range image is defined as [p] = {~tl(N,, Np) _< V~t },
(6)
where ct, 1 _<~ < 13 denotes the serial number of the DNP's, N= is the unit normal vector to the DNP ~, Np is the normal to the surface at point p and (N=,Np) is the angular difference between the two vectors N, and Np. Thus [p] actually lists the serial numbers of those neighborhood planes whose normals have angular deviations from Np (the normal to the surface at point p) less than the tolerance V~'t. If II[p] H= 1 (where [I'll represents the cardinality of a finite set), then [p] represents the DNP to which the point p is quantized. However if I[[P] [I > 1, confusion arises regarding to which of the DNP's contained in [p] the point p is actually quantized. In such cases the point p is quantized to a plane whose orientation is given by the vector sum of the unit normal vectors to all the DNP's contained in [p]. Such situations occur in the regions where more than one DNP's overlap. Thus though there are only 13 neighborhood planes (DNP's) but due to the overlapping of neighbor-
hood planes the number of quantized orientations which can actually be obtained is more than 13.
6. THE ALGORITHM
With this definition and the physical interpretation of neighborhood plane set (NPS) of a point p on a surface in a range image, a parallel algorithm for mesh connected computer (MCC) is now presented, which computes the NPS values of all the points in a range image. In a mesh connected computer the processors are arranged in a two-dimensional array of size n × n, where each processor is connected to its eight neighbors (except for edge processors which have fewer connections). The input image is an n x n array of depth values, stored one pixel per processor. Each processor computes the NPS value of its own pixel. We have shown in Example 1, the points along a vector t on a curve v = f(u), can be extracted by shifting the curve by vector t. Since two nonparallel intersecting vectors hi and h 2 define a plane in three dimensions, the same concept can be extended in three dimensions to segment piecewise smooth surfaces in range images into piecewise planar regions. Let hi and h Ebe two nonparallel intersecting vectors. Let R~ and R 2 be the sets of points (R~, R E~_Z 3, Z 3 is three-dimensional digital space) lying along the directions of h~ and hE, respectively. Then R~ c~ R 2 is the set of points lying on the plane containing both h I and h 2. If the tolerance in angular deviation given for computing DNP's in nil (Vet = 0 °) then shift along any two nonparallel vectors hi and h 2 would have been sufficient to ensure that all the extracted surface points lie on the plane containing h i and h 2. The following example shows that when V~t > 0 ° shifts in only two directions extract some extra points the normais at which deviate from the orientation of a DNP ct by an angle more than the tolerance in angular deviation V0t. Example 3. Consider the illustration in Fig. 7(a). The normal to a DNP ct is represented by N=. Surface points are extracted using shifts along two vectors h t and h2, both lying in DNP ct. With a shift along h x and V~t as tolerance in angular deviation, the set of surface points R 1lying within two planes A and B are extracted,
An SIMD algorithm for range image segmentation
26
,z
z
.--~ ....
///'P4. ++ I ((-" v~'t i -..~
----.~,
x
(i)
(b)
Fig. 7. (a) Portion of a hemispherical surface extracted by applying shift along two intersecting h I and both lying in DNP ~. (b) Expanded view of the extracted surface patch.
0
/
,T-T.R--~ I,~A/I
/ I
,
,4 . . . . I
,
I
13 ,/\,
I
I
13 71"~'-.l---f "71 ~
,
1/i4-t-~
1
L#
III
2
5
I
I
0/I '/
0
"
t,'.____
~"
I
0
0
//~i~
....
I
I
6
71
I '
i
- ,
3
0 __
T
/'~
'
I \ / i
2
6
o
,
r-~t--,--f
'.l/l~l,
~
(----'l~"-X
'
,
2/c,,,/'
,,~
¢--,~fdfl-( i IklV! I
.4
74
,,/"~
,,"
7
8
9
0
0
F - ; ~ , S ~ ,
! /I
,I . . . .
/ i1
211
0~-4,
1/=.j~"
~'~,~
"1~
1
,zf
-~
~
-
-
, \ i _Ji-1 t
,
z
i
;.. 10
11
!
i
2
Jr
_
<.',," .... 12
13
Fig. 8. Straight lines on different digital neighborhood planes passing through the point p.
h2
P.K. BlSWAS et al.
262
where both the planes A and B pass through the centre of curvature of the hemisphere, perpendicular to the plane containing h a and N,, and are inclined at an angle of Vq6 with N,. Similarly, shift along h2 extracts the set of surface points R 2 lying within two planes C and D, where both C and D pass through the centre of curvature of the hemisphere, perpendicular to the plane containing h 2 and N,, and make angle V~kt with N,. Then R 1 n R 2 gives the set of points lying on the portion of the surface enclosed by the planes A, B, C and D and is shown in Fig. 7(b) in expanded form. Consider a point p at a corner of the surface element. Clearly ( N ~ , N p ) > V~t. Thus shifts along two nonparallel vectors only include some points the angular deviations at which are more than the tolerance V~bt. [] It follows from Example 3, if R i is the set of points extracted with a shift along a vector h i (h i lying in the plane of ct), then the point p is quantized to D N P ~ if
direction i, 0 < i < k (k = 3 for D N P ' s 1-9 and k = 2 for D N P ' s 10-13) on D N P ~t, then the set of points lying along the direction of t,~ is given by I(rh7 (x) + zio) - r(x)[n~w I(r_ g(x) - ze) - r(x) ln', (9) where 6, represents the threshold for D N P ct (recall that 6, = 6/cos 0n). The set of points quantized to D N P ~t is given by (~
(l(rh;(x) + z e) - r ( x ) l ~"
O
w ](r _ hT(x) - z e) - r(x)]~),
(10)
where k is a constant and its value is different for different neighborhood planes (k = 3 for DNP's 1 9 and k = 2 for D N P ' s 10-13). Ediuation 10 gives the set of points quantized to D N P ct as a binary image B _~ Z 2 such that Vx~B, b(x) = 1 [-where b(x) represents the binary value of x] iff (x,r(x))EZ 3 is quantized to D N P ~t.
p~Ri Yh i lying on ct.
(7)
Since there are infinite n u m b e r of straight lines, all passing through the point p and lying in a plane parallel to the D N P ~, the above method ideally requires infinite number of shifts. Fortunately, in digital domain the number of shifts in each of the D N P ' s are finite. In a 3 x 3 x 3 digital neighborhood of a point p, there are four straight lines in each of the D N P ' s 1 9 (these are represented by directions 0, I, 2 and 3) and three straight lines in each of the D N P ' s 10-13 (represented by directions 0, 1 and 2), all passing through the point p. These straight lines in all the D N P ' s are shown in Fig. 8. Thus a D N P ~ can also be defined by a set of straight lines {PI 0 < i < k}; where k is a constant (k = 3 for 1 < c t < 9 and k = 2 for 10 < ct < 13). As we have already mentioned, a range image is modelled as a function r(x), where x e Z 2 and (x, r(x)) Z 3. The function r(x) represents the range value at point x e Z 2. Thus the information stored in a range image is essentially a three-dimensional surface map. A shift t in three-dimension is denoted as [h, z], where the component of the shift on the image plane is h and on the depth axis is z. Note that z can assume both positive and negative values. Such a shift in threedimension can be realized by applying a shift by a vector h to the image r(x) along the image plane and addition of the scalar quantity z. An image r(x) shifted by a vector t is represented by rt(x), and the image r(x) thresholded by 6 is represented by r(x) ~. Thus a picture r(x) shifted by a vector t can be represented by rt(x ) = (rh(x) + Z). Similarly a picture r(x) shifted by a vector - t can be represented as r_t(x ) =(r_h(X)--z ). Thus the set of surface points lying along the direction of t is given by
[(rh(x)+z)--r(x)16ul(r_h(X)--z)--r(x)l '~,
(8)
where [.[ represents the absolute value. If t~, = [hzo,z~o] represents the shift vector along the
7. S M O O T H I N G
In the above discussion it was assumed that the image is noise free. However in practice, due to the presence of noise in the range data, discreteness in the directions of shift vectors and the amount of threshold, discreteness in the orientations of D N P ' s and the error occuring from the tolerance given in the definition of the DNP's, the points belonging to a D N P are obtained in the form of clusters of white pixels. In between clusters of such white pixels, n u m b e r of black pixels may be found. Similarly some isolated or small patches of white pixels may also occur which do not actually belong to that DNP. Hence the binary image needs to be smoothed to remove the undesired white or black patches (if the size of such patches are not too large) to produce segments of definite boundary and uniform interior. The size of the patches to be smoothed depends upon the noise level and also on the size of the objects present in the image. If this size is too large then some of the valid planar segments may also be smoothed out. A patch size of 5 × 5 gives good result with the images used in the present work. Smoothing is a two stage process. In the first stage all the black pixels inside a cluster of white pixels are removed. The second stage removes all the undesired white pixels. This is done using morphological opening and closing operations. Equation (10) gives an algorithmic method to extract the set of points in Z 3 which can be quantized to a D N P ~. Ifa point p e Z 3 is quantized to more than one neighborhood planes, then the set of all such neighborhood planes determines the neighborhood plane set (NPS) of point p denoted by [p]. As pointed out in example 1, that this technique can not be used to extract line segments which are parallel to v-axis, by similar reason since the normals to the DNP's 1, 2, 4 and 5 are perpendicular to (i.e., these D N P ' s are parallel to) the depth axis (Z-axis), points quantized to
An SIMD algorithm for range image segmentation these planes can not be determined by this method. An algorithmic description of the technique is given below. A loorithm S E G M E N T Input: Range image I _ Z 2
Output: Array of NPS values N / / A t the end of the algorithm entry N(x) contains the NPS value of the point (x, r(X))EZ 3, VX~I. Initially N(x) contains N U L L for all x . / / / / A r r a y s rl, r2, r3, r4, r l and f2 are used to store intermediate results.// for e = 1 to 13 do / / E x c e p t e = 1, 2, 4 and 5. denotes the serial number of D N P . / / (1) Vx e I do in parallel r4(x ) , - 1 (2) for i = l to k do / / i denotes the direction of the shift vector on DNP ~.// (a) V x e l do in parallel r x ( x ) ~ r (x)+ z i, ' hl (b) Vx e l do in'parallel if Ir1(x) - r(x) l < 6, then ra (x) = 1 else fl (x) = 0 (c) V x e l do in parallel r ~ ( x ) ~ - r_ ~;(x) - x i,
(d) ¥ x e l do in parallel iflrz(x) - r(x)l < 6, then f2(x) = 1 else f2(x) = 0 (e) V x e l do in parallel r3(x)*--rx(x) v f2(x) (f) ¥ x e l do in parallel r4(x)~-r4(x) ^ r3(x) end for (3) V x e l do in parallel if r4(x) = 1 then N(x) ~- N(x) u {~} end for End S E G M E N T Theorem 1. Algorithm S E G M E N T takes constant time on a mesh connected computer. Proof. In algorithm SEGMENT, Step 1 takes constant time. Each of the Steps 2(a) through 2(f) takes constant time and each of them is repeated for (k + 1) times where the value ofk is bounded (2 or 3 depending upon the D N P being computed). Thus the time taken by Step 2 is also bounded and is constant. The time taken by Step 3 is also constant. Thus each of the Steps I - 3 takes constant time and since each of them is repeated nine times (DNP's 1, 2, 4 and 5 being excluded), the algorithm S E G M E N T also takes constant time. Q.E.D. Computation of NPS value alone does not guarantee proper segmentation, since two adjacent parallel planes separated by a jump boundary will also have the same NPS value, and thus will form a single segment. Hence some restriction must be imposed on the algorithm S E G M E N T to prevent such segment merging. For
263
every pixel 'x' both its NPS value [x] and actual range value r(x) are considered for this purpose. For two neighboring pixels 'x' and 'y' if [x] = [y] and It(x) r(y) l > 6, where 6 is some threshold (in this case J is taken as 5), then 'x' and 'y' are boundary pixels of two adjacent parallel segments separated by a jump boundary. The NPS values of both 'x' and 'y' are set to NULL. This shrinks the boundaries of such segments by one pixel on all sides, but avoids the problem of segment merging. S. RESULTS AND DISCUSSIONS
Experimental results are shown for a synthetic range image of a sphere (Fig. 9a) and three real range images referred to as 'Poly' (Fig. 10a), 'Quad' (Fig. 11a)and 'Ring' (Fig. 12a). The images of Poly and Quad are of size 256 x 256 and are obtained from PRIMULA range data base of Osaka University, Japan. The other image Ring is of size 128 x 128 and is obtained from Environmental Research Institute of Michigan (ERIM), USA. The image Sphere contains only quantization noise, whereas the SNR value in P01Y and Quad is around 20dB and in Ring the SNR value is around 30dB. The results are shown with a shift of 3. The value of threshold was calculated accordingly. For the synthetic range image of sphere, the sets of points quantized to different digital neighborhood planes are shown in Fig. 9(b)-(j). In each of these figures the region of points belonging to a D N P is indicated by the closed boundary inside an outer circle, where the outer circle represents the boundary of the hemisphere. The positions of these regions with respect to the centre of the hemisphere give an indication of the range of the surface normals which are quantized. to the corresponding DNP. For example, Fig. 9(b) indicates that all the surface points where the directions of the surface normals are closed to the positive direction of the Z-axis are quantized to the D N P 3. Closed boundaries in Fig. 9(k) represent the regions of points with the same NPS value. Each such region is marked with an integer. The NPS values for these marked regions are listed in Table 1. Though region 1 is background, the points in it assume the NPS value of 3. Points belonging to region 2 are the points at which the directions of the surface normals are nearly perpendicular to the Z-axis. These points do not assume any NPS value because the DNP's 1, 2, 4 and 5 whose normals are perpendicular to Z-axis cannot be realized using this technique. For the image Poly the segmentation result without smoothing is shown in Fig. 10(b). Fig. 10(c) and (d) illustrate the segmented output when smoothing is carried out over a neighborhood of 3 x 3 and 5 x 5, respectively. It can be seen that without the smoothing operation the segment boundaries are not clearly defined (Fig. 10b). The segmentation results with range images Quad and Ring are illustrated in Figs 1 l(b) and 12(b) where the smoothing is carried out over a neighborhood of 5x5.
264
P.K. BISWAS et al.
(a)
(b) DNP 3
(c) DNP 6
(d) DNP 7
(e) DNP 8
(f) DNP 9
(g) DNP 10
(h) DNP 11
(i) DNP 12
(j) DNP 13
Fig. 9. (a) Synthetic ramge image of sphere. (b)-(j) Sets of surface points on synthetic sphere quantized to (k) different DNP's.
An SIMD algorithm for range image segmentation
265
Table 1. Regions of uniform NPS value in a hemispherical surface
Fig. 9. (Continued.)
(a) rt _
~'-' ~ t
-
-!
Region no.
NPS value
Region no.
NPS value
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
Background Nil
19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35
{3} {3,9} {3,8} {9,12} {8,11} {3,7,9,12} {3,7,8, 11} {3,7,9} {3,7} {3,7,8} {7,9,12} {7,8,12} {12} ill} {12,7} {7,11} {7}
{6}
{10} {13} {6,10} {6,13} {3,6} {3,6,9,10} {3,6,8,13} {6,9, 10} {6,8,13} {9,10} {8,13} {9} {8} {3,6,9} {3,6,8}
(b)
¢
(c)
CO)
Fig. 10. (a) Range image poly. (b) Segmentation output without smoothing. (c)Segmentation output with smoothing over 3 x 3 neighborhood. (d) Segmentation output with smoothing over 5 x 5 neighborhood.
266
P . K . BISWAS et al.
.bC
~--i
I;
I
I
(
(a)
(b)
Fig. 11. (a) Range image quad. (b) Segmentation output with smoothing over 5 × 5 neighborhood.
• o . . . . . . . . .
. . . . . . . . . . . . . . . . . .
.
. . . . .
. . o
~¢~,~,,,v......st ~tv*,*0ta...,a~,.,_
f.(t" "t
;to I
t
i,"
Z'I ii
. . . . . . .
: P
oN
~,
a
1
I
p
.", IO
•
s,II
JI
@
!
..~,.w,.,._.,.~.,..,..,.:
(a)
. . . . . . . . .
(b)
Fig. 12. (a) Range image Ring. (b) Segmentation output with smoothing over 5 x 5 neighborhood.
9. CONCLUSION
In this paper the problem of segmentation of surfaces in range images is modelled as a quantization problem, where the orientations of the surface normals are quantized to the orientations of normals to some predefined planes (DNP's). Though the algorithm quantizes the surface normals it does not compute the normals at all the points on a surface. Rather the segmentation is achieved by simple shift, subtract and threshold operations. Hence this approach for surface segmentation in range images is computationally more efficient than the existing range image segmentation algorithms(3,6-10.12) which normally attempt either to fit the surface points to a surface function or to calculate a local property (such as normal, curvature etc.). The algorithm described in this paper is parallel and its time complexity is shown to be of 0(1) on a mesh connected computer. The algorithm presented in this paper produces excellent segmentation results with polyhedral ob-
jects (Fig. 10d) and reasonably good results with curved objects (Fig. l i b and 12b). The segmentation produced by our algorithm is comparable to the ones in Mukherjee e t al. (13) Morphological operations of opening and closing are conveniently incorporated in our technique to smooth the regions quantized to different DNP's, thus making our technique robust to noise. Smoothing takes care of the noise present in the range data and the noise arising out of the discreteness in shift vectors, discreteness in the threshold and the tolerance given in the DNP's. Excellent segmentation results show the power of our technique. Acknowledgement--The authors are grateful to the anonymous reviewer for his/her valuable comments on an earlier
version of this paper. REFERENCES
1. P.J. Besl and R.'C. Jain, Three-dimensional object recognition, A C M Comp. Surv. 17(1), 75-145 (1985).
An SIMD algorithm for range image segmentation
2. P.K. Biswas, J. Mukherjee, B.N. Chatterji and P.P. Chakraborty, Qualitative description of three-dimensional scenes, Int. J. Pattern Recognition and Artificial Intelligence 6(4), 65t-672 (1992). 3. D. L. Milgrim and C. M. Bjorklund, Range image processing planar surface extraction, Proc. 5th Int. Conf. Pattern Recognition, pp. 912-919, Miami, Florida (1-4 December 1980). 4. T. C. Henderson, Efficient 3-D object representation for industrial vision systems, IEEE Trans. Pattern Anal. Mach. InteU. 5(6), 609 617 (1983). 5. T. C. Henderson and B. Bhanu, Three-point seed method for the extraction of planar faces from range data, Proc. Workshop Industrial Applications of Machine Vision Research, pp. 181 186, Triangle Park, North Corolina (May 1982). 6. J. H. Han, R. A. Volz and T. H. Mudge, Range image segmentation and surface parameter extraction for 3-D object recognition of industrial parts, Proc. on 1987 IEEE Int. Conf. on Robotics and Automation, 1,380-386 (March 1987). 7. M. Herbert and J. Ponce, A new Method for segmenting 3-D scenes into primitives, Proc. Int. Conf. Pattern Recognition, pp. 836-838, Munich, West Germany, (21-22 October 1982).
267
8. I. K. Sethi and S. N. Jayaramamurthy, Surface classification using characteristic contours, Proc. 7th Int. Conf. Pattern Recognition, pp. 438-440, Montreal, P. Q., Canada (30 July-2 August 1984). 9. P. J. Besl and R. C. Jain, Segmentation through variableorder surface fitting, IEEE Trans. Pattern Anal. Much. Intell. 10(2), 167-192 (1988). 10. P. J. Besl and R. C. Jain, Invariant Surface Characteristics for 3-D object recognition in range images, Computer Vision Graphics and Image Processing, 33, 33-80 (1986). 11. M. Brady, J. Ponce, A. Yuille and H. Asada, Describing surfaces, Comput. Vision Graphics Image Process. 32, 1-28 (1985). 12. Russel, W. Taylor, M. Savini and A. P. Reeves, Fast segmentation of range imagery into planar regions, Computer Vision Graphics and Image Process. 45, 42-60 (1989). 13. Jayanta Mukherjee, P. P. Das and B. N. Chatterji, Segmentation of Range Images, Pattern Recognition 25(10), 1141-1156 (1992). 14. P. K. Biswas, S. S. Biswas, A. K. Roy and B. N. Chatterji, Segmentation of range images using morphological operators, Proc. IEEE Region 10 Intl. Conf. Energy, Computer, Communication and Control Systems, pp. 338343, New Delhi (28 30 August 1991).
About the Author--PRABIR KR. BISWAS completed his B.Tech(Hons), M.Tech and Ph.D. in electronics and electrical communication eng. from Indian Institute of Technology, Kharagpur, in the years 1985, 1988 and 1991, respectively. From 1985 to 1987 he was with Bharat Electronics Ltd., Ghaziabad, India. Currently he is an Assistant Professor in the Dept. of Electronics and Electrical Communication Eng., Indian Institute of Technology, Kharagpur. His areas of interest include pattern recognition, computer vision, parallel processing, distributed control and computer network.
About the Author--S. S. BISWAS completed his B.Tech (hons), M. Tech and Ph.D. in electronics and electrical communication engineering from Indian Institute of Technology, Kharagpur in the years, 1985, 1986 and 1993, respectively. Currently he is a Software Engineer in Tata Consultancy Services Ltd., India. His areas of interest include pattern recognition, computer vision and systems engineering.
About the Author--B. N. CHATTERJI completed his B.Tech(hons) and Ph.D. in electronics and electrical communication eng. from Indian Institute of Technology in the years 1965 and 1970. He was with Telerad Pvt. Ltd. Bombay in 1965 and Central Electronics Engineering Research Institute, Pilani in 1966. Since 1967 he has worked as a faculty member in the Dept. of Electronics and Electrical Communication Eng., Indian Institute of Technology, Kharagpur where he has been Professor since 1980 and was Head of the Department during 1987-1991. He is a fellow of several professional societies like the Institution of Engineers (India), the Institution of Electronics and Telecommunication Engineers (India) etc. Dr Chatterji has published extensively in journals and is the author of two text books. He is the recipient of Hari Om Ashram Prerit Dr Vikram Sarabhai Award for the year 1983 for his contributions in the area of Electronics and Telecommunications. His areas of interest are pattern recognition, image processing, computer vision, parallel processing, system identification and digital signal processing.