A Coarse-to-Fine Multiscale Mesh Representation and its Applications

10/08/2018 ∙ by Hao-Chiang Shao, et al. ∙ 0

We present a novel coarse-to-fine framework that derives a semi-regular multiscale mesh representation of an original input mesh via remeshing. Our approach differs from the conventional mesh wavelet transform strategy in two ways. First, based on a lazy wavelet framework, it can convert an input mesh into a multiresolution representation through a single remeshing procedure. By contrast, the conventional strategy requires two steps: remeshing and mesh wavelet transform. Second, the proposed method can conditionally convert input mesh models into ones sharing the same adjacency matrix, so it is able be invariant against the triangular tilings of the inputs. Our experiment results show that the proposed multiresolution representation method is efficient in various applications, such as 3D shape property analysis, mesh scalable coding and mesh morphing.



There are no comments yet.


page 8

page 9

page 10

page 11

page 12

This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

The mesh surface model, sometimes called the wireframe model and considered as a representation of 3D image data, is widely-used in applications, such as animation, game design, facial recognition

[1], scene reconstruction [2], medical modeling [3], and bio-image atlasing and visualization [4, 5, 6]. Traditionally, fundamental mesh processing methods, e.g., remeshing, simplification, compression, and editing, have been developed as hierarchical decomposable operators, which allow users to make a change at a coarser resolution and refine it at a finer resolution. Accordingly, multiresolution and wavelet analysis have generated considerable interest in the field of mesh surface representation.

Lounsbery et al. [7] and Eck et al. [8] extended Mallat’s multiresolution wavelet analysis method [9] to the mesh surface space of an arbitrary topological type. The extension decomposes a complicated surface into a sequence of surfaces, each of which approximates the input surface at a lower resolution, and the coarsest approximation is generally referred to as the base mesh [10, 11]. Mesh wavelet analysis requires that the input surface is a semi-regular mesh. However, in most cases, meshes are not semi-regular, irrespective of whether they are derived by 3D scanning devices or obtained from tiling 3D isosurfaces. Therefore, a remeshing step, called the remesher, which maps an arbitrary mesh to a semi-regular one, is required before wavelet analysis can be performed [12, 13, 11]. We call this conventional mesh surface analysis strategy forward wavelet approximation because the remesher is applied first, followed by wavelet analysis in a fine-to-coarse manner. In general, the remeshing process begins from a base mesh, which is derived from simplification of the original input mesh. The isosurface defined by the original input mesh may therefore be approximated based on various remeshing results, each of which corresponds to a specific base mesh. Consequently, the wavelet decomposition derived by forward approaches cannot be invariant against triangular tiling, i.e., the topological information111In this paper, we say that two meshes have the same topological structure/information if they have the same number of vertices, edges, and faces and the same adjacency matrix. of the isosurface defined by the original input mesh.

Instead of the forward approach, we can exploit a backward

process to synthesize the wavelet transform from a base mesh to the input mesh surface via subdivision, i.e., remeshing. Because the proposed backward approach operates in a coarse-to-fine manner, two major considerations need to be clarified in advance: (1) the way to select vertices for constructing the base mesh; and (2) where the new vertices should be interpolated at finer levels. Without imposing any constraint on the shape similarity of the base mesh, any decision about the first consideration would be correct. Users can manually select vertices from the input mesh to construct a base mesh; or they can derive a base mesh from the input by either a mesh simplification algorithm

[14, 15] or a base mesh optimization algorithm [16, 17]. As for the second consideration, the piercing procedure [12]

provides a feasible solution. Based on the piercing procedure, the problem can be solved by a subdivision scheme, which can derive displacement vectors to deliver newly interpolated vertices to where they are supposed to be on the isosurface defined by the original input mesh.

In this paper, we propose a backward wavelet remesher (BWR) that uses the backward wavelet process to generate a semi-regular approximation of the original input mesh from a user-specified input base mesh. BWR has two properties. First, it derives the wavelet information by only implementing the remeshing process, whereas the forward approach requires two steps: remeshing and mesh wavelet analysis. Second, because BWR works in a coarse-to-fine manner, it guarantees that the multiresolution representation of two homomorphic meshes will share the same topological structure if the vertices’ order and adjacency matrix of their base meshes are identical. Because of the second property, the multiresolution mesh representation derived by BWR is useful in a larger number of applications than compression and approximation. When two sequences of multiresolution mesh approximations have the same topological structure, one mesh can be morphed into another at any level by simply using vertex interpolation. It is also straightforward to make a 3D structural comparison between the two sequences of the multiresolution meshes. Our experiment results demonstrate the efficacy of the two properties, and show that the PSNR-performance of the proposed approach is comparable to that of forward approaches for scalable coding.

The remainder of this paper is organized as follows. We review related work In Section 2; describe the backward wavelet framework in Section 3; explain parameter determination in Section 4; consider applications of the proposed approach in Section 5; present our experiment results in Section 6; and discuss limitations of the approach and future work in Section 7. Section 8 contains our concluding remarks.

2 Related Work

The concept of the mesh wavelet transform presented by Lounsbery et al. in 1997 [7] motivated subsequent multiresolution mesh representation studies. The authors focused on the relationship between the mesh wavelet transform and the inverse of the subdivision process. Their method processes a semi-regular mesh in a fine-to-coarse manner, and the analysis matrix required to decompose the mesh surface is the inverse of the synthesis matrix used in the lazy wavelet framework [7, 11]:


In Eq. (1), and are the analysis matrices; and are the synthesis matrices; and is the subdivision matrix (described in Section 3) that characterizes the subdivision scheme.

Based on [7], Khodakovsky et al. designed the Progressive geometry compression (PGC) method [18], which applies the mesh wavelet transform on the semi-regular remeshed approximation derived by Mulitresolution adaptive parameterization of surfaces (MAPS) [10] for progressive coding. Because the connectivity of a remeshed approximation is encoded with the subdivision method, Khodakovsky et al. developed in PGC a hierarchical edge-tree structure (shown in Figure 4(b)) to link edges across different scales. With such a structure, it becomes possible to extend the conventional zero-tree coding method from 2D images to 3D meshes. Subsequently, Guskov et al. and Khodakovsky et al. designed Normal meshes (NM) approaches [12, 13] to overcome the inconvenience of compressing three-entry coefficient vectors generated by the mesh wavelet transform. Based on the coarsest mesh derived by MAPS and a concomitant mapping function that records how vertices on the original input mesh are projected on all coarser resolution meshes and approximated via barycentric coordinate system, NM modifies the vertices’ positions on the coarsest mesh to obtain its base mesh. From such a base mesh, hopefully, vertices on a finer level mesh can be represented as normal offsets of the previous coarser mesh. In other words, the NM algorithm generates a multiresolution mesh in which each level can be written as a normal offset from its coarser version. As a result, the wavelet information of such meshes can be represented as a scalar rather than a three-entry vector. There are other methods of multiresolution mesh representation. For example, Guskov et al. [19] and Botsch et al. [20] developed algorithms based on Burt and Adelson’s image pyramid decomposition [21]; while Dong et al. [22] used a tight wavelet frame representation for mesh surface denoising. However, these methods are not closely related to mesh wavelet decomposition, so we will not discuss them further.

Conventional multiresolution analysis methods operate in a fine-to-coarse manner. Consequently, the multiscale representations derived by such approaches vary with the triangular tilings of the remeshed approximation of the original input. Hence, it is difficult for forward mesh wavelet analysis to derive a multiscale comparison of two homomorphic isosurfaces when their triangular tilings comprise different numbers of vertices, edges, and faces. For example, to characterize local shape variations of the human hippocampus via the spherical wavelet transform [23], Nain et al. designed a very elaborate pre-processing routine for their source hippocampus mesh surface models [24]. To obtain a tiling-invariant multiresolution analysis result, Nain et al.’s approach uses conformal mapping to map the hippocampus models to a sphere and then remeshes the mapped models so that they share the same topological information. In the final step, the spherical wavelet transform is applied after registering the remeshing results. Accordingly, applications based on forward approaches have a bias toward multicale approximation issues, such as compression, progressive transmission and level-of-detail control of a single mesh model [11]. They are less suitable for warping, registration and structural comparison of the multiscale properties of two or more isosurfaces.

To resolve the shortcomings of forward wavelet approaches, we posit that mesh wavelet analysis should be synthesized in a coarse-to-fine manner. The concept of the proposed method is illustrated in Figure 1. Briefly, detailed information about each resolution can be recorded by the subdivision coefficients and the corresponding unit direction vectors . The latter should be adaptive to the local structure of the coarse resolution meshes. The advantages of the proposed BWR are twofold. First, from the given base meshes with the same connectivity, the subdivision procedure can guarantee that the finer level approximation meshes will share the same adjacency matrix. Because of this tiling-invariant property, BWR acts as a transformation procedure that can convert input meshes into a standard reference domain. Second, BWR derives the wavelet information directly from remeshing. By contrast, forward wavelet analysis is implemented in two steps: remeshing and mesh wavelet analysis. In addition, based on BWR, newly interpolated vertices on finer level meshes can be represented as a scalar denoting the length of the displacement vector, similar to that achieved by NM. The coding efficiency of BWR is therefore comparable to that of the NM-based method. In the next section, we review the frameworks of the conventional subdivision scheme and mesh wavelet transform, and then describe the proposed backward wavelet remesher (BWR) in detail.

Fig. 1: BWR representation of a sphere with an octahedron as the base mesh. denotes the midpoint subdivisions, and denotes the procedure of placing the midpoint vertices on the surface of the original mesh.

3 The Backward Wavelet Framework

Based on the lazy wavelet framework described in [7], the vertex interpolated in the subdivision processes can be represented as the sum of its initial position and a suitable displacement vector denoting the detailed information; that is:


This equation is the basis for our backward wavelet remesher (BWR) design. The detailed information about each resolution can be further represented as (as shown in Figure 1) if the displacement vector can be expressed as the product of a unit direction vector and a scalar coefficient , i.e., .

In the following subsections, we first review the matrix notations of the subdivision framework and conventional mesh wavelet transform and then describe the BWR method. This section provides the theoretical derivation to explain how we obtain a multiresolution mesh representation of an irregular input mesh by subdividing a user-specified base mesh.

3.1 Subdivision and Its Matrix Notations

Let be the matrix of the coordinates of the vertices in the approximation subspace . The subdivision process determines the vertices in the finer approximation subspace

by applying a linear transform

on the vertices in as follows:


If the number of vertices in is and that in is , then and are an matrix and an matrix respectively. Note that is the number of newly interpolated vertices in , and the -th row of denotes the coordinate vector of , i.e., the -th vertex in . Additionally, is an subdivision matrix that can be represented as two block submatrices, viz. (of size ) and (of size ). is used to map the vertices from to . When the subdivision is an interpolating scheme,

is an identity matrix; hence, all existing vertices in

will not be re-positioned in . In contrast, attempts to insert new vertices based on ’s vertices coordinates and connectivities. For example, the of the butterfly subdivision scheme is an sparse matrix in which each row only contains eight non-zero entries with the values , , or . The locations and values of non-zero entries in are related to the vertices that are exploited to derive the coordinate of a new vertex, as expressed in Eqs. (12)-(14) and illustrated in Figure 2.

3.2 Wavelet Transform Framework

In the wavelet transform framework, the coordinate matrix of the vertices in is derived according to


Here, are obtained by applying a subdivision process on ; and are the wavelet coordinate matrices [7] that contain the “detailed information” of the newly interpolated vertices.

Let , , and . Then, Eq. (4) can be re-written as follows:


where and denote the coordinates of vertices in and respectively; are the newly inserted vertices in ; and are the wavelet coordinate vectors. If we choose and , Eq. (6) can be expressed in a lazy wavelet framework and becomes




Eq. (8) shows that the mesh wavelet representation can be regarded as a subdivision process in which the vertices’ coordinates in are obtained from a linear combination of those in plus the “new” information in .

If we use an interpolating subdivision that lets , we obtain


Eq. (9) indicates that (1) the positions of the vertices in are unchanged in after the subdivision process because ; and (2) the positions of the newly inserted vertices in are obtained by the subdivision process and corresponding translation offset . As the coordinates of the interpolated vertices are determined by the matrices and , Eq. (9) leads to the proposed backward wavelet remesher for multiresolution mesh representation.

3.3 BWR Framework

Eq. (9) implies that if a subdivision procedure can send the newly interpolated vertices to where they are supposed to be located, can be positioned on the target mesh surface without the detailed information . That is, by letting , Eq. (9) becomes


Note that, as mentioned in [11], there is no limitation on the subdivision method in conventional multiresolution mesh representation approaches. Our only restriction here is the presumption that and (that is, the subdivision method itself must be an interpolating scheme), and the form of does not matter. Therefore, the goal of BWR is to develop a subdivision procedure that can ensure (1) are on the isosurface defined by the original input mesh ; and (2) will not fold patches and distort the geometry.

Although the subdivision strategy required by BWR is inexplicit at this point, we should first check what will happen if the butterfly subdivision method [25] is adopted. Let be the initial coordinate of the -th interpolated vertex in . Then, the lower part of Eq. (9), i.e., , can be rewritten in the vertex-wise form as follows:


As shown in Figure 2, the butterfly subdivision scheme yields


where is a scalar parameter; and are the eight points in the neighborhood of for the butterfly weighting average kernel. Eq. (12) can then be re-written as




Therefore, we have


The above equation states that the coordinate of the vertex is rooted at the mid-point of and and stretched in the direction .

Rewriting Eq. (10) according to Eq. (16), the coordinates of the vertices in can be represented as the following matrix:


In Eq. (17), is the mid-point subdivision kernel [26], defined as


where is an identity matrix, and is an sparse matrix with two nonzero entries and in each row. is an diagonal matrix whose -th diagonal element is , and is an matrix whose -th row is the direction vector, , corresponding to . By letting , Eq.(17) can be represented as a two-scale relationship of the wavelet transform as follows:


The matrix


can be interpreted as a wavelet coordinate matrix. Eq. (19) shows a one-step, coarse-to-fine procedure for BWR. The locations of the vertices in are kept the same in the finer subspace, a strategy that follows the interpolating subdivision scheme. Meanwhile, the newly interpolated vertices are derived from the vertices in , modified with the mid-point subdivision matrix and the detailed information in . As shown in Figure 1, the BWR construction is comprised of a sequence of midpoint subdivision procedures; and for each subdivision, the newly interpolated vertices are placed on the surface of the original mesh.

Eq. (19) implies that the problem of designing a proper , which is described in Eq. (10), can be converted to another problem of finding and . In the next section, we describe how to derive these two matrices. Additionally, because is evaluated according to the neighborhood of the newly interpolated vertices, the subdivision procedure of BWR is adaptive to local geometric information. Finally, we should explain that although we let in Eq.(19), this term should be preserved as per Eq. (17) for possible future user-interactive extensions and robustness.

Fig. 2: The eight neighboring vertices for butterfly subdivision interpolation of the vertex .
(a) (b) (c)
Fig. 3: Determining the weight along the direction vector . The bold solid line is the contour of the original mesh; the dashed line connecting and is a portion of the coarser mesh ; is located at the midpoint ; and the new vertex on is at the intersection of and the original mesh. If there are multiple intersections, the best candidate is always one of the two nearest intersections, i.e., the smallest and the largest . (a) There are two intersections in the positive direction of and the smallest is chosen. For (b) and (c), the intersections have one positive candidate (the smallest ) and one negative candidate (the largest ). The value of is determined by checking whether or points in the direction of . Here, is the vertex normal of the midpoint ; and and are the normal vectors of the triangular patches associated with and respectively.

4 Parameter Determination

BWR is comprised of three steps: selecting the vertices of the base mesh in the original mesh; determining the direction matrix ; and deriving the weighting matrix from the vertices in and the original mesh.

4.1 Selecting the Vertices of the Base Mesh in the Original Mesh

The vertices of the base mesh must be chosen from the original mesh because BWR is based on the interpolating subdivision scheme. The base mesh can be derived by using a mesh simplification method [14, 15, 10], a base mesh optimization algorithm [16, 17], or through user intervention. The base mesh and the original input mesh must have the same genus type.

4.2 Direction Matrix

To preserve the topological structure during synthesis, BWR requires that a newly interpolated vertex must intersect the original mesh in a way that preserves the neighborhood topology of the original mesh. This property is determined by the rows in the direction matrix .

There are two scenarios in which a new vertex derived by the vertices may not be located in the region enclosed by : 1) when the direction vector used to place the newly interpolated vertex lies on a flat or nearly flat area; and 2) when the neighborhood of the newly interpolated vertex covers sharp edges. In the first case, Eq. (14) states that if are coplanar or nearly coplanar, will lie on (or almost lie on) the same plane. This implies that the interpolated vertex might be moved along some traverse direction of the plane of . In the second case, the neighboring vertices on any side of the sharp edge may shift to an improper direction, such as the traverse direction of the local patch, and cause patch folding or patch overlapping.

In both cases, the solution is to modify the direction vector in order to preserve the neighborhood topology in a finer approximation subspace. We determine if the interpolated vertex is located in a flat area or in a neighborhood of an edge; then the direction of the vertex is replaced by the normal vector of according to the method in [11, 27].

4.3 Weighting Matrix

The parameters in the diagonal weighting matrix are chosen so that the newly interpolated vertices are incident on the original mesh. This is called the piercing procedure [12], and it can be achieved by setting the -th element of as so that


and belongs to the original mesh. The entry is determined by extending the direction vector

to intersect the original mesh. The steps of the procedure are as follows: 1) find the hyperplane that includes a triangular patch in the original mesh; 2) extend

to find the intersection point on the hyperplane; and 3) determine if the intersection point is inside the triangular patch. Let the triangular patch in the original mesh be . The hyperplane that contains with a normal vector can be represented as follows:


where , and “t” denotes matrix transpose. The point where intersects the hyperplane is at . Here, can be derived by replacing in Eq. (22) with . If the intersecting point is in the triangle , then can be represented by a barycentric coordinate system of , and ; that is,


where and .

The shape of the original mesh determines whether it can contain multiple valid intersections with triangles. As shown in Figure 3, the best candidate is always one of the two nearest intersections: one from the smallest and the other from the largest . The following procedure can be used to make the final decision. Let and be the normal vectors of two candidate triangular patches corresponding to the smallest and the largest respectively. In addition, let be the normal vector at the midpoint (see Eq. (21)). The triangular patch whose normal vector points in the direction of is chosen to determine the value of .

5 Applications

In this section, we discuss two practicable applications of the proposed method, namely, scalable coding and morphing.

5.1 Scalable Coding

BWR can derive a multiresolution approximation of an original input mesh and represent the location of every newly interpolated vertex as a scalar. Consequently, in terms of storage efficiency, BWR is as efficient as the Normal Meshes (NM) method 222In the Normal Meshes approach, as reported in [12], more than of the vertices can be represented as scalar and rest of them are still represented as vector. [12, 13] for scalable (progressive) coding applications.

Fig. 4: An example of mapping the vertices of a multi-resolution mesh to an image wavelet transform. The subdivision coefficient used to insert a new vertex in edge in (b) is placed in the pixel position, which is labeled , in (a). As shown in the mesh, the edge is subdivided into four edges: and ; and the corresponding coefficients are shown as children of in the image. In general, the children of are , and . The pixel labeled is set at the mean of the three values in its underside, dextral, and diagonal positions.

Scalable coding can be regarded as a special case without any user input at the locations of the newly interpolated vertices in . This corresponds to setting in Eq. (17

) as a zero matrix. Our mesh scalable coding framework regards the scalar matrices,

’s, as a series of wavelet transform coefficients [7] and exploits Khodakovsky et al’s edge-tree structure [18]. In addition, we map all the scalars (at different resolutions of the meshes) to the wavelet transform coefficients of an image. Then, we can use a common wavelet-based scalable codec for image compression to progressively transmit the vertices of our multiscale mesh representation. Progressive transmission can be performed by sending the coarsest mesh first, followed by the wavelet coefficients in a bitplane-by-bitplane manner for the refinement. Figure 4 shows a simple example to demonstrate the parent-child relationship between a vertex at scale and its children at scale . With the map of mesh vertices to image wavelet transforms, any wavelet-based scalable coding algorithm would be suitable for transmitting and encoding the information about the mesh vertices.

5.2 Morphing

Morphing is the process that gradually changes a source object through intermediate objects into a target object [28]. In the conventional mesh morphing process, an essential step is to map both the source and the target meshes into the intermediate one. After deriving a function , which can map the isosurface defined by mesh- to that defined by the intermediate mesh, the morphing result can be obtained according to Eq.(24):


where and are the source and the target meshes respectively, denotes the vertices’ coordinates, and is the inverse function of .

BWR can reduce the complexity of deriving the mapping function. Because BWR can remesh meshes of the same genus type into ones that share the same topological information, it is straightforward to construct a one-to-one mapping between vertices of two remeshed approximations. Given two mesh models, and , we can generate an intermediate mesh between them based on a one-to-one mapping function by


Hence, with the aid of BWR, the conventional mesh morphing problem can be simplified to a “mesh blending” problem. For cases where metamorphism among multiple meshes needs to be derived simultaneously333For example, in order to produce a standard average atlas model for 3D image atlasing, one may have to register, warp and average multiple 3D models simultaneously, as described in [6]. Such an average can be considered as an equal-weighted metamorphism among all models., the blended intermediate meshes can be obtained by


with and .

6 Experiment Results

In this section, we demonstrate the multiresolution approximations of BWR on some meshes and compare their performance with that of other remeshers. We also compare the coding performance of several multi-resolution mesh representation methods and analyze the morphing results. In the following experiments the mesh-to-mesh distance was measured by a well-known method, namely, METRO [29, 30]; and the bit-rate was measured by the bit-per-vertex (bpv).

(a) (b) (c) (d) (e)
(f) (g)
Fig. 5: Comparison of errors in meshes derived at different resolutions. (a)-(d) the coarsest mushroom body mesh plus meshes , , and derived by our method; (e) the reference mesh; (f) the forward root-mean-square and mean distances of the reconstructed meshes compared to those of the reference mesh; (g) the backward root-mean-square and mean distances of the reference mesh compared to those of the reconstructed meshes. The PM algorithm yielded the best performance, as shown in (f) and (g). However, for all methods, the average error per voxel was less than . Since the bounding box diagonal of the reference mesh was voxels, less than voxel difference indicates that the performances of all methods are comparable. The number of vertices in the meshes derived by PM is the same as that in the corresponding . Note that .
(a) (b) (c) (d) (e)
(f) (g) (h) (i) (j)
Fig. 6: The reconstruction results of Venus Head. From left to right, top to bottom: the base-domain octahedron , the subdivision results , and the reference surface model . The numbers of vertices and patches of are and respectively; and the reference contains vertices and patches.
Fig. 7: Synthesis results of ”Vebbit”. The coarser meshes are from Venus, and the subdivision coefficients are from Rabbit, that is, . From left to right, . The top row shows the , which serves as the beginning of the reconstruction shown in each column; and the bottom two rows show the synthesis results of different views.

6.1 Remesh and multi-resolution representation

First, we compared the remesh performance of BWR with that of PM and MAPS. The error between the original mesh and the reconstructed mesh was measured by METRO. Figures 5(a)-(d) show the mesh models of a mushroom body, and a neuropil with a pair of tri-axial structures in the Drosophila brain [6], reconstructed by our method at different resolutions. The original reference mesh shown in Figure 5(e) was reconstructed from a confocal image stack by Chen’s method [31]. The errors measured between the reference mesh and the reconstructed mesh of various mesh representations are compared in Figures 5(f) and (g). Our method’s performance was between those of MAPS and PM444PM (Progressive Mesh) method is not a remesher. However, since the vertex removal strategy of PM is to remove iteratively a vertex whose absence results in minimal change in current mesh geometry, a simplified mesh derived by PM is regarded to the best approximation of the source mesh at that resolution.. Moreover, Figures 5(f) and (g) show that the average error per voxel of the three methods was less than because the bounding box diagonal for the mushroom body model is voxels. Therefore, the performances of these methods are comparable. Finally, note that the -errors of the remeshed Venus Head and Rabbit models derived by BWR (shown in Figures 6 and 10) are and with respect to the bounding box diagonal lengths and respectively. The models’ PSNR values are thus dB and dB respectively. The remesh performance of BWR is the same as that of NM and MAPS.

Next, we consider the semi-regular approximation results of BWR. Figure 6 shows the multi-resolution approximation of the Venus Head model. Because BWR can start from a base mesh that has the same genus as the original mesh, we chose an octahedron as the base mesh. The experiment results indicate that different kinds of shape components and details appear in different scales. From to , the global shape of the input model is reconstructed, whereas facial features and detailed components are reconstructed when and respectively.

Fig. 8: Synthesis results of ”Ranus”. The coarser meshes are from Rabbit, and the subdivision coefficients are from Venus. In this example, respectively.
(a) Skull (b) (c)
(d) (e) (f)
Fig. 9: Visualization of subdivision coefficients derived from the Skull model. Sub-figures (b)-(f) were generated by MATLAB’s built-in function trisurf(Tri, X, Y, Z, C), where (Tri, X, Y, Z) defines the -th scale unit sphere that the subdivision coefficient matrix (i.e., the input argument C) was mapped to.

The following experiment results demonstrate that (1) the proposed method can achieve multiresolution representation; and (2) the behavior of is similar to that of conventional wavelet coefficients. Let denote the subdivision coefficients used to reconstruct from . We borrow the “direct sum” symbol to link the subdivision procedure and the coefficient matrix. Then, in Figure 7, each column represents a synthesis result of . The top row shows the ; and the middle and bottom rows show, respectively, the top views and front views of the same synthesis results. Figure 8 illustrates the results of synthesizing and , and the bottom row is the rear view. These experiment results demonstrate that the obtained coefficients can describe how the surface changes from a coarse scale to a fine scale. Figure 7 shows that if we give the subdivision coefficients of Rabbit to Venus Head’s base-domain earlier, the synthesized Vebbit would be more similar to Rabbit. This finding implies that the subdivision coefficients obtained in coarser resolutions represent global shape components, whereas those obtained in finer resolutions represent local details. Moreover, based on Figure 8, we find that the hair bud can be roughly decomposed into two components: shape (the left two columns) and textures (the right three columns). The smaller local variations of the texture components make the Rabbit model Godzilla-like from the rear view. This observation suggests that the proposed method could be applicable to mesh editing, provided that the user-specified region of interest and the corresponding subdivision coefficients are given.

(a) (b)
(c) (d)
(e) (f) (g) (h) (i)
(j) (k) (l) (m) (n)
Fig. 10: Comparison of the scalable coding performance on Venus Head and Rabbit. Sub-figures (a)-(b) are the bpv-PSNR curves of various methods. BWR’s performance is comparable to that of the state-of-the-art methods. Sub-figures (c)-(d) show the three-dimensional bpv-PSNR-resolution curves of BWR. The thick red curve denotes the recommended transmission strategy for the model. Sub-figures (e)-(n) show the BWR reconstruction results of Rabbit.
(a), 0.125 bpv (b), 0.5 bpv (c), 2.0 bpv (d), 0.125 bpv (e), 0.5 bpv (f), 2.0 bpv
Fig. 11: Comparison of the visual qualities of the Venus Head and Rabbit meshes derived by BWR at different bpv values. The meshes derived at bpv and bpv can be distinguished (the ears and chest in the Rabbit mesh, and the mouth and surface texture in the Venus Head mesh); however, meshes derived at bpv and bpv are almost indistinguishable.
(a) (b)
Fig. 12: Morphing sequences of the Skull, Rabbit, and Venus Head models. (a) Front view. (b) Back view.
Fig. 13: The remeshed results derived by BWR share the same topological structure, i.e., . In the figure, and are, respectively, the vertices’ coordinates and the connectivity information of each mesh model. Consequently, we posit that BWR acts as if it is a transformation procedure capable of converting input meshes into a standard reference space.

Finally, Figure 9 shows that BWR can capture local features in a multiscale fashion. In contrast to the experiment results presented in Figures 7 and 8, we mapped each to a unit sphere and then visualized it to highlight the local features of each scale. The Skull model in Figure 9(a) has five features: 1) eye socket areas, 2) nose area, 3) teeth, 4) jaw, and 5) “headus” and “snake-like pattern”. The visualization of in Figure 9(b) proves once again that subdivision coefficients of coarser scales represent global shape structures. The coefficients within the eye socket areas are all negative so that they can make these regions appear concave. Thereafter, sharpens the boundaries of the eye sockets and jaw. Positive subdivision coefficients in these zones lift up newly-interpolated vertices to raise the local surface. Meanwhile, characterizes the shapes of the teeth, generates the boundaries of the nose, and continues to sharpen the eye sockets and jaw areas. The two layers of the coefficients and represent very local details, such as the “headus” and “snake-like pattern” on the forehead and boundaries of teeth. We also find that the fringes in Figures 9(e)-(f) are thinner and narrower than those in Figures 9(c)-(d) because and represent finer structures on the mesh surface. To summarize, the experiment results presented in this subsection show that BWR can generate a multiscale approximation of the input reference mesh, as well as decompose the input reference mesh to capture local features in a multiscale manner. Such local and multiresolutional features can further enable users to evaluate structural commonness and variations of 3D surface models. Examples of deriving 3D shape properties via BWR can be found in [32].

6.2 Scalable Coding

We compared the coding performance of our method with that of existing methods on two benchmark models: Venus Head and Rabbit. Our implementation uses the SPIHT algorithm for scalable coding [33]. The source code for SPIHT can be downloaded from the official Matlab forum, Matlab Central [34]. To ensure a fair comparison, we did not implement the compared methods. Instead, we used the performance curves of those methods provided in [13] and compared them with the curves derived by our method, as shown in Figures 10(a) and (b). The measurements of our curves were derived from a sequence of mesh approximations ( and for Venus Head; and, and for Rabbit). The multi-resolution approximations of Rabbit from to are shown in the third and fourth rows of Figure 10.

Figures 10(c) and (d) show, respectively, the three-dimensional (bpv-PSNR-resolution) performance curves for the Venus Head and Rabbit meshes. On both curves, the PSNR first increases with increases in the bpv and resolution, and then forms a plateau when the bpv is slightly more than and the resolution is higher than . Thus, the Venus Head and Rabbit meshes can be decoded with a satisfactory mesh when the bpv and resolution are set at and respectively. Figure 11 shows that the visual quality of meshes reconstructed at bpv is almost indistinguishable from that derived at a higher bpv.

We observe that, at a low bit rate, our Venus Head curve in Figure 10(a) has a better PSNR than other curves, whereas the PSNR curve of our Rabbit model in Figure 10(b) is not as good as that of the others. This is because we did not impose any shape constraints while constructing the base mesh, so the surface-to-surface error between our , an octahedron, and the original input mesh is surely larger. Hence, it is reasonable that the PSNR curve of Rabbit is not as good as those of other methods in a low bit-rate environment due to an octahedral base mesh. This weakness can be corrected by compressing a higher resolution mesh. As indicated by Figures 10(c) and (d), the of the Venus Head model yields a better scalable coding performance at a low bit-rate, whereas the performance of models is worse, even at a high bit-rate. Note that the PSNR values of the two uncompressed models are about 86dB. This fact implies that when the bandwidth is limited, it is still possible to raise the PSNR effectively by first sending low bpv subdivision coefficients to increase the total number of vertices, and then refining the subdivision coefficients gradually555For example, in our experiment an model contains vertices and an contains vertices. To transmit an model in 0.25bpv (bit-per-vertex) environment requires about bits. This amount is almost equivalent to that required to transmit an model in 16.0 bpv.. This transmission strategy is illustrated by the thick red lines in Figures 10(c) and (d).

6.3 Morphing

Figures 12 illustrates the metamorphism of three remeshed models. The top apex of the triangle is Skull, the bottom-left apex is Venus Head, and the bottom-right apex is Rabbit. We used an octahedron as the base mesh for each input mesh. The morphing objects in the triangle were obtained as follows: (1) from each input mesh, six vertices were selected in the order of top, right, front, left, back, and bottom; (2) the six vertices were used to construct an octahedron as the base mesh for each input mesh; (3) the BWR algorithm was applied to the octahedrons to derive the remeshed results; (4) we let the mapping function mentioned in Eq.(25) be because the remeshed results share the same vertex order and the same topological information; and, (5) the vertices of the three remeshed models were blended as follows:


where , , are barycentric coefficients of the Skull-Venus-Rabbit triangle system shown in Figure 12.

This experiment shows that BWR can convert input meshes into a standard reference domain, as illustrated in Figure 13. By ordering the apexes, edges, and faces of the base mesh correctly, the user-specified pre-alignment process ensures that the connectivity information of the multiresolution approximation meshes is the same. Therefore, the mesh morphing problem as well as mesh registration and warping can be solved in a straightforward manner.

6.4 Time Complexity

The piercing procedure used to determine the subdivision coefficient matrix causes a computational bottleneck in BWR. In practice, the procedure can be computed rapidly if its implementation is optimized. For example, Guskov et al. reported that the NM method, comprised of a base mesh generating phase and a piercing procedure phase, can generate remeshing results containing vertices from a -vertex base mesh in minutes by using a Pentium III CPU [12]. However, because we adopted the most straightforward strategy, i.e., a full-search strategy, to implement the piercing procedure, the computational complexity of our program is theoretically , which is higher than that of Guskov et al.’s implementation. Here, “f” is the number of triangular faces of the input reference mesh, and “m” is the number of newly-interpolated vertices of the current scale. Table I and Figure 14 show the execution times of our experiments, each of which begins by reading the and files (.txt) and ends by writing the file (.txt) to the hard disk. The experiment results show that the computational complexity of our full-search-based implementation is , an upper bound for the piercing procedure.

In sum, the computational time cost of BWR can be reduced to several minutes if the implementation is optimized in the same way as the NM method. In addition, we should emphasize that BWR can produce multiscale mesh representation and remeshing results simultaneously, whereas the total time cost of the conventional mesh wavelet transform should include the computation times of remeshing and mesh wavelet analysis. In this paper, we do not consider the problem of how to design a data structure to facilitate the piercing procedure. We will address the problem in a future study.

Scale and Venus Head Rabbit Skull
# of new vertices
: 12 7.26 4.94 1.78
: 48 11.82 7.78 2.77
: 192 31.07 20.58 6.81
: 768 114.5 71.11 22.88
: 3,072 420.5 275.3 89.92
: 12,288 1649 1093 373.3
: 49,152 6804 4478 1670
: 196,608 29689 22079 9673
TABLE I: Time costs of the examples discussed in Section 6. All times are in seconds on a virtual Windows XP machine with 2GB RAM and one core shared by an Intel i5-3470 CPU. The numbers of vertices and faces of the input references Venus Head, Rabbit, and Skull models are , , and respectively.
Fig. 14: Computational time cost. Each data point in the figure includes the program execution time and I/O time.

7 Limitation and Future Work

The limitation of the proposed method is that the and the original input mesh must have a correct intersection. For models with conspicuously protruding structures, e.g., the legs of the horse model and the ears of the bunny model, it is difficult to generate the multiresolution approximation from a very simple base mesh, e.g., an octahedron. There are two possible ways to resolve this issue. The first is to derive a multiresolution representation by BWR from a more complicated base mesh, e.g. a simplified mesh derived by PM or MAPS, similar to the experiment results we described in Figure 5. The second way is to develop a more robust algorithm to derive and iteratively adjust for each interpolated vertex. The constraints on developing a robust algorithm for of the scale are as follows. (1) The of adjacent newly interpolated vertices should vary gradually and homogeneously; that is, the angle between the of any two adjacent vertices should be small so that triangular patches will not be folded. (2) The must be derived solely by ’s local information; otherwise, the decoder of the scalable transmission will impose an extra overhead, e.g., either an overhead on or the in Eq.(17), to reconstruct . The first constraint is the reason that we chose the vertex normal to replace the derived by the butterfly scheme in flat areas and regions covering shape edges. The second constraint is optional because such overhead information only influences the efficiency of scalable coding and will not alter the topological structure of our multiresolution mesh.

Given the above limitation and the discussion in Section 6.4, there are two aspects of BWR that could be improved in our future work. The first is to develop an algorithm that can guarantee the direction vectors of adjacent newly interpolated vertices will vary gradually and homogeneously. The second is to design a data structure to integrate BWR and the piercing procedure so that the time complexity can be reduced.

8 Concluding Remarks

We have proposed a backward coarse-to-fine method (called BWR) that can construct a multi-resolution approximation of a mesh surface. With the help of the original mesh, our method starts by using any base mesh of the same genus type as the original mesh and terminates in a semi-regular approximation of the original mesh. In addition, we have shown that the locations of new vertices in BWR can be derived as a simple closed form and represented as a scalar. Our experiment results show that this representation is effective for scalable coding and its performance is comparable to or better than that of some state-of-the-art methods. We have also shown that BWR can (1) generate a tiling-invariant multiscale mesh representation; and (2) act as if it is a transformation procedure capable of converting the original input mesh into a standard reference domain. As a result, morphing becomes a straightforward interpolation of corresponding vertices in the approximating meshes derived by BWR. We believe that the simplicity and flexibility of BWR makes it suitable for mesh editing operations and applicable to biomedical applications; for example, tracking deformations of a 3D beating heart model and deriving a 3D calibrated reference space for image atlasing.


  • [1] V. Blanz and T. Vetter, “Face recognition based on fitting a 3d morphable model,” IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 25, no. 9, pp. 1063–1074, 2003.
  • [2] D. Gallup, J.-M. Frahm, and M. Pollefeys, “Piecewise planar and non-planar stereo for urban scene reconstruction,” in

    2010 IEEE Conference on Computer Vision and Pattern Recognition (CVPR)

    .   IEEE, 2010, pp. 1418–1425.
  • [3] Y. Zheng, A. Barbu, B. Georgescu, M. Scheuering, and D. Comaniciu, “Four-chamber heart modeling and automatic segmentation for 3-d cardiac ct volumes using marginal space learning and steerable features,” IEEE Transactions on Medical Imaging, vol. 27, no. 11, pp. 1668–1681, 2008.
  • [4] A.-S. Chiang, C.-Y. Lin, C.-C. Chuang, H.-M. Chang, C.-H. Hsieh, C.-W. Yeh, C.-T. Shih, J.-J. Wu, G.-T. Wang, Y.-C. Chen et al., “Three-dimensional reconstruction of brain-wide wiring networks in Drosophila at single-cell resolution,” Current Biology, vol. 21, no. 1, pp. 1–11, 2011.
  • [5] G.-Y. Chen, C.-C. Wu, H.-C. Shao, H.-M. Chang, A.-S. Chiang, and Y.-C. Chen, “Retention of features on a mapped drosophila brain surface using a bézier-tube-based surface model averaging technique,” IEEE Transactions on Biomedical Engineering, vol. 59, no. 12, pp. 3314–3326, 2012.
  • [6] H.-C. Shao, C.-C. Wu, G.-Y. Chen, H.-M. Chang, A.-S. Chiang, and Y.-C. Chen, “Developing a stereotypical Drosophila brain atlas,” IEEE Transactions on Biomedical Engineering, vol. 61, no. 12, pp. 2848–2858, 2014.
  • [7] M. Lounsbery, T. D. DeRose, and J. Warren, “Multiresolution analysis for surfaces of arbitrary topological type,” ACM Transactions on Graphics, vol. 16, no. 1, pp. 34–73, 1997.
  • [8] M. Eck, T. DeRose, T. Duchamp, H. Hoppe, M. Lounsbery, and W. Stuetzle, “Multiresolution analysis of arbitrary meshes,” in Proceedings of the 22nd annual conference on Computer graphics and interactive techniques.   ACM, 1995, pp. 173–182.
  • [9] S. Mallat, “A theory for multiresolution signal decomposition: the wavelet representation,” IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 11, no. 7, pp. 674–693, 1989.
  • [10] A. W. Lee, W. Sweldens, P. Schröder, L. Cowsar, and D. Dobkin, “Maps: Multiresolution adaptive parameterization of surfaces,” in Proceedings of the 25th annual conference on Computer graphics and interactive techniques.   ACM, 1998, pp. 95–104.
  • [11] E. J. Stollnitz, T. D. DeRose, and D. H. Salesin, Wavelets for computer graphics: theory and applications.   Morgan Kaufmann, 1996.
  • [12] I. Guskov, K. Vidimče, W. Sweldens, and P. Schröder, “Normal meshes,” in Proceedings of the 27th annual conference on Computer graphics and interactive techniques.   ACM Press/Addison-Wesley Publishing Co., 2000, pp. 95–102.
  • [13] A. Khodakovsky and I. Guskov, “Compression of normal meshes,” in Geometric modeling for scientific visualization.   Springer, 2003, pp. 189–206.
  • [14] H. Hoppe, “Progressive meshes,” in Proceedings of the 23rd annual conference on Computer graphics and interactive techniques.   ACM, 1996, pp. 99–108.
  • [15] M. Garland and P. S. Heckbert, “Surface simplification using quadric error metrics,” in Proceedings of the 24th annual conference on Computer graphics and interactive techniques.   ACM Press/Addison-Wesley Publishing Co., 1997, pp. 209–216.
  • [16] M. Marinov and L. Kobbelt, “Optimization techniques for approximation with subdivision surfaces,” in Proceedings of the ninth ACM symposium on Solid modeling and applications.   Eurographics Association, 2004, pp. 113–122.
  • [17] M. Marinov and L. Kobbelt, “Automatic generation of structure preserving multiresolution models,” in Computer Graphics Forum, vol. 24, no. 3.   Wiley Online Library, 2005, pp. 479–486.
  • [18] A. Khodakovsky, P. Schröder, and W. Sweldens, “Progressive geometry compression,” in Proceedings of the 27th annual conference on Computer graphics and interactive techniques.   ACM Press/Addison-Wesley Publishing Co., 2000, pp. 271–278.
  • [19] I. Guskov, W. Sweldens, and P. Schröder, “Multiresolution signal processing for meshes,” in Proceedings of the 26th annual conference on Computer graphics and interactive techniques.   ACM Press/Addison-Wesley Publishing Co., 1999, pp. 325–334.
  • [20] M. Botsch and L. Kobbelt, “A remeshing approach to multiresolution modeling,” in Proceedings of the 2004 Eurographics/ACM SIGGRAPH symposium on Geometry processing.   ACM, 2004, pp. 185–192.
  • [21] P. J. Burt and E. H. Adelson, “The laplacian pyramid as a compact image code,” IEEE Trans. Commun., vol. 31, no. 4, pp. 532–540, 1983.
  • [22] B. Dong, Q. Jiang, C. Liu, and Z. Shen, “Multiscale representation of surfaces by tight wavelet frames with applications to denoising,” Applied and Computational Harmonic Analysis, 2015.
  • [23] P. Schröder and W. Sweldens, “Spherical wavelets: Efficiently representing functions on the sphere,” in Proceedings of the 22nd annual conference on Computer graphics and interactive techniques.   ACM, 1995, pp. 161–172.
  • [24] D. Nain, S. Haker, A. Bobick, and A. Tannenbaum, “Multiscale 3-d shape representation and segmentation using spherical wavelets,” IEEE Transactions on Medical Imaging, vol. 26, no. 4, pp. 598–618, 2007.
  • [25] N. Dyn, D. Levine, and J. A. Gregory, “A butterfly subdivision scheme for surface interpolation with tension control,” ACM transactions on Graphics (TOG), vol. 9, no. 2, pp. 160–169, 1990.
  • [26] J. Warren and H. Weimer, Subdivision methods for geometric design: A constructive approach.   Morgan Kaufmann, 2001.
  • [27] M. Meyer, M. Desbrun, P. Schröder, and A. H. Barr, “Discrete differential-geometry operators for triangulated 2-manifolds,” in Visualization and mathematics III.   Springer, 2003, pp. 35–57.
  • [28] A. W. Lee, D. Dobkin, W. Sweldens, and P. Schröder, “Multiresolution mesh morphing,” in Proceedings of the 26th annual conference on Computer graphics and interactive techniques.   ACM Press/Addison-Wesley Publishing Co., 1999, pp. 343–350.
  • [29] P. Cignoni, C. Rocchini, and R. Scopigno, “Metro: Measuring error on simplified surfaces,” in Computer Graphics Forum, vol. 17, no. 2.   Wiley Online Library, 1998, pp. 167–174.
  • [30] Metro, “
  • [31] Y.-C. Chen, Y.-C. Chen, A.-S. Chiang, and K.-S. Hsieh, “A reliable surface reconstruction system in biomedicine,” Computer methods and programs in biomedicine, vol. 86, no. 2, pp. 141–152, 2007.
  • [32] H.-C. Shao and W.-L. Hwang, “Deriving 3d shape properties by using backward wavelet remesher,” in Proceeding of The 5th IEEE Global Conference on Signal and Information Processing.   IEEE, 2017, pp. 191–195.
  • [33] A. Said and W. A. Pearlman, “A new, fast, and efficient image codec based on set partitioning in hierarchical trees,” IEEE Transactions on Circuits and Systems for Video Technology, vol. 6, no. 3, pp. 243–250, 1996.
  • [34] MatlabCentral, “