Fast Implementation of Level Set Method andIts Real-time Applications 
Yumi Iwashita, Ryo Kurazume,
Tokuo Tsuji and Tsutomu Hasegawa
Department of Intelligent Systems
Kyushu University
Fukuoka City, JAPANyumi,kurazume,tsuji,hasegawa@irvs.is.kyushu-u.ac.jp
Abstract – The level set method (LSM) has been widelyused for various applications such as motion tracking and
3D geometrical modeling. However, the calculation cost ofreinitialization and updating of an implicit function is considerably expensive as compared with conventional active
contour models such as ”Snakes”. To tackle this problem,
we propose an ef cient algorithm of the LSM named the Fast
Level Set Method(FLSM). This paper introduces some experiments based on the FLSM, including 2D real-time trackingof moving objects in video images, and 3D simultaneous motion capture system of multiple targets using stereo range images.
Keywords: Active contour, Snakes, Deformable surfaces,
Level Set, Motion capture, Tracking.

1 Introduction
The level set method(LSM) proposed by S. Osher and J.
A. Sethian[7] has attracted much attention as a method thatrealizes a topology free active contour modeling by detecting and tracking a boundary. Various applications have been
proposed in the past including motion tracking[5] and 3Dmodeling[8][9].
However, the calculation cost of reinitialization and updating of the implicit function is considerably expensive ascompared with conventional active contour models such as”Snakes”[2]. Thus, this has been considered as a heavy burden for applying the level set method in real-time applications such as motion tracking in video images. To overcome
this problem, this paper presents a new high speed execution technique of the LSM named the Fast Level Set Method(FLSM), which makes it possible to apply the LSM for realtime applications. This paper introduces two examples of the
FLSM applications: 2D real-time tracking of moving objects
in video images and 3D real-time tracking of multiple humanbodies using stereo range images.

­

 0-7803-8566-7/04/$20.00 c 2004 IEEE.

Kenji Hara
Department of Visual Communication Design
Kyushu University
Fukuoka City, JAPANhara@design.kyushu-u.ac.jp

2 The Level Set Method
The LSM utilizes an implicit function which is definedin a space one dimensional higher than that of where a contour of interest is described. This function , which is defined as a distance function from a current contour in general, is updated according to a next PDE (Partial Differential
Equation).

    (1)where,  is a local curvature of , and is a speed function.
The contour to be tracked is detected as the cells with a valueof zero of the implicit function (zero level set), that is, thecontour line of  . In the implementation of the LSM,
the space is uniformly split by cells, and Eq.(1) is solved iteratively using numerical schemes such as the upwind scheme.
To solve Eq.(1), the speed function  has to be determined at each cell for every update process of . The distribution of the , which is known as the extension velocity
field[1], is constructed as follows: i) at the current zero levelset cell, is calculated according to the intensity of the current image at first; ii) next, at each cell except the zero levelset cell, the speed function is copied from the nearest zerolevel set cell. However, finding the nearest zero level set cellneeds large calculation cost.
In addition, since integral errors are accumulated during
the calculation of updating, the reinitialization process is required in which the proper quantity (the distance from thecurrent zero level set cell) at each cell is re-calculated periodically after several update processes. The constructed fieldis called the distance field. However, this reinitialization process is time consuming, because the nearest zero level set has
to be determined at each cell. Therefore, this is also a major
obstacle for the high speed implementation of the LSM.
To overcome these problems, several techniques have been
proposed in the past, such as the Narrow Band Method, theFast Marching Method, SFA (Sparse Field Algorithm), andHERMES. The most popular and efficient method is the one

6302

proposed by Adalsteinsson and Sethian[1] in 1999, which
combines the Narrow Band Method and the Fast Marching
Method. The calculation cost of this method is  ¾  (n is the number of cells along a side of the voxel space in
3D space).

3 Fast Level Set Method
In this section, we propose a new efficient algorithm
named the Fast Level Set Method (FLSM) with the calculation cost of ¾ . The key idea of the proposed FLSM is
the use of a reference map (Figure 1(a)). This map indicatesclassification in which each cell is categorized according toa distance from a center cell. For example, the class   consists of cells which are located  away from the center cell.
Rr
10

9

10

8

5

4

5

8

10

5

2

1

2

5

10

9

4

1

0

1

4

9

10

5

2

1

2

5

10

8

5

4

5

8

10

9

10

F1

F1

F2
F3

F1

F1
F2

F1

F1

F3

F1

F1

F7

F4

F1

F1

F3

F1

F1

F6
F7

F1

F6

F1

F5

F7
F6

F2

F5

F6

F5
F2

F6

F5

F3

F5

F4
F1

F1

F1

F3
F7

F1

F7
F5

zero level set(a) Refence map(1) Rr = 10

constructed consequently(Figure 1(c)).
In this process, the speed function is registered in orderof a distance from the zero level set. Thus, by overwriting
the distance stored in the reference map at the same time,the distance field can also be created. Additional cost by thisreinitialization process is only memory access time, and total calculation cost is the almost same as the one without thereinitialization process. Therefore, the reinitialization process can be performed frequently, even up to every updateprocess, with trivial calculation cost.
Moreover, the method can be improved by using a divided
reference map(Figure 1(b)). For example, let us consider azero level set cell  . If a left cell  ¼ adjacent to the cell
 is also a zero level set cell, there must be no cells which
are nearer to  than  ¼ in the left side region of the cell (A,B, and H). Therefore, by skipping the registration process
in these regions, it is possible to reduce the number of useless overwriting and execute the construction process of theextension velocity field more efficiently.
The comparison of the calculation cost in 3D (    )is shown in Table 1. The cost of the conventional level set
method using the fast marching method[1] is  ¾  .
On the other hand, the calculation cost of our method is¾ . Moreover, the additional cost of the reinitialization
process is trivial in our method.

(2) Rr = 10

Table 1: Comparison of calculation costs. Æ   is thewidth of the narrow band.
B C D

F2

B B C D D

F3

B B B0 C0 D0 D D

F1A A A0 Z E 0 E E
H H H0 G0 F 0 F F

F1H H G F F

F4

H G F

F1

F1

F1

F3

F1

F5

F7
F2

F3

F1

F5

F7
F6

F2

F5

F6

F1

F6
F3

F1

F1

(3) Rr = 9

F5

F2

F2

F1

F5

F5

F6

F2

F2

F2

F1

F5

F5

F5

F1

F3

F2

F2

F1

F5

F5

F6

F1

F3

F3

F2

F1

F5

F6

F6

F1

F4

F3

F3

F1

F6

F6

F7

F7

F4

F4

F3

F1

F6

F7

F7

F4

F4

F4

F1

F7

F7

F7

LSM
Construction ofvelocity field
Reinitialize
Updating process
Detection ofzero level set(4) Final result(b) Divided reference map (c) Procedure of constructing extension velocity field

Figure 1: Reference map and the construction process of theextension velocity field.

The extension velocity field is constructed efficiently using the reference map. Here, we assume that the speed function at the zero level set cell has already been determined.
At first, one zero level set cell is chosen. Then, using theclass Æ´Æ·½µ in the reference map, all cells which are located Æ Æ   away from this zero level set cell are selected, and the speed function at the zero level set cell is registered to these cells tentatively. This procedure is repeated
for all the zero level set cells. Next, the same procedure is
performed using the next class  Æ´Æ·½µ ½ . This process isrepeated until all the classes in the reference map are investigated. At the end of the registration processes for all theclasses, the speed function of the nearest zero level set cellis registered at each cell, and the extension velocity field is



FLSM

Æ¾   Æ¿ ¾ Æ¾  Æ¾ Æ¾ Æ¾ Æ¾ Figure 2 shows the 3D surface construction process of awired basket. The range data consists of 5530 points. Comparison of calculation time is shown in Figure 3 for two conventional level set methods, that is, the no extension and extension velocity fields[1], and the proposed FLSM. We performed several experiments for various 3D models with complex shapes. These results suggest that the FLSM is executed
two or three times faster that the methods proposed so far.

4 2D real-time tracking of moving objects in video images
We developed a real-time contour detection and tracking
system of moving objects in video images using the FLSM.
In this system, the size of the image is    pixels,
and the speed of image acquisition is  frames/second. The

6303

Figure 4 shows the results of detecting and tracking multiple objects in real scene. This shows that boundary detection
and real-time tracking of moving objects which move up to pixels/second are successfully realized. Though the detected contours are separated at first, they are merged to asingle contour when one object occludes the other. Whenthe objects split, the contours separate again.
initial surfacet = 0.11 s (4 updating)t = 2.93 s

t = 3.70 s

t = 4.53 s

t = 4.76 s

t = 3.45 s (16 updating)t = 0.73 s (9 updating)Figure 2: The process of 3D surface reconstruction.

400
350

Time [s.]

300

345.7

FLSM
LSM ( No extension velocity field )LSM ( Extension velocity field )

258.3

250

Figure 4: Simultaneous tracking of multiple objects.

200
150

90.5

100
50

1.6 3.4 4.5

12.3

28.8 41.5

5 3D real-time tracking of multiple
human bodies

0
64x64x64

128x128x128
Resolution [grids]

256x256x256

Figure 3: The comparison of the processing time (PentiumIV, GHz).
process of the FLSM is calculated at  frames/second using
Pentium IV  GHz processor.
The speed function is determined as  (2)where  is the local curvature of the implicit function  , and and  are constant values. is determinedby intensity     and its subtraction of the background image             as

         
The maximum value of the speed function,
1 update in the experiments.

(3), was cell per

We also implemented the proposed method for the
3D tracking of human bodies. The multi-view stereo
techniques[6] and the visual hull[3][4] are well-known 3D
all-round reconstruction methods of a real object using multiple cameras. However, when two or more objects are in
the scene at once, the reconstruction process of these objects
separately is quite difficult because of the mutual occlusion.
In addition, for the visual hull technique, to extract a silhouette of the object robustly is quite laborious, and reconstruction of a non-convex object is also impossible.
This paper introduces a prototype model of a 3D motion
tracking system of multiple targets using the FLSM with twostereo cameras. Characteristics of this system are as follows: i) By maintaining shapes of objects obtained previously, closed and smooth surfaces are given even if parts ofthe depth images cannot be captured temporarily. Therefore,
though the mutual occlusion occurs when multiple objects
are in the scene, we can reconstruct each object separately.
Thus, the proposed system is considered to be robust for mutual occlusion; ii) by using the FLSM, complex models with
arbitrary topology can be reconstructed; iii) by restricting the
maximum value of the speed function and updating of an
implicit function frequently, the proposed system is robust

6304

Stereocamera 2

detected
contour

2

2
2.3m

1.2m

Centervoxel

Objects

1

1

Stereocamera 1(b)(a)
stopping
regionOccludingregions

Objects

stopping
regiond2

zero level set A

D

Objects

Next, we apply the FLSM for the region-S. Here, we define voxels of the region-S which are adjacent to voxels voted
OUT as stopping region, and the rest of the region-S as inner
region(Fig.5(c)). A closed and nonintersecting initial surface
including the voxel space is set at first, and it is shrunk andsplit until the front reaches to the stopping region.
If more than two objects are in the scene, mutual occlusion might occur between them(Fig.5(b)). In addition, apart of the body may be occluded by his hands or legs. In
these cases, a connected 3D model shown as a dotted line
in Fig.5(b) is reconstructed along to the occluding boundary.
To avoid the production of connected objects and to reconstruct each models separately, a quite small value is stored
to the speed function of the cells which are considered to belocated in the occluding regions. To judge whether a cell isin the occluding regions or not, we use a significantly simple method in this approaches. Here,  and  are distancesfrom inside zero level set cells at previous time to the outside
zero level set cells as shown in Fig.5(d).  is the sum of and  , and  is the distance radio defined as the following
equation;

 

d1

The procedure for determining

inner region

inner region(c)   


(d)

(4)

 is as follows:

1. When zero level set voxels are in the stopping region, ;

Figure 5: Occluding boundary.

2. When zero level set voxels are in the inner region and  ;(a)    , thenfor spike noise within stereo range data; iv) the calibration
process of mutual camera position is quite easy by utilizing
iterative registration algorithms such as ICP algorithm forobtained 3D points.

(b)    , then

3. When zero level set voxels are in the inner region and  ,   ½ ;5.1 Tracking of multiple targets using the
FLSM
Figure 5(a) shows the camera and workspace configuration of our system. Two stereo cameras(Point Grey Research, from the cenBumblebee) are placed at a distance ofter of the workspace. The workspace is divided into small
voxels (hereafter denoted the voxel space). After capturing
depth images using stereo cameras, multiple depth images
are combined in the voxel space.
Each voxel in the voxel space can hold one of two status,IN(inner human body) or OUT(outer human body). At first,
IN is assigned to all voxels. Then, each pixel (x,y) of a depth
image is backprojected onto the voxel space, and OUT is
voted until the distance from the camera is equal or smaller
than the depth value of (x,y). This process is repeated for all
depth images, and a region entirely composed of voted IN is
extracted, named a region-S .

 ¼
 ½

4. When zero level set voxels are not in the stopping andinner region,   ¾ .
Here, ¼  ½  ¾ , and we used    ,    inexperiments. By this procedure, the position of the previous
zero level set cells are maintained and two separate region
can be extracted even if mutual occlusion occurs.

5.2 Experiments
The experiments of detecting and tracking multiple humanbodies using the proposed system are carried out. In the experiment, the size of the depth images is    pixels,
and the FLSM calculation is done by the Pentium IV PC, GHz.
Figure 6 shows an example of the 3D modeling of a humanbody. We examined two kinds of resolutions of the voxel
space, that is,      voxels and      voxels.

6305

Stereo
Camera 2

StereoCamera 2

Box
Box
1

1
2

Box
1
2

StereoCamera 1(a) Object(b) Voxel data
(200x200x200)

2

StereoCamera 1(c) Texture mapping(a) Objects

Figure 6: 3D reconstruction of human body.

The processing cycle is  Hz(     voxels), and  Hz(      voxels), respectively.
Figure 7 shows an example of the 3D motion tracking of
two human bodies with a box. The resolutions of the voxel
space is  voxels. Though the detected surfaces are
separated at first, they are merged to a single surface when
passing the box. After the box is passed, the surfaces separate again.
Figure 8 shows the examples when 2 people walk acrossa scene. Though the mutual occlusion between them occurs
when one of them is occluded by the other, separate shapes
are detected and kept tracking. These examples suggest thatthe proposed system is capable of detecting and tracking
multiple human bodies in a scene individually even if a partof them is occluded.

6 Conclusions
We have proposed a new high speed execution technique of the level set method named the Fast Level Set
Method(FLSM). The efficiency of the FLSM was verified
through 2D real-time tracking of moving objects in video
images and 3D real-time tracking of multiple human bodies in stereo range images. Our future work includes theimplementation of the 3D real-time tracking system on a PCcluster, and the development of the high speed and high quality 3D real-time tracking system.

Acknowledgment
This research was partly supported by ”Ministry of Public Management, Home Affairs, Posts and Telecommunications” and Grants-in-Aid for Scientific Research(16700185).

(b) Top view(c) Cross view

Figure 7: 3D reconstruction of multiple human bodies witha box.

References[1] D. Adalsteinsson and J. A. Sethian: The fast construction of extension velocities in level set methods,
J.
Computational Physics, 148:2–22, 1999.
[2] M. Kass, A. Witkin and D. Terzopoulos: Snakes, active
contour models,Int. J. Computer Vision, 1(4):321–
331, 1988.
[3] A. Laurentini: The visual hull concept for silhouettebased image understanding,IEEE Trans. on PAMI,
16(2):150–162, 1994.
[4] W. N. Martin and J. K. Aggarwal: Volumetric description of objects from multiple views,
IEEE Trans. on
PAMI, 5(2):150–158, 1983
[5] N. Paragios and R. Deriche: Geodesic active contotraurs
and level sets for detection and tracking moving objects,
IEEE Trans. on PAMI, 22:266–280, 2000
[6] H. Saito, S. Baba, M. Kimura, S. Vedula and T.
Kanade: Appearance-based virtual view generation of
temporally-varying events from multi-camera images inthe 3d Room, Proc. 3DIM99,516–525, 1999.
[7] J. A. Sethian: Level Set Method and Fast Marching
Methods, second ed.
Cambridge University Press,
UK,1999

6306

2
1

2

2

2

1

1

2

1

1

StereoCamera 2
2

2

2

1

1

1

2

2

1

1

StereoCamera 1
StereoCamera 2
2

2

2
2
1

StereoCamera 1

(a)

1

1

(b)

(c)

(d)

Figure 8: 3D reconstruction of multiple human bodies.

[8] S. Yui, K. Hara, H. Zha and T. Hasegawa: A fast narrow
band method and its application in topology-adaptive 3D modeling, Proc. ICPR02, 4:122–125, Aug, 2002
[9] H. Zhao, S. Osher and R. Fedkiw: Fast surface reconstruction using the level set method,
Proc. 1stIEEE Workshop on Variational and Level Set Methods
in Computer Vision, 2001.

6307

2

1

1

(e)