Global constrained parameterization of triangulated surfaces

[0002] The present invention pertains to the field of three dimensional model generation, and more particularly to a method for reducing distortion of a feature when the feature is mapped onto a contoured surface.

[0003] Texture mapping is widely used to improve the visual richness of 3D surfaces in computer generated images. Each 3D surface is put in correspondence with a planar image through a function called a mapping. Such a mapping assigns a pair of coordinates (u,v) referring to a pixel of the planar image to each point of a surface. Thus, for instance, the latitude and longitude can define a trivial mapping of a sphere. This technique was described by E. Catmull in A Subdivision Algorithm for Computer Display of Curved Surfaces, PhD thesis, Dept. of Computer Sciences, University of Utah, December 1974, which is hereby incorporated by reference. This technique was first applied to bicubic patches using a recursive subdivision algorithm. Unfortunately, these methods often produce highly distorted textures in the resulting images, similar to FIG. 10 (which is not prior art, but is only referred to for illustrative purposes).

[0004] First attempts to minimize these distortions were described by R. Darwyn Peachey in Solid Texturing of Complex Surfaces, SIGGRAPH85 Conference Proceedings, volume 19, pages 287-296, ACM, July 1985, and E. Bier and K. Sloan in Two-Part Texture Mapping, IEEE Computer Graphics and Applications, pages 40-53, September 1986, both of which are incorporated herein by reference. These attempts involved separating the process into two steps. The texture pattern is first applied to a simple intermediate surface such as a box or a cylinder for which texture mapping is trivial. Then, this intermediate surface is projected on the target object. The choice of the intermediate surface and its orientation, together with the projection method selected, dramatically affect the results, and great deal of user interaction was therefore required.

[0005] Marcel Samek, C. Slean, and H. Weghorst in Texture Mapping and Distortions in Digital Graphics, The Visual Computer, 2(5):313-320, September 1986, which is hereby incorporated by reference, taught that assigning texture coordinates to any surface is equivalent to flattening it. That method consists of unfolding a polygonal surface from a user selected seed. A similar technique was dislosed by C. Bennis, J. M. V'ezien, and G.Igl'esias in Piecewise Surface Flattening for Non-Distorted Texture Mapping, SIGGRAPH91 Conference Proceedings, volume 25, pages 237-246, ACM, July 1991, which is hereby incorporated by reference. This latter method means a parametric surface may be unfolded by allowing cuts to appear on the mapped texture when the discrepancy of the geodesic curvature goes beyond a given threshold.

[0006] Minimizing the distortions induced by texture mapping can be also realized using optimization techniques. In the method described by S. D. Ma and H. Lin in Optimal Texture Mapping, EUROGRAPHICS'88, pages 421-428, September 1988, which is incorporated herein by reference, a mapping of any surface is constructed by starting from a grid of points sampled on the surface. The grid is then iteratively optimized by minimizing a global distortion criterion. V. Krishnamurthy and M. Levoy in Fitting Smooth Surfaces to Dense Polygon Meshes, SIGGRAPH96 Conference Proceedings, pages 313-324, ACM, August 1996, and incorporated herein by reference, disclosed a similar approach for converting a triangulated mesh into a set of B-Spline surfaces. It is also possible to construct a mapping by assigning (u,v) coordinates to the vertices of the mesh. This naturally leads to the use of harmonic maps, as described by M. Eck, T. DeRose, T. Duchamp, H. Hoppes, M. Lounsbery and W. Stuetzle in Multiresolution Analysis of Arbitrary Meshes, SIGGRAPH95 Conference Proceedings, pages 173-182, ACM, August 1995, which is hereby incorporated herein by reference. This method consists of minimizing a metric dispersion criterion. Unfortunately, this does not always preserve angles accurately. Another approach disclosed by M. S. Floater in Parametrization and Smooth Approximation of Surface Triangulations, Computer Aided Geometric Design, 14(3):231-250, April 1997, and incorporated herein by reference, generalizes the barycentric mapping method described by W. T. Tutte in Convex Representation of Graphs, Proc. London Math. Soc., volume 10, 1960, also incorporated herein by reference. The (u,v) texture coordinates are found to be the solution of a linear system, where each (u,v) point is a convex combination of its neighbors. M. S. Floater (referenced above) taught a method for choosing the coefficients of these convex combinations to mimic the chord length parametrization for curves. These global methods give good results for most simple surfaces, but suffer from several limitations when applied to more complex surfaces. Further, since the criterion to be minimized is integral to the optimization algorithm, it is often difficult to take into account user defined information. For instance, as most surfaces are not developable, distortions will still remain, and the user may want to specify the distribution of these distortions.

[0007] Additionally, it is desirable to have a mapping function that handles “cuts” in the surface, such as for example are found when mapping a geological surface, which may represent a faulted stratigraphic grid. One technique for handling this is described by J. Bloomenthal in Modeling the Mighty Maple, SIGGRAPH85 Conference Proceedings, volume 19, pages 305-311, ACM, July 1985. Bloomenthal uses several distinct “patches” on the surface and makes the edges of the patches match to provide continuity of a texture over the discontinuity (“cut”). This is a computationally intense method, and it is desirable to provide a simpler method for dealing with discontinuities along a line in the surface.

[0008] Prior interactive mapping methods, where the parametrization is partially or completely defined by the user, are described by H. K. Pedersen in Decorating Implicit Surfaces, SIGGRAPH95 Conference Proceedings, pages 291-300, ACM, 1995, P. Litwinowicz and G. Miller in Efficient Techniques for Interactive Texture Placement, SIGGRAPH94 Conference Proceedings, pages 119-122, ACM, July 1994, and J. Maillot, H. Yahia, and A. Verroust in Interactive Texture Mapping, SIGGRAPH93 Conference Proceedings, volume 27, ACM, 1993, all of which are hereby incorporated herein by reference.

[0009] Therefore, what is needed in the art is a method for texture mapping onto a surface which reduces distortion of the texture, and which preferably allows for user input to distribute any remaining distortion to preferred regions. Additionally, it is desirable to have a texture mapping method which is capable of handling discontinuities in the surface.

[0010] The invention discloses a novel global optimization method, and more particularly, a texture mapping method which reduces distortion encountered by prior texture mapping methods. Contrary to the prior art, the method disclosed herein is based on a modular approach enabling the way the mapping is constructed to be customized. For instance, it is possible to tune the perpendicularity and homogeneous spacing of isoparametric curves over the entire surface, thus specifying the surface zones where distortions should be minimized in order of preference. The disclosed method also allows the user to make the mapping respect a set of user specified isoparametric curves. Moreover, the mapping can be made continuous through cuts, hence allowing the mapping of a texture on a complex cut surface in a single step. The method can be extended relatively easily to honor other kinds of linear constraints. All these constraints allow the method to take into account user specified information while being much more automatic than prior interactive mapping methods, where the parametrization is partially or completely defined by the user.

[0011] The disclosed method can easily take into account various additional information. For example, it is possible with this method to specify the zones where distortions should be minimized in order of preference, to make a set of isoparametric passes through user specified curves, and to “sew” “cuts” of a surface in texture space. Moreover, the method can be extended by defining new constraints, once these constraints are expressed as linear (or linearizable) relations.

[0012] The method can be easily implemented, since it only requires an efficient representation of triangulated meshes, which is provided by most C.A.D. packages. Thus, the algorithm has been integrated as a basic algorithm into a widely used C.A.D. software dedicated to geology, and several methods other than these associated with texture mapping have been developed based on this algorithm, such as unfolding geological layers and performing computations in texture space.

[0013] Given a surface represented by a set of adjacent triangles, i.e., a “triangulated surface”, one embodiment of the invention comprises a method for building a global parameterization (u(x), v(x)) honoring constraints in a least square sense. By way of example, the method can be used to honor linear constraints by specifying that the isoparametric curves are regularly spaced on the surface. Alternately, the constraint can comprise specifying that the isoparametric curves C(u₀) are orthogonal to the isoparametric curves C(v₀). Another constraint that can be used in the method comprises specifying that an isoparametric curves C(u₀) or C(v₀) are coincident with a given curve which can be drawn on the surface. Yet another constraint that can be used in the method, which is particulalry useful in a “cut” surface (i.e., one comprising discontinuities, particularly curvilinear discontinuities), comprises making the parameterization continuous through the “cut”.

[0014] The invention can be described as a method for mapping a feature onto a simulated surface having contoured portions, the surface being definable by a plurality of generally adjacent triangulations, the surface being further definable by a plurality of intersecting isoparametric curves linking identified values of the feature to identified points on the surface, the isoparametric curves characterized in that when the surface is flat, the isoparametric curves intersect one another in a perpendicular manner, the method comprising maintaining the perpendicularity of the intersections of at least some of the isoparametric curves on the contoured portions of the surface. This embodiment of the invention may be described as maintaining the orthogonality of the isoparametric curves.

[0015] The invention may be further described as a method for mapping a feature onto a simulated surface having contoured portions, the surface being definable by a plurality of generally adjacent triangulations, the surface being further definable by a plurality of intersecting isoparametric curves linking identified values of the feature to identified points on the surface, the isoparametric curves characterized in that when the surface is flat, the isoparametric curves are spaced a predetermined distance apart, the method comprising maintaining the predetermined spaced distance between at least some of the isoparametric curves on the contoured portions of the surface. This embodiment of the invention may be described as maintaining the homogeneous spacing of the isoparametric curves.

[0016] The invention may yet further be described as a method for mapping a feature onto a simulated surface, the surface being definable by a plurality of generally adjacent triangulations, the surface being further definable by a plurality of intersecting isoparametric curves linking identified values of the feature to identified points on the surface, the surface characterized by a defined secondary isoparametric curve projected onto the surface, the method comprising aligning at least one of the plurality of isoparametric curves to be coincident with the defined secondary isoparametric curve. This embodiment of the invention may be described as aligning details of the texture to details of the surface by specifying an isoparametric curve, and may be used to locally constrain a mapping (versus a global constraint).

[0017] The invention may yet further be described as a method for mapping a feature onto a simulated surface having contoured portions, the surface being definable by a plurality of generally adjacent triangulations, each triangulation be defined by vertices terminating at points on the surface, the surface being further definable by a plurality of intersecting isoparametric curves linking identified values of the feature to identified points on the surface, the surface being characterized by a linear discontinuity such that a point on the surface along the discontinuity may be identified with two vertices of non-adjacent triangulations, the method comprising maintaining the identified value of the feature associated with a point on the surface along the discontinuity at the two both vertices. This embodiment of the invention may be described as mapping two borders of a cut on the surface to the same isoparametric curve, or “sewing” the cut in (u,v) domain space.

[0018] Preferred embodiments of the invention are described below with reference to the following accompanying drawings.

[0019] FIG. 1 is a diagram showing mapping Φ from a surface S of R³to D⊂R².

[0020] FIG. 2 is a diagram showing mapping the function Φ interpolated over a triangle.

[0021] FIG. 3A is a representation of a three dimensional, contoured surface showing isoparametric curves obtained without non-distorted constraints.

[0022] FIG. 3B is a representation of a three dimensional, contoured surface showing isoparametric curves obtained with non-distorted constraints.

[0023] FIG. 4 is a diagram showing the local orthonormal basis (X, Y) of a triangle T=(α₀, α_l, α₂).

[0024] FIG. 5 is a diagram showing a constant gradient across the common edge of two triangles T and _{{haeck over (T)}}.

[0025] FIG. 6 is a diagram showing the aligning details of the texture to details of the three dimensional model by specifying an isoparametric curve.

[0026] FIG. 7 is a diagram showing extrapolating a mapping onto a three dimensional, contoured surface from four user specified isoparametric curves.

[0027] FIG. 8 is a diagram showing how two borders of a three dimensional, contoured cut surface can be connected in texture-space.

[0028] FIG. 9 is a diagram showing the representation of a three dimensional surface (a face) by a triangulated mesh.

[0029] FIG. 10 is a diagram showing a texture mapped over the triangulated mesh surface of FIG. 9 wherein the isoparametric contour were plotted without the use of non-distortion constraints.

[0030] FIG. 11 is a diagram showing the representation of a three dimensional surface (a face) wherein orthogonality and homogeneity non-distortion constraints were used on the isoparametric contour lines.

[0031] FIG. 12 is a diagram showing a texture mapped over the triangulated mesh surface of FIG. 11.

[0032] FIG. 13 is a diagram showing a texture mapped over a three dimensional, cut contoured surface such as might be used to represent a faulted stratigraphic geological layer.

[0033] FIG. 14 is a detail of the upper left corner of the three dimensional surface of FIG. 13, shown as represented by a simplified mesh.

[0034] FIG. 15 is a diagram of a three dimensional model of a salt dome, showing a texture map where the orthogonality of the isoparametric curves is respected, but the sizes of the squares differ a great deal over the contoured surface of the salt dome.

[0035] FIG. 16 is variation of the diagram of a three dimensional model of a salt dome shown in FIGS. 15 and 17, with an average solution between that shown for those two figures, and where the same weighting for the non-distortion criteria has been used for the two constraints (orthogonality and homogeneity) but with finer segments mapped over the contoured surface than in FIG. 15.

[0036] FIG. 17 is variation of the diagram of a three dimensional model of a salt dome shown in FIG. 15 and 16, but the squares have approximately the same size while the isoparametric curves are far from orthogonal as compared with FIG. 15.

[0037] FIG. 18 is a diagram of a texture mapped three dimensional, contoured surface with an isoparametric curve and non-aligned texture.

[0038] FIG. 19 is a variation of the texture mapped three dimensional surface of FIG. 18 with the texture matched to the isoparametric curve.

[0039] By way of overview, we will first describe the notion of mapping a function defined on a triangulated mesh, as well as a new method for constructing such mapping functions based on an iterative optimization algorithm. We will then describe how these mappings can be optimized in order to minimize the distortions.

[0040] In the following description, we wish to be clear about terms. Therefore, we will now provide the following definitions:

[0041] By “mapping” we mean any operation which continuously associates a point (u,v) of a parametric 2-dimensional space to a point x(u,v) located on a 3-dimensional surface.

[0042] By “global parameterization” of a surface embedded in a 3 dimensional space, we mean any operation which continuously associates any point x on the surface to a point (u(x),v(x)) of a parametric space. Note that this is the inverse of “mapping”.

[0043] By saying that a given surface is said “parameterized” we mean that there is a given global parameterization associating any point x of this surface to a pair of parameters (u(x),v(x)) of a parametric domain. It should be noticed that several parameterizations can be associated to a given surface.

[0044] Given a parameterized surface and a value u₀in the parametric domain, by “isoparametric curve C(u₀)”, we mean the set of points located on the surface and generated by x(u₀,v) when v varies in the parametric domain.

[0045] Given a parameterized surface and a value v₀in the parametric domain, by “isoparametric curve C(v₀)”, we mean the set of points located on the surface and generated by x(u,v₀) when u varies in the parametric domain.

[0046] The Mapping Function Φ and Discrete Mapping φ

[0047] Turning to FIG. 1, the concept of a mapping function, or “mapping”, Φ, from a surface S of in three a dimensional region R³to a domain on a two dimensional region, D⊂R², is shown. As shown, given an open surface S of R³, a mapping Φ is a one-to-one transform that maps the surface S to a subset D of R², and is given by the following formula (Formula 1): $\left(x,y,z\right)\in S\to \Phi \left(x,y,z\right)=\left[\begin{array}{c}{\Phi }^u\left(x,y,z\right)\\ {\Phi }^v\left(x,y,z\right)\end{array}\right]$

[0048] wherein the following definitions are given:

[0049] D is called the parametric (u,v) domain.

[0050] As Φ is, by definition, a one-to-one function, it has an inverse function x=Φ⁻¹, called a parametrization of the surface, such that the following formula (Formula 2) applies: $\left(u,v\right)\in 𝒟\to x\left(u,v\right)={\Phi }^{-1}\left(u,v\right)=\left[\begin{array}{c}x\left(u,v\right)\\ y\left(u,v\right)\\ z\left(u,v\right)\end{array}\right]$

[0051] If a surface has a parametrization x defined, the inverse Φ=x⁻¹ of this parametrization provides a mapping function. A method of applying this technique to cubic splines is disclosed in A Subdivision Algorithm for Computer Display of Curved Surfaces (E. Catmull), incorporated above, in which a recursive subdivision scheme is described, making it possible to avoid inverting of the parametrization directly. In the following description, the surface S is provided with a triangulation G={Ω, T}, where Ω is the set of the vertices of the triangulation, and T the set of the triangles of G, defined as vertex triplets. For the sake of simplicity, Ω will be identified with the interval [1 . . . M] of integers, where M=|106 | denotes the number of vertices of the triangulation. The geometric location at a vertex αεΩ is denoted as p(α) in the following description.

[0052] For surfaces relevant to the invention, we define the value of Φ at the vertices Ω of the triangulation G. This information can be stored as a set of (u_i, v_i) values, where 1≦i≦M. The selection of the values for (u_i, v_i) is explained below. This defines a discrete function φ:Ω→R²such that ∀α_iεΩ, φ(α_i)={φ^u(α_i), φ^v(α_i)}=(u_i, v_i). As shown in FIG. 2, a mapping function Φ can be then defined as the linear interpolation of φ over each triangle T=(α_i,α_j,α_k, of T. For each point p in T, Φ is given the following formula (Formula 3): $\begin{array}{c}\{\begin{array}{cccccc}\Phi \left(p\right)& =& & \left(1-\lambda -\mu \right)& ·& \varphi \left({\alpha }_i\right)\\ & & +& \lambda & ·& \varphi \left({\alpha }_j\right)\\ & & +& \mu & ·& \varphi \left({\alpha }_k\right)\end{array}\end{array}$

[0053] where:

[0054] λ and μ are the local barycentric coordinates at the point p in T

[0055] φ(α_i)=(u_i, v_i); φ(α_j)=(u_j, v_j); φ(α_k)=(u_k, v_k)

[0056] Discrete Smooth Interpolation

[0057] Given a triangulation G={Ω, T}, we preferably to assign (u, v) coordinates to each vertex αεΩ. It is known that the following two sufficient conditions construct a mapping:

[0058] 1. The image of the border of the surface through φ in the parametric (u, v) domain is a convex polygon; and

[0059] 2. Each internal node is a convex combination of its neighbors. It is also known that these two conditions are sufficient to define mappings, but are not necessary to define mappings. For example, as demonstrated below, the first condition (concerning the image of the border) can be replaced by a less restrictive condition.

[0060] More formally, the second condition (regarding the internal nodes) can be written in the form of the following equation (Equation 1): $\begin{array}{cc}\forall k\in \Omega ,\sum _{\alpha \in N\left(k\right)}{v}^{\alpha }\left(k\right)·\varphi \left(\alpha \right)=0& \left(1\right)\end{array}$

[0061] where:

[0062] N(k) denotes the set of nodes directly connected to k, including k, and

[0063] the v^α(k) are given coefficients such that the following equation set (Equation 2) applies: $\begin{array}{cc}\{\begin{array}{cc}{v}^{\alpha }\left(k\right)>0& \forall \alpha \in N\left(k\right)-\left\{k\right\}\\ v^k\left(k\right)=-\sum _{\underset{\alpha \ne k}{\alpha \in N\left(k\right)}}{v}^{\alpha }\left(k\right)\ne 0& \forall k\in \Omega \end{array}& \left(2\right)\end{array}$

[0064] Once boundary nodes have been mapped to a convex polygon in parametric domain space, (u, v) coordinates are assigned to the internal nodes of the triangulation. Instead of finding φ by directly solving Equation 1, the method described herein comprises minimizing a global criterion in a least square sense, honoring at the same time a set of linear constraints, as will be shown below. The algorithm is preferably based on the Discrete Smooth Interpolation algorithm (D.S.I.), described in Discrete Smooth Interpolation in Geometric Modeling, ACM-Transactions on Graphics, 8(2):121-144, 1989, and Discrete Smooth Interpolation, Computer Aided Design Journal, 24(4):263-270, 1992, both by J. L. Mallet, and both incorporated herein by reference. Those publications also describe the notions of generalized roughness, linear constraints, and the iterative D.S.I. algorithm, discussed below. The criterion minimized by the D.S.I. method is called the roughness R (“roughness”), and is defined in Equation 3 below: $\begin{array}{cc}R\left(\varphi \right)=\sum _{k\in \Omega }\sum _{v\in \left\{u,v\right\}}{\left\{\sum _{\alpha \in N\left(k\right)}{v}^{\alpha }\left(k\right)·{\varphi }^{\nu }\left(k\right)\right\}}^2& \left(3\right)\end{array}$

[0065] The minimum of this functional R(φ) is reached if ∂R(φ)/∂φ^υ(α)=0 for each αεΩ and for each υε{u, v}, where υ denotes one of the two components of φ. This yields the following equation (Equation 4): $\begin{array}{cc}{\varphi }^{\nu }\left(\alpha \right)=-\frac{G^{\nu }\left(\alpha |\varphi \right)}{g^{\nu }\left(\alpha \right)}\text{where:}|\begin{array}{ccc}{G}^{\nu }\left(\alpha |\varphi \right)& =& \sum _{k\in N\left(\alpha \right)}\left\{v^{\alpha }\left(k\right)·\sum _{\underset{\beta \ne \alpha }{\beta \in N\left(k\right)}}{v}^{\beta }\left(k\right)·{\varphi }^{\nu }\left(\beta \right)\right\}\\ g^{\nu }\left(\alpha \right)& =& \sum _{\alpha \in N\left(\alpha \right)}{\left\{v^{\alpha }\left(k\right)\right\}}^2\end{array}& \left(4\right)\end{array}$

[0066] The following algorithm (Algorithm 1) computes iteratively the assignments of (u, v) coefficients minimizing the roughness given in Equation 3. It is known that Algorithm 1 converges to a unique solution as soon as at least one node α has its value φ(α) fixed, and provided that the chosen v^α(k) coefficients honor Equation 2. As described further below, this method can be enhanced using D.S.I. constraints.

[0067] Algorithm 1: $\mathrm{let}I\mathrm{be}\mathrm{the}\mathrm{set}\mathrm{of}\mathrm{nodes}\mathrm{where}\varphi \mathrm{is}\mathrm{unknown}$ $\mathrm{let}{\varphi }_{\left[0\right]}\mathrm{be}a\mathrm{given}\mathrm{initial}\mathrm{approximated}\mathrm{solution}$ $\mathrm{while}\left(\mathrm{more}\mathrm{iterations}\mathrm{are}\mathrm{needed}\right)\{\mathrm{for\_all}\left(\alpha \in I\right)\{\mathrm{for\_all}\left(\nu \in \left\{u,v\right\}\right)\{{\varphi }^{\nu }\left(\alpha \right):=-\frac{G^{\nu }\left(\alpha \right)}{g^{\nu }\left(\alpha \right)}\}\}\}$

[0068] Where the v^α(k) coefficients are concerned, several choices are available. One possible choice is referred to as the shape preserving weighting, and is described in Parametrization and Smooth Approximation of Surface Triangulations (M. S. Floater, referenced and incorporated above). This type of weighting ensures that the location of a vertex in parametric space relative to its neighbors mimics the local geometry around the vertex being considered. By separating the criteria minimizing the distortions from those which ensure that a valid mapping is constructed, we can obtain a finer control on the way the surface is parametrized. For this reason, the simple harmonic weighting defined by the following formula (Formula 4) is used for the {v^α(k)}: $v^{\alpha }\left(k\right)=\{\begin{array}{cc}1& \mathrm{if}\alpha \in N\left(k\right)-\left\{k\right\}\\ -\mathrm{degree}\left(k\right)& \mathrm{if}\alpha =k\end{array}$

[0069] where degree(k) denotes the number of neighbors of k. Less preferably, one of the previously mentioned weightings such as the shape preserving or gaussian weightings can be used instead, since they both satisfy Equation 2, but as described below, by using linear constraints, the same effect can be obtained with higher flexibility.

[0070] NON-DISTORTED MAPPING

[0071] We will now describe the criterion to be minimized in order to construct a texture mapping having reduced distortion, which will be referred to further herein as a “non-distorted texture mapping”, although it will be appreciated that some distortion may still be present, although the distortion will generally be much reduced over mappings performed by prior methods. In an overview of the method, this criterion preserves the perpendicularity and constant spacing of the isoparametric curves traced on the surface. In other words, the gradients of u and v should be perpendicular one to another and constant all over the surface. This is clearly demonstrated in FIGS. 3A and 3B. FIG. 3A shows a surface 10 having isoparametric curves (shown by representative lines 11) obtained without non-distortion constraints, while FIG. 3B shows a surface 10 having isoparametric curves (shown by representative lines 12) obtained with non-distortion constraints. FIG. 3B clearly shows more even spacing of the grid lines and the cells, and the perpendicularity of the grid lines, while FIG. 3B shows irregular spacing and a high frequency of non-perpendicular grid lines. To accomplish such non-distortion, we define the gradient of a function interpolated over a triangulated mesh from the vertices of the triangulation, as described more fully below. Modification of Algorithm 1 to account for the non-distortion criterion will also be more fully described below.

[0072] Data in FIG. 3 shows the effect of the constraints described in this section as applied while parametrizing a triangulated mesh. The isoparametric curves obtained when applying Algorithm 1 are shown in FIG. 3A, whereas the constraints described further below give the result shown in FIG. 3B.

[0073] Gradient of a Discrete Function φ Interpolated over a Triangulation G

[0074] As shown in FIG. 4 (and reproduced in the text following this paragraph), each triangle T=(α₀,α₁,α₂) of T can be provided with a local orthonormal basis (p(α₀), X, Y). The function φ^ν_T(X,Y) denotes the linear interpolation of φ over the triangle T, where νε{u, v} represents one of the two components of φ and where (X,Y) are the local coordinates in the orthonormal basis (p(α₀), X, Y) of T. The following three formulaes, which are representative of vector quantities and correlate to FIG. 5, shall collectively be known as “Formula 5”. $\begin{array}{c}X=\frac{p\left({\alpha }_1\right)-p\left({\alpha }_0\right)}{\left|\right|p\left({\alpha }_1\right)-p\left({\alpha }_0\right)\left|\right|}\\ N=\frac{X\times \left\{p\left({\alpha }_2\right)-p\left({\alpha }_0\right)\right\}}{\left|\right|X\times \left\{p\left({\alpha }_2\right)-p\left({\alpha }_0\right)\right\}\left|\right|}\end{array}$

[0075] In this ortho normal basis, it is possible to confirm that the gradient of φ^ν_T(X, Y) is constant over T and is a linear combination of the values of (φ^ν_Tat the three vertices of the triangle T. The six coefficients D_x(α_j) and D_y(α_j) given in Equation 5 below are solely dependent on the geometry of the triangle T. $\begin{array}{cc}|\begin{array}{c}\frac{\partial {\varphi }_T^v}{\partial X}=\sum _{j=0}^2{D}_X\left({\alpha }_j\right)·{\varphi }^v\left({\alpha }_j\right)\\ \frac{\partial {\varphi }_T^v}{\partial Y}=\sum _{j=0}^2{D}_Y\left({\alpha }_j\right)·{\varphi }^v\left({\alpha }_j\right)\end{array}\mathrm{where}:\{\begin{array}{c}\begin{array}{c}{D}_X\left({\alpha }_0\right)=\left(y_1-y_2\right)/d\\ D_X\left({\alpha }_1\right)=\left(y_2-y_0\right)/d\\ D_X\left({\alpha }_2\right)=\left(y_0-y_1\right)/d\\ D_Y\left({\alpha }_0\right)=\left(x_2-x_1\right)/d\\ D_Y\left({\alpha }_1\right)=\left(x_0-x_2\right)/d\\ D_Y\left({\alpha }_2\right)=\left(x_1-x_0\right)/d\end{array}\\ d=\left(x_1-x_0\right)·\left(y_2-y_0\right)-\left(x_2-x_0\right)·\left(y_1-y_0\right)\\ \begin{array}{c}{y}_j=\left(p\left({\alpha }_j\right)-p\left({\alpha }_0\right)\right)·Y\\ x_j=\left(p\left({\alpha }_j\right)-p\left({\alpha }_0\right)\right)·X\end{array}|\forall \mathrm{j\epsilon }\left\{0,1,2\right\}\end{array}& \left(5\right)\end{array}$

[0076] Using the definition given in Equation 5 of the gradient of φ, we can write the equations corresponding to the orthogonality and homogeneous spacing of the isoparametric curves. The orthogonality of the iso-u and iso-v curves in a triangle T is given by the following Equation 6 as: $\begin{array}{cc}\left[\begin{array}{cc}\frac{\partial {\varphi }_T^u}{\partial X}& \frac{\partial {\varphi }_T^u}{\partial Y}\end{array}\right]·\left[\begin{array}{c}\frac{\partial {\varphi }_T^v}{\partial X}\\ \frac{\partial {\varphi }_T^v}{\partial Y}\end{array}\right]=0& \left(6\right)\end{array}$

[0077] If we consider that φ^uis fixed and that φ^vis to be determined, replacing in Equation 6 the gradient of φ^vwith its expression given in Equation 5 yields the following equation (Equation 7), which linearly combines the values of φ^vat the three vertices (α₀,α₁,α₂) of T. The equation to be used when φ^uis interpolated can be obtained by exchanging u and v in Equation 7. $\begin{array}{cc}\sum _{\mathrm{j\epsilon }\left\{0,1,2\right\}}\left\{{\varphi }^v\left({\alpha }_j\right)·\left(\frac{\partial {\varphi }_T^u}{\partial X}·D_X\left({\alpha }_j\right)+\frac{\partial {\varphi }_T^u}{\partial Y}·D_Y\left({\alpha }_j\right)\right)\right\}=0& \left(7\right)\end{array}$

[0078] The remaining condition on φ concerns the homogeneous spacing of the isoparametric curves. In other words, the gradient must not vary from one triangle to another. This requires that a common basis for two adjacent triangles T and {tilde over (T)} be defined, as shown in FIG. 5 wherein the gradient 15 is indicated as being the same for both T and {tilde over (T)}. The same expressions as introduced in FIG. 4 and Formula 5 apply for FIG. 5. The vector X is shared by the two bases, and {tilde over (Y)} is such that Y and {tilde over (Y)} would become colinear if the pair of triangles (T, {tilde over (T)}) were unfolded along their common edge [α₀,α₁].

[0079] The homogeneous spacing of the isoparametric curves is verified if, and only if, for each edge of T the equation below (Equation 8) is verified: $\begin{array}{cc}\{\begin{array}{c}\frac{\partial {\varphi }_T^u}{\partial X}=\frac{\partial {\varphi }_{\tilde{T}}^u}{\partial X};\frac{\partial {\varphi }_T^u}{\partial Y}=-\frac{\partial {\varphi }_{\tilde{T}}^u}{\partial \tilde{Y}}\\ \frac{\partial {\varphi }_T^v}{\partial X}=\frac{\partial {\varphi }_{\tilde{T}}^v}{\partial X};\frac{\partial {\varphi }_T^v}{\partial Y}=-\frac{\partial {\varphi }_{\tilde{T}}^v}{\partial \tilde{Y}}\end{array}& \left(8\right)\end{array}$

[0080] By replacing in Equation 8 the gradients of φ^uand φ^vby their expressions in T and {tilde over (T)}, the following four linear equations concerning the two components X and Y of the gradients of φ^uand φ^uare obtained (collectively, Equation 9). The term δ_Wtakes into account the fact that Y and {tilde over (Y)} point in an opposite direction. $\begin{array}{cc}\begin{array}{c}\forall v\epsilon \left\{u,v\right\},\\ \forall W\epsilon \left\{X,Y\right\},\end{array}|\begin{array}{c}\begin{array}{c}\begin{array}{c}{\varphi }^v\left({\alpha }_0\right)·\left\{D_W\left({\alpha }_0\right)+{\delta }_W·{\tilde{D}}_W\left({\alpha }_0\right)\right\}+\\ {\varphi }^v\left({\alpha }_1\right)·\left\{D_W\left({\alpha }_1\right)+{\delta }_W·{\tilde{D}}_W\left({\alpha }_1\right)\right\}+\end{array}\\ {\varphi }^v\left({\alpha }_2\right)·\left\{D_W\left({\alpha }_2\right)\right\}+\end{array}\\ {\varphi }^v\left({\tilde{\alpha }}_1\right)·\left\{{\delta }_W·{\tilde{D}}_W\left({\tilde{\alpha }}_2\right)\right\}=0\end{array}& \left(9\right)\\ \mathrm{where}{\delta }_W=\{\begin{array}{cc}-1& \mathrm{if}W=X\\ +1& \mathrm{if}W=Y\end{array}& \end{array}$

[0081] Honoring Linear Constraints

[0082] As described above, D.S.I. can be used to construct a mapping of a triangulated mesh. We next want to take into account the two principle criteria minimizing the distortions of the mapping, namely the perpendicularity and homogeneity criteria previously discussed. These two criteria can be written as a set of linear equations. As it is not possible to honor these constraints for a non-developable surface, they will be respected in a least square sense, thus minimizing the distortions. The general form of such a constraint is given in Equation 10 below: $\begin{array}{cc}\sum _{\alpha \epsilon \Omega }\left\{A_{c^v}\left(\alpha \right)·{\varphi }^v\left(\alpha \right)\right\}=b_{c^v}& \left(10\right)\end{array}$

[0083] where the values A_cν (α) and the scalar b_cν are constant given coefficients defining the constraint c.

[0084] Equation 7, corresponding to the perpendicularity of the isoparametric curves in the triangle T=(α₀,α₁,α₂), yields two constraints, c_T^uand c_T^v, to be honored when interpolating φ^uand φ^vrespectively. The expression of c_T^vis given below in Equation 11 following. The expression of the twin constraint c_T^ucan be obtained by permuting u and v in this Equation 11. $\begin{array}{cc}|\begin{array}{c}\begin{array}{c}\forall \mathrm{j\epsilon }\left\{0,1,2\right\},A_{c_T^v}\left({\alpha }_j\right)=\frac{\partial {\varphi }_T^u}{\partial X}·D_X\left({\alpha }_j\right)+\frac{\partial {\varphi }_T^u}{\partial Y}·D_Y\left({\alpha }_j\right)\\ \forall \mathrm{\alpha \epsilon }\left\{{\alpha }_0,{\alpha }_1,{\alpha }_2\right\},A_{c_T^v}\left(\alpha \right)=0\end{array}\\ b_{c_T^v}=0\end{array}& \left(11\right)\end{array}$

[0085] The homogeneity criterion specified by Equation 9 can be expressed by the following four constraints c^uX_E, c^uY_E, c^vX_Eand c^vY_Eyielded by Equation 12 below, to be taken into account at each edge E=(α₀,α₁) of the triangulation G. The vertices α₂and {tilde over (α)}₂denote the two remaining vertices of the two triangles T and {tilde over (T)} sharing the edge E. $\begin{array}{cc}|\begin{array}{c}{A}_{c_{{\rm E}}^{{\nu }_W}}\left({\alpha }_0\right)=\left\{D_W\left({\alpha }_0\right)+{\delta }_W·\widetilde{D_W}\left({\alpha }_0\right)\right\}\\ A_{c_{{\rm E}}^{{\nu }_W}}\left({\alpha }_1\right)=\left\{D_W\left({\alpha }_1\right)+{\delta }_W·\widetilde{D_W}\left({\alpha }_1\right)\right\}\\ A_{c_{{\rm E}}^{{\nu }_W}}\left({\alpha }_2\right)=D_W\left({\alpha }_2\right)\\ A_{c_{{\rm E}}^{{\nu }_W}}\left(\widetilde{{\alpha }_2}\right)={\delta }_W·\widetilde{D_W}\left({\alpha }_2\right)\\ A_{c_{{\rm E}}^{{\nu }_W}}\left(\alpha \right)=0\forall \alpha \notin \left\{{\alpha }_0,{\alpha }_1,{\alpha }_2,\widetilde{{\alpha }_2}\right\}\\ b_{c_{{\rm E}}^{{\nu }_W}}=0\end{array}\text{where:}\nu \in \left\{u,v\right\};W\in \left\{X,Y\right\};{\delta }_W=\{\begin{array}{c}-1\mathrm{if}W=X\\ +1\mathrm{if}W=Y\end{array}& \left(12\right)\end{array}$

[0086] The roughness criterion which D.S.I. minimizes can be generalized in order to honor a set C of linear constraints in a least square sense. Here, the set C of constraints is given by Equation 13 (following), where E denotes the set of the edges of the triangulation G(Ω,T). $\begin{array}{cc}C=\left(\bigcup _{{\rm T}\in \tau }\left\{c_{{\rm T}}^u,c_{{\rm T}}^v\right\}\right)\bigcup \left(\bigcup _{{\rm E}\in ϵ}\left\{c_{{\rm E}}^{u_X},c_{{\rm E}}^{u_Y},c_{{\rm E}}^{v_X},c_{{\rm E}}^{v_Y}\right\}\right)& \left(13\right)\end{array}$

[0087] The generalized roughness R*(φ), taking into account the degree of violation of the constraints C, is given by Equation 14 below. In addition to the equation of the roughness given in Equation 3, several terms correspond to the linear constraints, as described further below. $\begin{array}{cc}{R}^{*}\left(\varphi \right)=R\left(\varphi \right)+\phi ·\sum _{c\in C}{\varpi }_c·{\left\{\left(\sum _{\nu \in \left\{u,v\right\}}\sum _{\alpha \in \Omega }{A}_c^{\nu }\left(\alpha \right)·{\varphi }^{\nu }\left(\alpha \right)\right)-b_c\right\}}^2& \left(14\right)\end{array}$

[0088] In Equation 14, the term R(φ) is the roughness (see Equation 3), and the second term represents the degree of violation of the linear constraints. Each constraint c is weighted by a given _c>0 coefficient, allowing to tune the relative importance of the constraints. For example, it is possible to make the mapping respect the perpendicularity rather than the homogeneity. Moreover, since each triangle T and edge E has an individual constraint defined, as well as an individual associated _ccoefficient, it is possible to select the surface zones where the distortions are to be minimized in order of preference. The remaining coefficient φε[0, +∞] is a given parameter called the fitting factor and representing the importance of the constraints relative to the roughness.

[0089] The functional R*(φ) is a quadratic form, whose minimum is reached if ∂R*(φ)/∂φ^ν(α)=0 for each νε{u, v} and for each αεΩ. This yields the following equations (Equation 15-17), which solution minimizes R*(φ): $\begin{array}{cc}{\varphi }^{\nu }\left(\alpha \right)=-\frac{G^{\nu }(\alpha |\varphi )+\phi ·{\Gamma }^{\nu }(\alpha |\varphi )}{g^{\nu }\left(\alpha \right)+\phi ·{\gamma }^{\nu }\left(\alpha \right)}& \left(15\right)\\ |\begin{array}{c}{\Gamma }^{\nu }(\alpha |\varphi )=\sum _{c\in C}{\varpi }_c·{\Gamma }_c^{\nu }(\alpha |\varphi )\\ {\gamma }^{\nu }\left(\alpha \right)=\sum _{c\in C}{\varpi }_c·{\gamma }_c^{\nu }\left(\alpha \right)\end{array}\text{with:}& \left(16\right)\\ \{\begin{array}{c}{\gamma }_c^{\nu }\left(\alpha \right)=({A_c^{\nu }\left(\alpha \right)}^2\\ {\Gamma }_c^{\nu }(\alpha |\varphi )=A_c^{\nu }\left(\alpha \right)·\left\{\sum _{\beta \ne \alpha }{A}_c^{\nu }\left(\beta \right)·{\varphi }^{\nu }\left(\beta \right)-b_c\right\}\end{array}& \left(17\right)\end{array}$

[0090] The orthogonality constraint can be modified in the iterative D.S.I. algorithm. The two internal loops iterating on the components of φ and on the nodes of Ω respectively have been inverted. At each iteration, φ^uis interpolated while φ^vis considered to be constant, then the roles of φ^uand φ^vare permuted. The resulting algorithm given below (Algorithm 2) assigns (u, v) coordinates to the vertices of the triangulation while respecting the specified set of constraints. $\mathrm{let}I\mathrm{be}\mathrm{the}\mathrm{set}\mathrm{of}\mathrm{nodes}\mathrm{where}\varphi \mathrm{is}\mathrm{unknown}$ $\mathrm{let}{\varphi }_{\left[0\right]}\mathrm{be}a\mathrm{given}\mathrm{initial}\mathrm{approximated}\mathrm{solution}$ $\mathrm{while}\left(\mathrm{more}\mathrm{iterations}\mathrm{are}\mathrm{needed}\right)\{\mathrm{for\_all}\left(\nu \in \left\{u,v\right\}\right)\{\mathrm{for\_all}\left(\alpha \in I\right)\{{\varphi }^{\nu }\left(\alpha \right):=-\frac{G^{\nu }\left(\alpha \right)+{\Gamma }^{\bigvee }(\alpha |\phi )}{g^{\nu }\left(\alpha \right)+{\gamma }^{\nu }\left(\alpha \right)}\}\}\}$

[0091] The constraints defined and described so far herein provide a user of the method with a global control on the mapping function. Even if the orthogonality and perpendicularity constraints can be weighted locally to specify the zones where distortions are preferably to be minimized, this may be not sufficient for some applications, where a more precise set of local constraints is required. For example, it may be desirable to align some details of textures with details of models, which can be achieved by specifying isoparametric curves. Moreover, the model to be texture mapped can present cuts, and the user may want to define a single mapping function for a cut model instead of “sewing together” several patches. This desire can be fulfilled by making the mapping continuous through cuts, as will be described further below.

[0092] Specifying an Isoparametric Curve

[0093] As shown in FIG. 6, which shows the aligning details of the texture to the details of the model by specifying an isoparametric curve, we consider that we have a given polygonal curve L={P₀, . . . P_m} associated with a given value u₀of the parameter u. We will now describe the constraints to be honored for making the isoparametric of the mapping defined by (u=u₀) correspond to the projection of L on the surface S. Each point p_iof L yields a constraint c_piensuring that the isoparametric curve u=u₀of the mapping φ passes near the projection p′_iof p_ion S.

[0094] The triangle T=(α₀,α₁,α₂) is the triangle of S of FIG. 6 that contains p′_i, and (λ₀,λ₁,λ₂) are the barycentric coordinates of p′_iin T. The linear relation to be honored is given in Equation 18 below: $\begin{array}{cc}|\begin{array}{c}\sum _{j\in \left\{0,1,2\right\}}{\lambda }_j·{\varphi }^u\left({\alpha }_j\right)=u_0\\ \mathrm{where}:\{\begin{array}{c}\sum _{j\in \left\{0,1,2\right\}}{\lambda }_j·p^u\left({\alpha }_j\right)=p_i^{\prime }\\ \sum _{j\in \left\{0,1,2\right\}}{\lambda }_j=1\end{array}\end{array}& \left(18\right)\end{array}$

[0095] Equation 19 below gives the expression of the constraint c_piin the form of Equation 10. Such a constraint per point p_iis added to the set C to be honored by D.S.I., described above. $\begin{array}{cc}\{\begin{array}{cc}{A}_{c_{p_i}}\left({\alpha }_j\right)={\lambda }_j& \forall j\in \left\{0,1,2\right\}\\ A_{c_{p_i}}\left({\alpha }_j\right)=0& \forall \alpha \notin \left\{{\alpha }_0,{\alpha }_1,{\alpha }_2\right\}\\ A_{c_{p_i}}\left({\alpha }_j\right)=u_0& \end{array}& \left(19\right)\end{array}$

[0096] As described above, the two following sufficient conditions ensure that a discrete function φ defines a mapping:

[0097] The image of the border of the surface through φ in the parametric (u, v) domain is a convex polygon; and

[0098] Each internal node is a convex combination of its neighbors.

[0099] Introducing the constraints to ensure that the isoparametric curves are orthogonal, with homogeneous spacing techniques the first condition can be replaced by a less restrictive one. As shown in FIG. 7, it is then sufficient for a user to specify four arcs of isoparametric curves {u₁,u₂, v₁, v₂} using the constraint previously introduced to extrapolate a mapping on the surface S, as indicated by the dashed exemplary lines 16. Thus, by using the algorithm not only as an interpolator, but also as an extrapolator, it is possible to construct mappings for surfaces having complex shaped borders 18 by leaving φ unspecified on the border.

[0100] Constructing a Mapping for a Cut Surface

[0101] Let us now consider that the surface has linear discontinuities (“cuts”), and that we want the mapping function Φ to be continuous through these cuts. FIGS. 13 and 14 show good examples of such “cuts”, while FIG. 8 shows the detail of such a cut, separating the surface S into two surfaces S₁and S₂along lines 18′ and 18″. In order to accommodate the discontinuity, the surface will be considered as a single patch (as it was before being “cut”). The set of constraints described below allows us to assign (u, v) coordinates to the vertices of the triangulation in such a way that the two borders of a cut are mapped to the same curve by the interpolated Φ mapping function. In other words, the cuts are “sewn” in (u, v) domain space (texture space).

[0102] As shown in FIG. 8, twin set of points {q_i, i=0 . . . n} and {{tilde over (q)}_i, i=0 . . . n} are sampled on the twin borders of the cut. We describe now how to make the mapping match at each pair (q_i, {tilde over (q)}_i) of points. More precisely, we want to respect the following conditions: $\begin{array}{cc}\forall \nu \in \left\{u,v\right\}|\begin{array}{cc}\left(1\right)& {\varphi }_T^{\nu }\left(q_i\right)={\varphi }_T^{\nu }\left({\tilde{q}}_i\right)\\ \left(2\right)& \mathrm{grad}{\varphi }_T^{\nu }=\mathrm{grad}{\varphi }_{\tilde{T}}^{\nu }\end{array}& \left(20\right)\end{array}$

[0103] where T and {tilde over (T)} denote the triangles containing q_iand {tilde over (q)}_irespectively. The gradient grad φ^ν_Tis computed as described in Section 3 (see Equation 5), using the basis shown in FIG. 4. Using the methods introduced in the previous two sections, one can translate these two conditions into the D.S.I. constraints c^ν_{qi,{tilde over (q)}i}and c^νW_{qi,{tilde over (q)}i}, given below in the equations 21 and 22 respectively. $\begin{array}{cc}|\begin{array}{cc}{A}_{c_{q_i,\widetilde{q_i}}^{\nu }}\left({\alpha }_j\right)={\lambda }_j& \forall j\in \left\{0,1\right\}\\ A_{c_{q_i,\widetilde{q_i}}^{\nu }}\left(\widetilde{{\alpha }_j}\right)=-\widetilde{{\lambda }_j}& \forall j\in \left\{0,1\right\}\\ A_{c_{q_i,\widetilde{q_i}}^{\nu }}\left({\alpha }_j\right)=0& \forall \alpha \notin \left\{{\alpha }_0,{\alpha }_1,\widetilde{{\alpha }_0},\widetilde{{\alpha }_1}\right\}\\ b_{c_{q_i,\widetilde{q_i}}^{\nu }}\left({\alpha }_j\right)=0& \end{array}& \left(21\right)\end{array}$

[0104] where λ_j,jε0,1and {tilde over (λ)}_j,jε0,1denote the barycentric coordinates of q_iin [p(α₀), p(α₁)] and {tilde over (q)}_iin [p({tilde over (α)}₀), p{tilde over (α)}₁)] respectively.

[0105] The four constraints c^uX_{qi,{tilde over (q)}i}, c^uY_{qi,{tilde over (q)}i}, c^νX_{qi,{tilde over (q)}i}and c^νY_{qi,{tilde over (q)}i}yielded by Equation 22 below ensure a constant gradient of the mapping through the cut. In other words, an isoparametric curve points in the same direction in the two corresponding triangles T and {tilde over (T)}. $\begin{array}{cc}|\begin{array}{cc}{A}_{c_{\mathrm{qi},\tilde{q}i}^{v_W}}\left({\alpha }_j\right)=D_j& \forall \mathrm{j\epsilon }\left\{0,1,2\right\}\\ A_{c_{\mathrm{qi},\tilde{q}i}^{v_W}}\left({\tilde{\alpha }}_j\right)={\delta }_W·{\tilde{D}}_j& \forall \mathrm{j\epsilon }\left\{0,1,2\right\}\\ A_{c_{\mathrm{qi},\tilde{q}i}^{v_W}}\left(\alpha \right)=0& \forall \mathrm{\alpha \epsilon }\left\{{\alpha }_0,{\alpha }_1,{\alpha }_2,{\tilde{\alpha }}_0,{\tilde{\alpha }}_1,{\tilde{\alpha }}_2\right\}\\ b_{c_{\mathrm{qi},\tilde{q}i}^{v_W}}\left(\alpha \right)=0& \end{array}\mathrm{where}:& \left(22\right)\\ v\epsilon \left\{u,v\right\};W\epsilon \left\{X,Y\right\};{\delta }_W=\{\begin{array}{cc}-1& \mathrm{if}W=X\\ +1& \mathrm{if}W=Y\end{array}& \end{array}$

[0106] In addition to the preferred embodiment described above wherein the invention comprises a method for building a global parameterization (u(x), v(x)) honoring constraints in a least square sense, the invention further includes a computer configured to perform the steps of the method. The computer comprises a memory, such as a disk drive, configured to store a series of computer executable instructions for performing the method, and a processor for processing the computer executable instructions. The computer can further comprise an output device for outputting the results after the method has been performed. The computer is configured to store in the memory a numerical representation of a two-dimensional feature (surface), as well as a three dimensional surface onto which the feature is to be mapped.

[0107] In yet another embodiment of the invention, the invention comprises a computer readable medium such as a magnetic disk, tape, or diskette onto which are stored a series of computer executable instructions for performing the steps of the method.

[0108] Yet another embodiment of the invention comprises a representation of a three dimensional surface mapped with a two-dimensional parametric feature (surface) using the method of the invention. Such may be for example a computer plot on paper or a screen display, or a file saved on a computer readable medium, which represents the surface.

[0109] One can see in FIGS. 9-12 the results of the method applied to a triangulated mesh 21 representing a face 20 (see FIG. 9). The effect of the orthogonality and homogeneity constraints can be brought to the fore by comparing FIG. 10 (no constraint used) and FIG. 12 (orthogonality and homogeneity enforced), where a checker pattern is mapped to the mesh. The isoparametric curves corresponding to FIG. 12 are shown in FIG. 11, where one can see that the iso-u curves (exemplary lines 22, generally horizontal) are perpendicular to the iso-v curves (exemplary lines 24, generally vertical) at the intersections of the two curves. For all these FIGS. 9-12, the constraints ensuring the continuity of the mapping through cuts have been specified at the mouth 26 and the eyes 28 of the model. This model has 3000 triangles (FIG. 9), and has been parametrized after 100 iterations in approximatively one minute using a Sun R4000 workstation computer.

[0110] Industrial utility; applications

[0111] The methods described herein can be applied to problems other than those associated with texture mapping. For example, in the realm of geosciences, the methods described herein have been implemented into a geology oriented C.A.D. software program. Other exemplary applications include:

[0112] Unfolding surfaces representing the boundaries of geological layers while preserving the volume of the layers;

[0113] Generating grids suitable for finite elements analysis;

[0114] Beautifying triangulated meshes by remeshing in (u, v) domain space;

[0115] Constructing Spline surfaces from triangulated meshes; and

[0116] Performing computations such as geostatistical simulations in (u, v) domain space.

[0117] Not only do these applications highly desire that mappings present non-distorting properties, which is addressed by the method disclosed herein, but in addition, these applications can benefit from the ability of the disclosed method to take into account additional information expressed in the form of linear constraints.

[0118] The method applied to geological data is demonstrated in FIGS. 13-19. In FIG. 13, one can see a mapping of a complex surface 30 having cuts 32, 33 and 34, corresponding to boundaries of a geological layer presenting faults. In FIG. 14, the isoparametric curves (exemplary lines 34 and 36) of the mapping are displayed, and one can see that the mapping is continuous through the cuts 32 and 33 of the surface 30. In FIGS. 15-17, a surface 40 representing a dome of salt is parametrized. For this kind of surfaces which are far from developable, distortions will still remain, and one can choose a compromise between the orthogonality and the homogeneity of the mapping by tuning the weightings _cof the two constraints. In FIG. 15, the orthogonality is respected, but the sizes of the squares differ a great deal, whereas in FIG. 17 the squares have approximatively the same size while the isoparametric curves are far from orthogonal. An average solution is shown in FIG. 16, where the same weighting has been used for the two constraints. One can see in FIG. 18 a mapped surface 50 with an isoparametric curve 52 specified. As shown in FIG. 19, the texture has been aligned to this curve 52. Such constitutes making the isoparametric curves C(u0) or C(v0) coincident with the curve 52. This is yet another method for honoring a linear constraint to reduce distortions in a mapped feature on a three-dimensional triangulated surface.

[0119] In compliance with the statute, the invention has been described in language more or less specific as to structural and methodical features. It is to be understood, however, that the invention is not limited to the specific features shown and described, since the means herein disclosed comprise preferred forms of putting the invention into effect. The invention is, therefore, claimed in any of its forms or modifications within the proper scope of the appended claims appropriately interpreted in accordance with the doctrine of equivalents.