Viewpoint-invariant image matching and generation of three-dimensional models from two-dimensional imagery

[0001] This application claims priority to and the benefits of U.S. Provisional Applications Ser. Nos. 60/452,429, 60/452,430 and 60/452,431 filed on Mar. 6, 2003 (the entire disclosures of which are hereby incorporated by reference).

[0002] The present invention relates to object modeling and matching systems, and more particularly to the generation of a three-dimensional model of a target object from two- and three-dimensional input.

[0003] In many situations, it is useful to construct a three-dimensional (3D) model of an object when only a partial description of the object is available. In a typical situation, one or more two-dimensional (2D) images of the 3D object may be available, perhaps photographs taken from different viewpoints. A common method of creating a 3D model of a multi-featured object is to start with a base 3D model which describes a generic or typical example of the type of object being modeled, and then to add texture to the model using one or more 2D images of the object. For example, if the multi-featured object is a human face, a 3D “avatar” (i.e., an electronic 3D graphical or pictorial representation) would be generated by using a pre-existing, standard 3D model of a human face, and mapping onto the model a texture from one or more 2D images of the face. See U.S. Pat. No. 6,532,011 B1 to Francini et al., and U.S. Pat. No. 6,434,278 B1 to Hashimoto. The main problem with this approach is that the 3D geometry is not highly defined or tuned for the actual target object which is being generated.

[0004] A common variant of the above approach is to use a set of 3D base models and select the one that most resembles the target object before performing the texture mapping step. Alternatively, a single parameterized base model is used, and the parameters of the model are adjusted to best approximate to the target. See U.S. Pat. No. 6,556,196 B1 to Blanz et al. These methods serve to refine the geometry to make it fit the target, at least to some extent. However, for any target object with a reasonable range of intrinsic variability, the geometry of the model will still not be well tuned to the target. This lack of geometric fit will detract from the verisimilitude of the 3D model to the target object.

[0005] Conventional techniques typically also require that the 2D images being used for texturing the model be acquired from known viewpoints relative to the 3D object being modeled. This usually limits the use of such approaches to situations where the model is being generated in a controlled environment in which the target object can be photographed. Alternatively, resort may be had to human intervention to align 2D images to the 3D model to be generated. See U.S. Patent Publication No. 2002/0012454 to Liu et al. This manual step places a severe limit on the speed with which a 3D model can be generated from 2D imagery.

[0006] Accordingly, a need exists for an automated approach that systematically makes use of available 2D source data for a 3D object to synthesize an optimal 3D model of the object.

[0007] The present invention provides an automated method and system for generating an optimal 3D model of a target multifeatured object when only partial source data describing the object is available. The partial source data often consists of one or more 2D projections of the target object or an obscuration of a single projection, but may also include 3D data, such as from a 3D camera or scanner. The invention uses a set of reference 3D representations that span, to the extent practicable, the variations of the class of objects to which the target object belongs. The invention may automatically identify feature items common to the source data and to the reference representations, and establish correspondences between them. For example, if the target object is a face, the system may identify points at the extremities of the eyes and mouth, or the nose profile, and establish correspondences between such features in the source data and in the reference representations. Manual identification and matching of feature items can also be incorporated if desired. Next, all possible positions (i.e., orientations and translations) for each 3D reference representation are searched to identify the position and reference representation combination whose projection most closely matches the source data. The closeness of match is determined by a measure such as the minimum mean-squared error (MMSE) between the feature items in the projection of the 3D representation, and the corresponding feature items in the source projection. A comparison is performed in 3D between the estimated deprojected positions of the feature items from the 2D source projection and the corresponding feature items of the 3D representation. The closest-fitting 3D reference representation may then be deformed to optimize the correspondence with the source projection. Each point in the mesh which defines the geometry of the 3D representation is free to move during the deformation. The search for the best-fitting position (i.e., orientation and translation) is repeated using the deformed 3D representation, and the deformation and search may be repeated iteratively until convergence occurs or terminated at any time.

[0008] Thus the geometry of the 3D model is tailored to the target object in two ways. First, when more than one reference representation is available, the selection of the best-fitting reference representation from a set of references enables the optimal coarse-grain choice to be made. Second, deformation enables fine scale tuning in which errors introduced by inaccurate choice of viewpoint are progressively reduced by iteration. The invention requires no information about the viewpoint from which the 2D source projection was captured, because a search is performed over all possible viewpoints, and the viewpoint is taken to be that which corresponds to the closest fit between the projected 3D representation and the 2D source data.

[0009] In a first aspect, the invention comprises a method of comparing at least one source 2D projection of a source multifeatured object to a reference library of 3D reference objects. In accordance with the method, a plurality of reference 3D representations of generically similar multifeatured objects is provided, and a viewpoint-invariant search of the reference 3D representations is performed to locate the reference 3D representation having a 2D projection most resembling the source projection(s). In some embodiments, resemblance is determined by a degree of alignment between feature items in the 3D representation and corresponding feature items in the source 2D projection(s). Each reference 3D representation may be searched over a range of possible 2D projections of the 3D representation without actually generating any projections. The search over a range of possible 2D projections may comprise computing a rigid motion of the reference 3D representation optimally consistent with a viewpoint of the source multifeatured object in at least one of the 2D projections. The rigid motions may comprise pitch, roll, yaw, and translation in three dimensions. Automatic camera calibration may be performed by estimation of camera parameters, such as aspect ratio and field of view, from image landmarks.

[0010] In some embodiments, the optimum rigid motion may be determined by estimating a conditional mean pose or geometric registration as it relates to feature items comprising points, curves, surfaces, and subvolumes in a 3D coordinate space associated with the reference 3D representation such that the feature items are projectionally consistent with feature items in source 2D projection(s). MMSE estimates between the conditional mean estimate of the projected feature items and corresponding feature items of the reference 3D representation are generated. The rigid motion may be constrained by known 3D position information associated with the source 2D projection(s).

[0011] In some embodiments, the feature items may include curves as well as points which are extracted from the source projection using dynamic programming. Further, areas as well as surfaces and or subvolumes may be used as features generated via isocontouring (such as via the Marching Cubes algorithm) or automated segmentation algorithms. The feature items used in the matching process may be found automatically by using correspondences between the 2D source projection(s) and projected imagery of at least one reference 3D object.

[0012] The invention may further comprise the step of creating a 3D representation of the source 2D projection(s) by deforming the located (i.e., best-fitting) reference 3D representation so as to resemble the source multifeatured object. In one embodiment, the deformation is a large deformation diffeomorphism, which serves to preserve the geometry and topology of the reference 3D representation. The deformation step may deform the located 3D representation so that feature items in the source 2D projection(s) align with corresponding features in the located reference 3D representation. The deformation step may occur with or without rigid motions and may include affine motions. Further, the deformation step may be constrained by at least one of known 3D position information associated with the source 2D projection(s), and 3D data of the source object. The deformation may be performed using a closed form expression.

[0013] In a second aspect, the invention comprises a system for comparing at least one source 2D projection of a source multifeatured object to a reference library of 3D reference objects. The system comprises a database comprising a plurality of reference 3D representations of generically similar multifeatured objects and an analyzer for performing a viewpoint-invariant search of the reference 3D representations to locate the reference 3D representation having a 2D projection most resembling the source projection(s). In some embodiments, the analyzer determines resemblance by a degree of alignment between feature items in the 3D representation and corresponding feature items in the source 2D projection(s). The analyzer may search each reference 3D representation over a range of possible 2D projections of the 3D representation without actually generating any projections. In some embodiments, the analyzer searches over a range of possible 2D projections by computing a rigid motion of the reference 3D representation optimally consistent with a viewpoint of the source multifeatured object in at least one of the 2D projections. The rigid motions may comprise pitch, roll, yaw, and translation in three dimensions. The analyzer may be configured to perform automatic camera calibration by estimating camera parameters, such as aspect ratio and field of view, from image landmarks.

[0014] In some embodiments, the analyzer is configured to determine the optimum rigid motion by estimating a conditional mean of feature items comprising points, curves, surfaces, and subvolumes in a 3D coordinate space associated with the reference 3D representation such that the feature items are projectionally consistent with feature items in the source 2D projection(s). The analyzer is further configured to generate MMSE estimates between the conditional mean estimate of the projected feature items and corresponding feature items of the reference 3D representation. The rigid motion may be constrained by known 3D position information associated with the source 2D projection(s).

[0015] In some embodiments, the analyzer is configured to extract feature items from the source projection using dynamic programming. In further embodiments, the analyzer may be configured to find feature items used in the matching process automatically by using correspondences between source imagery and projected imagery of at least one reference 3D object.

[0016] The invention may further comprise a deformation module for creating a 3D representation of the at least one source 2D projection by deforming the located (i.e., best-fitting) reference 3D representation so as to resemble the source multifeatured object. In one embodiment, the deformation module deforms the located reference 3D representation using large deformation diffeomorphism, which serves to preserve the geometry and topology of the reference 3D representation. The deformation module may deform the located 3D representation so that feature items in the source 2D projection(s) align with corresponding features in the located reference 3D representation. The deformation module may or may not use rigid motions and may use affine motions. Further, the deformation module may be constrained by at least one of known 3D position information associated with the source 2D projection(s), and 3D data of the source object. The deformation module may operate in accordance with a closed form expression.

[0017] In a third aspect, the invention comprises a method of comparing a source 3D object to at least one reference 3D object. The method involves creating 2D representations of the source object and the reference object(s) and using projective geometry to characterize a correspondence between the source 3D object and a reference 3D object. For example, the correspondence may be characterized by a particular viewpoint for the 2D representation of the 3D source object.

[0018] In a fourth aspect, the invention comprises a system for comparing a source 3D object to at least one reference 3D object. The system comprises a projection module for creating 2D representations of the source object and the reference object(s) and an analyzer which uses projective geometry to characterize a correspondence between the source 3D object and a reference 3D object.

[0019] In a fifth aspect, the above described methods and systems are used for the case when the 3D object is a face and the reference 3D representations are avatars.

[0020] In a sixth aspect, the invention comprises a method for creating a 3D representation from at least one source 2D projection of a source multifeatured object. In accordance with the method, at least one reference 3D representation of a generically similar object is provided, one of the provided representation(s) is located, and a 3D representation of the source 2D projection(s) is created by deforming the located reference representation in accordance with the source 2D projection(s) so as to resemble the source multifeatured object. In some embodiments, the source 2D projection(s) is used to locate the reference representation. In further embodiments, the set of reference representations includes more than one member, and the reference most resembling the source 2D projection(s) is located by performing a viewpoint-invariant search of the set of reference representations, without necessarily actually generating any projections. The search may include computing a rigid motion of the reference representation optimally consistent with a viewpoint of the source multifeatured object in at least one of the source projections.

[0021] In a preferred embodiment, a 3D representation of the source projection(s) is created by deforming the located reference representation so as to resemble the source multifeatured object. The deformation may be a large deformation diffeomorphism. In some embodiments, the deformation deforms the located reference so that feature items in the source projection(s) align with corresponding feature items in the located 3D reference representation. In some embodiments, the deformation is performed in real time.

[0022] In a seventh aspect, the invention comprises a system for creating a 3D representation from at least one source 2D projection of a source multifeatured object. The system includes a database of at least one reference 3D representation of a generically similar object, and an analyzer for locating one of the provided representation(s). The system further includes a deformation module for creating a 3D representation of the source 2D projection(s) by deforming the located reference representation in accordance with the source 2D projection(s) so as to resemble the source multifeatured object. In some embodiments, the analyzer uses the source 2D projection(s) to locate the reference representation. In further embodiments, the set of reference representations includes more than one member, and the analyzer locates the reference most resembling the source 2D projection(s) by performing a viewpoint-invariant search of the set of reference representations, without necessarily actually generating any projections. The search may include computing a rigid motion of the reference representation optimally consistent with a viewpoint of the source multifeatured object in at least one of the source projections.

[0023] In a preferred embodiment, the deformation module creates a 3D representation of the source projection(s) by deforming the located reference representation so as to resemble the source multifeatured object. The deformation may be a large deformation diffeomorphism. In some embodiments, the deformation module deforms the located reference so that feature items in the source projection(s) align with corresponding feature items in the located 3D reference representation. In some embodiments, the deformation module operates in real time.

[0031] FIG. 1 illustrates the basic operation of the invention in the case where the 3D target multifeatured object is a face and the set of reference 3D representations are avatars. The matching process starts with a set of reference 3D avatars which represent, to the extent practicable, the range of different types of heads to be matched. For example, the avatars may include faces of men and women, faces with varying quantities and types of hair, faces of different ages, and faces representing different races. Typically the reference set includes numerous (e.g., several hundred or more) avatars, though the invention works with as few as one reference object, and as many for which storage space and compute time are available. In some situations, only a single reference avatar will be used. This case may arise, for example, when the best-fitting avatar has been selected manually, or by some other means, or when only one reference avatar is available. In FIG. 1, the source data of the target face is illustrated as a single 2D photograph of the target 3D face taken from an unknown viewpoint. First, selected features to be used for matching are identified in the source photograph. The features may be points, such as the extremities of the mouth and eyes, or curves representing a profile, an eyebrow, or other distinctive curve, or subareas such as an eyebrow, nose or cheek. The corresponding features are identified in the reference avatars. The selection of features in the target photograph, and the identification of corresponding features in the reference avatars may be done automatically according to the invention. Next, a viewpoint-invariant search is conducted in which each 3D avatar is notionally subjected to all possible rigid motions, and the features projected into 2D. The positions of the projected avatar features are compared to the feature positions in the target photograph. The avatar for which a particular rigid motion provides the closest fit between projected features and those of the source photograph is selected as the best reference avatar. FIG. 1 illustrates the best-fitting reference avatar to be the middle one.

[0032] Next, the best reference avatar is deformed to match the target photograph more closely. First the features of the photograph are reverse projected to the coordinates of the best reference avatar in the orientation and position corresponding to the best match. The mesh points of the avatar are then deformed in 3D to minimize the distances between the reverse-projected features of the photograph and the corresponding avatar features. The avatar resulting from this deformation will be a closer approximation to the target 3D face. The rigid motion search and deformation steps may be repeated iteratively, e.g., until the quality of fit no longer improves appreciably. The resulting 3D model is the optimal match to the target face.

[0033] The invention can be used effectively even when the source imagery includes only a part of the target face, or when the target face is partially obscured, such as, for example, by sunglasses or facial hair. The approach of the invention is suitable for any multifeatured object, such as faces, animals, plants, or buildings. For ease of explanation, however, the ensuing description will focus on faces as an exemplary (and non-limiting) application.

[0034] FIGS. 2A, 2B, and 2C show the components of a representative avatar. In one embodiment of the invention, the geometry of the avatar is represented by a mesh of points in 3D which are the vertices of set of triangular polygons approximating the surface of the avatar. FIG. 2A illustrates a head-on view 202 and a side view 204 of the triangular polygon representation. In one representation, each vertex is given a color value, and each triangular face may be colored according to an average of the color values assigned to its vertices. The color values are determined from a 2D texture map 206, illustrated in FIG. 2B, which may be derived from a photograph. In FIG. 2C, the final avatar with texture is illustrated in a head-on view 208 and side view 210. The avatar is associated with a coordinate system which is fixed to it, and is indexed by three angular degrees of freedom (pitch, roll, and yaw), and three translational degrees of freedom of the rigid body center in three-space. In addition, individual features of the avatar, such as the chin, teeth and eyes may have their own local coordinates (e.g., chin axis) which form part of the avatar description. The present invention may be equally applied to avatars for which a different data representation is used. For example, texture values may be represented as RGB values, or using other color representations, such as HSL. The data representing the avatar vertices and the relationships among the vertices may vary. For example, the mesh points may be connected to form non-triangular polygons representing the avatar surface.

[0035] The invention may include a conventional rendering engine for generating 2D imagery from a 3D avatar. The rendering engine may be implemented in OpenGL, or in any other 3D rendering system, and allows for the rapid projection of a 3D avatar into a 2D image plane representing a camera view of the 3D avatar. The rendering engine may also include the specification of the avatar lighting, allowing for the generation of 2D projections corresponding to varying illumination of the avatar. Lighting corresponding to a varying number of light sources of varying colors, intensities, and positions may be generated.

[0036] The feature items in the 2D source projection which are used for matching are selected by hand or via automated methods. These items may be points or curves. When the source projection includes a front view, suitable points may be inflection points at the lips, points on the eyebrows, points at the extremities of the eyes, or extremities of nostrils, and suitable curves may include an eyebrow or lip. When the source projection includes a side view, the feature points corresponding to the profile are used and may include the tip of the nose or chin. Suitable feature item curves may include distinct parts of the profile, such as nose, forehead, or chin.

[0037] When the feature items are determined manually, a user interface is provided which allows the user to identify feature points individually or to mark groups of points delineated by a spline curve, or to select a set of points forming a line.

[0038] The automated detection of feature items on the 2D source projection is performed by searching for specific features of a face, such as eyeballs, nostrils, and lips. As understood by those of ordinary skill in the art, the approach may use Bayesian classifiers and decision trees in which hierarchical detection probes are built from training data generated from actual avatars. The detection probes are desirably stored at multiple pixel scales so that the specific parameters, such as for orientation of a feature, are only computing on finer scales if the larger-scale probes yield a positive detection. The feature detection probes may be generated from image databases representing large numbers of individuals who have had their features demarcated and segregated so that the detection probes become specifically tuned to these features. The automated feature detection approach may use pattern classification, Bayes nets, neural networks, or other known techniques for determining the location of features in facial images.

[0039] The automated detection of curve features in the source projection may use dynamic programming approaches to generate curves from a series of points so as to reduce the amount of computation required to identify an optimal curve and maximize a sequentially additive cost function. Such a cost function represents a sequence of features such as the contrast of the profile against background, or the darkness of an eyebrow, or the crease between lips. A path of N points can be thought of as consisting of a starting node x₀ and a set of vectors v₀, v₁, . . . , v_N−1 connecting neighboring nodes. The nodes comprising this path are defined as

[0040] $x_i=\sum _{j=0}^{i-1}{v}_j+x_0.$
Rather than searching over all paths of length N, dynamic programming may be used to generate maximum (or minimum) cost paths. This reduces the complexity of the algorithm from K^N to NK² where N is the length of a path and K is the total number of nodes, as dynamic programming takes advantage of the fact that the cost is sequentially additive, allowing a host of sub-optimal paths to be ignored. Dynamic programming techniques and systems are well-characterized in the art and can be applied as discussed herein without undue experimentation.

[0041] Next, the 3D rotation and translation from the avatar coordinates to the source projection is determined. This corresponds to finding the viewpoint from which the source projection was captured. In preferred embodiments, this is achieved by calculating the position of the avatar in 3D space that best matches the set of selected feature items in the 2D source projection. Generally, these feature items will be points, curves, or subareas and the source projection will be a photograph on which the position of these items can be measured, either manually or automatically. The position calculation may be based on the computation of the conditional mean estimate of the reverse projection positions in 3D of the 2D feature items, followed by the computation of MMSE estimates for the rotation and translation parameters in 3D, given the estimates of the 3D positions of the feature items. Since position in 3D space is a vector parameter, the MMSE estimate for translation position is closed form; when substituted back into the squared error function, it gives an explicit function in terms of only the rotations. Since the rotations are not vector parameters, they may be calculated using non-linear gradient descent through the tangent space of the group or via local representation using the angular velocities of the skew-symmetric matrices.

[0042] In addition to or in lieu of the least squares or weighted least squares techniques described herein, the distance metrics used to measure the quality of fit between the reverse projections of feature items from the source imagery and corresponding items in the 3D avatar may be, for example, Poisson or other distance metrics which may or may not satisfy the triangle inequality.

[0043] If feature items measured in 3D are available, such as from actual 3D source data from 3D cameras or scanners, the feature item matching may be performed directly, without the intermediate step of calculating the conditional mean estimate of the deprojected 2D features. The cost function used for positioning the 3D avatar can be minimized using algorithms such as closed form quadratic optimization, iterative Newton descent or gradient methods.

[0044] The 3D positioning technique is first considered without deformation of the reference avatar. In the following, a 3D reference avatar is referred to as a CAD (computer-aided design) model, or by the symbol CAD. The set of x_j=(x_j,y_j,z_j), j=1, . . . , N features is defined on the CAD model. The projective geometry mapping is defined as either positive or negative z, i.e., projection occurs along the z axis. In all the projective geometry

[0045] $p_j=\left(\frac{{\alpha }_1{x}_j}{-z_j},\frac{{\alpha }_2{y}_j}{-z_j}\right)$
(for negative z-axis projection), or

[0046] $p_j=\left(\frac{{\alpha }_1{x}_j}{z_j},\frac{{\alpha }_2{y}_j}{z_j}\right)$
(for positive z-axis projection) is the projected position of the point x_j where α is the projection angle. Let the rigid transformation be of the form A=O,b:xOx+b centered around x_c=0. For positive (i.e., z>0) mapping and n=0,

[0047] $p_i=\left(\frac{{\alpha }_1{x}_i}{z_i},\frac{{\alpha }_2{y}_i}{z_i}\right),i=1,\dots ,N,$
where n is the cotangent of the projective angle. The following data structures are defined throughout:

[0048] $\begin{array}{cc}{P}_i=\left(\begin{array}{c}\begin{array}{c}\frac{p_{i1}}{{\alpha }_1}\\ \frac{p_{i2}}{{\alpha }_2}\end{array}\\ 1\end{array}\right),Q_i=\left(I-\frac{{P_i\left(P_i\right)}^t}{{P_i}^2}\right),\stackrel{\_}{Q}=\sum _{i=1}^N{Q}_i,X_Q=\sum _{i=1}^N{Q}_i{X}_i& \left(\mathrm{Equation}1\right)\end{array}$

[0049] $\begin{array}{cc}{X}_j=\left(\begin{array}{ccccccccc}{x}_j^1& x_j^2& x_j^3& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& x_j^1& x_j^2& x_j^3& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& x_j^1& x_j^2& x_j^3\end{array}\right),& \left(\mathrm{Equation}2\right)\end{array}$
with (·)^t matrix transpose, and the identity matrix

[0050] $I=\left(\begin{array}{ccc}1& 0& 0\\ 0& 1& 0\\ 0& 0& 1\end{array}\right).$
For negative (i.e., z<0) mapping,

[0051] $p_j=\left(\frac{{\alpha }_1{x}_j}{-z_j},\frac{{\alpha }_2{y}_j}{-z_j}\right),$
and the change

[0052] $P_i={\left(-\frac{p_{i1}}{{\alpha }_1},-\frac{p_{i2}}{{\alpha }_2},1\right)}^t$
is made.

[0053] The basis vectors Z₁, Z₂, Z₃ at the tangent to the 3×3 rotation element O are defined as:
Z₁=1₁O^old=[o₂₁, o₂₂, o₂₃, −o₁₁, −o₁₂, −o₃₁, 0, 0, 0]^t  (Equation 3)
Z₂=1₂O^old=[o₃₁, o₃₂, o₃₃, 0, 0, 0, −o₁₁, −o₁₂, −o₁₃]^t  (Equation 4)
Z₃=1₃O_old=[0, 0, 0, o₃₁, o₃₂, o₃₃, −o₂₁, −o₂₂, −o₂₃]^t  (Equation 5)

[0054] $\begin{array}{cc}\mathrm{where}{1}_1=\left(\begin{array}{ccc}0& 1& 0\\ -1& 0& 0\\ 0& 0& 0\end{array}\right),1_2=\left(\begin{array}{ccc}0& 0& 1\\ 0& 0& 0\\ -1& 0& 0\end{array}\right),1_3=\left(\begin{array}{ccc}0& 0& 0\\ 0& 0& 1\\ 0& -1& 0\end{array}\right).& \left(\mathrm{Equation}6\right)\end{array}$

[0055] The reverse projection of feature points from the 2D projection may now be performed. Given the feature points p_j=(p_j1,p_j2), j=1, 2, . . . in the image plane, the minimum norm-estimates for {circumflex over (z)}_i are given by,

[0056] $\hat{O},\hat{b}\mathrm{as}{\hat{z}}_i=\frac{<O{x}_i+b,P_i>}{{P_i}^2},$
and the MMSE Ô,{circumflex over (b)} satisfies

[0057] $\begin{array}{cc}\underset{z,O,b}{\min }\sum _{i=1}^N{O{x}_i+b-z_i{P}_i}^2=\underset{O,b}{\min }\sum _{i=1}^N{\left(O{x}_i+b\right)}^t{Q}_i\left(O{x}_i+b\right).& \left(\mathrm{Equation}7\right)\end{array}$

[0058] During the process of matching a source image to reference avatars, there may be uncertainty in the determined points x, implying that cost matching is performed with a covariance variability structure built in to the formula. In this case, the norm has within it a 3×3 matrix which represents this variability in the norm.

[0059] The optimum rotation and translation may next be estimated from feature points. Given the projective points p_j, j=1, 2, . . . , the rigid transformation has the form O,b: xOx+b (centered around center x_c=0). Then for positive (z>0) mapping and n=0,

[0060] $p_j=\left(\frac{{\alpha }_1{x}_j}{z_j},\frac{{\alpha }_2{y}_j}{z_j}\right),$

[0061] $\begin{array}{cc}\underset{z,O,b}{\min }\sum _{i=1}^N{O{x}_i+b-z_i{P}_i}^2=\underset{O,b}{\min }\sum _{i=1}^N{\left(O{x}_i+b\right)}^t{Q}_i\left(O{x}_i+b\right).& \left(\mathrm{Equation}8\right)\end{array}$

[0062] The optimum translation/rotation solutions are preferably generated as follows. Compute the 3×9 matrix M_i=X_i− Q⁻¹X_Q and evaluate exhaustively the cost function choosing minimum Ô and computing the translation {circumflex over (b)}=−(Σ_i=1^NQ_i)⁻¹Σ_i=1^NQ_iÔx_i with minimum Ô attained, for example, via brute force search over the orthogonal group which may for example be parameterized by pitch, roll or yaw or running the gradient search algorithm to convergence as follows.

[0063] $\begin{array}{cc}\mathrm{Brute}\mathrm{Force}\text{:}\hat{O}=\arg \underset{O}{\min }{O}^t\left(\sum _{i=1}^N{M}_i^t{Q}_i{M}_i\right)O;& \left(\mathrm{Equation}9\right)\\ \begin{array}{c}\mathrm{Gradient}\text{:}{O}^{\mathrm{new}}=e^{\sum _{i=1}^3{\alpha }_i^{\mathrm{new}}{1}_i}{O}^{\mathrm{old}},{\alpha }_j^{\mathrm{new}}\\ =<2\left(\sum _{i=1}^N{M}_i^t{Q}_i{M}_i\right)O^{\mathrm{old}},Z_j>,j=1,2,3\\ \mathrm{with}<f,g>=\sum _{i=1}^3{f}_i{g}_i.\end{array}& \left(\mathrm{Equation}10\right)\end{array}$

[0064] In a typical application some information about the position of the object in 3D space is known. For example, in a system which takes a succession of photographs of a moving source object, such as in a tracking system, the position from a previous image may be available. The invention may incorporate this information into the matching process as follows. Given a sequence of points p_i, i=1, . . . , N and a rigid transformation of the form A=O,b:xOx+b (centered around x_c=0), the MMSE of rotation and translation Ô,{circumflex over (b)} satisfies:

[0065] $\begin{array}{cc}\underset{z,O,b}{\min }\sum _{i=1}^N{O{x}_i+b-z_i{P}_i}^2+{\left(b-\mu \right)}^t{\Sigma }^{-1}\left(b-\mu \right)=\underset{O,b}{\min }\sum _{i=1}^N{\left(O{x}_i+b\right)}^t{Q}_i\left(O{x}_i+b\right)+{\left(b-\mu \right)}^t{\Sigma }^{-1}\left(b-\mu \right).& \left(\mathrm{Equation}11\right)\end{array}$
The 3×9 matrix M_i and a 3×1 column vector N are computed:
=M_i=X_i− Q_Σ⁻¹X_Q, N= Q_Σ⁻¹X_Q, Q_Σ=( Q+Σ⁻¹), ψ= Q_Σ⁻¹Σ_μ, φ= Q_Σ⁻¹Σ_μ−μ, Σ_μ=Σ⁻¹μ.  (Equation 12)
The translation is then determined {circumflex over (b)}=− Q_Σ⁻¹X_QÔ+ Q_Σ⁻¹Σ_μ at the minimum Ô obtained by exhaustive search or gradient algorithm run until convergence:

[0066] $\begin{array}{cc}\begin{array}{c}\mathrm{Brute}\mathrm{Force}\text{:}\hat{O}=\arg \underset{O}{\min }{O}^t\left(\sum _{i=1}^N{M}_i^t{Q}_i{M}_i+N^t{\Sigma }^{-1}N\right)O+\\ 2O^t\left(\sum _{i=1}^N{M}_i^t{Q}_i\psi -N{\Sigma }^{-1}\varphi \right)\end{array}& \left(\mathrm{Equation}13\right)\\ \begin{array}{c}\mathrm{Gradient}\text{:}{O}^{\mathrm{new}}=e^{\sum _{i=1}^3{\alpha }_i^{\mathrm{new}}{1}_i}{O}^{\mathrm{old}},{\alpha }_j^{\mathrm{new}}\\ =<2\left(\sum _{i=1}^N{M}_i^t{Q}_i{M}_i+N^t{\Sigma }^{-1}N\right)O^{\mathrm{old}}+\\ 2\sum _{i=1}^N{M}_i^t{Q}_i\psi -N^t{\Sigma }^{-1}\varphi ,Z_j>\end{array}& \left(\mathrm{Equation}14\right)\end{array}$
with the projection onto the basis vectors Z₁, Z₂, Z₃ of equations 3-5 defined at the tangent to O^old in the exponential representation where α^new are the directional derivatives of cost.

[0067] The rotation/translation data may be indexed in many different ways. For example, to index according to the rotation around the center of the object, rather than fixed external world coordinates, the coordinates are just reparameterized by defining {tilde over (x)}←x−x_c. All of the techniques described herein remain the same.

[0068] The preferred 3D algorithm for rigid motion efficiently changes states of geometric pose for comparison to the measured imagery. The preferred 3D algorithm for diffeomorphic transformation of geometry matches the geometry to target 2D image features. It should be understood, however, that other methods of performing the comparison of the 3D representation to the source imagery may be used, including those that do not make use of specific image features.

[0069] Once the rigid motion (i.e., rotation and translation) that results in the best fit between 2D source imagery and a selected 3D avatar is determined, the 3D avatar may be deformed in order to improve its correspondence with the source imagery. The allowed deformations are generally limited to diffeomorphisms of the original avatar. This serves to preserve the avatar topology, guaranteeing that the result of the deformation will be a face. The deformations may also enforce topological constraints, such as the symmetry of the geometry. This constraint is especially useful in situations where parts of the source object are obscured, and the full geometry is inferred from partial source information.

[0070] FIGS. 3A and 3B illustrate the effect of avatar deformation on the matching of the avatar to the source imagery. In FIG. 3A, feature points are shown as black crosses on the source image 302. An example is the feature point at the left extremity of the left eyebrow 304. The projections of the corresponding feature points belonging to the best-matching reference avatar with optimal rigid motion prior to deformation are shown as white crosses. It can be seen that the projected point corresponding to the left extremity of the left eyebrow 306 is noticeably displaced from its counterpart 304. In FIG. 3B, the same source image 302 is shown with feature points again indicated by black crosses. This time, the best-fitting avatar feature points shown as white crosses are now projected after deformation. The correspondence between source feature points and avatar feature points is markedly improved, as shown, for example, by the improved proximity of the projected left eyebrow feature point 308 to its source counterpart 304.

[0071] The 3D avatar diffeomorphism calculation starts with the initial conditions for placement of the avatar determined by the feature item detection and computation of the best-fitting rigid motion and the original geometry of the avatar. It then proceeds by allowing all of the points on the avatar to move independently according to a predefined formula to minimize the distance between the deformed avatar points in 3D and the conditional mean estimates of the 2D landmark points reverse projected to the 3D coordinates. Once this diffeomorphism is calculated, the 3D landmark rigid motion algorithm is applied again to the source projections and feature items to find the best guess of the camera positions given this newly transformed avatar with its new vertex positions. Subsequently, a new diffeomorphism is generated, and this process is continued until it converges. Alternatively, iteration may not be used, with the rigid motion calculation being performed only a single time, and just one diffeomorphism transformation applied. In the case where camera orientations (i.e., the viewpoint of the measured source projections) are known precisely, these can be used as fixed inputs to the calculation, with no rigid transformation required. When the measured sets of feature items are in 3D, such as from a cyber scan or 3D camera observations of the candidate head, the avatar is transformed onto the candidate sets of points without any intermediate generation of the candidate points in 3D space via the conditional mean algorithm for generating 3D points from 2D sets of points.

[0072] The diffeomorphic deformation of an avatar proceeds as follows. Given the set x_j=(x_j,y_j,z_j), j=1, . . . , N feature items defined on the CAD model, with the projective geometry mapping with

[0073] ${\alpha }_1=\frac{2n}{w},{\alpha }_2=\frac{2n}{h},$
where n is cotangent of the angle, and w, h are aspect ratio widths and heights,

[0074] $\left(x,y,z\right)\mapsto p\left(x,y,z\right)=\left(\frac{{\alpha }_{1}x}{z},\frac{{\alpha }_{2}y}{z}\right),$
with observations of the feature items through the projective geometry

[0075] $p_j=\left(\frac{{\alpha }_1{x}_j}{z_j},\frac{{\alpha }_2{y}_j}{z_j}\right).$
The goal is to construct the deformation of the CAD model x→x+u(x), xεCAD with unknown camera rigid motions corresponding to the measured projective image feature items. The projective points for each orientation v=1, . . . , V, and smoothing matrices

[0076] $\begin{array}{cc}\begin{array}{c}{P}_i^{\left(v\right)}=\left(\begin{array}{c}\begin{array}{c}\frac{p_{\mathrm{i1}}^{\left(v\right)}}{{\alpha }_1}\\ \frac{p_{\mathrm{i2}}^{\left(v\right)}}{{\alpha }_2}\end{array}\\ 1\end{array}\right),Q_i^v=I-\frac{P_i^{\left(v\right)}{P}_i^{\left(v\right)t}}{{P_i^{\left(v\right)}}^2},\\ {\stackrel{\_}{Q}}_i=\sum _{v=1}^V{O}^{\left(v\right)t}{Q}_i^{\left(v\right)}{O}^{\left(v\right)},{\stackrel{\_}{Q}}_{Oi}=\sum _{v=1}^V{O}^{\left(v\right)t}{Q}_i^{\left(v\right)},\end{array}& \left(\mathrm{Equation}15\right)\\ \stackrel{\_}{Q}=\underset{\underset{3N\times 3N}{︸}}{\left(\begin{array}{cccc}{\stackrel{\_}{Q}}_1& 0& \cdots & 0\\ 0& {\stackrel{\_}{Q}}_2& 0& \cdots \\ ⋮& ⋮& ⋮& ⋮\\ 0& 0& \cdots & {\stackrel{\_}{Q}}_N\end{array}\right)},{\stackrel{\_}{Q}}_o=\underset{\underset{3N\times 3N}{︸}}{\left(\begin{array}{cccc}{\stackrel{\_}{Q}}_{O_1}& 0& \cdots & 0\\ 0& {\stackrel{\_}{Q}}_{O_2}& 0& \cdots \\ ⋮& ⋮& ⋮& ⋮\\ 0& 0& \cdots & {\stackrel{\_}{Q}}_{O_N}\end{array}\right)},& \left(\mathrm{Equation}16\right)\\ K=\underset{\underset{3N\times 3N}{︸}}{\left(\begin{array}{cccc}{K}_{11}& K_{12}\dots & K_{1N}& \\ K_{21}& K_{21}& \cdots & K_{2N}\\ ⋮& ⋮& ⋮& ⋮\\ K_{\mathrm{N1}}& K_{\mathrm{N1}}& \cdots & K_{\mathrm{NN}}\end{array}\right)},K_{ij}=\left(\begin{array}{ccc}{K}_{ij}\left(1,1\right)& K_{ij}\left(1,2\right)& K_{ij}\left(1,3\right)\\ K_{ij}\left(2,1\right)& K_{ij}\left(2,2\right)& K_{ij}\left(2,3\right)\\ K_{ij}\left(3,1\right)& K_{ij}\left(3,2\right)& K_{ij}\left(3,3\right)\end{array}\right),& \left(\mathrm{Equation}17\right)\end{array}$
are constructed, where for example K_ij=diag(e^{−α∥xi−xj∥},e^{−α∥xi−xj∥},e^{−α∥xi−xj∥}) corresponds to the square root inverse Laplacian operator L=diag(−∇²+c).

[0077] In one embodiment, the avatar may be deformed with small deformations only and no rigid motions. For this embodiment, it is assumed that the measured feature items are all points from a single camera viewing which generated the projected source image in which the feature points were measured. The goal is to construct the deformation of the CAD model constructing the mapping x→x+u(X), xεCAD:

[0078] $\begin{array}{cc}\underset{u,z_n}{\min }{\mathrm{Lu}}^2+\sum _{n=1}^N{\left(x_n+u\left(x_n\right)\right)-z_n{P}_n}^2=\underset{u}{\min }{\mathrm{Lu}}^2+\sum _{n=1}^N{\left(x_n+u\left(x_n\right)\right)}^t{Q}_n\left(x_n+u\left(x_n\right)\right).& \left(\mathrm{Equation}18\right)\end{array}$

[0079] First, the transformation of the model x→x+u(x) with u(x)=Σ_n=1^NK(x_n,x)β_n and where

[0080] $\begin{array}{cc}\left(\begin{array}{c}{\beta }_1\\ ⋮\\ {\beta }_N\end{array}\right)=K^{-1}\left(-\frac{1}{{\sigma }^2}{\left(K^{-1}+\frac{1}{{\sigma }^2}\stackrel{\_}{Q}\right)}^{-1}\left(\stackrel{\_}{Q}\left(\begin{array}{c}{x}_1\\ ⋮\\ x_N\end{array}\right)\right)\right)& \left(\mathrm{Equation}19\right)\end{array}$
is computed. Next, rigid motions are added and the following equation solved for the optimizer:

[0081] $\begin{array}{cc}\underset{u}{\min }{\mathrm{Lu}}^2+\sum _{n=1}^N\sum _{v=1}^V(O^{\left(v\right)}\left(x_n+u\left(x_n\right)\right)+{\hspace{1em}{b}^{\left(v\right)})}^t\hspace{1em}{Q}_n^{\left(v\right)}\left(O^{\left(v\right)}\left(x_n+u\left(x_n\right)\right)+b^{\left(v\right)}\right).& \left(\mathrm{Equation}20\right)\end{array}$
The transformation of the model using small deformation x→x+u(x) is computed, where

[0082] $\begin{array}{cc}u\left(x\right)=\sum _{n=1}^{N}K\left(x_n,x\right){\beta }_n\mathrm{and}\begin{array}{c}\left(\begin{array}{c}{\beta }_1\left(44\right)\\ ⋮\\ {\beta }_N\end{array}\right)=K^{-1}(-\frac{1}{{\sigma }^2}{\left(K^{-1}+\frac{1}{{\sigma }^2}\stackrel{\_}{Q}\right)}^{-1}(\stackrel{\_}{Q}\left(\begin{array}{c}{x}_1\\ ⋮\\ x_N\end{array}\right)+\\ \left(\begin{array}{c}\sum _{v=1}^V{O}^{\left(v\right)t}{Q}_1^{\left(v\right)}\left(b^{\left(v\right)}\right)\\ ⋮\\ \sum _{v=1}^V{O}^{\left(v\right)t}{Q}_N^{\left(v\right)}\left(b^{\left(v\right)}\right)\end{array}\right))),\end{array}& \left(\mathrm{Equation}21\right)\end{array}$

[0083] In another embodiment, diffeomorphic deformations with no rigid motions of the avatar are applied. In the case that the change in shape of the face is extensive the large deformation φ:xφ(x) satisfying φ=φ₁, φ_t=∫₀¹v_s(φ_s(x))ds+x, xεCAD is generated. The deformation of the CAD model constructing the mapping xφ(x), xεCAD is constructed:

[0084] $\begin{array}{cc}\underset{v_1,t\in \left[0,1\right],z_n}{\min }{\int }_0^1{{\mathrm{Lv}}_t}^{2}dt+\sum _{n=1}^N{\varphi \left(x_n\right)-z_n{P}_n}^2={\min }_{v_1,t\in \left[0,1\right]}{\int }_0^1{{\mathrm{Lv}}_t}^{2}dt+\sum _{n=1}^N{\varphi \left(x_n\right)}^t{Q}_n\varphi \left(x_n\right).& \left(\mathrm{Equation}22\right)\end{array}$
Using the initialization v^new=0, φ^new(x)=x, xεCAD, mappings are repeatedly generated from the new vector field by running an iterative algorithm until convergence:

[0085] $\begin{array}{cc}{v}_t^{\mathrm{new}}(·)=\sum _{n=1}^{N}K\left({\varphi }_t^{\mathrm{new}}\left(x_n\right),·\right)D_{{\varphi }_t^{\mathrm{new}}\left(x_n\right)}{\varphi }_{t,1}^{\mathrm{new}}{Q}_n{\varphi }^{\mathrm{new}}\left(x_n\right),& \left(\mathrm{Equation}23\right)\\ {\varphi }^{\mathrm{new}}\left(x\right)={\int }_0^1{v}_t^{\mathrm{new}}\left({\varphi }_t^{\mathrm{new}}\left(x\right)\right)dx+x,& \left(\mathrm{Equation}24\right)\end{array}$
where

[0086] $D_{{\varphi }_t\left(x_n\right)}{\varphi }_{t,1}=\left(\frac{\partial {\varphi }_1\circ {\varphi }_t^{-1}}{\partial y_j}{❘}_{{\varphi }_t\left(x_n\right)}\right).$
The addition rigid motions to the large deformation.
xφ(x), xεCAD is accomplished as follows:

[0087] $\begin{array}{cc}\underset{v_t:{\stackrel{.}{\varphi }}_t=v_t\left({\varphi }_t\right),t\in \left[0,1\right],z_n}{\min }{\int }_0^1{{\mathrm{Lv}}_t}^{2}dt+\sum _{n=1}^N\sum _{v=1}^V{O^{\left(v\right)}\varphi \left(x_n\right)+b^{\left(v\right)}-z_n^{\left(v\right)}{P}_n^{\left(v\right)}}^2=\underset{v_t,t\in \left[0,1\right]}{\min }{\int }_0^1{{\mathrm{Lv}}_t}^{2}dt=\sum _{n=1}^N\sum _{v=1}^V{\left(O^{\left(v\right)}\varphi \left(x_n\right)+b^{\left(v\right)}\right)}^t{Q}_n^{\left(v\right)}\left(O^{\left(v\right)}\varphi \left(x_n\right)+b^{\left(v\right)}\right).& \left(\mathrm{Equation}25\right)\end{array}$
Using the initialization v^new=0, φ^new(x)=x, xεCAD, a mapping is generated from the new vector field by running an iterative algorithm until convergence:

[0088] $\begin{array}{cc}{v}_t^{\mathrm{new}}(·)=\sum _{n=1}^{N}K\left({\varphi }_t^{\mathrm{new}}\left(x_n\right),·\right)D_{{\varphi }_t^{\mathrm{new}}\left(x_n\right)}{\varphi }_{t,1}^{\mathrm{tnew}}\left(\sum _{v=1}^V{O}^{\left(v\right)t}{Q}_n^{\left(v\right)}\left(O^{\left(v\right)}{\varphi }^{\mathrm{new}}\left(x_n\right)+b^{\left(v\right)}\right)\right);& \left(\mathrm{Equation}26\right)\\ {\varphi }^{\mathrm{new}}\left(x\right)={\int }_0^1{v}_t^{\mathrm{new}}\left({\varphi }_t^{\mathrm{new}}\left(x\right)\right)dx+x,\mathrm{where}{D}_{{\varphi }_t\left(x_n\right)}{\varphi }_{t,1}=\left(\frac{\partial {\varphi }_1\circ {\varphi }_t^{-1}}{\partial y_j}{❘}_{{\varphi }_t\left(x_n\right)}\right).& \left(\mathrm{Equation}27\right)\end{array}$

[0089] In a further embodiment, the deformation may be performed in real time for the case when the rigid motions (i.e., the rotation/translation) which bring the avatar into correspondence with the one or more source 2D projection are not known. A similar approach to the one above is used, with the addition of an estimation of the rigid motions using the techniques described herein. The initialization u^new=0 is used. Rigid motions are calculated using the rotation/translation techniques above to register the CAD model xx+u^new(x) to each photograph, generating rigid motions O^(v)new,b^(v)new, v=1, 2, . . . O^(v)new,b^(v)new are fixed from the previous step, and the deformation of the CAD model xx+u^new(x) or large deformation x→φ(x) are computed using the above techniques to solve the real-time small deformation or large deformation problem:

[0090] $\begin{array}{cc}\left(\mathrm{small}\right)\underset{u}{\min }{\mathrm{Lu}}^2+\sum _{n=1}^N\sum _{v=1}^V{\left(O^{\left(v\right)}\left(x_n+u\left(x_n\right)\right)+b^{\left(v\right)}\right)}^t{Q}_n^{\left(v\right)}\left(O^{\left(v\right)}\left(x_n+u\left(x_n\right)\right)+b^{\left(v\right)}\right).& \left(\mathrm{Equation}28\right)\\ \left(\mathrm{large}\right)\underset{v_t,t\in \left[0,1\right]}{\min }{\int }_0^1{{\mathrm{Lv}}_t}^{2}dt+\sum _{n=1}^N\sum _{v=1}^V({O^{\left(v\right)}\left(\varphi \left(x_n\right)+b^{\left(v\right)}\right)}^t{Q}_n^{\left(v\right)}\left(O^{\left(v\right)}\left(\varphi \left(x_n\right)\right)+b^{\left(v\right)}\right).& \left(\mathrm{Equation}29\right)\end{array}$

[0091] In another embodiment, the avatar is deformed in real-time using diffeomorphic deformations. The solution to the real-time deformation algorithm generates a deformation which may be used as an initial condition for the solution of the diffeomorphic deformation. Real-time diffeomorphic deformation is accomplished by incorporating the real-time deformations solution as an initial condition and then performing a small number (in the region of 1 to 10) iterations of the diffeomorphic deformation calculation.

[0092] The deformation may include affine motions. For the affine motion A:x→Ax where A is the 3×3 generalized linear matrix so that

[0093] $\begin{array}{cc}\underset{A,z_n}{\min }\sum _{n=1}^N\sum _{v=1}^V{O^{\left(v\right)}{\mathrm{Ax}}_n+b^{\left(v\right)}-z_n^{\left(v\right)}{P}_n^{\left(v\right)}}^2=\underset{A}{\min }\sum _{n=1}^N\sum _{v=1}^V{\left(O^{\left(v\right)}{\mathrm{Ax}}_n+b^{\left(v\right)}\right)}^t{Q}_n^{\left(v\right)}\left(O^{\left(v\right)}{\mathrm{Ax}}_n+b^{\left(v\right)}\right),& \left(\mathrm{Equation}30\right)\end{array}$
the least-squares estimator Â:x→Âx is computed:

[0094] $\begin{array}{cc}\hat{A}=-{\left(\sum _{n=1}^N{X}_n^t{\stackrel{\_}{Q}}_n{X}_n\right)}^{-1}\left(\sum _{n=1}^N{X}_n^t\sum _{v=1}^V{O}^{\left(v\right)t}{Q}_n^{\left(v\right)}{b}^{\left(v\right)}\right).& \left(\mathrm{Equation}31\right)\end{array}$

[0095] In many cases, both feature items in the projective imagery as well as the imagery itself can be used to drive the deformation of the avatar. Augmentation of source data to incorporate source imagery may improve the quality of the fit between the deformed avatar and the target face. To implement this, one more term is added to the deformation techniques. Let I be the measured imagery, which in general includes multiple measured images I^(v), v=1, 2, . . . corresponding to an indexed sequence of pixels indexed by pε[0,1]², with the projection mapping points

[0096] $x=\left(x,y,z\right)\in {\mathrm{IR}}^3\mapsto p\left(x\right)=\left(p_1\left(x\right)=\frac{{\alpha }_{1}x}{z},p_2\left(x\right)=\frac{{\alpha }_{2}y}{z}\right).$
For the discrete setting of pixels in the source image plane with color (R,G,B) template, the observed projective Π(p) is an (R,G,B) vector and the projective matrix becomes

[0097] $P_x=\left(\begin{array}{ccc}\frac{{\alpha }_1}{z+n}& 0& 0\\ 0& \frac{{\alpha }_2}{z+n}& 0\end{array}\right),$
operating on points (x, y, z)εIR³ according to the projective matrix

[0098] $P_x:\left(x,y,z\right)\mapsto \left(p_1\left(x,y,z\right),p_2\left(x,y,z\right)\right)=\left(\begin{array}{ccc}\frac{{\alpha }_1}{z+n}& 0& 0\\ 0& \frac{{\alpha }_2}{z+n}& 0\end{array}\right)\left(\begin{array}{c}x\\ y\\ z\end{array}\right),$
the point x(p) being the revealed point which is not occluded (closest point to the projection on the ray) on the 3D CAD model which projects to the point p in the image plane. Next, the projected template matrices resulting from finite differences on the (R,G,B) components at the projective coordinate p of the template value are required. The norm is interpreted component wise:

[0099] $\begin{array}{cc}{\nabla }^t\prod \left(p\right)=\left(\begin{array}{cc}\frac{\partial }{\partial p_1}\prod {\left(p\right)}^r& \frac{\partial }{\partial p_2}\prod {\left(p\right)}^r\\ \frac{\partial }{\partial p_1}\prod {\left(p\right)}^g& \frac{\partial }{\partial p_2}\prod {\left(p\right)}^g\\ \frac{\partial }{\partial p_1}\prod {\left(p\right)}^b& \frac{\partial }{\partial p_2}\prod {\left(p\right)}^b\end{array}\right),\stackrel{.}{I}\left(p\right)\stackrel{.}{=}I\left(p\right)-\prod \left(p\right),& \left(\mathrm{Equation}32\right)\\ {\tilde{\nabla }}^t\prod \left(p\right)={\nabla }^t\prod \left(p\right)P_{x\left(p\right)}=\left(\begin{array}{ccc}\frac{\partial \prod {\left(p\right)}^r}{\partial p_1}\frac{{\alpha }_1}{z\left(p\right)+n}& \frac{\partial \prod {\left(p\right)}^r}{\partial p_2}\frac{{\alpha }_2}{z\left(p\right)+n}& 0\\ \frac{\partial \prod {\left(p\right)}^g}{\partial p_1}\frac{{\alpha }_1}{z\left(p\right)+n}& \frac{\partial \prod {\left(p\right)}^g}{\partial p_2}\frac{{\alpha }_2}{z\left(p\right)+n}& 0\\ \frac{\partial \prod {\left(p\right)}^b}{\partial p_1}\frac{{\alpha }_1}{z\left(p\right)+n}& \frac{\partial \prod {\left(p\right)}^b}{\partial p_2}\frac{{\alpha }_2}{z\left(p\right)+n}& 0\end{array}\right),& \left(\mathrm{Equation}33\right)\end{array}$
with matrix norm ∥A−B∥²=|A^r−B^r|²+|A^g−B^g|²+|A^b−B^b|².

[0100] Associated with each image is a translation/rotation assumed already known from the previous rigid motion calculation techniques. The following assumes there is one 2D image, with O,b identity, and let any of the movements be represented as x→x+u(x). Then u(x)=Ox−x is a rotational motion, u(x)=b is a constant velocity, u(x)=Σ_ie_iE_i(x) is a constrained motion to a basis function such as “chin rotation,” “eyebrow lift,” etc., and the general motion u is given by:

[0101] $\begin{array}{cc}\underset{u}{\min }\sum _{p\in {\left[0,1\right]}^2}{I\left(p\right)-\prod \left(p\right)-{\nabla }^t\prod \left(p\right)P_{x\left(p\right)}u\left(x\left(p\right)\right)}^2=\underset{u}{\min }\sum _{p\in {\left[0,1\right]}^2}\stackrel{.}{I}\left(p\right)-{\tilde{\nabla }}^t\prod \left(p\right){u(x\left(p\right)}^2& \left(\mathrm{Equation}34\right)\\ =\underset{u}{\min }-2\sum _{p\in {\left[0,1\right]}^2}{u\left(x\left(p\right)\right)}^t\stackrel{.}{I}\left(p\right)\tilde{\nabla }\prod \left(p\right)\sum _{p\in {\left[0,1\right]}^2}{u\left(x\left(p\right)\right)}^t\tilde{\nabla }\Pi \left(p\right){\tilde{\nabla }}^t\prod \left(p\right)u\left(x\left(p\right)\right).& \left(\mathrm{Equation}35\right)\end{array}$
This is linear in u so closed-form expressions exist for each of the forms of u, for example, for the unconstrained general spline motion,
û(x(p))=({tilde over (∇)}Π(p) ∇^tΠ(p))⁻¹İ(p){tilde over (∇)}Π(p).  (Equation 36)
This approach can be incorporated into the other embodiments of the present invention for the various possible deformations described herein.

[0102] In the situation where large numbers of feature points, curves, and subvolumes are to be automatically generated from the source projection(s) and 3D data (if any), image matching is performed directly on the source imagery or on the fundamental 3D volumes into which the source object can be divided. For the case where the avatar is generated from 2D projective photographic images, the measured target projective image has labeled points, curves, and subregions generated by diffeomorphic image matching. Defining a template projective exemplar face with all of the labeled submanifolds from the avatar, the projective exemplar can be transformed bijectively via the diffeomorphisms onto the target candidate photograph, thereby automatically labeling the target photographs into its constituent submanifolds. Given these submanifolds, the avatars can then be matched or transformed into the labeled photographs. Accordingly, in the image plane a deformation φ: xφ(x) satisfying φ=∫₀¹v_t(φ_t(x))dt+x, xεR² is generated. The template and target images I₀, I₁ are transformed to satisfy

[0103] $\begin{array}{cc}\underset{v_t,t\in \left[0,1\right],z_n}{\min }{\int }_0^1{{\mathrm{Lv}}_t}^{2}dt+\sum _{n=1}^N{I_0\circ {\varphi }^{-1}-I_1}^2.& \left(\mathrm{Equation}37\right)\end{array}$
The given diffeomorphism is applied to labeled points, curves, and areas in the template I₀ thereby labeling those points, curves, and areas in the target photograph.

[0104] When the target source data are in 3D, the diffeomorphisms are used to define bijective correspondence in the 3D background space, and the matching is performed in the volume rather than in the image plane.

[0105] The following techniques may be used to select avatar models automatically. Given a collection of avatar models {CAD^a, a=1, 2, . . . }, and a set of measured photographs of the face of one individual, the task is to select the avatar model which is most representative of the individual face being analyzed. Let the avatar models be specified via points, curves, surfaces, or subvolumes. Assume for example N feature initial and target points, x_n^a, x_nεIR^d, n=1, . . . , N, with x_n=x_n^a+u(x_n^a), one for each avatar a=1, 2, . . . .

[0106] In one embodiment, the avatar is deformed with small deformations only and no rigid motions. For this embodiment, it is assumed that the measured feature items are all points from a single camera viewing which generated the projected source image in which the feature points were measured. The matching x_n^a x_n^a+u(x_n^a), n=1, . . . , N is constructed and the CAD^a model which is of smallest metric distance is selected. The optimum selected avatar model is the one closest to the candidate in the metric. Any of a variety of distance functions may be used in the selection of the avatar, including the large deformation metric from the diffeomorphism mapping technique described above, the real-time metric described above, the Euclidean metric, and the similitude metric of Kendall. The technique is described herein using the real-time metric. When there is no rigid motion, then the CAD model is selected to minimize the metric based on one or several sets of features from photographs, here described for one photo:

[0107] $\begin{array}{cc}\begin{array}{c}{\mathrm{CAD}}^{\stackrel{\_}{a}}=\arg \underset{{\mathrm{CAD}}^a,a=1,2,\dots }{\min }\underset{u,z_n}{\min }{\mathrm{Lu}}^2+\\ \sum _{n=1}^N{\left(x_n^a+u\left(x_n^a\right)\right)-z_n{P}_n}^2\\ =\arg \underset{{\mathrm{CAD}}^a,a=1,2,\dots }{\min }\underset{u}{\min }{\mathrm{Lu}}^2+\\ \sum _{n=1}^N{\left(x_n^a+u\left(x_n^a\right)\right)}^t{Q}_n\left(x_n^a+u\left(x_n^a\right)\right).\end{array}& \left(\mathrm{Equation}38\right)\end{array}$
In other embodiments both metrics including unknown or known rigid motions, large deformation metrics, or affine motions can be used such as described in equations 28 and 29, respectively.

[0108] For selecting avatars given 3D information such as features of points, curves, surfaces and or subvolumes in the 3D volume, then the metric is selected which minimizes the distance between the measurements and the family of candidate CAD models. First, the K matrix defining the quadratic form metric measuring the distance is computed:

[0109] $\begin{array}{cc}K=\left(\begin{array}{cccc}K\left(x_1^a,x_1^a\right)& K\left(x_1^a,x_2^a\right)& & K\left(x_1^a,x_N^a\right))\\ ⋮& ⋮& ⋮& ⋮\\ K\left(x_N^a,x_1^a\right)& K\left(x_N^a,x_2^a\right)& & K\left(x_N^a,x_N^a\right)\end{array}\right);& \left(\mathrm{Equation}39\right)\end{array}$
K may be, for example, K(x_i,x_j)=diag(e^{−a∥xi−xj∥}). Next, the metric between the CAD models and the candidate photographic feature points is computed according to x_n^ax_n=Ax_n^a+b+u(x_n^a), n=1, . . . , N and CAD^a of small distance is selected. Exact matching or inexact matching (σ=0 or inexact σ≠0) can be used:

[0110] $\begin{array}{cc}{\mathrm{exactCAD}}^{\hat{a}}=\arg \underset{{\mathrm{CAD}}^a,a=1,2,\dots }{\min }\underset{A,b}{\min }\sum _{\mathrm{ij}=1}^N{\left({\mathrm{Ax}}_i+b-x_i^{\prime }\right)}^t{\left(K^{-1}\right)}_{\mathrm{ij}}\left({\mathrm{Ax}}_j+b-x_j^{\prime }\right);& \left(\mathrm{Equation}40\right)\\ {\mathrm{inexact}}^{\hat{a}}=\arg \underset{{\mathrm{CAD}}^a,a=1,2,\dots }{\min }\underset{A,b}{\min }\sum _{\mathrm{ij}=1}^N{\left({\mathrm{Ax}}_i+b-x_i^{\prime }\right)}^t{\left({\left(K+{\sigma }^{2}I\right)}^{-1}\right)}_{\mathrm{ij}}\left({\mathrm{Ax}}_j+b-x_j^{\prime }\right).& \left(\mathrm{Equation}41\right)\end{array}$
The minimum norm is determined by the error between the CAD model feature points and the photographic feature points.

[0111] The present invention may also be used to match to articulated target objects. The diffeomorphism and real-time mapping techniques carry the template 3D representations bijectively onto the target models, carrying all of the information in the template. The template models are labeled with different regions corresponding to different components of the articulated object. For example, in the case of a face, the articulated regions may include teeth, eyebrows, eyelids, eyes, and jaw. Each of these subcomponents can be articulated during motion of the model according to an articulation model specifying allowed modes of articulation. The mapping techniques carry these triangulated subcomponents onto the targets, thereby labeling them with their subcomponents automatically. The resulting selected CAD model therefore has its constituent parts automatically labeled, thereby allowing each avatar to be articulated during motion sequences.

[0112] In the case when direct 3D measurement of the source object is

[0113] available x_n, y_nεR³, n=1, . . . , N from points, curves, surfaces, or subvolumes, the techniques for determining the rotation/translation correspondences are unchanged. However, since the matching terms involve direct measurements in the volume, there is no need for the intermediate step to determine the dependence on the unknown z-depth via the MMSE technique.

[0114] Accordingly, the best matching rigid motion corresponds to:

[0115] $\begin{array}{cc}\underset{O,b}{\min }\sum _{n=1}^N{{\mathrm{Ox}}_n+b-y_n}^2.& \left(\mathrm{Equation}42\right)\end{array}$

[0116] The real-time deformation corresponds to:

[0117] $\begin{array}{cc}\underset{u,z_n}{\min }{\mathrm{Lu}}^2+\sum _{n=1}^N{\left(x_n+u\left(x_n\right)\right)-y_n}^2.& \left(\mathrm{Equation}43\right)\end{array}$

[0118] The diffeomorphism deformation corresponds to with φ=∫₀¹v_t(φ_t(x))dt+x, xε³:

[0119] $\begin{array}{cc}\underset{v_t,t\in \left[0,1\right],z_n}{\min }{\int }_0^1{{\mathrm{Lv}}_t}^{2}dt+\sum _{n=1}^N{\varphi \left(x_n\right)-y_n}^2.& \left(\mathrm{Equation}44\right)\end{array}$

[0120] The techniques described herein also allow for the automated calibration camera parameters, such as the aspect ratio and field of view. The set of x_j=(x_j,y_j,z_j), j=1, . . . , N features are defined on the CAD model. The positive depth projective geometry mapping with

[0121] ${\alpha }_1=\frac{1}{{\gamma }_1},{\alpha }_2=\frac{1}{{\gamma }_2}$
is defined, according to

[0122] $\left(x,y,z\right)\mapsto p\left(x,y,z\right)=\left(\frac{x}{{\gamma }_{1}z},\frac{y}{{\gamma }_{2}z}\right),$
zε[0,∞), n>0. Given are observations of some features through the projective geometry

[0123] $p_j=\left(\frac{x_j}{{\gamma }_2{z}_j},\frac{y_j}{z_j}\right).$

[0124] The calibration of the camera is determined under the assumption that there is no transformation (affine or other) of the avatar. The z value is parameterized by incorporating the n frustrum distance so that all depth coordinates are the above coordinates plus frustrum depth. Videos can show different aspect ratios

[0125] $\mathrm{AR}=\frac{{\alpha }_1}{{\alpha }_2}$
and fields-of-view

[0126] $\mathrm{FOV}=2{\tan }^{-1}\frac{1}{{\alpha }_1}.$
The technique estimates the aspect ratios γ₁,γ₂ from measured points P_i=(γ₁p_i1,γ₂p_i2,1)^t, i=1, . . . , N:

[0127] $\begin{array}{cc}\underset{O,b{\gamma }_1,{\gamma }_2,z_i,i=1,\dots }{\min }\sum _{i=1}^N{{\mathrm{Ox}}_i+b-z_i{P}_i}^2=\underset{O,b{\gamma }_1,{\gamma }_2}{\min }\sum _{i=1}^N{\mathrm{Ox}}_i+b^2-\sum _{i=1}^N\frac{{{\left({\mathrm{Ox}}_i+b\right)}^t{P}_i}^2}{{P_i}^2}.& \left(\mathrm{Equation}45\right)\end{array}$
Using the initialization γ₁^new=γ₂^new=1, the calculation is run to convergence.

[0128] In the first step, the data terms

[0129] $P_i=\left(\begin{array}{c}{\gamma }_1^{\mathrm{old}}{p}_{\mathrm{i1}}\\ {\gamma }_2^{\mathrm{old}}{p}_{\mathrm{i2}}\\ 1\end{array}\right)$
are solved, and the following rotations/translations are computed:

[0130] $\begin{array}{cc}\begin{array}{c}{Q}_i=\left(I-\frac{{P_i\left(P_i\right)}^t}{{P_i}^2}\right),\stackrel{\_}{Q}\\ =\sum _{i=1}^N{Q}_i,X_Q=\sum _{i=1}^N{Q}_i{X}_i,M_i=X_i-{\stackrel{\_}{Q}}^{-1}{X}_Q;\end{array}& \left(\mathrm{Equation}46\right)\\ \begin{array}{c}\hat{O}=\arg \underset{O}{\min }{O}^t\left(\sum _{i=1}^N{M}_i^t{Q}_i{M}_i\right)O,\hat{b}\\ =-{\left(\sum _{i=1}^N{Q}_i\right)}^{-1}\sum _{i=1}^N{Q}_i\hat{O}{x}_i,{\psi }_i=\left(\begin{array}{c}{\left(\hat{O}{x}_i+\hat{b}\right)}_1{p}_{\mathrm{i1}}\\ {\left(\hat{O}{x}_i+\hat{b}\right)}_2{p}_{\mathrm{i2}}\\ {\left(\hat{O}{x}_i+\hat{b}\right)}_3\end{array}\right).\end{array}& \left(\mathrm{Equation}47\right)\end{array}$