Cristian Sminchisescu
email: cristian.sminchisescu@inrialpes.frdnm@graphics.cis.upenn.edusven@cs.toronto.edu

Dimitris Metaxas

Sven Dickinson

Improving the Scope of Deformable Model Shape and Motion Estimation

Keywords: Deformable models, geometric constraints, object tracking, model grouping, bundle adjustement

Previous approaches to deformable model shape estimation and tracking have assumed a fixed class of shapes representation (e.g., deformable superquadrics), initialized prior to tracking. Since the shape coverage of the model is fixed, such approaches do not directly accommodate incremental representation discovery during tracking. As a result, model shape coverage is decoupled from tracking, thereby limiting both processes in terms of scope and robustness. We present a novel deformable model framework that accommodates the incremental incorporation during tracking of new geometric primitives (lines, in addition to points) that are not explicitly captured in the initial deformable model but that are moving consistently with its image motion. As these new features are detected via consistency checks, they are added to the model, providing incremental soft constraints on the estimation of its rigid parameters. The consistency checks are based on trilinear relationships between geometric primitives. Consequently, we not only increase both model scope and, ultimately, its higher-level shape coverage, but improve tracking robustness and accuracy, by directly employing the new features in both forward prediction and reconstruction. Our new formulation is a step towards automating model shape estimation and tracking, since it requires significantly reduced initial model hand-crafting. We demonstrate our approach on two separate image-based tracking domains, each involving complex 3D object shape and motion.

Introduction

Approach

Deformable models offer a powerful, data-driven framework for recovering the shape and motion of an object through numerical integration of a dynamical system that encodes the fitting error. The process is a physics-based interpretation To appear IEEE CVPR 2001 of an optimization process, in which the initial model shape and pose are supplied by an initialization process that aligns the model with the object of interest in the initial image. The model subsequently moves under the action of forces derived from image measurements which, in turn, are converted (through a Jacobian matrix) into generalized forces that act on the model's rigid and non-rigid degrees of freedom.

The quest for robustness and accuracy during model estimation naturally leads to the search for geometrically and photometrically derived image forces (cues), for model parameterizations, and for ways to combine them in a consistent manner. Although previous approaches have investigated the use of particular combinations of cues in the context of deformable models, they assumed a fixed and known initial model shape representation.

In this paper, we present a novel two-step technique in which new model structure is revealed during tracking using geometric consistency checks. These checks combine heterogeneous information, involving an initial model's rigid parameters as well as independently tracked line features in the image. Our goal is to enhance the higher level model shape representation with additional deformable shapes, and enrich its basic discretization power through the incremental addition of line features.

In our approach, line features (not belonging to the initial deformable model) are tracked in an image sequence. Those features whose motion is consistent with the model are reconstructed, integrated into the initial (high-level) model representation, and subsequently used to increase tracking robustness. The technique consequently attempts to bridge the gap between model-based, top-down estimation techniques and classical bottom-up, feature-based reconstruction techniques, by relaxing some constraints on each side. Complete model representation is no longer fully known a-priori, while knowledge of feature motion is acquired during reconstruction, once their identity as being part of the model is established.O u r final model is a mixture of lines, points and parameterized shapes, and generalizes the shape coverage of previouslyused deformable models. We show how a deformable model framework naturally accommodates rigid and non-rigid highlevel shape representations and their associated constraints, and various model discretizations such as points and lines, by combining the simplicity of linear reconstruction methods with the refinement and bias removal of non-linear estimation methods.

We present experiments involving complex objects and motion in monocular video sequences and show that the proposed method is able to recover model structure accurately and efficiently. We also quantitatively show how the novel recovered structure significantly improves the accuracy and speed of the deformable models estimation process, providing important additional constraints especially during difficult to track, towards camera motion.

Relation to Previous Work

In the area of physics-based deformable models, various formulations have been proposed ([25, 13, 17, 16]). As powerful as these techniques are, they typically assume the model representation is fixed and known a-priori, sometimes imposing a heavy burden on the model initialization/recovery process ( [START_REF] Dickinson | Shape recovery using distributed aspect matching[END_REF][START_REF] Dickinson | Integrating qualitative and quantitative shape recovery[END_REF]). Furthermore, a representational "gap" exists between the coarse, parametric shapes used to model the objects and the point-based discretizations used to bind them to the image; no obvious method exists to bridge this gap through the inclusion of other basic geometric primitives, e.g., lines. On the other hand, it has been recognized ( [START_REF] Seitz | Implicit Scene Reconstruction from Probability Density Functions[END_REF]) that the rather different parameterizations of various features encountered in different approaches or applications usually lead to difficult problems when attempting to integrate those features within a particular representation or optimization procedure.

In the area of constraint integration, particularly in a deformable model setting, contours and stereo ( [START_REF] Terzopoulos | Constraints on deformable models: Recovering 3-D shape and non-rigid motion[END_REF]), shading and stereo ( [START_REF] Fua | Taking Advantage of Image-Based and Geometry-based constraints to recover 3-D surfaces[END_REF]), contours and optical flow ( [START_REF] Decarlo | The integration of optical flow and deformable models with applications to human face shape and motion estimation[END_REF][START_REF] Decarlo | Combining Information in Deformable Models Using Hard Constraints[END_REF]), and shading ( [START_REF] Samaras | Incorporating illumination constraints in deformable models[END_REF]) have been introduced, with good results. Beyond the particular choice of sources of information to use, the approaches differ in the way they fuse their component information. Some combine the information in a symmetric manner and weight them statistically (soft constraints). Others favor a particular hierarchical constraint satisfaction order with an exact policy, such that inconsistent contributions to the solution from constraints further down in the hierarchy are pruned away by constraints higher-up (hard constraints).

Finally, in the area of feature-based rigid motion and structure estimation, various approaches based on the type of correspondences (2-D to 2-D, 2-D to 3-D, or 3-D to 3-D) as well as the type of features (lines, points, or corners), have also been proposed (see [10] for a review). Several algorithms, commonly based on Kalman filtering, are also known ( [START_REF] Crowley | Measurements and integration of 3-D structures by tracking edge lines[END_REF][START_REF] Vieville | Feed-forward recovery of motion and structure from a sequence of 2-D line matches[END_REF]). Approaches to simultaneous structure and motion estimation, based on the inversion of the so-called forward model, have been presented. Many of these methods, such as [START_REF] Szeliski | Recovery 3-D shape and motion from image streams using non-linear least-squares[END_REF][START_REF] Azarbayejani | Recursive Estimation of motion, structure and focal length[END_REF][START_REF] Taylor | Structure and Motion from Line Segments in Multiple Images[END_REF], could be related to the classical relative orientation paper ( [START_REF] Horn | Relative Orientation[END_REF]).

The approach we present here is based on a deformable model formulation similar to ( [START_REF] Terzopoulos | Dynamic 3D models with local and global deformations:Deformable superquadrics[END_REF]), in which we employ an articulated chain of deformable superquadrics with global deformations, such as tapering bending, shearing and twisting. For broader modeling coverage, additional local deformation functions could also be used. As in the standard formulation, we employ local contour-based forces derived from image gradients for basic model motion estimation. We diverge from previous approaches in that we assume our initial model is incomplete, and attempt to recover additional model structure through the use of geometric consistency checks embedded in a tracking framework. Specifically, we integrate constraints in the form of new line features discovered in the image as moving consistently with the model, with constraints in the form of point features derived from a discretization of the model. The line feature consistency checks we use are based on the conditions present in separate, bottom-up structure and motion estimation under 2-D to 2-D line correspondences in the Euclidean calibrated case ([10, 15, 28]). A later body of work for calibrated, uncalibrated, and projective reconstruction cases ([22, 21, 8] addresses trilinear constraints between points and lines in 3 views, and presents linear, bottom-up, and model-free rigid reconstruction methods (see [START_REF] Hartley | Multiple View Geometry in Computer Vision[END_REF] for a comprehensive review).

However, unlike these approaches, we assume an incomplete deformable and Euclidean model and we don't attempt to solve for the rigid parameters. On the contrary, by knowing the estimated rigid parameters of the model, we only attempt to check whether the consistency conditions are indeed satisfied for arbitrarily independent tracked lines in at least three frames. Once consistent lines have been identified, they are reconstructed, and the model representation is enhanced to include them. Their prediction in subsequent frames under the action of the Euclidean group is derived, while their structure is continuously re-estimated together with the deformable structure and motion of the other model parts. Their contribution to the solution for the model's estimated rigid parameters is linearized by means of their corresponding Jacobian matrix and integrated as an additional soft constraint (see section 3.5).

Deformable Model Formulation

In the next two sections, we shall briefly review the deformable model formulation in terms of its underlying geometry and dynamics (see [START_REF] Metaxas | Physics-Based Deformable Models[END_REF] for details).

Model Geometry

The reference shape of the model is defined over a domain Ω as p = T(q d , u),whereT defines a global deformation based on the parameters, q d and u ∈ Ω. Following [START_REF] Solina | Recovery of Parametric models from range images: the case of superquadrics with local and global deformations[END_REF], we employ a deformable superquadric ellipsoid, having global tapering and bending deformations, as a reference shape. The posi-tion of a point on the model is expressed with respect to a world coordinate system Φ as x = D(q r , p),w h e r eD represents a rigid displacement parameterized by q r . Finally, both rigid and non-rigid parameters are assembled in a vector, q =(q r , q d ), sometimes referred to as the model's generalized coordinates.

Dynamics and Generalized Forces

The velocity of points on the model is given by:

ẋ(u)=L(q, u) q (1)
where L is the model Jacobian matrix [START_REF] Terzopoulos | Dynamic 3D models with local and global deformations:Deformable superquadrics[END_REF]. By omitting the mass term (since we want a model free of inertia), the Lagrangian dynamics equations become [START_REF] Terzopoulos | Dynamic 3D models with local and global deformations:Deformable superquadrics[END_REF]:

Dq + Kq = f q , f q = L ⊤ fdu (2)
where D = γL ⊤ L,andK = diag(k si ), with k si being the stiffness associated with the global parameter i.

When the estimation of the model is based on image contours, the above formulation has to be extended to include: 1) the perspective projection Jacobian, L p , computed in terms of the model Jacobian and the perspective transform Jacobian evaluated at points on the model discretization, and 2) the definition of local forces f p derived from image potentials (see [START_REF] Metaxas | Physics-Based Deformable Models[END_REF] for derivations). Consequently, the generalized forces with respect to edges are defined as:

f q e = L ⊤ p f p du (3)

Line Feature Formulation

The formulation of the tracking process in the previous sections is based on a physics-based deformable framework, where the model moves under the action of local forces derived from image potentials. In this section, we extend that formulation by integrating new line features into the modelline features that are initially not part of the model, but for which evidence exists in the image. We begin by tracking a minimal model in a sequence of images and attempt, incrementally over time, to: 1) identify line features moving consistently with the model, and 2) augment the model with these consistent features to improve its tracking. At each iteration, the model not only increases in scope but covers more of the image features, resulting in more robust tracking of the object. Our approach consequently attempts to relate two processes at the image and model levels, respectively. At the image level, we employ purely image-based techniques to track lines in a sequence of frames. We call such lines imagetracked lines (ITL). At the model level, we use contour based gradient forces to estimate the model (rigid and non-rigid) parameters. First, we decide whether an ITL represents a line belonging to the object (but not present in its model) by means of two geometric consistency checks, derived from tracking the ITL in at least three successive frames. The lines that pass this test are called consistent image-tracked lines (CITL). Second, for CITLs detected in the previous step, we attempt to robustly recover their structure in terms of an underlying parameterization in a model-centered frame. Third, we predict how a CITL will appear in a subsequent image based on the current estimate of the model's rigid motion. We call this prediction a model-predicted line (MPL). Finally, we use the error between a CITL and a MPL to define additional image alignment forces. They are subsequently linearized by means of their corresponding Jacobian matrices and combined into the model estimation procedure as soft constraints. The initially reconstructed line representation is subsequently re-estimated, together with the entire model structure and motion, within the deformable model framework.

Line parameterization

In the following formulation, we denote a line in 3-D by lower-case letters, l i (i =1 ..n), and its corresponding projected line (or segment) in the image plane by capital letters, L i (i =1 ..n). A 3-D line is parameterized by a unit vector v, representing one of two possible directions on the 3-D line l, and a vector d terminating on l and perpendicular to it (see Fig. 1a). This line representation forms a 6-dimensional parameter space with 4 degrees of freedom. Consequently, the underlying relation between the line primitives can be identified as a 4-dimensional manifold embedded in the abstract parameter space, and any line can be identified with (actually two) points on this manifold ( [START_REF] Kanantani | Statistical optimization for geometric computation: Theory and Practice[END_REF][START_REF] Taylor | Structure and Motion from Line Segments in Multiple Images[END_REF]). The 3-D line, l, and the optical center of the camera determine a plane (line's interpretation plane) with normal, N = (N x ,N y ,N z ). This plane intersects the image plane, defined by the equation z = f (f being the camera focal length) at line L. The equation of line L in the image plane can be consequently written as:

N x X + N y Y + N z f =0 (4)
The above relation, used in an inverse manner, yields an immediate representation for the normal of the plane containing a 3-D line in space by knowing the image plane equation of its projected (observed) line, or segment. In vector terms, given a plane, P =( N x ,N y ,N z , 0) = (N, 0),a n d any point belonging to the interpretation plane of a given line, X =( X, Y, Z, 1), we can write the plane equation as P ⊤ X = 0.

Model-Based Consistency Checks

The standard formulations for recovering motion and structure using 2-D to 2-D line correspondences (e.g., [10,[START_REF] Yen | Determining 3-D Motion and Structure of a Rigid Body Using Straight line Correspondences[END_REF][START_REF] Mitiche | Line-based computation of structure and motion using angular invariance[END_REF]) rely on three frames and at least six line correspondences (although no formal proof is yet available, [10]) for uniquely recovering the structure and motion of a rigid object (within a scale factor for translation and structure parameters). Furthermore, for the two view case, the resulting system of equations is underconstrained, so one would expect infinitely many solutions.

Consider now the motion of a 3-D line, l, in three successive frames (l i ,i =1 , 2, 3, with direction support v i ,i = 1, 2, 3) and it's corresponding image projected lines (L i ,i = 1, 2, 3). The motion between frames 1 and 2 is described by the translation and rotation, t 12 and R 12 , and for frames 1 and 3,b yt 13 and R 13 , respectively. The corresponding normals for the interpretation planes, P 1 ,P 2 ,P 3 , determined by the line L and the center of projection in the three frames, are N 1 , N 2 , N 3 , respectively (see Fig. 1b).

The following two relations can be derived either geometrically or algebraically, from the quantities presented above ([10, 15, 28]):

N 1 • (R -1 12 N 2 × R -1 13 N 3 )=0 (5)
-t 12 • (R 12 N 1 )= N 2 × R 12 N 1 N 2 × R -1 23 N 3 • R -1 23 t 23 • R -1 23 N 3 (6)
In this model-based formulation, we are not interested in solving for the rotation and translation but, given a model with known motion and some independent ITLs, checking whether those lines are moving consistently with our model (CITLs). More specifically, given an ITL in three frames (that is, knowing N 1 , N 2 and N 3 ) as well as the motion of the model (that is, R 12 , t 12 and R 13 , t 13 ), we use equations ( 5) and ( 6) to check whether the motion of the line is consistent with the model's motion. If so, we hypothesize that the line is part of the object and therefore should be added to the model. One can prove that the above two relations are necessary and sufficient conditions for consistency. In practice, the above two relations are never satisfied exactly, due to slight errors in the image tracked lines and errors in the estimated model rigid parameters, so a threshold must be chosen, based on the covariance of both model rigid parameters estimation and the image tracked lines, respectively.

Robust Unbiased Model-Based Structure Recovery

Once a moving image line has been assigned to the model through the above consistency check, the next step is to recover its structure, i.e., the vector pair (v, d) in a modelcentered coordinate system. In order to increase the robustness of the recovery process, we can use as many line correspondences in as many frames (at least two) as are available. The process can be formulated as follows: all interpretation planes for the line correspondences in a camera frame are transformed to a common, model-centered coordinate frame. Each line, l i , having the interpretation plane, P i , is subject to the displacement, D -1 c D -1 i ,where:

D c = R c t c 01
is the displacement corresponding to the camera, and D i is the displacement corresponding to the model (in the world coordinate system) in frame i. Then, the equation of the plane in the object-centered frame is:

P i ⊤ • D -1 c D -1 i
• X = 0.B y stacking together the equations for corresponding lines, we obtain:

A • X =         P ⊤ 1 • D -1 c D -1 1 P ⊤ 2 • D -1 c D -1 1 ... ... ... P ⊤ k • D -1 c D -1 k         • X = 0
Since we are looking for the line intersection of all the above planes, the above [k x 4] matrix A should have rank 2. Any point p on the intersecting line can be written as a linear combination of the singular vectors corresponding to the 2 smallest singular values of the matrix A:

p = a • X s1 + b • X s2 (7)
The corresponding line parameterization can be subsequently recovered as:

v = X s1 -X s2 d =(I - v • v ⊤ v 2 ) • X s1 (8)
As mentioned before, a 3D line has 4 intrinsic degrees of freedom while a projected image line has just 2. Measurements are collected in 3 frames, so this will determine 3x2-4=2 independent relations.

The stability of the reconstruction can be verified in terms of the ratio of the 2nd and 3rd singular values of A (remember that in the noise free case the last two singular values should be zero), being satisfactory when this ratio is high. The above linear method, although robust, might be prone to bias in the initial line parameter estimates, due to fixed rigid displacements. However, we are working in a deformable model framework, in which both model structure and model motion are jointly estimated. Consequently, any such initial linear bias is automatically removed during estimation in subsequent frames.

Forward Model Line Prediction

Once consistent lines (CITL) have been identified, recovered, and effectively added to the model, we use these new lines to improve the tracking of the enhanced model by imposing further constraints on its alignment with the data. The approach we follow is to define alignment forces between a CITL in the image and the projection (MPL) of its corresponding (new) model line, as predicted under the forward motion of the model.

v 2 = R 12 v 1 (9)
d 2 =(R 12 d 1 + t 12 ) -v 2 ((R 12 d 1 + t 12 ) • v 2 ) (10)
N 2 = v 2 × d 2 v 2 × d 2 (11)
As mentioned in the introductory section, N 2 fully identifies the MPL in the second frame. The Jacobian matrix associated with the line parameterization is derived analytically (it is omitted here, due to lack of space). Given a line representation, u l =(v ⊤ , d ⊤ ) ⊤ ,andthe rigid parameterization of the model q r , the Jacobian matrix of the line parameters with respect to the model parameters is a [6x7] matrix:

L u l = ∂u l ∂q r (12)

Image Forces and Soft Constraints

Given an MPL, identified by equation ( 4), and a CITL, identified by its endpoints P 1 =( X 1 ,Y 1 ) and P 2 =( X 2 ,Y 2 ), we define 2-D forces to compensate for their alignment error. Since the CITL and MPL each define an interpretation plane, we define the line alignment forces in terms of the alignment of their corresponding interpretation planes (see Fig. 2b). Given N i , the normal of the CITL interpretation plane, and N m , the normal of the MPL interpretation plane, the corresponding alignment force can be written as:

f l (u l i )=k N i -N m (13)
Notice that a final transformation maps a line representation to an observable "normal" on which the line alignment is actually performed. The corresponding transformation is given by equation [START_REF] Horn | Relative Orientation[END_REF], and involves the computation of another [3x6] Jacobian matrix:

L N = ∂N ∂u l (14)
The Jacobian L l = ∂N ∂qr corresponding to the mapping between the line representation in the model frame and the observation (the normal of the corresponding line interpretation plane) is computed via the chain rule in terms of Jacobians [START_REF] Kanantani | Statistical optimization for geometric computation: Theory and Practice[END_REF] and [START_REF] Metaxas | Physics-Based Deformable Models[END_REF].

The normal alignment forces map the observation error corresponding to the line, l i , to the model's parameter space by means of the corresponding line Jacobian, L l (u l i ).A saresult, the combined generalized forces (due to image gradient and line alignment) acting on the model's degrees of freedom sum over both point and line forces mapped by their corresponding Jacobians:

f q = i L p (u i ) ⊤ f p (u i )+ lj L l (u l j ) ⊤ f l (u l j ) (15)

Experiments

The experiments consist of two sequences, each containing of 4 seconds of video (200 frames recorded at 50 fps) of a moving bike and of a space robotics end-effector grapple fixture. Both sequences involve significant translational and rotational motion in the camera frame. Part of the bike frame is modeled and tracked by means of a model consisting of 3 parts (Fig. 3a), while the grapple fixture consists of only a square parallelepiped (Fig. 4a). Newly recovered models Figure 4: Initial and reconstructed model, with new superquadrics parts fitted are displayed in Fig. 3b and Fig. 4b, in which additional deformable shapes have been fitted based on 3-D reconstructed model lines, their covariance, and their corresponding image contours (an algorithm for such higher-level shape-recovery and grouping is beyond the scope of this paper, and will be described in a subsequent paper).

The models have been manually initialized in the first frame of the sequence, and are displayed in subsequent frames aligned (overlaid and rendered in flat-shaded grey for the bike and brown wireframe for the grapple fixture, respectively) with the part of the object they model (for tracking results, see Fig. 6). Additional lines, not part of the initial models, are tracked through the sequence by means of an independent gradient-based line tracker, involving interest pointbased tracking, followed by line fitting and re-normalization in each frame. Lines that are determined to be CITLs are displayed in green. The model reconstructed/predicted lines are displayed in yellow. The consistency check is evaluated, using a threshold τ =0.05, based on the lines present in frames [START_REF] Solina | Recovery of Parametric models from range images: the case of superquadrics with local and global deformations[END_REF]40,60) for the bike sequence. For the space shuttle sequences, lines are incrementally tracked and recovered (some are not visible initially) and the consistency checks are performed several times for different lines, in the frames [START_REF] Solina | Recovery of Parametric models from range images: the case of superquadrics with local and global deformations[END_REF]40,60), (120, 140, 160) and (160,175,190).

It is important that the estimation of the models' rigid motion is accurate as it affects both the validity of consistency checks and the quality of the line reconstruction. In practice, since the models are initialized manually, there is always uncertainty associated with their initialization. For this reason, we do not attempt to check and reconstruct lines immediately after model initialization, but rather allow a slight delay (around 20 frames in our experiments), such that the model locks onto the data.

Any initial bias in the linear reconstruction is subsequently eliminated through the re-estimation of the line representation jointly with all the rigid and non-rigid model parameters in the deformable model framework.

The least-squares-based reconstruction for the bike sequence is based on 12 frames within the interval, 20-60, while for the grapple fixture, we use 12 frames in 20-60, 120-160, 160-190, respectively. In our experiments, this provided us with sufficiently important motion and sufficient additional lines for accurate, constrained reconstruction. The stability of the reconstruction is checked, as explained in §3.3 by the ratio of the 2nd and 3rd singular values associated with the A matrix. We thus employ a principled criteria for selecting sufficient lines and corresponding informative displacements for the reconstruction. In both sequences, the frames prior to new line reconstruction are tracked using the prior model along with contour-based image forces, while the rest of the sequence (once the CITLs have been reconstructed) is tracked using the enhanced model (consisting of both the initial model and the reconstructed lines, and using both contour and alignment forces between CITLs and MPLs).

Although no lines are reconstructed until frame 60 in the bike sequence, their re-projection is displayed over the entire sequence (note that since they are reconstructed in a modelcentered coordinate frame, they can be transferred backwards, since the inter-frame motion of the model has already been estimated). For the grapple fixture sequence, we avoided plotting the lines over the entire sequence in order to show how they are incrementally tracked and recovered. Notice that for both sequences, the reconstructed lines' projections are correct.

The tracking results emphasize increased stability, as the line features are reconstructed and integrated into the tracking process as additional cues, especially during difficult to track towards camera motion of the grapple fixture. The number of iterations in the integration of motion equation ( 2) decreases with the integration of additional line cues (see Fig. 5 for plots). The average per-node error decreased from 0.9 (initially) to 0.4 pixels (frame 60) in the bike tracking sequence, and from 1.2 pixels (initially) to 0.8 (frame 60), 0.6 (frame 160) and 0.2 (frame 190) in the grapple fixture sequence, respectively, reflecting improved accuracy due to the integration of line soft constraints in the deformable model estimation process.

Conclusions

We have presented a framework for incremental model acquisition and tracking using deformable models. We have relaxed the constraint that the representation of the model has to be fully known a-priori, and have enhanced its basic discretization structure in terms of points and lines while preserving higher-level shape information (in terms of parametric shapes). This has allowed us to formulate geometric constraints for feature consistency, reconstruction, and tracking via cue integration in a flexible manner through the integration of top-down, unbiased defomable model-based estimation techniques and bottom-up, linear feature reconstruction techniques. Our initial experimental results are promising and show that the proposed framework is able to deal with objects with complex shape and motion, and to track and incrementally recover structure with improved accuracy. Subsequent work will include more rigorous quantitative evaluation of the reconstruction and tracking improvement results, as well as grouping and high-level model abstraction.

Figure 1 :

 1 Figure 1: Line Representation and Motion

Figure 2 :

 2 Figure 2: Lines transfer and alignment

Figure 3 :

 3 Figure 3: Initial and Reconstructed Model Tracking

Figure 5 :

 5 Figure 5: Number of iterations per frame (Bezier interpolation) decreases due to the addition of line constraints.

  frame 0, no reconstructed lines (h) frame 40, 3 reconstructed lines (i) frame 80, 3 reconstructed lines (j) frame 120, 4 reconstructed lines (k) frame 160, 7 reconstructed lines (l) frame 190, 10 reconstructed lines

Figure 6 :

 6 Figure 6: Model Tracking (the bike model in grey flat-shaded and the grapple fixture model in brown wireframe) with CITLs (green) and MPLs (yellow)

Acknowledgments

The authors would like to thank Bill Triggs for many helpful discussions and guidance, Doug DeCarlo for explanations and discussions on deformable models, and Alex Telea for implementation assistance. The authors also give special thanks to Piotr Jasiobedzki, MacDonald Dettwiler Space and Advanced Robotics Ltd., for kindly providing the grapple fixture sequence used in the paper.


Cristian Sminchisescu
email: cristian.sminchisescu@inrialpes.fr

Alexandru Telea

A Framework for Generic State Estimation in Computer Vision Applications

Experimenting and building integrated, operational systems in computational vision poses both theoretical and practical challenges, involving methodologies from control theory, statistics, optimization, computer graphics, and interaction. Consequently, a control and communication structure is needed to model typical computer vision applications and a flexible architecture is necessary to combine the above mentioned methodologies in an effective implementation. In this paper, we propose a three-layer computer vision framework that offers: a) an application model able to cover a large class of vision applications; b) an architecture that maps this model to modular, flexible and extensible components by means of object-oriented and dataflow mechanisms; and c) a concrete software implementation of the above that allows construction of interactive vision applications. We illustrate how a variety of vision techniques and approaches can be modeled by the proposed framework and we present several complex, application oriented, experimental results.

Introduction

Experimenting and building systems in computational vision poses several major challenges. First, such systems involve methodologies from various areas, such as object modeling, optimization, control theory, statistics and computer graphics. Secondly, a control and communication structure is necessary to build complete classes of vision application in terms of the above methodologies. Finally, in order to satisfy the different demands of particular vision applications, one needs a software architecture where they can be built in a complete, extensible, flexible and interactive manner.

These are some reasons for which relatively few generic-purpose vision software systems exist. Most such systems are monolithic software architectures built around specialized linear algebra or numerical optimization modules such as Netlib [START_REF]NETLIB REPOSITORY[END_REF] and LAPACK [START_REF] Anderson | LAPACK User's Guide[END_REF], image processing tools such as SUSAN [START_REF] Smith | SUSAN -a new approach to low level image processing[END_REF] and Intel's Image Processing Library [START_REF][END_REF], or basic visualization tools such as OpenGL [START_REF] Woo | The OpenGL Programming Guide[END_REF] and Geomview [START_REF] Phillips | MUNZNER Geomview -An Interactive Geometry Viewer[END_REF]. Few such systems, if any, have a high-level architecture able to integrate the many aspects To appear ICCV-ICVS 2001 of a computer vision problem (computations, control, visualization, user interaction). Moreover, their design often lacks the simplicity, extensibility and completeness required to build various vision applications out of generic, reusable components. For instance, Intel's Open Source Computer Vision Library [START_REF][END_REF] offers various functionalities, such as camera calibration, image thresholding, image-based gesture recognition, but does not structure them in a modular, extensible, and customizable way. Target Jr [19], intended object-oriented, employs a one-way communication to its Netlib-based Fortran optimization routines. Since this control is hard-wired in its classes, Target Jr is not usable for applications that require user input or data monitoring during the optimization process. Finally, almost no computer vision software we know allows adding interactive manipulation and visualization of its data structures to the computational code in a simple, yet generic manner. Visual interactivity is important in several respects, such as the setting of correspondences between 3D vision models and 2D images, and for monitoring and control of the time evolution of vision applications.

Summarizing, many vision systems lack a generic control model. Secondly, they provide heterogenous functionalities that cannot be easily adapted to new application contexts.

In order to address these problems, we propose a generic high-level application model that covers a large class of vision approaches, a software architecture of this generic model that combines the flexibility, extensibility, and interactivity requirements, and an efficient and effective implementation of this architecture. The proposed application model is based on the idea that a large range of vision applications share the concept of generic optimal state estimation. Such applications involve:

-a model characterized by its state -a generative transformation which predicts discretized model features in the observation space, based on a current state configuration -an association of predicted and extracted features in the observation space to evaluate a configuration cost -a control strategy that updates the model state such that the evaluation cost meets an optimality criterium This application model is powerful enough to cover many techniques such as deformable models/dynamical systems, optimization based methods and sampling based methods, as well as combinations of them (mixed methods). The application model is flexible and extensible as it does not make any assumptions about the ways the model is represented or discretized, the type and number of extracted features (cues), how the association is done, and the underlying strategy used to estimate and evolve the model state. From the proposed model, we derive an architecture that provides desired extensibility, flexibility and modularity requirements in terms of object-oriented and dataflow mechanisms, and finally, we propose an implementation of the above architecture which allows building complex vision applications with visual parameter monitoring and interactive control. The rest of the paper describes the proposed vision framework, as follows. Section 2 describes the generic vision application model we propose. Section 3 describes an architecture that combines the model's genericity with the flexibility, extensibility, and modularity requirements, and presents a concrete implementation of the architecture. Section 4 presents several vision applications constructed in the proposed framework.

Application Model

The proposed application model exploits the observation that many vision applications are based on the generic state estimation concept. This concept involves optimal estimates of a model's parameters based on a (possibly temporal) sequence of observations. In detail, our application model involves the following elements (see also Fig. 1 which depicts these elements and the data streams between them):

1. a "representational" discretization in a spatial domain. 2. a composite (generally non-linear) generalized transformation (T), parameterized in terms of the current configuration (typically an instance of the state) which generates predictions in the observation space for points in the discretized domain. This item and the previous one form the model representation. 3. a sequence of (possibly temporal) observations or extracted features. 4. a way to associate predictions to features to evaluate a configuration cost or higher order operators associated with it (their computation typically needs equivalent quantities of T). 5. a strategy to evolve the model state based on the evaluation of configuration costs or higher-order operators associated to it, such as its gradient or Hessian, in order to match an optimality criterium (see Fig. 1). 6. several user interaction and visualization components responsible for the application building and scene exploration, including interactive model-data couplings, 2D and 3D renderings, as well as classical graphics user interfaces (GUI).

The application model elements outlined above are detailed in the following sections.

Configuration Evaluator

The basic task of the configuration evaluator is to compute the cost of a configuration or higher order operators (gradient or Hessian) associated with it. Their evaluation is typically performed in terms of equivalent quantities computed for the model/data generalized transform (T).

Model Representation

The model is spatially discretized in a domain . For any point ¡ £¢ , we can compute a prediction in the observation space ¤ ¦¥ ¨ § © ¡ .

The Jacobian matrix of the generalized transformation makes the connection between differential quantities in parameter and observation spaces:

¤ ¥ § !¥ #" (1)
The process of model estimation involves a data association problem between individual model feature predictions $ &% and one or more observations that we shall gener- ically denote ' $ % (with additional subscripts if these are several). We refer to ( )$ % © 0¤ 1 2¥ ' $ % 43 $ % © 5¤ 6 as the feature prediction error.

Feature Extraction, Data Association and Error Evaluation

These components extract the various types of information used in vision applications. Typical datasets include 2D image data, namely edge/contour/silhouette information, optical flow information, or 3D range data obtained, for instance, from a multi-camera system.

A subsequent data associator or matching stage establishes correspondences between model predictions and data features. The matching is either explicit such as in the case of an optical flow module or a 3D model to range features, or implicit, when every predicted model feature has already a computed cost on a potential surface (see below and Fig. 3). In the case of explicit matched features a separate error evaluation stage computes the feature prediction error, based on a certain error distribution. Common error distributions include both well-known non-robust Gaussian ones and robustified ones that model the total (inlier plus outlier) distribution for the observation, e.g.:

where ¥ ¤£ control the outlier threshold and robust influence, respectively, ¡ is the squared feature prediction error, ¡ ¥ ( )$ &% % ( )$ % , and % is a symmetric positive definite weighting matrix associated with feature ! . For implicit feature matching, let us consider the case of edges, where a potential surface can be build such that it has gaps in the places in the images corresponding to edge features or steep changes in intensity. The image potential is computed for every image frame from the image intensity " © 0¤ $# as follows:

% © 5¤ &# ¥ 3 (' 0) 21 © ¢3 54 76 " © 5¤ &# ) (3)
where £ determines the width of the Gaussian function 3 4 , 6 denotes the convolution operator, and ' determines the potential surface steepness. The potential is generating a 2D force field, given by: 8

% @9 (A B DC © 5¤ $# ¥ 3 E1 % © 0¤ $# (4)

State Representation and Control

The state representation and control components constitute the core of the application model. Various state representations correspond to various control strategies. Typical ones include unimodal, Gaussian state representation and multiple-hypothesis (samplebased) representations. Generic control strategies include:

-continuous ones (deformable models, continuous optimization) which generally assume unimodal state representation and evaluation of the cost function and its higher order operators (gradient, Hessian, etc.);

-discrete ones which only involve cost function evaluations and sampling methods for focusing the search effort; -mixed strategies that involve a combination of the first two.

Regardless of the state representation and control employed, there is a unique interface between these components and the configuration evaluator (Fig. 1).

Dynamical Systems and Deformable Models

A deformable model estimates the state of the system by numerical integrating a synthetic dynamical system that encodes the fitting error. The Lagrangian equations governing its evolution can be written as:

¢¡ ¤£ ¦¥ ¤£ ¨ § !¥ 8 © 8 © ¥ " © (5)
where ¥ " " , ¥ ¥ " " , and § ¥ ! "! © $# &% % , with and being tuning parameters, # '% % being the stiffness associated with the parameter ! , and 8 © are general- ized "forces", acting on the state parameters and © is a distance error in the observation space (typically an " )( norm).

Continuous Optimization Based Methods

Optimization based methods perform a non-linear estimation in terms of the generalized transformation gradient and Hessian. Popular optimizers include second order damped Newton trust region methods that choose a descent direction by solving the regularized system [START_REF] Fletcher | Practical Methods of Optimization[END_REF]:

© 10 2£ 43 5 ¥ 3 ! (6)
where is a symmetric positive-definite matrix and 3 is a dynamically chosen weight- ing factor. For robust error distribution specific gradient and Hessian approximations have to be derived in terms of the robustifiers [START_REF] Sminchisescu | A Robust Multiple Hypothesis Approach to Monocular Human Motion Tracking[END_REF]. For least squares problems, ! ¥ )" © and 0 76 8 " " .

Sampling Methods Sampling methods, usually known in vision under the generic name of CONDENSATION [START_REF] Isard | CONDENSATION -conditional density propagation for visual tracking[END_REF], are discrete methods that propagate the entire parameter distribution in time as a set of hypotheses, or samples, with their associated probability. In each frame, the entire distribution is resampled (i.e. recomputed and reweighted) based on new image observations. As these methods do not traditionally have a continuous component, evaluating the distribution in our scheme requires only the evaluation of a configuration, but no other associated higher-order operators. The computational burden lies in the strategies for sampling the parameter distribution in order to locate typical sets, i.e. areas where most of the probability mass is concentrated. Known strategies include importance, partitioned (Gibbs) or annealing based sampling methods (see [START_REF] Mackay | An Introduction to Monte Carlo Methods[END_REF]). We developed an uniform interface such that various discrete sampling methods can be parameterized by the search strategy.

Mixed Continuous/Discrete Methods Mixed continuous/discrete methods have been proposed recently [START_REF] Sminchisescu | A Robust Multiple Hypothesis Approach to Monocular Human Motion Tracking[END_REF] in an effort to combine generic robustness properties of samplebased techniques with the local informed, accuracy and speed properties of continuous ones. The state is represented effectively as a set of hypotheses. In each frame, each hypothesis is subject to a continuous optimization followed by a hypothesis generation based on the uncertain directions of the continuous estimate. Finally, the current set of hypotheses is pruned and propagated to the next time step.

Interaction/Visualization

The interaction and visualization modules are responsible for scene exploration and manipulation, parameter monitoring and control, and interactive application building. While further details are provided in Section 3, we point out that the dual, simultaneous type of control, user driven vs. application driven, imposes particular demands on the system design, namely consistency and control problems: user input has to be consistently integrated and propagated into the application data structures (even during application driven execution), while the structuring of individual application components has to allow external interaction during the execution of their different operations. We address these problems by means of automatically enforcing dataflow consistency and by using object-orientation in order to design components with open, non-encapsulated, control.

Architecture

The above sections present a generic model that can accommodate a large class of vision applications. These applications share the conceptual structure described in Fig. 1 in terms of component functionalities and intercommunication. To make this application model viable, one must produce a software implementation of it that complies with the requirements of modularity, simple application construction and extension, and flexible control specification, discussed in the previous sections.

We have achieved the above combination by choosing a specific software architecture for implementing the discussed conceptual model. This architecture is based on the combination of two design principles: object-orientation and dataflow, as follows.

Object Orientation

Object orientation (OO) provides a sound, complete framework for modeling a system as a set of related software modules, by what is called a class hierarchy [START_REF] Gamma | Design Patterns: Elements of Reusable Object-Oriented Software[END_REF]. Specific mechanisms such as subclassing allow the incremental construction of application functionality by specializing a few basic concepts, as well as an effective reuse of the written code. In our case, we have implemented our vision application model as a C++ class library that specializes a few base classes that correspond to the generic concepts shown in Fig. 1. The fixed characteristics of the model, such as the data interfaces between the generic components in Fig. 1, reside the library's base classes. Subclasses add specific functionality to the basic components, in an orthogonal way. For example, different control strategies or model parameterizations, error evaluations, or feature detection techniques can be added easily, independently on each other, and without modifying the basic software architecture.

Dataflow

Object orientation effectively models the static, structural relations between the framework's components. To model the dynamic, control relations, we added the dataflow concept to the C++ classes. We structure a vision application as a network of classes that have data inputs, outputs, and an update operation. Application execution is driven by the network structure: once a class's input change, the class reads the new input, updates its output, and triggers the execution of the other classes connected to its output, thus enforcing the dataflow consistency. A dataflow architecture allows constructing complex control strategies in a simple, yet flexible way, by connecting together the desired components in the desired network. Keeping the vision computational code inside the components and the control strategy outside them, in the network structure, has two advantages: the vision components are simple to write, extend, and understand, and they are directly reusable in applications having different control scenarios.

Implementation

Writing a vision application in the above setup implies constructing and updating the desired dataflow network. We have made this process flexible, by integrating our vision C++ library in the VISSION dataflow application environment [START_REF] Telea | VISSION: An Object Oriented Dataflow System for Simulation and Visualization[END_REF]. In VISSION, the dataflow networks are constructed interactively by assembling iconic representations of C++ classes in a GUI network editor (see Fig. 5 b). Figure 4 a,b show two such networks for the vision applications discussed in the next section. The graphical icons are actual subclasses of the vision components shown in Fig. 1.

Besides providing a simple, intuitive way to construct the application, this solution offers several other advantages as an implementation for our architecture. First, VISSION provides graphics user interfaces (GUIs) automatically for all the classes of a network, thus allowing for parameter changes and monitoring. Figure 5 a shows such an interface for the 3D range data tracking application discussed in the next section. Secondly, once a parameter is changed, the traversal and update of the dataflow network is performed automatically, thus enforcing the dataflow consistency. Thirdly, VISSION can dynamically load different C++ class libraries. Consequently, we integrated our vision library with the Open Inventor library which provides several direct manipulation tools and 2D and 3D viewers by which monitoring the time evolution of a vision experiment and setting up model-data correspondences can be done in a simple, interactive way, as described in the next section.

There exist several dataflow application environments similar to VISSION, such as AVS [START_REF] Upson | The Application Visualization System: A Computational Environment for Scientific Visualization[END_REF] or Khoros [START_REF] Research | Khoros Professional Student Edition[END_REF]. From an end user point of view, these environments offer the same visual programming, dataflow, and user interface facilities as VISSION. However, integrating our vission C++ library in such environments would pose several architectural problems. First, these environments are primarily designed to integrate C or FOR-TRAN code. Integrating C++ libraries, such as our vision library or the Open Inventor library, is a complex task. Secondly, these environments assume that the code to be integrated is developed based on an application programmer interface (API) provided by the environment. Such an API usually contains data types and functions via which the user code and the environment communicate. In our case, however, we wish to keep our vision library independent on a specific environment API. In this way, our vision library can be deployed in different environments, such as custom turnkey applications, with little or no modification.

Applications

In this section we discuss applications involving shape estimation and tracking both on 2D monocular image data and 3D range data. These applications that we have implemented in our framework are typical examples for the state estimation approach based on temporal observations.

In order to model a large class of objects or structures, we have implemented a hierarchical multi-object parametric representation (see Fig. 2), by subclassing the model representational components (Fig. 1). More precisely, any node ¢¡ ¤£ ¦¥ correspond- ing to one of the objects discretization domain can be transformed into a 3D point § ¡ ©¨ § ¡ , and subsequently into an image prediction ¡ ¨ ¡ , by means of a composite non-linear transformation: ¡ !#" $¡ %'& ( § ¡ )10 2 3 54 6 7¡ 64 68 9 @ A4 ¡ B [START_REF] Gamma | Design Patterns: Elements of Reusable Object-Oriented Software[END_REF] where 8 represents a sequence of parametric deformations which construct the corre- sponding part in a self-centered reference frame, 0 represents a chain of transforma- tions that position the corresponding part in the hierarchy, and & represents the perspec- tive projection of the camera, in case we work with 2D data. The entire transformation " (see Subsection 2.1) is assembled in terms of the individual transforms " !¡ .

2D Image-Based Temporal Estimation

In this subsection, we illustrate the flexibility of the presented framework by a temporal estimation application based on 2D image sequences. The examples we show involve various types of model representations, various types of feature extractors and data associators as well as different control strategies and interaction modules.

Human Body Tracking

We use human body motion sequences that consist of an arm and an entire human body motion (see Figs. 9 and8). They both involve complex motions with significant self-occlusion, in cluttered backgrounds, including subjects wearing loose clothing. In order to accommodate the complexity of such an application, we have derived a novel mixed-type control method, combining robust (Eqn.2) continuous optimization and sampling along the uncertain directions of the continuous estimate (see Section 2.2).

For model representation, we use the hierarchical scheme mentioned above and superquadric ellipsoids with tapering and bending deformations as basic representational primitives. The estimation involves the articulations degrees of freedom, namely 6 for the arm and 30 for the human model.

We employ both edge and intensity features as cues during estimation. The optical flow has been extracted based on a robust multi-scale computation and the edges have been weighted based on the motion boundary map extracted from that computation (see Fig. 3). Furthermore, we used a multiple association scheme for the assignment of model contour predictions to image observations. Consequently, the Feature Extractor and Data Associator modules have been implemented by chaining together the corresponding Extractors/Associators for contours and optical flow, weighted accordingly (see Fig. 4 a).

A 3D model to 2D image Manipulator interaction module has been employed for specifying correspondences between the model joints and the subject joints in the image, in order to initialize the model in the first frame of the sequence. Subsequent estimation, both for the model initial pose and proportions as well as articulation parameters during the next frames, is fully automated based on image cues. Finally, the model geometry and deforming forces are visualized by the GeomMapper and the Hedge-hogMapper modules respectively, within an interactive Viewer module.

Incremental Model Acquisition and Tracking

The bicycle sequence in Figure 10 uses a mixed model representational structure and a control strategy based on deformable models for the state estimates (Section 2.2). Only contours are used as image cues initially, based on a single association scheme using implicit matching in a computed image potential (Section 2.1). User interaction was needed during model initialization as for the human sequence.

The tracking started with a minimal hierarchical model, a frame consisting of 3 parametric shapes. New model components, with different discretization and parameterizations, are discovered and recovered during tracking based on geometric consistency checks (see [START_REF] Sminchisescu | DICKINSON Incremental Model-Based Estimation Using Geometric Consistency Constraints[END_REF] for details). Particularly, lines moving consistently with the bicycle frame, but not part of the initial parameterization, are identified, reconstructed, This application needs flexibility at all levels of the vision system: we need various representations and discretizations for model components, various feature extractors and data associators corresponding to each component prediction/observation couple, and a way to integrate each component contribution as a cue in the model control. For instance, in a deformbale model control, Jacobians corresponding to each component representational mapping need to be evaluated.

3D Flow Feature Tracking

In this section, we present a 3D feature tracking and reduced modeling application based on information extracted from fluid dynamics simulations. We employ a model representation based on a superquadric ellipsoid with tapering, bending, shearing and twisting parametric transformations. Next, we employ a control strategy based on deformable models (see Section 2.2). The application's network is shown in Fig. 4 a.

Conclusions

We have presented a layered framework for the modeling, design and implementation of integrated computer vision applications. Our proposal includes: (1) an application model in the form of generic state estimation, involving a control and communication structure able to cover a large class of vision paradigms; (2) an architecture able to support flexibility, extensibility and interactivity demands based on object-oriented and We have subsequently presented how various vision technique classes, like the ones based on dynamical systems/deformable models, continuous optimization, random sampling, as well as combination of them, can be naturally integrated in the proposed framework. Finally, we have presented several complex integrated vision applications which have been easily built within our framework.

Fig. 1 .

 1 Fig. 1. Vision application model

Fig. 2 .

 2 Fig. 2. Various hierarchical models used during tracking

Fig. 3 .

 3 Fig. 3. Image processing operators: original image (a), edge detection (b), motion boundaries (c), robust flow in horizontal (d) and vertical (e) directions % © ¢¡ ¤£ 4 ¥ ¦¥ 4© ¨ § 3 ©

Fig. 4 .

 4 Fig. 4. Tracking networks for human motion (a) and 3D range data (b)

Fig. 5 .

 5 Fig. 5. Data viewer and model manipulator interface (a). A vision application in the visual dataflow network editor (b)

Fig. 6 .Fig. 7 .

 67 Fig. 6. 3D range data tracking examples

Fig. 8 .Fig. 9 .

 89 Fig. 8. Frames during human tracking sequence

Dataset courtesy of N. Zabusky, D.Silver and X. Wang[START_REF] Silver | WANG Tracking and Visualizing Turbulent 3D Features[END_REF]


C Sminchisescu
email: cristian.sminchisescu@inria.fr

A Telea

HUMAN POSE ESTIMATION FROM SILHOUETTES A CONSISTENT APPROACH USING DISTANCE LEVEL SETS

Keywords: human tracking, model-based estimation, constrained optimization, level set methods, fast marching methods

We present a novel similarity measure (likelihood) for estimating three-dimensional human pose from image silhouettes in model-based vision applications. One of the challenges in such approaches is the construction of a model-to-image likelihood that truly reflects the good configurations of the problem. This is hard, commonly due to the violation of consistency principle resulting in the introduction of spurious, unrelated peaks/minima that make the search for model localization difficult. We introduce an entirely continuous formulation which enforces model estimation consistency by means of an attraction/explanation silhouette-based term pair. We subsequently show how the proposed method provides significant consolidation and improved attraction zone around the desired likelihood configurations and elimination of some of the spurious ones. Finally, we present a skeleton-based smoothing method for the image silhouettes that stabilizes and accelerates the search process.

INTRODUCTION AND PREVIOUS WORK

Human pose estimation from images is an active area of computer vision research with many potential applications ranging from computer interfaces to motion capture for character animation, biometrics or intelligent surveillance. One promising approach, called model based [Smin01b, Deut00, Heap98, Smin01a, Gavr96, Breg98, Kakad96, Rehg95], relies on a 3D articulated volumetric model of the human body to constrain the localization process in one or several images. The goal in human pose estimation applications is to estimate the model's articulation and possibly structural parameters such that the projection of the 3D geometrical model closely fits a human in one or several images. Typically, model localization is a multi-dimensional expensive search process in the model parameter space for good cost configurations defined in terms of maxima of a likelihood, or minima of an energy function. Such costs are defined in terms of the association of model predictions with extracted image features. The search process produces a parameter configuration which brings the 3D model close to the tracked 2D image in the metric of the predefined likelihood model. The above problem is hard since likelihood cost surfaces are typically multi-peaked, due to factors like multiple scence objects, ambiguous feature assignments, occlusions, and depth uncertainties.

Search strategies for locating good peaks in the model parameter space based on local and global search methods, possibly in temporal sequences, have received significant attention [Smin01b, Deut00, Heap98, Smin01a, Gavr96, Breg98] and are not addressed here. However, the dificulty and intinsically ill-posed nature of such search problems raise two complementary questions about the design of the cost surface whoose minima are to be found: what are good image features which will readily qualify for likelihood terms for sampling and continuous evaluation ? how to define such terms to limit the number of spurious minima in parameter space and render the search more efficient and effective.

Likelihood models defined in terms of edges [START_REF] Deutscher | Articulated Body Motion Capture by Annealed Particle Filtering[END_REF][START_REF] Sminchisescu | A Robust Multiple Hypothesis Approach to Monocular Human Motion Tracking[END_REF][START_REF] Kakadiaris | Model-Based Estimation of 3D Human Motion with Occlusion Prediction Based on Active Multi-Viewpoint Selection[END_REF], silhouettes [START_REF] Deutscher | Articulated Body Motion Capture by Annealed Particle Filtering[END_REF][START_REF] Sminchisescu | A Robust Multiple Hypothesis Approach to Monocular Human Motion Tracking[END_REF][START_REF] Heap | Wormholes in Shape Space: Tracking through discontinuities changes in shape[END_REF] or intensities [START_REF] Sminchisescu | A Robust Multiple Hypothesis Approach to Monocular Human Motion Tracking[END_REF][START_REF] Sidenbladh | Stochastic Tracking of 3D Human Figures Using 2D Image Motion[END_REF][START_REF] Rehg | Model-Based Tracking of Self Occluding Articulated Objects[END_REF] are the most common. While image intensities seem to be good cues for various types of optical-flow based local search, they are not invariant to lighting changes, and typically rely on low inter-frame intensity variations and motion. It is consequently difficult to sample configurations out of the region where such photometric model is valid. Edges and/or silhouettes have therefore been more used in approaches that employ, at least partially, some form of parameter-space sampling [START_REF] Deutscher | Articulated Body Motion Capture by Annealed Particle Filtering[END_REF][START_REF] Sminchisescu | A Robust Multiple Hypothesis Approach to Monocular Human Motion Tracking[END_REF][START_REF] Heap | Wormholes in Shape Space: Tracking through discontinuities changes in shape[END_REF].

Deutscher [START_REF] Deutscher | Articulated Body Motion Capture by Annealed Particle Filtering[END_REF] uses a silhouette based term for his cost function design in a multi-camera setting. However, this term peaks if the model is inside the silhouette without demanding that the silhouette area is fully explained (see Sec. 4.1). Consequently, an entire family of undesired configurations situated inside the silhouette will generate good costs under this likelihood model. Moreover, the term is purely discrete, not suitable for continuous estimation. The situation is alleviated by the use of the additional cues and sensor-fusion from multiple cameras with good results. Delamarre [START_REF] Delamarre | 3D Articulated Models and Multi-View Tracking with Silhouettes[END_REF] uses silhouette contours in a multi-camera setting and computes assignments using a form of ICP (Iterative Closest Point) algorithm and knowledge of normal contour directions. The method is local and not necessarily enforces globally consistent assignments, but again relies on fusing information from many camera to ensure consistency. Brand [START_REF] Brand | Shadow Puppetry, Proc. ICCV[END_REF] and Rosales [START_REF] Rosales | Inferring Body Pose without Tracking Body Parts[END_REF] use silhouettes to infer temporal and static human poses. However, their motivation is slightly different in using silhouettes as inputs to a system which directly learns 3D to 2D mappings. Summarizing, many likelihood terms used in modelbased vision applications have the undesirable property that they not only peak around the desired model configurations, which correspond to subject localization in the image, but also in totally unrelated, false configurations. This poses huge burdens on any search algorithm, as the number of spurious minima could grow unbounded and therefore discriminating them from "good peaks" can only be done via temporal processing. Consequently, any finite samples/hypothesis estimator has a great chance to miss significant, true minima.

In practice, extracting pose from silhouette using single images remains an under-constrained problem with potential multiple solutions. A more global search method, multiple cameras, temporal disambiguation and/or additional features have thus to be used in conjunction with the local method we propose in this work, to robustify the search for good cost configurations [START_REF] Sminchisescu | Covariance-Scaled Sampling for Monocular 3D Body Tracking[END_REF][START_REF] Deutscher | Articulated Body Motion Capture by Annealed Particle Filtering[END_REF][START_REF] Heap | Wormholes in Shape Space: Tracking through discontinuities changes in shape[END_REF][START_REF] Sminchisescu | A Robust Multiple Hypothesis Approach to Monocular Human Motion Tracking[END_REF][START_REF] Gavrila | 3-D Model Based Tracking of Humans in Action:a Multiview Approach[END_REF]. In this paper, we assume a reasonable initialisation and restrict our attention to the design of likelihoods with larger basin of attraction zones and globally consistent responses around the desirable cost minima. We achieve this by means of an entirely continuous formulation and a new likelihood term for silhouettes in model-based applications. The proposed term allows a globally consistent response for the subject localization in the image by means of a pair of attraction/explanation components that a) push the geometric model inside the subject's silhouette and b) demand that the area associated with the silhouette is entirely explained by the model. We subsequently show how this proposal significantly improves the pose estimation results compared to previously used similarity measures.

In Section 2, we describe the human body model we employ. Section 3 outlines the search process for optimal configurations. Section 4 introduces our new likelihood terms and details its two components. Section 5 presents a new technique for smoothing the image-acquired silhouettes that stabilizes and accelerates the search process. Finally, Section 6 concludes the paper and proposes directions for future work.

HUMAN MODEL

MODEL DESCRIPTION

Our human body model (Fig. 1) consists of kinematic 'skeletons' of articulated joints controlled by angular joint parameters ¡ £¢ , covered by 'flesh' built from superquadric ellipsoids with additional tapering and bending parameters [START_REF] Barr | Global and local deformations of solid primitives[END_REF]. A typical model has around 30 joint parameters, plus 8 internal proportion parameters ¡ ¥¤ encoding the positions of the hip, clavicle and skull tip joints, plus 9 deformable shape parameters for each body part, gathered into a vector ¡ ¥¦ . The state of a complete model is thus given as a single parameter vector ¡ ¨ § © ¡ ¢ ¡ ¥¦ ¡ ¤ . We note, however, that only joint parameters are typically estimated during object localization and tracking, the other parameters remaining fixed.

Although this model is far from photo-realistic, it suffices for a high-level interpretation and realistic occlusion prediction. Moreover, it offers a good tradeoff between computational complexity and coverage in typical motion tracking applications.

MODEL TO IMAGE FITTING

The model is used in the human pose estimation application as follows (see also Fig. 2 for an overview of the application pipeline).

£ ¤ § ¤£ ¤ © ¡ using com- posite nonlinear transformations £ ¤ © ¡ § ¦¥ ¤ © ¡ § § © ¡ ¤ © ¡ § § © ©¨© ¡ ¥¢ ¡ ¤ © ¡ ¥¦ ¤
, where is a sequence of parametric deformations that construct the corresponding part in its own reference frame, ïs a chain of rigid transformations that map it through the kinematic chain to its 3D position, and § is the perspective projection.

During parameter estimation (see Sec. 3), predictionto-image matching cost metrics are evaluated for predicted image feature £ ¤ , and the results are summed to produce the image contribution to the overall parameter space cost function. For certain likelihood terms like edge based ones, predictions £ ¤ are associ- ated with nearby image features £ ¤ . The cost is then a function of the prediction errors £ ¤ © ¡ § £ ¤ £ ¤ © ¡ .

For other likelihood terms (like the silhouette attraction term we employ here), a potential surface is built for the current image, and the prediction is only evaluated at a certain location on this surface.

PARAMETER ESTIMATION

We aim towards a probabilistic interpretation and optimal estimates of the model parameters by maximizing the total probability according to Bayes rule:

¡ © ¡ ! " £ $# ¡ © £ % ¡ ¡ © ¡ § '& )( 10 32 © ©4 ¢ 65 4 87 @9 ¡ © ¡ (1)
where 4 ¢ and 4 87 are the new silhouette likelihood terms we propose, defining similiarity criteria between the model projection and the image silhouette to be defined in the next section, and ¡ © ¡ is a prior on model parameters. The prior encodes static knowledge on humans, such as anatomical joint angle limits for parameters or non-penetration constraints on the body parts (see [START_REF] Sminchisescu | Covariance-Scaled Sampling for Monocular 3D Body Tracking[END_REF][START_REF] Sminchisescu | A Robust Multiple Hypothesis Approach to Monocular Human Motion Tracking[END_REF] for details).

In a maximum a-posteriori estimate (MAP) approach, we spatially discretize the continuous formulation in Eqn. 1, and attempt to minimize the negative loglikelihood, or 'energy', for the total posterior probability. The energy is expressed as the following cost

function: A © ¡ § CB ED F © ¡ © £ % ¡ ¡ © ¡ § CB ED F ¡ © £ G ¡ H IB ED F ¡ © ¡ § 4 ¢ 5 4 7 5 A © ¡
where A © ¡ is the negative log of the model prior. In the following, we shall concentrate on the behavior and properties of the negative log-likelihood 4 ¢ P5 4 87 V is a symmetric positive-definite stabilization matrix (ofter set to identity)

T is a dynamically chosen weighting factor `a Yc is a matrix containing joint angle limits constraints acting as effective priors, defining an addmissible subspace to search for model parameters (see [START_REF] Sminchisescu | Covariance-Scaled Sampling for Monocular 3D Body Tracking[END_REF][START_REF] Sminchisescu | A Robust Multiple Hypothesis Approach to Monocular Human Motion Tracking[END_REF] for details).

The parameter T controls the descent type: T ¡ £¢ leads to a gradient descent, while T ¤ h leads to a Newton-Raphson step. The optimization routine automatically decides over the type and size of the optimal step within the admissible trust radius (see [START_REF] Fletcher | Practical Methods of Optimization[END_REF][START_REF] Triggs | Bundle Adjustment -A Modern Synthesis[END_REF] for details).

OBSERVATION LIKELIHOOD

Whether continuous or discrete, the search process depends critically on the observation likelihood component of the parameter-space cost function. Besides smoothness properties, necessary for the stability of the local continuous descent search, the likelihood should be designed to limit the number of spurious local minima in parameter space. We propose a new likelihood term, based on two components: the first component maximizes the modelimage silhouette area overlap.

the second component pushes the model inside the image silhouette.

The above pair of cost terms produces a global and consistent response. In other words, this term enforces the model to remain within the image silhouette, but also demands that the image silhouette is entirely explained, i.e. that all silhouette parts contribute to the cost function that drives the fitting process. In the following, we detail the two cost components.

SILHOUETTE-MODEL AREA OVERLAP TERM

This term maximizes the model-image area overlap.

The area of the predicted model can be computed from the model's projected triangulation by summing over all visible triangles ¥ §¦ © (triangles having all the vertices

© ¡ ¤ ¤ ¤ visible). ¢ § ! #" %$ '& ! ¤ © ¡ ¤ )( © ¤ 0 ( ¤ 0 ¡ ( (2)
where 1 describes the modulo operation, and the computation assumes the triangle vertices are sorted in counter-clockwise order to preserve positive area sign. In subsequent derivations we drop the modulo notation for simplicity.

Let 32

be the area of the target silhouette. The area alignment cost, i.e. the difference between the model and image silhouette areas, is:

4 ¢ § 45 76 ¡ © ! #" %$ 8& ¢ 92 ¡ (3)
The gradient and Hessian for the area-based cost-term can subsequently be derived (by dropping the scaling term):

Q ¢ § ¤@ 4 ¢ @ ¡ § © ! #" %$ '& ¢ 32 ! #" %$ '& A ¢ A ¡ (4)
where:

A ¢ A ¡ § ! ¤ B A ¡ ¤ A ¡ C © ¤ B0 ¤ B0 ¡ (5) 5 ! ¤ ¡ ¥¤ © A ¤ B0 A ¡ A ¤ 0 ¡ A ¡ C (6) R ¢ § D@ ¡ 4 ¢ @ ¡ ¡ § ! #" %$ & A ¢ A ¡ C A ¢ A ¡ (7) 5 © ! #" %$ '& ¢ ¤ ! #" %$ '& A ¡ ¢ A ¡ ¡ (8)
One should notice that the individual partial derivatives E GF 'H E GF and E GI H E GF represent the columns of the indi- vidual Jacobian matrix evaluated at the corresponding prediction for the mesh node P , £ ¤ © ¡ § © ¡ ¤ Q ¤ .

In practice, computing node visibility and area differences is rather fast, as we use the frame and z buffers to this end.

SILHOUETTE ATTRACTION TERM

This second term pushes the model inside the image silhouette. Adding over all projected model nodes P , this term writes:

4 7 § 45 R6 ¡ ! ¤ 4 87 H (9)
where 4 87 H is the distance from a predicted model point £ ¤ © ¡ to a given silhouette 92 . We estimate 4 7 H by computing the distance transform of the silhouette

32

and evaluating it in the points P :

4 7 H © £ ¤ © ¡ 2 § © £ ¤ © ¡ (10)
We use a level-set based approach to quickly and robustly estimate , as follows. We initialize to zero on 32 , i.e. regard 92 as the zero level set of the function . Next, we compute by solving the Eikonal equation [Seth99]: § 4

(11) for all points outside 92

. The solution of equation 11 has the property that its isolines, or level sets, are at equal distance from each other in the 2D space (Fig. 3). Consequently, is a good approximation of the distance transform . Equation 11 is efficiently solved by using the fast marching method (FMM), introduced by Sethian in [START_REF] Sethian | A Fast Marching Level Set Method for Monotonically Advancing Fronts[END_REF]. We briefly outline here the FMM. A detailed description of the FMM, up to the implementation details we have ourselves used, is given in [START_REF] Sethian | A Fast Marching Level Set Method for Monotonically Advancing Fronts[END_REF][START_REF] Sethian | Level Set Methods and Fast Marching Methods[END_REF]. First, is initialized to zero in all points on the silhouette 2

. Next, the solution is built outwards starting from the smallest known value. This is done by evolving a so-called narrow band of pixels, initially identical to 2 , in normal direction to 2 , with unit constant speed. As the narrow band evolves, it takes the shape of the consecutive, equidistant level sets, or isolines, of the function (Fig. 3 d).

Using the FMM to compute the distance has sev- eral advantages. First, the function obtained is con- tinuous over the 2D plane, which is important as we need to evaluate its first and second order derivatives, as explained below. Secondly, the FMM performs robustly even for noisy silhouettes 92

. This is essential for practical applications, as the silhouettes extracted from real images have many disconnected, spurious pixels (Fig. 3 b is a typical example). Thirdly, the FMM is very efficient, as it needs ¡ © £¢ ¥¤ §¦ ©¨Q opera- tions, where ¢ is the number of image pixels and is the average number of pixels in the narrowband, of the same order as the number of pixels on the silhouette's contour. is computed in real time for h h ¡ pixel images on an SGI O2 R5000 machine. Finally, implementing the FMM is straightforward, as described in [START_REF] Sethian | A Fast Marching Level Set Method for Monotonically Advancing Fronts[END_REF]. Overall, we believe that using the FMM to compute is a more efficient and effective method than e.g. chamfer based methods widely used in vision and imaging applications.

The gradient and Hessian of the corresponding silhouette attraction term are computed from the modelimage Jacobian, as follows:

Q 7 § ! ¤ @ © £ ¤ © ¡ @ ¡ § ! ¤ " %$ C ¤ A A £ ¤ (12) R 7 § ! ¤ @ ¡ @ ¡ ¡ ! ¤ " %$ C ¤ A ¡ A £ ¡ ¤ ¤ (13)
Figure 4 shows the effect of the silhouette attraction and area overlap terms for two images taken from a longer tracking sequence. The figure shows the initial images (a,e), the initial model configuration (b,f), and the fitting results obtained when using only the silhouette attraction term (c,g) and finally both the silhouette attraction and the area overlap terms (d,h). One can notice that the silhouette attraction term does not suffice for a good fit. Indeed, any parameter configuration which places the model inside the image silhouette can be potentially chosen. Adding the area overlap term stabilizes the estimation and drives it towards relatively satisfactory results. Moreover, the cost term has the desired properties of a wide attraction zone. This makes it a good candidate in tracking applications where recovery from tracking failures is highly desirable.

SILHOUETTE SMOOTHING

The gradient Q and Hessian R introduced in the pre- vious sections are at the core of the optimization process that fits the model to the observed image features. The stability of the optimization is influenced by the behavior of Q and R : if the silhouette data are noisy, then the cost terms 4 7 and 4 ¢ , and their derivatives Q and R , are not smooth functions. In such cases, the optimization process might fail or take too long to converge or might fit the model erroneously to the image silhouette.

We alleviate this problem by performing a smoothing on the silhouettes acquired from the image data. The smoothing aims to produce silhouettes from the image data that can be easier approximated by our human body models than the original 'raw' silhouettes. The process runs as follows (see Fig. 6 for an overview). In second next step, we compute the skeleton of the silhouette, as follows. We apply the FMM algorithm inwards on the raw silhouette and compute the distance map ¢¡ of all the points inside the silhouette to its boundary (Fig. 6 b). The silhouette skeleton is then computed as being those points of the evolving narrow band that meet other similar points due to the band's evolution under normal speed. In other words, the skeleton points are those points where the narrow band collapses onto itself during its evolution driven by the FMM algorithm. We identify these points using a technique similar to the ones described in [START_REF] Siddiqi | The Hamilton-Jacobi Skeleton[END_REF][START_REF] Ogniewicz | O: Hierarchic Voronoi Skeletons[END_REF].

In the third step, the obtained skeleton (Fig. 6 c) is pruned of its small, less significant branches by retaining only its points that originate from points on the initial narrow band situated at a distance larger than a given threshold [Ogni95a, Ogni95b, Sidd99]. The above pruning scheme is based on two observations: a) every skeleton point is generated by the collapsing of a compact segment of the original boundary [START_REF] Siddiqi | The Hamilton-Jacobi Skeleton[END_REF][START_REF] Kimmel | Skeletonization vis Distance Maps and Level Sets[END_REF], and b) the importance of a skeleton point can be measured by the length of the boundary segment out of which it originates [START_REF] Ogniewicz | Automatic Medial Axis Pruning by Mapping Characteristics of Boundaries Evolving under the Euclidean Geometric Heat Flow onto Voronoi Skeletons[END_REF][START_REF] Ogniewicz | O: Hierarchic Voronoi Skeletons[END_REF].

In the last step, we 'inflate' the pruned skeleton to ob- To do this, we execute again the FMM algorithm outwards from the skeleton, as follows. We initialize the narrow band to the skeleton points and the function to the value of at those points, where is the distance from the skele- ton to the silhouette, computed in the previous step. We stop the FMM execution when points of the outwards evolving narrow band reach a value of zero.

At that moment, the inflated skeleton matches the initial silhouette (Fig. 6 d). However, due to its pruning, most of the noise of the initial raw silhouette has been removed, as seen in the examples in Fig. 5.

Since the FMM algorithm performs in real time, as noted in Sec. 4.2, the whole skeleton-based smoothing process takes less than a second for our typical images. By adjusting the skeleton pruning threshold, we obtain different smoothing levels. Smoother silhouettes, produced by a higher threshold, lead in practice to a more stable and sensibly faster convergence of the model parameter estimation. Moreover, pruned skeletons typically lead, due to the properties of the Eikonal equation used in the reconstruction, to silhouettes having rounded edges. These shapes are easier appproximated by the superquadric shapes used in our human body model than the raw, arbitrarily shaped silhouettes. However, if the skeletons are pruned too much, the smoothed silhouettes might miss important image cues, such as the orientation of a limb. Conversely, less smoothed silhouettes are closer to the observed data, thus more accurate, but, as mentioned, may lead to numerically unstable derivative estimations. Currently we estimate, by trial and error, a good value for the pruning threshold for a given application configuration (camera parameters, lighting, raw silhouette extraction parameters, optimization method parameters, etc). This value works well for the various images we have tried it on. However, a better strategy we plan to investigate is to use an adaptively optimal threshold for each image.

CONCLUSIONS

We have presented a method to build more consistent likelihood terms for silhouettes, and applied it for human pose estimation in a model based context. Aiming to build cost surfaces whose minima accurately reflect the good configurations in the problem, we define a novel likelihood model composed of an attraction term and an area overlap term which ensures consistent model localization in the image with improved attraction zones. Secondly, we propose a smoothing method for the silhouettes extracted from the image data that stabilizes the optimization process used for pose estimation. Both the likelihood attraction term and silhouette smoothing method are based on distance functions extracted using level-set techniques for evolving boundaries under constant speed in the normal direction. In particular, the fast marching method allows us to calculate distance transforms, skeletons, and to reconstruct silhouettes from their skeletons in a simple to implement and efficient way.

Our future work aims at employing silhouette skeletons, extracted with level set methods, directly as likelihood terms for human pose estimation applications.

Together with this, we aim to develop an automatic procedure of setting the prunung threshold for the skeleton-based smoothing we employ on the imageextracted silhouettes.

Figure 1 :

 1 Figure 1: Human model: flat shaded (a,b) and discretization (c,d)

Figure 2 :

 2 Figure 2: Human pose estimation application pipeline

.

  Various search methods attempt to identify the minima of the function A , by either local continous descent, stochastic search, parameter space subdivision or combinations of them [Smin01b, Deut00, Heap98, Smin01a, Gavr96, Breg98]. All these methods require the evaluation of A . Continuous methods require supplementary evaluations of the first order gradient Q and sometimes the second order Hessian R of A . In this paper, we use a second order local continuous method, where a descent direction is chosen by solving the regularized subproblem [Flet87]: © R S5 UT WV YX ¡ § Q subject to ba Yc ed ¡ gf 'h where:

Figure 3 :

 3 Figure 3: Distance transform computation: original image (a), silhouette (b), distance plot (c) and distance level sets (d)

Figure 4 :

 4 Figure 4: Model estimation based on various silhouette terms original images (a,e), initial models (b,f), silhouette attraction term only (c,g), silhouette attraction and area overlap terms (d,h)

Figure 5 :

 5 Figure 5: Examples of raw silhouettes, skeletons, and smoothed silhouettes

Figure 6 :

 6 Figure 6: Skeleton-based silhouette smoothing pipeline