UFDC Home  Search all Groups  UF Institutional Repository  UF Institutional Repository  UF Theses & Dissertations  Vendor Digitized Files   Help 
Material Information
Subjects
Notes
Record Information

Full Text 
A DYNAMIC FRAMEWORK FOR SUBDIVISION SURFACES By CHHANDOMAY MANDAL A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA @Copyright 1998 by Chhandomay Mandal To My Parents, and My Brother ACKNOWLEDGMENTS I am deeply grateful to my advisor Dr. Baba C. Vemuri, whose guidance, en couragement and support made the completion of this work possible. I thank him for his academic as well as nonacademic advice and assistance during my study at the University of Florida, and for his impressive example of integrity and academic responsibility. I wish to express sincere appreciation to my coadvisor, Dr. Hong Qin, for introducing me to the wonders of computer graphics, for patiently guiding me all the way till the end, and for assisting me in various ways. I would also like to thank Dr. Sartaj Sahni, Dr. Paul A. Fishwick and Dr. Douglas A. Cenzer for serving in my supervisory committee and advising me on various aspects of this dissertation. My dissertation research was supported in part by the NSF grant ECS9210648 and the NIH grant RO1LM05944 of Dr. Vemuri, and by the NSF CAREER award CCR9702103 and DMI9700129 of Dr. Qin. I wish to acknowledge Dr. Tim McIn erney, Dr. Gregoire Malandain, Dr. Hughes Hoppe, Dr. Kari Pulli, Dr. Dimitry Goldgof, Dr. Christina Leonard and Dr. ShangHong Lai for helping me at various stages of the work by providing various data sets, software and figures. I would like to thank the faculty, staff and students of the computer and in formation science and engineering department, who were always there to help me. In my years with the computational vision, graphics and medical imaging group, I have enjoyed working with a set of bright students, including Yanlin Guo, Yi Cao, Li Chen, Li Wang, Shuangying Huang, Chinar Kapoor, Fengting Chen, Arun Srini vasan, Jundong Liu and Jun Ye. Special thanks goes to my long time roommates, Raja Chatterjee and Kingshuk Majumdar, for helping me in various ways during my stay at Gainesville. Finally, I would like to thank my parents and my brother for their constant en couragement and support, for cheering me up during the hard days and for providing an ever increasing incentive to finish my graduate study. I would have been nowhere near my goals without them. TABLE OF CONTENTS ACKNOWLEDGMENTS ........... .................. iv LIST OF FIGURES .. ..... ............... .. .... ix ABSTRACT ............................... ...... xiii CHAPTERS 1 INTRODUCTION ..................... ........ 1 1.1 Problem Statement ............................ 2 1.2 Proposed Solution ............ ..... ........... 3 1.3 Contributions ............................... 5 1.4 Outline of Dissertation .......................... 8 2 BACKGROUND ................................ 11 2.1 Subdivision Surfaces ........................... 11 2.2 Physicsbased Deformable Surface Models ............... 17 2.3 Shape Modeling Using Physicsbased Subdivisionsurface Model 20 2.4 Shape and Motion Estimation Using Physicsbased Subdivisionsurface Model ........... ....................... 22 2.5 Unified Framework for Shape Recovery and Shape Modeling ..... 24 3 DYNAMIC CATMULLCLARK SUBDIVISION SURFACES ....... 27 3.1 Overview of the CatmullClark Subdivision Scheme ......... 27 3.2 Formulation ................. ..... ......... 30 3.2.1 Assigning Patches to Regular Faces .............. 31 3.2.2 Assigning Patches to Irregular Faces .............. 33 3.2.3 Kinematics ............. ................ 36 3.2.4 Dynamics ................... ...... 44 3.2.5 Multilevel Dynamics ............ ........ .... 46 3.3 Finite Element Implementation ... .... .. 49 3.3.1 Data Structures .................. .. ..... .. 49 3.3.2 Normal Elements .................. .... .. 51 3.3.3 Special Elements ........................ .. 51 3.3.4 Force Application .................. ....... ..52 3.3.5 Discrete Dynamic Equation .................. ..54 3.3.6 Model Subdivision .......... ........ 55 3.4 Generalization of the Approach ...... .. 56 4 DYNAMIC BUTTERFLY SUBDIVISION SURFACES ... 57 4.1 Overview of the (Modified) Butterfly Subdivision Scheme ...... 57 4.2 Formulation ............. 61 4.2.1 Local Parameterization ... ..62 4.2.2 Dynamics ................ 72 4.2.3 Multilevel Dynamics ...................... .. 74 4.3 Finite Element Implementation ..... ... 76 4.3.1 Elemental Mass and Damping matrices .... 77 4.3.2 Elemental Stiffness Matrix ..... .. 79 4.3.3 Other Implementation Issues ... ..... 81 4.4 Generalization of the Approach ... ..... 82 5 UNIFIED DYNAMIC FRAMEWORK FOR SUBDIVISION SURFACES 83 5.1 Overview of the Unified Approach .. .. 83 5.2 Unified Dynamic Framework for CatmullClark Subdivision Surfaces 85 5.2.1 Local Parameterization ... ... 88 5.2.2 Dynamics ........ . 92 5.2.3 Finite Element Implementation ..... 92 5.3 Unified Dynamic Framework for Loop Subdivision Surfaces 95 5.3.1 Local Parameterization ............ .... .. .. 97 5.3.2 Dynamics ............. ..... .... ......... 100 5.3.3 Finite Element Implementation .... 100 5.4 A General Outline of the Framework for Approximating Subdivision Schemes ............... .................. 101 5.5 A General Outline of the Framework for Interpolatory Subdivision Schemes ............ ............ 102 6 MULTIRESOLUTION DYNAMICS .................. .. ..103 6.1 Overview of Multiresolution Analysis and Wavelets .... 109 6.2 Multiresolution Analysis for Surfaces of Arbitrary Topology 114 6.2.1 Nested Spaces through Subdivision ... 116 6.2.2 Inner Product ................ ........ 118 6.2.3 Biorthogonal Surface Wavelets on Arbitrary Manifold 119 6.3 Multiresolution Representation ... ... .. 125 6.4 Dynamics ................................ .. 128 6.5 Implementation Details .................. ....... .. 129 7 SYSTEM ARCHITECTURE ........... ............... 133 7.1 Topological Information Processing Module ........... 135 7.2 Subdivision Module .................. .. ....... .. 135 7.3 Finite Element Analysis Module ... 135 7.4 Force Synthesis Module .... 136 7.5 Update Engine ................. 136 7.6 Display Module ............ ...... 137 8 APPLICATIONS ............... 138 8.1 Geometric Modeling .... . ..... 138 8.2 Shape Recovery from Range and Volume Data ... 142 8.3 Nonrigid Motion Estimation. ....... 148 8.4 Multiresolution Editing .................. ..... .. 150 8.5 Natural Terrain Modeling .................. ... ..151 9 CONCLUSIONS AND FUTURE WORK. ...... ..... 155 9.1 Conclusions ............... ....... 155 9.2 Future Directions .......... . 156 9.2.1 Automatic Change of Topology .... 156 9.2.2 Local Refinement ........ ...... 157 9.2.3 Automatic Model Parameter Selection .... 157 9.2.4 Constraint Imposition ..................... .. 157 9.2.5 Recovery of Sharp Features ... ..... 158 9.2.6 Automatic Model Initialization 158 9.2.7 Improved Synthesized Force Fields .. 158 REFERENCES .................. ... 159 BIOGRAPHICAL SKETCH .................. .......... .. 168 LIST OF FIGURES 2.1 Refinement of an initial control mesh to obtain the limit surface. (Courtesy : H. Hoppe) .................. ....... 12 3.1 CatmullClark subdivision : (a) initial mesh, (b) mesh obtained after one step of CatmullClark subdivision, and (c) mesh obtained after another subdivision step. .................. ...... ..28 3.2 A rectangular mesh and its limit surface consisting of 4 bicubic surface patches. .................................. 32 3.3 A mesh with an extraordinary point of degree 3 and its limit surface. 34 3.4 Local subdivision around the extraordinary point and the limit surface. 35 3.5 Local subdivision around the extraordinary point and the correspond ing patches in the limit surface from different levels of subdivision. .. 43 3.6 The data structure used for dynamic subdivision surface implementation. 50 3.7 Model subdivision to increase the degrees of freedom : (a) evolution of the initial mesh, and (b) the corresponding limit surface evolution perceived by the user. ...... ........ 55 4.1 (a) The control polygon with triangular faces; (b) mesh obtained after one subdivision step using butterfly subdivision rules. ... 58 4.2 (a) The contributing vertices in the jth level for the vertex in the j+1 th level corresponding to the edge between v and v2; (b) the weighing factors for different vertices. .. ...... 58 4.3 (a) The weighing factors of contributing vertex positions for an edge connecting two vertices of degree 6; (b) the corresponding case when one vertex is of degree n and the other is of degree 6. ... 60 4.4 The smoothing effect of the subdivision process on the triangles of the initial mesh. ................... ............... 62 4.5 Tracking a point x through various levels of subdivision : (a) initial mesh, (b) the selected section (enclosed by dotted lines) of the mesh in (a), after one subdivision step, (c) the selected section of the mesh in (b), after another subdivision step. 64 4.6 Different set of control points at a subdivided level obtained by apply ing different subdivision matrices on a given set of control points in a coarser mesh. ............. .... 67 68 4.7 (a) An initial mesh, and (b) the corresponding limit surface. The do mains of the shaded elements in the limit surface are the corresponding triangular faces in the initial mesh. The encircled vertices in (a) are the degrees of freedom for the corresponding element. ... 75 5.1 A control mesh with an extraordinary vertex of degree 5 and the corre sponding limit surface : (a) using the concepts developed in Chapter 3, where the limit surface consists of quadrilateral normal elements and a pentagonal special element; (b) using the unified concept developed in this chapter, where the limit surface consists of one single type of quadrilateral finite element. ................... 86 5.2 In CatmullClark subdivision, each nonboundary quadrilateral face in the control mesh has a corresponding quadrilateral patch in the limit surface : (a) control mesh, (b) limit surface. ... 87 5.3 (a) The marked 16 control vertices define the shaded quadrilateral patch associated with the shaded regular face in the control mesh; (b) the marked 14 control vertices define the shaded quadrilateral patch associated with the shaded irregular face in the control mesh.. 91 5.4 (a) The control polygon with triangular faces; (b) mesh obtained after one subdivision step using Loop's subdivision rules. ... 96 5.5 An initial mesh and the corresponding limit surface obtained using Loop's subdivision rules. The domains of the shaded triangular patches in the limit surface are the corresponding triangular faces in the initial mesh. The encircled vertices are the control vertices for the corre sponding triangular patch in the limit surface. .... 98 5.6 Each triangular patch in the limit surface can be associated with a nonboundary triangular face in the initial mesh, which in turn can be parameterized over a triangle with vertices at (0, 0), (1, 0) and (0, 1). 99 6.1 Schematic block diagram of the multilevel dynamics approach. 104 6.2 Multilevel dynamics vs. multiresolution dynamics. ... 106 6.3 Representation of the degrees of freedom in multilevel dynamics and multiresolution dynamics approach. ..... 108 6.4 The filter bank. .............. 112 6.5 Subdivision refinement of a tetrahedron. .. ...... 119 6.6 Wavelet construction. ................ ...... ..123 7.1 System architecture. ......................... .. 134 8.1 (a), (b), (c) and (d) : Initial shapes (obtained applying CatmullClark subdivision rules on control meshes); (e), (f), (g) and (h) : the corre sponding modified shapes after interactive force application. 139 8.2 (a), (b), (c) and (d) : Initial shapes (obtained applying modified but terfly subdivision rules on control meshes); (e), (f), (g) and (h) : the corresponding modified shapes after interactive sculpting via force ap plication ................... .. .............. 141 8.3 (a) Range data of a bulb along with the initialized model, (b) an intermediate stage of evolution, and (c) the fitted dynamic Catmull Clark subdivision surface model. ..143 8.4 (a) Range data of an anvil along with the initialized model, (b) an intermediate stage of evolution, and (c) the fitted dynamic Catmull Clark subdivision surface model. ..... 144 8.5 (a) Range data of a head along with the initialized model, (d) the fitted dynamic butterfly subdivision model, and (c) visualization of the shape from another view point. . 144 8.6 (a) A slice from a brain MRI, (b) initialized model inside the region of interest superimposed on the slice, (c) the fitted model superim posed on the slice, and (d) a 3D view of the dynamic CatmullClark subdivision surface model fitted to the cerebellum. .. ... 146 8.7 (a) Data points identifying the boundary of the caudate nucleus on a MRI slice of human brain, (b) data points (from all slices) in 3D along with the initialized model, and (c) the fitted dynamic butterfly subdivision model. .......................... .. 146 8.8 Snapshots from the animation of canine heart motion over a cardiac cycle using the dynamic butterfly subdivision model. ... 149 8.9 (a) The initialized model along with the data set; (b) the fitted model with two subdivisions on the initial mesh, along with attached springs for editing. The model after editing (c) at lower resolution, (d) at the same resolution of the fitted model, and (e) at higher resolution. 150 8.10 (a) Discrete elevation data set (4096 points), (b) fitted dynamic but terfly subdivision surface model with 841 vertices (without noise ad dition), and (c) fitted dynamic subdivision surface model with 841 vertices when Gaussian noise is added. 153 8.11 Synthesized natural terrain of different roughness using the dynamic butterfly subdivision surface model with 841 vertices from a data set of 3099 elevation values. ........ 154 Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy A DYNAMIC FRAMEWORK FOR SUBDIVISION SURFACES By Chhandomay Mandal December 1998 Chairman: Dr. Baba C. Vemuri Cochairman : Dr. Hong Qin Major Department: Computer and Information Science and Engineering Subdivision surfaces are extensively used to model smooth shapes of arbitrary topology. Recursive subdivision on a userdefined initial control mesh generates a visually pleasing smooth surface in the limit. However, users have to carefully select the initial mesh and/or manipulate the control vertex positions at different levels of subdivision hierarchy to satisfy the functional and aesthetic requirements in the smooth limit surface. This modeling drawback results from the lack of direct manip ulation tools for the limit surface. In this dissertation, techniques from physicsbased modeling are integrated with geometric subdivision methodology, and a dynamic framework is presented for direct manipulation of the smooth limit surface generated by the subdivision schemes using physicsbased "force" tools. In a typical subdivision scheme, the user starts with an initial control mesh which is refined recursively using a fixed set of subdivision rules, and a smooth sur face is produced in the limit. Most often this limit surface does not have an analytic expression, and hence poses challenging problems in incorporating mass and damp ing distribution functions, internal deformation energy, forces, and other physical quantities required to develop a physicsbased subdivision surface model. In this dissertation, local parameterization techniques suitable for embedding the geometric subdivision surface model in a physicsbased modeling framework have been devel oped. Specific local parameterization techniques have been fully developed for the CatmullClark, modified butterfly and the Loop subdivision schemes. Techniques for assigning material properties to geometric subdivision surfaces are derived, and a motion equation for the dynamic model has been formulated using Lagrangian dynamics. Furthermore, advantages of the physicsbased deformable models are in corporated into the conventional subdivision schemes, and a dynamic hierarchical control of this model is introduced. Finally, a multiresolution representation of the control mesh is developed and a unified approach for deriving subdivision surface based finite elements is presented. The proposed dynamic framework enhances the applicability of the subdivision surfaces in modeling applications. It is also useful for hierarchical shape recovery from large range and volume data sets, as well as for nonrigid motion tracking from a temporal sequence of data sets. Multiresolution representation of the initial mesh controlling the smooth limit surface enables global and local editing of the shape as desired by the modeler. This dynamic framework has also been used for synthesizing natural terrain models from sparse elevation data. CHAPTER 1 INTRODUCTION Generating smooth surfaces of arbitrary topology poses a grand challenge for the computer graphics, computeraided geometric design and scientific visualization researchers. The existing techniques for modeling smooth complex shapes can be broadly classified into two distinct categories namely, (1) explicit patching using parametric surfaces and (2) subdivision surfaces. The explicit patching technique involves partitioning the model surface into a collection of individual parametric surface patches. Adjacent surface patches are then explicitly stitched together using continuity constraints. This explicit stitching process is very complicated for modeling smooth surfaces of arbitrary topology due to the continuity constraints which need to be satisfied along the patch boundaries. Subdivision surfaces are simple procedural models which offer an alternate rep resentation for the smooth surfaces of arbitrary topology. A typical recursive subdi vision scheme produces a visually pleasing smooth surface in the limit by repeated application of a fixed set of refinement rules on an userdefined initial control mesh. The initial control mesh is a simple polygonal mesh of the same topological type as that of the smooth surface to be modeled. At each step of subdivision, a finer polygonal mesh with more vertices and faces is constructed from the previous one via the refinement process, and the smooth surface is obtained in the limit. However, a few subdivision steps on the initial mesh generally suffice to approximate the smooth surface for all practical purposes. Various sets of rules lead to different subdivision schemes which mainly differ in the smoothness property of the resulting limit surface and/or the type of initial mesh (i.e. triangular, quadrilateral etc.) chosen. 1.1 Problem Statement To model a smooth surface of arbitrary topology using subdivision surfaces, first a polygonal mesh of same topology needs to be chosen as the initial mesh. This ini tial mesh, also known as the control mesh, is refined via recursive subdivision using a fixed set of rules, and a smooth surface of the desired topology is obtained in the limit. However, users have to carefully select the initial mesh and/or manipulate the control vertex positions at different levels of subdivision hierarchy in order to satisfy the functional and aesthetic requirements on the smooth limit surface. For example, to obtain a desired effect on the smooth limit surface, it might be necessary to reposition a handful of vertices in the mesh obtained after one subdivision step, or a large number of vertices might need to be moved in the mesh produced after three subdivision steps! This process is not intuitive and at best laborious. Despite the presence of a variety of subdivision schemes in the computer graphics and geo metric modeling literature, there is no direct and natural way of manipulating the limit surface. The current stateoftheart only permits the modeler to interactively obtain the desired effects on the smooth surface by kinematically manipulating the vertex positions at various levels of the subdivision hierarchy. In this dissertation, the challenging problem of directly manipulating the smooth limit surface at arbitrary locations/areas is addressed and a novel solution is presented where the modeler can not only manipulate the smooth limit surface directly but also control the extent of the manipulation effect (i.e., global or local manipulation) on the limit surface. Subdivision surfaces have also been used for recovering shapes from a given set of points in 3D. However, most of the existing subdivision surfacebased shape recovery techniques resort to complex algorithms to derive a mesh for the underlying shape, and then mesh optimization techniques are used to obtain a compact representation of the same. Nevertheless, this process yields a control mesh of the smooth subdi vision surface representing the underlying shape that typically uses a large number of degrees of freedom (control vertices) for representation. In this dissertation, an efficient hierarchical method of recovering shapes from range and volume data sets is proposed where the control mesh of the smooth limit surface will use very few degrees of freedom for representation. 1.2 Proposed Solution Physicsbased modeling techniques offer a potential solution to the problem of directly manipulating the smooth limit surface generated by the recursive sub division procedure. In the physicsbased modeling paradigm, a deformable model is derived by assigning mass, damping, stiffness and other physical properties to a surface model. The model is deformed by applying synthesized forces and this de formation is governed by physical laws. Now, if the purely geometric subdivision surface models can be embedded in a physicsbased modeling paradigm, then the smooth limit surface can be deformed directly by using physicsbased "force" tools. However, this procedurebased surface model obtained through the subdivision pro cess does not have a closedform analytic formulation in general, and hence poses challenging problems to incorporate mass and damping distribution functions, inter nal deformation energy, forces, and other physical quantities required to develop a physicsbased model. Techniques of locally parameterizing the smooth limit surface generated by various subdivision schemes are proposed in this dissertation. Once a suitable local parameterization scheme is developed for a specific subdivision scheme, a dynamic framework is provided for the corresponding subdivision scheme where the modeler can directly manipulate the smooth limit surface by using synthesized forces. At the same time, a dynamic framework for the wavelets derived using subdi vision schemes will assist in adopting a multiresolution representation of the control mesh defining the smooth limit surface, and the modeler can control the extent of manipulation effect by choosing a desired level of editing. Synthesized force applica tion at a lower resolution will yield a global effect whereas manipulations at a finer resolution will have localized effects on the limit surface. The motion of this physics based deformable subdivision surface model is governed by a secondorder differential equation, which is solved numerically using the finite element method. New types of finite elements for the chosen subdivision scheme are also presented for representing the smooth limit surface. The proposed dynamic framework for the subdivision surfaces provides an effi cient solution to the shape recovery problem as well. A simple subdivision surface model with very few vertices in the control mesh can be initialized (positioned) fully inside a given set of points in 3D. The initialized model will be deformed by applying forces synthesized from the given data points. When an equilibrium is obtained, the number of vertices in the control mesh can be increased via a subdivision step on the current control mesh thereby increasing the degrees of freedom for model representa tion, and a new equilibrium with a better fit to the given data set can be obtained. This process can be repeated till a prescribed error in fit is achieved. Similar approach can be taken for shape recovery from volume data sets as well where a different type of synthesized force needs to be specified. The hierarchical shape recovery process ensures a compact representation of the recovered smooth limit surface using very few degrees of freedom. 1.3 Contributions In this dissertation, techniques from physicsbased modeling are integrated with geometric subdivision methodology to present a scheme for directly manipulating the smooth limit surface generated by the subdivision process. As a result, unlike the existing geometric solutions that only allow the operations on control vertices, the proposed methodology and algorithms permit the user to physically modify the shape of subdivision surfaces at desired locations via application of forces. This gives the user a "virtual" clay/playdough modeling environment. The proposed model can be edited directly in a hierarchical fashion using synthesized forces. Also, this physicsbased subdivision surface model efficiently recovers shapes as well as non rigid motions from large range and volumetric data sets. Note that this dissertation neither proposes a new subdivision technique nor provides a different interpretation of any existing subdivision technique, but integrates the advantages of subdivision surfacebased and physicsbased modeling techniques to solve important theoretical and practical problems. Although the principles of physicsbased modeling are well understood by the graphics experts and modeling researchers, this dissertation will greatly advance the state of the art in physicsbased shape modeling due to the contributions listed below. Local parameterization techniques for the smooth limit surface generated by various subdivision schemes are systematically derived in a hierarchical frame work, and subsequently the initial control polyhedron can be viewed as the parametric domain of the physicsbased smooth limit surface. The smooth dynamic subdivision surface in the limit is treated as a "real" physical object represented by a set of novel finite elements. The basis (shape) functions of these new variety of finite elements are derived using the subdivision schemes. The proposed finite element methods will prove to be useful not only in the areas of computer graphics and geometric design, but also in engineering analysis as well. The subdivision techniques are used to create a surface model that incorporates mass and damping distribution functions, internal deformation energy, forces, and other physical quantities. The motion equations are also systematically derived for this dynamic subdivision surface model whose degrees of freedom are the collection of initial userspecified control vertices. Therefore, the advantages of both the physicsbased modeling philosophy and the geometric subdivision schemes are incorporated within a single unified framework. Users will be able to manipulate this physicsbased model in an arbitrary region, and the model will respond naturally (just like a realworld object would) to this force application. This shape deformation is quantitatively characterized by the timevarying displacement of control points that uniquely define the geometry of the limit surface. The dynamic framework for wavelets derived using subdivision schemes en ables a multiresolution representation of the control mesh defining the smooth limit surface. This provides additional flexibility to the physicsbased modeling framework since the modeler is free to choose the desired editing level. If a lower resolution editing level is chosen, the synthesized force application on the smooth limit surface will have global effects, whereas editing with forcebased tools at a finer resolution will yield a localized effect. Algorithms and procedures are developed which approximate the proposed new finite elements using a collection of linear and/or bilinear finite elements subject to the implicit geometric constraints enforced by the subdivision rules. This hierarchicallystructured approximation can satisfy any userspecified error tol erance. The proposed dynamic framework enhances the applicability of subdivision sur face models in various application areas. It provides a direct and intuitive way of manipulating shapes in a hierarchical fashion for geometric modeling applications. It has also been successfully used for efficient and hierarchical shape recovery from range and volume data sets as well as for tracking a shape of interest from a time sequence of range or volume data sets. Finally, the dynamic framework of subdivision surfaces is combined with a variant of a fractal surface synthesis technique to present a novel natural terrain modeling method. 1.4 Outline of Dissertation The rest of the dissertation is organized as follows : Chapter 2 contains an overview of the subdivision surfaces along with a review of the related literature. The motivation for embedding the subdivision surface models in a physicsbased modeling framework is discussed, and the proposed model is com pared with the existing physicsbased models. The advantages of a unified framework for shape modeling and shape recovery are also pointed out in this chapter. Chapter 3 provides a dynamic framework for the CatmullClark subdivision scheme, which is one of the most popular subdivision techniques for modeling compli cated objects of arbitrary genus. Analytic formulation of the limit surface generated by the CatmullClark subdivision scheme is derived and the "physical" quantities required to develop the dynamic model are introduced [57, 60, 76]. The governing dynamic differential equation is derived using Lagrangian mechanics and is imple mented using a finite element technique. In Chapter 4, a dynamic framework is provided for the butterfly subdivision scheme, another popular subdivision technique for modeling smooth surfaces of arbi trary topology. A local parameterization scheme for the butterfly subdivision surface is derived in a hierarchical style. The physicsbased butterfly subdivision model is formulated as a set of novel finite elements which are optimally approximated by a collection of standard finite elements subjected to implicit geometric constraints [58, 61]. Chapter 5 presents an unified approach for providing a dynamic framework for subdivision surfaces in general. In particular, it has been shown that the limit surface obtained using any subdivision scheme can be viewed as a collection of either quadri lateral or triangular finite elements whose basis (shape) functions can be derived using the chosen subdivision scheme [59]. In Chapter 6, a dynamic framework is presented for the subdivision surface based wavelet schemes. This dynamic framework enables a multiresolution represen tation of the control mesh defining the smooth limit surface and consequently, the modeler can control the extent of the effect of manipulation on the limit surface by choosing a proper editing level. In Chapter 7, a system that integrates the implementation of the concepts pro posed in earlier chapters is presented. Several modules comprise the entire system, and the functionality of these modules are discussed in this chapter. 10 The proposed modeling framework is used in various application areas and the results are presented in Chapter 8. The proposed dynamic model has been success fully used in geometric modeling, shape recovery and nonrigid motion estimation applications. A novel technique for natural terrain modeling is also presented where the dynamic framework is combined with a variant of a fractal surface synthesis technique. Finally, conclusions are drawn in Chapter 9 where future directions of research are also pointed out. CHAPTER 2 BACKGROUND In this chapter, a brief overview of the subdivision surfaces is presented followed by a review of the previous work done in the area of subdivision surfaces. This is followed by an overview of the physicsbased models along with a review of the re lated literature. The motivating factors for embedding geometric subdivision surface models in a physicsbased modeling paradigm are presented in the next section. The advantages of the proposed dynamic models over the existing physicsbased models for shape recovery are discussed and finally advantages of a unified framework for shape modeling and shape recovery are also presented. 2.1 Subdivision Surfaces The input to any subdivision scheme is an initial mesh (also known as control mesh), MO = (Vo, Fo), which is a collection of vertices V and a collection of faces Fo. The subdivision surface S(MO) associated with the initial mesh MO is defined as the limit of the recursive application of the refinement R as shown in Fig.2.1. The refinement R, when applied on a mesh Mk = (VI, Fk), creates a refined mesh Mk+ l (Vk+l, Fk+1) where the vertices in Vk+' are computed as affine combinations of the vertices in Vk and the faces in Fk+1 are obtained by splitting each face in Fk into a fixed number of subfaces. (b) Meh M = R(M) (c) Mesh M = R(R(M)) (d) Limit snlrfax S(M) = R(M) Figure 2.1. Refinement of an initial control mesh to obtain the limit surface. (Cour tesy : H. Hoppe) It may be noted that the vertices introduced through the subdivision process have a fixed degree in general (4 in case of quadrilateral meshes and 6 in case of tri angular meshes). The vertices which do not have this degree are called extraordinary vertices. The number of extraordinary vertices are very few as these vertices are not introduced via subdivision refinement in general. For example, in the CatmullClark subdivision scheme (defined on quadrilateral meshes) the number of extraordinary vertices does not change after the first subdivision on the initial mesh, whereas in the modified butterfly subdivision scheme (defined on triangular meshes) the subdivision process never introduces an extraordinary vertex. The limit surface defined by the subdivision process is most often C (first derivative continuous) or C' (second deriva tive continuous) depending on the subdivision rules expect at very few extraordinary points corresponding to the extraordinary vertices in the mesh. Specific subdivision schemes along with the properties of their limit surfaces will be discussed in later chapters. A lot of research have been carried out on subdivision surfaces, which are mainly focused either on the development of a new subdivision technique or on analyzing the smoothness properties of the limit surface generated by a subdivision scheme. However, in this dissertation the focus is entirely different, namely, embedding the subdivision surfaces in a physicsbased modeling paradigm to enhance the applicabil ity of these surface models. In the rest of this section, these "traditional" techniques of subdivision surfaces are reviewed. Chaikin [14] first introduced the concept of subdivision to the computer graphics community for generating a smooth curve from a given control polygon. During the last two decades, a wide variety of subdivision schemes for modeling smooth surfaces of arbitrary topology have been derived following Chaikin's pioneering work on curve generation. In general, these subdivision schemes can be categorized into two distinct classes namely, (1) approximating subdivision techniques and (2) interpolatory sub division techniques. The limit surface obtained using an approximating subdivision technique only approximates the initial mesh, whereas the limit surface interpolates the initial mesh in case of interpolatory subdivision techniques. Among the approximating schemes, the techniques of Doo and Sabin [27, 83] generalize the idea of obtaining uniform biquadratic patches from a rectangular con trol mesh. This scheme is an example of vertex subdivision scheme where a vertex surrounded by k faces is split into k subvertices, one for each face. Doo and Sabin's subdivision technique can be applied to any mesh of arbitrary topology, and the re sulting smooth limit surface would be biquadratic Bsplines, except at extraordinary points corresponding to the extraordinary vertices in the control mesh. Catmull and Clark [10] developed a method for recursively generating a smooth surface from a polyhedral mesh of arbitrary topology. The CatmullClark subdivision surface, defined by an arbitrary initial mesh, can be reduced to a set of standard bicubic Bspline patches except at a finite number of degenerate points. This scheme is an example of face subdivision scheme, where a ksided face is split into k sub faces. The details of this subdivision scheme are discussed in Section 3.1. Halstead et al. [38] proposed an algorithm to construct a CatmullClark subdivision surface that interpolates the vertices of a mesh of arbitrary topology. Loop [51] presented a subdivision scheme based on the generalization of quartic triangular Bsplines for triangular meshes. Loop's subdivision rules are presented in Section 5.3. Hoppe et al. [40] extended Loop's work to produce piecewise smooth surfaces with discontinuities in selected areas of the limit surface. They introduced new subdivision rules that allow for sharp features such as creases and corners in the limit surface. Peters and Reif [73] proposed a simple subdivision scheme for smoothing poly hedra. Their refinement rules yield a C1 surface that has a piecewise quadratic pa rameterization except at isolated extraordinary points. Most recently, nonuniform DooSabin and CatmullClark surfaces that generalize nonuniform tensor product Bspline surfaces to arbitrary topologies were introduced by Sederberg et al. [88]. All the schemes mentioned above generalize recursive subdivision schemes for generating limit surfaces with a known parameterization. Various issues involved with character animation using these approximating subdivision schemes were discussed at length by DeRose et al. [25]. The most wellknown interpolationbased subdivision scheme is the "butterfly" algorithm proposed by Dyn et al. [30]. Butterfly subdivision method, like other sub division schemes, makes use of a small number of neighboring vertices for subdivision. It requires simple data structures and is extremely easy to implement. However, it needs a topologically regular setting of the initial (control) mesh in order to obtain a smooth C' limit surface. A variant of this scheme with better smoothness prop erties can be found in Dyn et al. [29]. Zorin et al. [109] developed an improved interpolatory subdivision scheme (which we call the modified butterfly scheme) that retains the simplicity of the butterfly scheme and results in much smoother surfaces even from irregular initial meshes. These interpolatory subdivision schemes have extensive applications in wavelets on manifolds, multiresolution decomposition of polyhedral surfaces and multiresolution editing. A variational approach for interpolatory refinement has been proposed by Kobbelt [43, 44] and by Kobbelt and Schrider [46]. In this approach, the vertex positions in the refined mesh at each subdivision step are obtained by solving an optimization problem. Therefore, these schemes are global, i.e., every new vertex position depends on all the vertex positions of the coarser level mesh. The local refinement property which makes the subdivision schemes attractive for implementation in the computer graphics applications is not retained in the variational approach. The derivation of various mathematical properties of the smooth limit surface generated by the subdivision algorithms is rather complex. Doo and Sabin [28] first analyzed the smoothness behavior of the limit surface using the Fourier transform and an eigenanalysis of the subdivision matrix. Ball and Storry [3, 4] and Reif [80] further extended Doo and Sabin's prior work on continuity properties of subdivision surfaces by deriving various necessary and sufficient conditions on smoothness for different subdivision schemes. Specific subdivision schemes were analyzed by Schweitzer [87], Habib and Warren [37], Peters and Reif [74] and Zorin [108]. Most recently, Stam [90] presented a method for exact evaluation of CatmullClark subdivision surfaces at arbitrary parameter values. It may be noted that the focus of this dissertation is not on deriving a new subdivision technique or analyzing a subdivision technique, but on deriving a deformable surface model using these subdivision schemes. 2.2 Physicsbased Deformable Surface Models Shape models can be broadly categorized into two types lumped parameter models and distributed parameter models. A lumped parameter model uses small number of parameters to represent model geometry, whereas a distributed parameter model uses large number of degrees of freedom for model representation. An example of a lumped parameter surface model is a superquadric. Spline surfaces serve as a good example for distributed parameter model. Deformable surfaces are distributed parameter models which use large number of degrees of freedom for representing the model geometry. A large variety of shapes can be modeled using this type of model, but handling a large number of degrees of freedom can be cumbersome. However, the large degrees of freedom of a deformable model are embedded in a physicsbased framework to allow only a "physically mean ingful" model behavior. Various types of energies are assigned to the model using the degrees of freedom, and the model "deforms" to find an equilibrium state with min imal energy. The motion equation is derived using Lagrangian dynamics [36], where various energies associated with model gives rise to internal and external forces. The equilibrium state of the model is a model position where the internal deformation force becomes equal to the externally applied force. These physicsbased deformable models are useful for modeling where the modeler can deform a surface by applying synthesized forces, and for data fitting where external forces are synthesized from a given data set such that the model approximates the given data at equilibrium. A deformable surface model is typically parameterized over the domain [0, 1]2. Let s(u, v, p) be a deformable surface model, where 0 < u, v < 1 and p is a collection of the degrees of freedom p,, i = 1, 2,..., n associated with the model (assuming the model has n degrees of freedom). The degrees of freedom pi(t) is a set of generalized coordinates which are functions of time and are assembled into the vector p. Let fi(t) be the generalized force represented by the vector fp and acting on p,. Let T be the kinetic energy, F be the dissipation energy and U be the potential energy of the deformable surface model. Then, the Lagrangian equation of motion for the model can be expressed as d T 8T OF 9U T + + f (2.1) dt( p, op, 8p, op, Let p(u, v) be the mass density function of the surface. Then the kinetic energy of the surface is T = J s sdudv = ipMp, (2.2) where M is called the mass matrix. Similarly, let y(u, v) be the damping density function of the surface. Then the dissipation energy is F = Tsi f dudv = TDp, (2.3) 2I 2 where D is called the damping matrix. The potential energy of the model can be defined using the thinplateundertension energy model [96], and is given by T f OT Ns s 9s'T s 82sT 2's U = ( + a22 + i1u1 82sT 02s _2sT 2s 1 +12 a + 22 )dudv pTKp, (2.4) u9uv Yaby a2v 82v 2 where aii(u, v) and Aj(u, v) are elasticity functions which control tension and rigidity respectively of the deformable surface model, and K is known as the stiffness matrix. The discretized equation of motion can be derived using the expressions of the kinetic, dissipation and potential energy in Eqn.2.1, which is given by Mp + Dp + Kp = fp, (2.5) where fp is the generalized force vector. The freeform deformable models discussed above were first introduced to com puter graphics and visualization in Terzopoulos et al. [99] and further developed by Terzopoulos and Fleischer [97], Pentland and Williams [72], Metaxas and Terzopou los [65] and Vemuri and Radisavljevic [104]. Celniker and Gossard [11] developed a system for interactive freeform design based on the finite element optimization of en ergy functionals proposed in Terzopoulos and Fleischer [97]. Bloor and Wilson [8, 9], Celniker and Welch [12] and Welch and Witkin [105] proposed deformable Bspline curves and surfaces which can be designed by imposing the shape criteria via the minimization of the energy functionals subject to hard or soft geometric constraints through Lagrange multipliers or penalty methods. Qin and Terzopoulos [77, 78, 100] developed dynamic NURBS (DNURBS) which are very sophisticated models suit able for representing a wide variety of freeform as well as standard analytic shapes. The DNURBS have the advantage of interactive and direct manipulation of NURBS curves and surfaces, resulting in physically meaningful hence intuitively predictable motion and shape variation. Deformable models are also widely used for shape recovery, segmentation, mo tion tracking and other computer vision and medical imaging applications. A detailed survey of deformable models used in these techniques can be found in McInerney and Terzopoulos [64] and the references therein. Some specific existing deformable models used for shape recovery and nonrigid motion estimation will be reviewed in Section 2.4. 2.3 Shape Modeling Using Physicsbased Subdivisionsurface Model Recursive subdivision surfaces are powerful for representing smooth geometric shapes of arbitrary topology. However, they constitute a purely geometric represen tation, and furthermore, conventional geometric modeling with subdivision surfaces may be difficult for representing highly complicated objects. For example, modelers are faced with the tedium of indirect shape modification and refinement through time consuming operations on a large number of (most often irregular) control vertices when using typical subdivision surfacebased modeling schemes. Despite the advent of advanced 3D graphics interaction tools, these indirect geometric operations remain nonintuitive and laborious in general. In addition, it may not be enough to obtain the most "fair" surface that interpolates a set of (ordered or unorganized) data points. A certain number of local features such as bulges or inflections may be strongly de sired while requiring geometric objects to satisfy global smoothness constraints in geometric modeling and computer graphics applications. In contrast, physicsbased modeling provides a superior approach to shape modeling that can overcome most of the limitations associated with traditional geometric modeling approaches. Freeform deformable models governed by the laws of continuum mechanics are of particular interest in this context. These dynamic models respond to externally applied forces in a very intuitive manner. The dynamic formulation marries the model geometry with time, mass, damping and constraints via a force balance equation. Dynamic models produce smooth, natural motions which are easy to control. In addition, they facilitate interaction especially direct manipulation of complex geometries. Fur thermore, the equilibrium state of the model is characterized by a minimum of the deformation energy of the model subject to the imposed constraints. The deformation energy functionals can be formulated to satisfy local and global modeling criteria, and geometric constraints relevant to shape design can also be imposed. The dynamic approach subsumes all of the aforementioned modeling capabilities in a formulation which grounds everything in realworld physical behavior. A severe limitation of the existing deformable models, including DNURBS, is that they are defined on a rectangular parametric domain. Hence, it can be very difficult to model surfaces of arbitrary genus using these models. In a recent work, DeRose et al. [25] assigned material properties to control meshes at various subdivi sion levels in order to simulate cloth dynamics using subdivision surfaces. Note that, they assign physical properties on the control meshes at various levels of subdivision and not on the limit surface itself, and hence can not solve the modeling goal we are trying to achieve. In this dissertation, a dynamic framework is presented for subdivision surfaces which combines the benefits of subdivision surfaces for modeling arbitrary topology as well as that of dynamic splines for interactive shape manip ulation by applying synthesized forces. The proposed dynamic framework presents the modeler a formal mechanism of direct and intuitive manipulation of the smooth limit surface, as if they were seamlessly sculpting a piece of realworld "clay." A dynamic framework for subdivision surfacebased wavelets adds flexibility in the pro posed modeling paradigm. It enables a multiresolution representation of the evolving control mesh defining the smooth limit surface, and the modeler can control the extent of manipulation effect by choosing a proper editing level. The formulation and implementation details of these dynamic frameworks are discussed in subsequent chapters. 2.4 Shape and Motion Estimation Using Physicsbased Subdivisionsurface Model The dynamic subdivision surface model has been developed primarily for mod eling arbitrary (known) topology where modelers can directly manipulate the limit surface by applying synthesized forces in an intuitive fashion. However, another im portant application of the dynamic subdivision surfaces is in nonrigid shape and motion reconstruction/recovery. Accurate shape recovery requires distributed pa rameter models which typically possess a large number of degrees of freedom. On the other hand, efficient shape representation imposes the requirement of geometry compression, i.e., models with fewer degrees of freedom. These requirements are conflicting and numerous researchers have been seeking to strike a balance between these contradicting requirements [5, 21, 47, 49, 62, 89, 100, 104]. Another important criterion in the design of shape models is that the initialization of the model during the shape recovery process should not be restricted to globally parameterized input meshes since it may be infeasible to globally parameterize shapes of arbitrary topol ogy. A physicsbased model best satisfying the above mentioned criteria is an ideal candidate for a solution to the shape recovery problem for obvious reasons. Deformable models which come in many varieties, have been used to solve this problem in the physicsbased modeling paradigm. These models involve the use of either fixed size [21, 66, 71, 98, 104] or adaptive size [15, 41, 47, 62, 95, 101] grids. The models with fixed grid size generally use less number of degrees of freedom for representation at the cost of accuracy of the recovered shape. On the other hand, models using adaptive grids generally need large number of degrees of freedom to recover the shapes accurately. Note that the shapes being recovered from the image data are smooth in most of the medical applications, i.e. the anatomical structures even with considerable amount of details have more or less a C1 surface. Under these circumstances, the finite element approaches as in Cohen and Cohen [21] and McInerney and Terzopoulos [62] need a large number of degrees of freedom for deriving a smooth and accurate representation. In addition, they can not represent shapes with known arbitrary topology. Moreover, almost all of these schemes require a globally parameterized mesh as their input which may be infeasible at times. The proposed model solves the shape recovery problem very effectively as it can recover shapes from large range and volume data sets using very few degrees of freedom (control vertices) for its representation and can cope with any arbitrary input mesh, not necessarily parameterized. The initialized model deforms under the influence of synthesized forces to fit the data set by minimizing its energy. Once the approximate shape is recovered, the model is further subdivided automatically and a better approximation to the input data set is achieved using more degrees of freedom. The process of subdivision after achieving an approximate fit is continued till a pre scribed error criteria for fitting the data points is achieved. The proposed dynamic subdivision surface models have also been successfully used in motion tracking and visualization of a dynamically deforming shape from a time sequence of volumetric data sets. 2.5 Unified Framework for Shape Recovery and Shape Modeling Currently, shape recovery and shape modeling are viewed as two distinct areas in computer graphics and geometric modeling literature. However, there are potential benefits if these two can be combined in an unified framework. For example, the modeler starts from scratch to build a specific model in a typical geometric modeling scenario. First, a rough shape is modeled and then it is finetuned by manipulating control vertex positions to obtain the desired effects. This turns out a cumbersome process in general. On the other hand, shape recovery using stateoftheart methods yield large polygonal meshes which are very difficult to manipulate, especially for global changes in shape. Most often the resulting meshes from a shape recovery application are not directly amenable to multiresolution analysis. Computationally expensive remeshing techniques are needed to convert these meshes into a specific type of meshes on which multiresolution analysis can be performed. This specific type of mesh is known as mesh with "subdivisionconnectivity", implying a topologically equivalent mesh with the same connectivity as of the given mesh can be obtained by recursive subdivision of a very simple known initial mesh. The proposed dynamic framework combines shape recovery and shape modeling in an unified framework where the modeler can scan in 3D points of a prototype model, recover the shape using the proposed dynamic subdivision surface model, and edit at any desired resolution using physicsbased force tools. Thus, the modeler is re lieved of the burden of both building an initial model and editing a cumbersome huge polygonal mesh. The shape recovery process starts with a smooth subdivision surface model which has a simple initial mesh. This model is deformed using synthesized forces from the given data points and is automatically refined using some predefined error in fit criteria. The initial mesh of the final recovered smooth shape has subdivi sion connectivity inbuilt, as it is obtained by subdivision refinement, and therefore no remeshing is needed for multiresolution analysis. The dynamic framework for subdivision surfacebased wavelets makes multiresolution editing using physicsbased 26 force tools easy to perform. These advantages of shape modeling and shape recovery in a unified framework will be further illustrated in later chapters. CHAPTER 3 DYNAMIC CATMULLCLARK SUBDIVISION SURFACES Subdivision surfaces, as mentioned earlier, can be broadly classified into two categories approximating and interpolatory subdivision schemes. The approxi mating subdivision schemes typically generalize recursive subdivision techniques for generating limit surfaces with known parameterizations. In this chapter, a dynamic framework will be presented for CatmullClark subdivision scheme, a popular approx imating subdivision technique. The CatmullClark subdivision scheme generalizes the idea of obtaining a bicubic Bspline surface by recursive subdivision of a rectangu lar initial mesh. Before discussing the dynamic framework, first the CatmullClark subdivision scheme is briefly reviewed. 3.1 Overview of the CatmullClark Subdivision Scheme CatmullClark subdivision scheme, like any other subdivision scheme, starts with a userdefined mesh of arbitrary topology. It refines the initial mesh by adding new vertices, edges and faces with each step of subdivision following a fixed set of subdivision rules. In the limit, a recursively refined polyhedral mesh will converge to a smooth surface. The subdivision rules are as follows: (b) (c) Figure 3.1. CatmullClark subdivision : (a) initial mesh, (b) mesh obtained after one step of CatmullClark subdivision, and (c) mesh obtained after another subdivision step. For each face, a new face vertex is introduced which is the average of all the old vertices defining the face. For each (nonboundary) edge, a new edge vertex is introduced which is the average of the following four points: two old vertices defining the edge and two new face vertices of the faces adjacent to the edge. For each (nonboundary) vertex V, a new vertex is introduced whose position is E + 2E + n3, where F is the average of the new face vertices of all faces adjacent to the old vertex V, E is the average of the midpoints of all the edges incident on the old vertex V and n is the number of the edges incident on the vertex. New edges are formed by connecting each new face vertex to the new edge vertices of the edges defining the old face and by connecting each new vertex corresponding to an old vertex to the new edge vertices of all old edges incident on the old vertex. New faces are defined as faces enclosed by new edges. An example of CatmullClark subdivision on an initial mesh is shown in Fig.3.1. The most important property of the CatmullClark subdivision surfaces is that a smooth surface can be generated from any control mesh of arbitrary topology. There fore, this subdivision scheme is extremely valuable for modeling various complicated geometric objects of arbitrary topology. CatmullClark subdivision surfaces include standard bicubic Bspline surfaces as their special case (i.e., the limit surface is a bicubic Bspline surface for a rectangular mesh with all nonboundary vertices of degree 4). In addition, the aforementioned subdivision rules generalize the recur sive bicubic Bspline patch subdivision algorithm. For nonrectangular meshes, the limit surface converges to a bicubic Bspline surface except at a finite number of ex traordinary points. These extraordinary points correspond to extraordinary vertices (vertices whose degree is not equal to 4) in the mesh. Note that, after the first sub division, all faces are quadrilaterals, hence all the new vertices created subsequently will have four incident edges. The number of extraordinary points on the limit sur face is a constant, and is equal to the number of extraordinary vertices in the refined mesh obtained after applying one step of the CatmullClark subdivision on the initial mesh. The limit surface is curvaturecontinuous everywhere except at extraordinary vertices, where only tangent plane continuity is achieved. In spite of the popularity of CatmullClark subdivision surfaces for representing complex geometric shapes of arbitrary topology, these subdivision surfaces may not be easily parameterizable and deriving a closedform analytic formulation of the limit surface can be very difficult. These deficiencies preclude their immediate pointwise manipulation and hence may restrain the applicability of these schemes. In the rest of the chapter, a dynamic framework is developed for the CatmullClark subdivision surfaces which offers a closedform analytic formulation and allows users to manipulate the model directly and intuitively. It may be noted that recently a scheme was proposed by Stam [90] for directly evaluating the CatmullClark subdivision surfaces at arbitrary parameter values, but the work presented here on dynamic CatmullClark subdivision surfaces [57, 60, 76] was published prior to the publication of the abovementioned work. 3.2 Formulation In this section, a systematic formulation of the new dynamic model based on CatmullClark subdivision surfaces is presented. This dynamic model treats the smooth limit surface as a function of its initial mesh. However, the control vertex po sitions need to be updated continually at every level of subdivision in order to develop the dynamic framework. Note that, all the vertices introduced through subdivision are obtained as an affine combination of control vertex positions in the initial mesh. Therefore, the dynamic behavior of the limit surface can be controlled by formulat ing the dynamic model on the initial mesh itself, the only exception being the case when the initial mesh has nonrectangular faces. This problem can be circumvented by taking the mesh obtained through one step of subdivision as the initial mesh. It may be noted that an additional subdivision step may be needed in some cases to isolate the extraordinary points, and the resulting mesh is treated as the initial mesh. A typical example of the above mentioned scenario is when the initial mesh is a tetrahedron where two subdivision steps are needed, and the dynamic framework can be developed by treating the mesh obtained after two subdivision steps as the initial mesh. To define the limit surface using the vertices of the initial mesh, the enumeration of the bicubic patches in the limit surface is necessary. In the next two sections, various schemes are presented for assigning the bicubic patches to various faces of the initial mesh. 3.2.1 Assigning Patches to Regular Faces In Fig.3.2, a rectangular control mesh is shown along with the bicubic Bspline surface (4 patches) in the limit after an infinite number of subdivision steps. Note that, each of the bicubic patches in the limit surface is defined by a rectangular face with each vertex of degree four, thereby accounting for 16 control points (from its 8 connected neighborhood) needed to define a bicubic surface patch in the limit. Therefore, for each rectangular face in the initial mesh with a degree of 4 at each vertex, the corresponding bicubic surface patch can be assigned to it in a straight forward way. In Fig.3.2, the surface patches S1,S2, S3 and S4 are assigned to face F1, F2, F3 and F4 respectively. The 16 control points for the patch S1, corresponding to face F1, are highlighted in Fig.3.2. Note that, the initial control mesh can be viewed as the parametric domain of the limit surface. Therefore, face Fi can be thought of as the portion of the parametric domain over which the patch Si is defined, i.e., has nonzero values. Nevertheless, each rectangular face (e.g. Fi) can be parametrically 32 UJ J   [  " ~ _ Figure 3.2. A rectangular mesh and its limit surface consisting of 4 bicubic surface patches. defined over [0, 1]2, and hence, all bicubic Bspline patches defined by 16 control points are locally parameterized over [0, 1]2. 3.2.2 Assigning Patches to Irregular Faces In Fig.3.3, a mesh containing an extraordinary point of degree 3 and its limit sur face are shown. The faces Fo, FI, ..., F are assigned to bicubic patches So, Si,..., Ss respectively (as they all have vertices of degree 4) following the aforementioned scheme. However, the central smooth surface enclosed by the patches So, S1,..., Ss consists of infinite number of bicubic patches converging to a point in the limit. A re cursive way of enumerating these bicubic patches and assigning them to various faces at different levels need to be developed in order to formulate the dynamic framework for CatmullClark subdivision surface model. The idea of enumerating the bicubic patches corresponding to faces having an extraordinary vertex is shown in Fig.3.4 where a local subdivision of the mesh con sisting of faces Fo, F,..., F. PF P P2 (and not the other boundary faces) of Fig.3.3 is carried out. Topologically, the resulting local subdivision mesh (shown as dotted mesh) is exactly the same as the mesh in Fig.3.3 and hence exactly the same number of bicubic patches can be assigned to its faces with vertices of degree 4 as is evident from Fig.3.4 (the new faces and the corresponding patches are marked by "p" and "n" respectively). This process of local subdivision and assignment of bicubic patches around an extraordinary point can be carried out recursively and in the limit, the enclosed patch corresponding to faces sharing the extraordinary point will converge to a point. However, there is no need to carry out an infinite number of subdivision 34 S 4 So S S S t   pi F 3n Figure 3.3. A mesh with an extraordinary point of degree 3 and its limit surface. selected mesh for S local subdivision   \7 Figure 3.4. Local subdivision around the extraordinary point and the limit surface. steps. This description is for formulation purposes only and the exact implementation will be detailed in a later section. 3.2.3 Kinematics In this section, the mathematics for the kinematics of the limit surface is devel oped via illustrative examples and then the generalized formulas are presented. It may be noted that a single bicubic Bspline patch is obtained as the limiting process of the CatmullClark subdivision algorithm applied to an initial 4 by 4 rectangular control mesh. Let sp(u, v), where (u, v) E [0, 1]2, denote this bicubic Bspline patch which can be expressed analytically as 3 3 sp(u, v) = (x(u, v), y(u, v), z(u, v)) = E E di,j B,4(u)Bj,4(v), (3.1) i=0 j=0 where dj represents a 3dimensional position vector at the (i,j)th control point location and Bi,4(u) and B,4 (v) are the cubic Bspline basis functions. The subscript p on s denotes the patch under consideration. Expressing Eqn.3.1 in a generalized coordinate system we have sp = Jq,, (3.2) where Jp is a transformation matrix storing the basis functions of a bicubic Bspline patch, and is of size (3, 48). Vector q, is the concatenation of all control points defining a Bspline patch in 3D. Note that in the concatenation of the control points, each control point has an (x, y, z) component. For example, the (x, y, z) components of the control point (i,j) correspond to positions 3k,3k + 1,3k + 2 where, k = 4i + j respectively in the vector qp. We can express the entries of Jp explicitly in the following way: Jp(0, k) = Jp(1, k + 1) = Jp(2, k + 2) = Bi,4(u)Bj,4(v) and J,(0, k + 1) = J,(0, k + 2) = Jp(1, k) = Jp(1, k + 2) = Jp(2, k) = Jp(2, k + 1) = 0. Limit surface with many bicubic patches from a rectangular initial mesh Now a limit surface consisting of many bicubic surface patches obtained after applying an infinite number of subdivision steps to a rectangular initial mesh is considered. For example, let the limit surface of Fig.3.2 be sm which can be written as 1 1 1 1 sm(u, v) = sm,(2u, 2v)+s,.(2(u ),2v)+s, (2(u),2(v2))+sa(2u,2(v )), (3.3) where s, (2u, 2v) sm(u, v) for 0 < u, v < and 0 otherwise. Similarly, s2,,Sm3 and s., are also equal to sm(u, v) for an appropriate range of values of u, v and 0 outside. It may be noted that SMn,, ,2,, sm,, sn correspond to patches S1, S, S3, S4 respectively in Fig.3.2. Eqn.3.3 in generalized coordinates becomes 4 Sm = Jlq + J2q2 + J3q3 + J44 = Z Jjqi, (3.4) i=1 where Jis are the transformation matrices of size (3, 48) and qis are the (x,y,z) com ponent concatenation of a subset of the control points of sm defining sm,, i = 1, 2,3 and 4. A more general expression for Sm is 4 sm = JAq + J2A2qm + JaA3qm + J4A4qm = JiAfqm = Jmqm. (3.5) i=l Where, q, is the 75component vector of 3D positions of the 25 vertex control mesh defining the limit surface si. Matrices Ai, 1 < i < 4, are of size (48, 75) each row consisting of a single nonzero entry (= 1) and the (3,75)sized matrix Jm = E=1 JA,. Limit surface with many bicubic patches from an arbitrary initial mesh The stage is now set to define the limit surface s using the vertices of initial mesh M for any arbitrary topology, assuming all faces are quadrilateral and no face contains more than one extraordinary point as its vertex (i.e., extraordinary points are isolated). As mentioned earlier, if these assumptions are not satisfied, one or two steps of global subdivision may be required and the resulting mesh can be treated as the initial mesh. Let the number of vertices in the initial mesh M be a, and let I of these be the extraordinary vertices. The number of faces in the initial mesh are assumed to be b, and k of these faces have vertices with degree 4 (henceforth termed a "normal face") and each of the remaining (b k) faces have one of the I extraordinary vertices (henceforth termed a "special face"). Let p be the 3a = N dimensional vector containing the control vertex positions in 3D. Using the formulations in Section 3.2.1 and Section 3.2.2, the smooth limit surface can be expressed as k s = En,+ s, (3.6) i=1 j=1 where n, is a single bicubic patch assigned to each of the normal faces and sj is a collection of infinite number of bicubic patches corresponding to each of the extraor dinary points. As s is a surface defined over the faces of the initial control mesh, each n, can be viewed as a bicubic patch defined over the corresponding regular face in the initial control mesh. Similarly, each sj can be defined over the corresponding irregu lar faces in the initial mesh (refer to Fig.3.2 Fig.3.4 for the detailed description on parametric domains of subdivision surfaces). For the simplicity of the following math ematical derivations on the dynamic model, from now on the parametric domain will not be explicitly provided in the formulation. Employing the same approach taken before to derive Eqn.3.5, it can be shown that k k k Sni = E(n')(P" ) = ( n J.)("A,))p = ("J)p, (3.7) I=1 i=1 i=l where "Ji,n pi and "nA are the equivalent of Ji, p, in Eqn.3.4 and A, in Eqn.3.5 respec tively. The presuperscript n is used to indicate that these mathematical quantities describe bicubic patch in the limit surface corresponding to normal faces. The following notational convention will be used for describing various math ematical quantities used in the derivation of the expression for a collection of infinite number of bicubic patches around an extraordinary vertex. The presuperscript s is used to represent a collection of bicubic patches around an extraordinary vertex, the subscript j is used to indicate the jth extraordinary point, the postsuperscript represents the exponent of a mathematical quantity and the level indicator (to rep resent various levels of subdivision in the local control mesh around an extraordinary vertex) is depicted via a subscript on the curly braces. The expression for sj is derived using the recursive nature of local subdivision around an extraordinary vertex as shown in Section 3.2.2. First, sj can be expressed as si = {'Jj}l{Spi + {s,}A, (3.8) where the first term of Eqn.3.8 is the generalized coordinate representation of the bicubic Bspline patches corresponding to the normal faces of the new local subdivi sion mesh obtained after one subdivision step on the local control mesh (similar to those patches marked n in Fig.3.4). {sj}, represents the rest of the infinite bicubic Bspline patches surrounding the extraordinary point (similar to the central patch enclosed by patches marked n in Fig.3.4). The vertices in the newly obtained lo cal subdivision mesh {'pjl} can be expressed as a linear combination of a subset of the vertices of the initial mesh M (which will contribute to the local subdivi sion) following the subdivision rules. We can name this subset of initial control vertices {'pj}o. Furthermore, there exists a matrix {'Bi,} of size (3c,3d), such that ({Bj}l('pj}o = {'Pj}i where {'pj}1 and {'pj}o are vectors of dimension 3c and 3d respectively. Applying the idea of recursive local subdivision again on {sj}l, sj can be further expanded as s = {*Jj3}1{'B3j},'pj} + {SJj03 'Bj}2{%jp} + {si}1. (3.9) In the above derivation, {'fj}1 is a vector of dimension 3d, comprising of a subset of the vertices defining the 3c dimensional vector {'pj}i. Note that, {'fij} has the same structure as {%pj}o, therefore, there exists a (3d,3d) matrix {*Cj}1 such that {'Cj),{'pj}, = {ij},. Each subdivision of a local mesh with d vertices creates a new local mesh with c vertices which contributes a fixed number of bicubic Bspline patches. Proceeding one step further, it can be shown that sj = {*Jl}1{'Bj,}'p}o+ {'Jj,}2{Bj}2{"'C}1{'pj}o +{'Jj}3{Bj}3{'Cj}1{'pj}o + {sj}3. (3.10) Because of the intrinsic property of the local recursive subdivision around the ex traordinary point, we have {'Jj} { SJj}2 = ... = {S'Jj} = ... = {'Jj}o In addition, the subdivision rules remain the same throughout the refinement process, we also have {1Bj}1 = {'Bj}2 = ... = {'Bj} = ... = {'Bj}m. So, the above equations can be further simplified leading to s = {SJj}1{Bj}{'pj+}o+ {'Jj}{'SBj}] {C,}{'pj}o +{"'Jj}1{'Bj}( {'Cj}){'pj}o +. = {*JJ}{iB,} ,( {SC)})P"{p "o, (3.11) =0 Vector sj can be written as sj = ('J, i'p,1, (3.12) where 'Jj = {'Jj,}l{'Bj}l(= 0o {'Cj,}) and 'pj = {'pj}0. The ideaoflocal recursive subdivision around an extraordinary point is illustrated in Fig.3.5. Note that, each vertex position in the subdivided mesh is obtained by an affine combination of some vertices in the previous level and hence any row of ({Cj}1 sums to 1. The largest eigenvalue of such a matrix is 1 and it can be shown that the corresponding infinite series is convergent following a similar approach as in Halstead et al. [38]. The rest of the derivation leading to an expression for s is relatively straight forward. Using the same approach used to derive the Eqn.3.7, it can be shown that sj = E(,'j)(p)) = (('Jj)('Aj))p = ('J)p. (3.13) j=1 j=1 j=1 From Eqn.3.6, Eqn.3.7 and Eqn.3.13, s = ("J)p + ('J)p. (3.14) Let J = ("J) + ('J), and hence s = Jp. (3.15) {SB} i J1 0 Figure 3.5. Local subdivision around the extraordinary point and the corresponding patches in the limit surface from different levels of subdivision. f<'>;?B^S>^ (Sj. WU ^^^ ^l ^ps''^^ * patches in the limit surface from different levels of subdivision. 3.2.4 Dynamics In the previous section, an expression is derived for the smooth limit surface obtained via infinite subdivision steps. This expression treats the smooth limit surface as a function of the control vertex positions in the initial mesh. Now the vertex positions in the initial mesh defining the smooth limit surface s are treated as a function of time in order to embed the CatmullClark subdivision model in a dynamic framework. The velocity of this surface model can be expressed as s(u,v,p)= Jp, (3.16) where an overstruck dot denotes a time derivative. The physics of the dynamic subdivision surface model is based on the workenergy version of Lagrangian dynamics [36] and is formulated in an analogous way to the dynamic framework presented in Section 2.2. Let ju(u, v) be the mass density function of the surface. Then the kinetic energy of the surface is T = Jf sdudv = p Mp, (3.17) where (using Eqn.3.16) M= f pJJdudv (3.18) is the mass matrix of size (N, N). Similarly, let y(u, v) be the damping density function of the surface. Then the dissipation energy is F = J fy sdudv = j1 Dp, (3.19) where D= f / 7JTjdudv (3.20) is the damping matrix of size (N, N). The potential energy of the smooth limit surface can be expressed as U = pTKp, (3.21) where K is the stiffness matrix of size (N, N). A thinplateundertension energy model [96] is used to compute the elastic potential energy of the dynamic subdivision surface. The corresponding expression for the stiffness matrix K is K = J(anJ J + a22J J +011 .TJ. + J12J J. + 223J J,)dudv. (3.22) where the subscripts on J denote the parametric partial derivatives. The aii(u, v) and i., ,,. v)s are elasticity functions controlling local tension and rigidity in the two parametric coordinate directions. Note that the thinplate energy expression might diverge at extraordinary points on the limit surface for CatmullClark subdivision scheme as shown in Halstead et al. [38]. Several methods have been suggested in Halstead et al. [38] to overcome the problem of the divergent series. However, in this dissertation, the problem is circumvented by setting the corresponding rigidity coef ficients to be zero at these points. Therefore, the thinplate energy at extraordinary points is zero. The effect is negligible to the overall thinplate energy as very few extraordinary points are present in the smooth limit surface. Using the expressions for the kinetic, dissipation and potential energy in Eqn.2.1, the same motion equation as in Eqn.2.5 can be obtained. The generalized force vector fp, which can be obtained through the principle of virtual work [36], is expressed as f, = f Jf(u, v, t)dudv (3.23) Different types of forces can be applied on the smooth limit surface, and the limit surface would evolve over time according to Eqn.2.5 to obtain an equilibrium position characterized by a minimum of the total model energy. 3.2.5 Multilevel Dynamics At this point, a dynamic framework is developed for the CatmullClark subdi vision scheme where the smooth limit surface evolves over time in response to the applied forces. The entire process can be described as follows. Given an initial mesh, a smooth surface is obtained in the limit. Users can directly apply synthesized forces to this smooth limit surface to enforce various functional and aesthetic constraints. This direct manipulation is then transferred back as virtual forces acting on the initial mesh through a transformation matrix (Eqn.3.23), and the initial mesh (as well as the underlying smooth limit surface) deforms continuously over time until an equilib rium position is obtained. The deformation of the surface in response to the applied forces is governed by the motion equation (Eqn.2.5). Within the physicsbased mod eling framework, the limit surface evolves as a consequence of the evolution of the initial mesh. One can apply various types of forces on the limit surface to obtain a desired effect, but the possible level of details appearing in a shape that can be obtained through evolution is constrained by the number of control vertices in the initial mesh. It might be necessary to increase the number of control vertices in the initial mesh in order to obtain a detailed shape through this evolution process. The number of control vertices defining the same smooth limit surface can be increased by simply replacing the initial mesh with a mesh obtained after one subdi vision step applied to the initial mesh. This new mesh has more number of vertices but defines the same limit surface. Therefore, the dynamic CatmullClark surface model can be subdivided globally to increase the number of vertices (control points) of the model. For example, after one step of global subdivision, the initial degrees of freedom p (refer to Eqn.3.15 and Eqn.3.16) in the dynamic system will be replaced by a larger number of degrees of freedom q, where q = Ap. A is a global subdivi sion matrix of size (M, N) whose entries are uniquely determined by CatmullClark subdivision rules (see Section 3.1 for the details about the rules). Thus, p, expressed as a function of q, can be written as p = (ATA)IATq = Atq, (3.24) where At = (ATA)'AT. Therefore, Eqn.3.15 and Eqn.3.16 can be rewritten as s = (JAt)q, (3.25) and s(u, v, q) = (JAt)q, (3.26) respectively. Now the equation of motion needs to be derived for this new subdivided model involving a larger number of control vertices namely q. The mass, damping and stiffness matrices need to be recomputed for this "finer" level. The structure of the motion equation as given by Eqn.2.5 remains unchanged, but the dimensionality and the entries of M, D, K, p and fp change correspondingly in this newly obtained subdivided level. In particular the motion equation, explicitly expressed as a function of q, can be written as Mq + Dq + Kq = fq, (3.27) where Mq f f p(At)TJTJAtdudv and the derivation of D,, K, and f, follow suit. It may be noted that further subdivision, if necessary, can be carried out in a similar fashion. Therefore, multilevel dynamics is achieved through recursive subdi vision on the initial set of control vertices. Users can interactively choose the level of detail representation of the dynamic model as appropriate for their modeling and de sign requirements. Alternatively, the system can automatically determine the level of subdivision most suitable for an application depending on some applicationspecific criteria. 3.3 Finite Element Implementation The evolution of the generalized coordinates for our new dynamic surface model can be determined by the secondorder differential equation as given by Eqn.2.5. An analytical solution of the governing differential equation can not be derived in general. However, an efficient numerical implementation can be obtained using the finite ele ment method [42]. The limit surface of the dynamic CatmullClark subdivision model is a collection of bicubic patch elements. We use two types of finite elements for this purpose normal elements (bicubic patches assigned to the normal faces of the initial mesh) and special elements (collection of infinite number of bicubic patches assigned to each extraordinary vertex of the initial mesh). The shape functions for the nor mal element are the basis functions of the bicubic surface patch, whereas the shape functions for the special element are the collection of basis functions corresponding to the bicubic patches in the special element. In the current implementation, the M, D and K matrices for each individual normal and special elements are calculated and they can be assembled into the global M, D and K matrices that appear in the corresponding discrete equation of motion. In practice, we never assemble the global matrices explicitly in the interest of time performance. The detailed implementation is explained in the following sections. 3.3.1 Data Structures The limit surface of the dynamic CatmullClark subdivision model is a collection of bicubic patches, and hence requires us to keep track of the control polygons defining such patches. We can get the control polygons for the normal elements in the limit Slist of faces list of edges Subdiv object S e list of vertices list of normal elements S subdivision aro list of faces list of facdges lubd of d  ofbdi ob Jim ofnorml elemCls i st of ononml element Figure 3.6. The data structure used for dynamic subdivision surface implementation. surface from the initial control mesh itself. However, we need to locally subdivide the initial control mesh around the extraordinary vertices to obtain the control polygons of the bicubic patches in the special elements on the limit surface. A subdivision surface defined by a control mesh at any level is designed as a class which has a pointer to its parent mesh, a set of pointers to its offspring meshes (arising out of local subdivision around the extraordinary vertices at that level), a list of faces, edges, vertices and normal elements Face, edge, vertex and normal elements are, in turn, classes which store all the connectivity and other information needed to either enumerate all the patches or locally subdivide around an extraordinary vertex in that level. The implementation takes the initial mesh as the base subdivision surface object (with its parent pointer set to NULL) and locally subdivides the initial mesh upto a userdefined maximum level around each extraordinary vertex to create offspring objects at different levels (Fig.3.6). At this point, let's take a closer look at the normal and special element data structures and computation of the corresponding local M, D and K matrices. 3.3.2 Normal Elements Each normal element is a bicubic surface patch and hence is defined by 16 vertices (from the 8connected neighborhood of the corresponding normal face) in the initial control mesh. Each normal element keeps a set of pointers to those vertices of the initial mesh which act as control points for the given element. For a normal element, the mass, damping and stiffness matrices are of size (16,16) and can be computed exactly by carrying out the necessary integration analytically. The matrix J in Eqn.3.18, 3.20 and 3.22 needs to be replaced by J, (of Eqn.3.2) for computation of the local M, D and K matrices respectively of the corresponding normal element. 3.3.3 Special Elements Each special element consists of an infinite number of bicubic patches in the limit. We have already described a recursive enumeration of the bicubic patches of a special element in Section 3.2.2. Let us now consider an arbitrary bicubic patch of the special element in some level j. The mass matrix M, of this patch can be written as M, = n' Mps (3.28) where Mp is the normal element mass matrix (scaled by a factor of y to take into account of the area shrinkage in bicubic patches at higher level of subdivision) and f, is the transformation matrix of the control points of that arbitrary patch from the corresponding control points in the initial mesh. The damping and stiffness matrices for the given bicubic patch can be derived in a similar fashion. These mass, damping and stiffness matrices from various levels of (local) subdivision can then be assembled to form the mass, damping and stiffness matrices of the special element. It may be noted that the stiffness energy due to plate terms diverges at the extraordinary points on the limit surface. We solve the problem using a slightly different approach than the one used in Halstead et al. [38]. When the area of the bicubic patch obtained via local subdivision of the initial mesh around an extraordinary vertex becomes smaller than the display resolution, the contribution from such a bicubic patch is ignored in computing the physical matrices of the corresponding special element. The number of extraordinary points in the limit surface is very few, and the above mentioned approximation is found to work well in practice. 3.3.4 Force Application The force f(u, v, t) in Eqn.3.23 represents the net effect of all externally applied forces. The current implementation supports spring, inflation as well as imagebased forces. However, other types of forces like repulsion forces, gravitational forces etc. can easily be implemented. To apply spring forces, a spring of stiffness k can be connected from a point do to a point xo on the limit surface (or, to the jth level approximation mesh), the net applied spring force being f(x, t) = s k(d0 s(x, t))6(x xo)dx, (3.29) where 6 is the unit impulse function implying f(xo, t) = kido s(xo, t) and vanishes elsewhere on the surface. However, the 6 function can be replaced with a smooth kernel to spread the force over a greater portion on the surface. The spring forces can be applied interactively using the computer mouse or the points from which forces need to be applied can be read in from a file. To recover shapes from 3D image data, we synthesize imagebased forces. A 3D edge detection is performed on a volume data set using the 3D MongaDeriche(MD) operator [69] to produce a 3D potential field P(x, y, z), which we use as an external potential for the model. The force distribution is then computed as f(x,y,z)= A VP(x, z) (3.30) 1 VP(lx, y, z) II' where A controls the strength of the force. The applied force on each vertex at the j th approximation level is computed by trilinear interpolation for evaluating Eqn.3.23 in Cartesian coordinates. More sophisticated imagebased forces which incorporate regionbased information such as gradients of a thresholded fuzzy voxel classification can also be used to yield better and more accurate shape recovery. It may be noted that we can apply spring forces in addition with the imagebased forces by placing points on the boundary of the region of interest in each slices of the 3D volume (MR, CT etc.) image data. 3.3.5 Discrete Dynamic Equation The differential equation given by Eqn.2.5 is integrated through time by dis cretizing the time derivative of p over time steps At. The state of the dynamic subdivision surface at time t + At is integrated using prior states at time t and t At. An implicit time integration method is used in the current implementation where discrete derivatives of p are calculated using ( + At) p(t + At) 2p(t) + p(t At) p(t+ At) =At2 (3.31) and S p(t + At) p(t At) p(t + At) = (3.32) 2At The elemental mass, damping and stiffness matrices can be assembled to get the global mass, damping and stiffness matrix for the smooth subdivision surface model. However, we do not assemble these global sparse matrices explicitly for efficiency reasons. For the time varying stiffness matrix, we recompute K at each time step. Using Eqn.2.5, Eqn.3.31 and Eqn.3.32, the discrete equation of motion is obtained as (2M + DAt + 2At2K)p(t + At) = 2At2(t + At) + (DAt 2M)p(t At) + 4Mp(t). (3.33) A evoslton evLion subdivson same limit surface at equilium with more patches eoion evolio (a) (b) Figure 3.7. Model subdivision to increase the degrees of freedom : (a) evolution of the initial mesh, and (b) the corresponding limit surface evolution perceived by the user. This linear system of equations is solved iteratively between each time step using the conjugate gradient method [34, 75]. 3.3.6 Model Subdivision The initialized model grows dynamically according to the equation of motion (Eqn.2.5). The degrees of freedom of the initialized model is equal to the number of control vertices in the initial mesh as mentioned earlier. When an equilibrium is achieved for the model, the number of control vertices can be increased by replacing the original initial mesh by a new initial mesh obtained through one step of Catmull Clark subdivision. This increases the number of degrees of freedom to represent the same (deformed) smooth limit surface and a new equilibrium position for the model can be obtained. This process is depicted schematically in Fig.3.7. Model subdivision might be needed to obtain a very localized effect on a smooth limit surface. For a shape recovery application, one may start with a very simple initial model, and when an approximate shape is recovered, the degrees of freedom can be increased to obtain a new equilibrium position for the model with a better fit to the given data set. The error of fit criteria for the discrete data is based on distance between the data points and the points on the limit surface where the corresponding springs are attached. In the context of imagebased forces, if the model energy does not change between successive iterations indicating an equilibrium for the given resolution, the degrees of freedom for the model can be increased by the abovementioned replacement scheme until the model energy is sufficiently small and the change in model energy between successive iterations becomes less than a prespecified tolerance. 3.4 Generalization of the Approach The proposed approach can be generalized for other approximating subdivision schemes. However, a more general approach is presented in Chapter 5 for deriving the dynamic framework. Dynamic CatmullClark subdivision surface model is reformu lated in Section 5.2 using this general approach. A dynamic framework for another popular approximating subdivision scheme namely, Loop's subdivision scheme, is also discussed in Chapter 5. CHAPTER 4 DYNAMIC BUTTERFLY SUBDIVISION SURFACES In the previous chapter, a dynamic framework has been presented for an approx imating subdivision scheme namely, CatmullClark subdivision scheme. In this chap ter, a dynamic framework is developed for (modified) butterfly subdivision scheme, which is an interpolatory subdivision technique. First, a brief overview of the (mod ified) butterfly subdivision scheme is presented. Next, a local parameterization tech nique for the limit surface obtained via (modified) butterfly subdivision is discussed. This parameterization scheme is then used to derive the dynamic model. The imple mentation details are also discussed. 4.1 Overview of the (Modified) Butterfly Subdivision Scheme The butterfly subdivision scheme [30], like many other subdivision schemes used in geometric design literature/applications, starts with an initial triangular mesh which is also known as the control mesh. The vertices of the control mesh are known as the control points. In each step of subdivision, the initial (control) mesh is refined through the transformation of each triangular face into a patch with four smaller triangular faces. After one step of refinement, the new mesh in the finer level retains the vertices of each triangular face in the previous level and hence, interpolates the coarser mesh in the previous level. In addition, every edge in each triangular face (a) Figure 4.1. (a) The control polygon with triangular faces; one subdivision step using butterfly subdivision rules. 6 v w V6 V V 1 5 4 12 3 Vj V2j VI w (b) (b) mesh obtained after w Figure 4.2. (a) The contributing vertices in the jth level for the vertex in the j+1 th level corresponding to the edge between v{ and v2; (b) the weighing factors for different vertices. I I " :. . is spilt by adding a new vertex whose position is obtained by an affine combination of the neighboring vertex positions in the coarser level. For instance, the mesh in Fig.4.1(b) is obtained by subdividing the initial mesh shown in Fig.4.1(a) once. It may be noted that all the newly introduced vertices corresponding to the edges in the original mesh have degree 6, whereas the position and degree of the original vertices do not change in the subdivided mesh. In the original butterfly scheme, the new vertices corresponding to the edges in the previous level are obtained using an eightpoint stencil as shown in Fig.4.2(a). To show how this stencil is used, a vertex to be introduced and the contributing vertices from the stencil are highlighted in Fig.4.1(a). The name of the scheme originated from the "butterfly"like configuration of the contributing vertices. The weighing factors for different contributing vertex positions are shown in Fig.4.2(b). The vertex e~1 in the j + 1th level of subdivision, corresponding to the edge connecting vertices v{ and v2 at level j, is obtained by ejf' = 0.5(vi + vN) + 2w(v~ + vJ) w(v3 + v) + v' + vj), (4.1) where 0 < w < 1, and vj denotes the position of the ith vertex at the jth level. The butterfly subdivision scheme produces a smooth C' surface in the limit except at the extraordinary points corresponding to the extraordinary vertices (ver tices with degree not equal to 6) in the initial mesh [109]. All the vertices introduced through subdivision have degree 6, and therefore, the number of extraordinary points in the smooth limit surface equals the number of extraordinary vertices in the initial (a) (b) Figure 4.3. (a) The weighing factors of contributing vertex positions for an edge connecting two vertices of degree 6; (b) the corresponding case when one vertex is of degree n and the other is of degree 6. mesh. Recently, this scheme has been modified by Zorin et al. [109] to obtain better smoothness properties at the extraordinary points. This modified scheme is known as modified butterfly subdivision technique. In Zorin et al. [109], all the edges have been categorized into three classes : (i) edges connecting two vertices of degree 6 (a 10 point stencil, as shown in Fig.4.3(a), is used to obtain the new vertex positions corresponding to these edges), (ii) edges connecting a vertex of degree 6 and a vertex of degree n 0 6 (the corresponding stencil to obtain new vertex position is shown in Fig.4.3(b), where q = .75 is the weight associated with the vertex of degree n 0 6, and si = (0.25 + cos(2ri/n) + 0.5cos(47ri/n))/n, i = 0, 1,..., n 1, are the weights associated with the vertices of degree 6), and (iii) edges connecting two vertices of degree n 5 6. The last case can not occur except in the initial mesh as the newly introduced vertices are of degree 6, and the new vertex position in this last case is obtained by averaging the positions obtained through the use of stencil shown in Fig.4.3(b) at each of those two extraordinary vertices. 4.2 Formulation In this section, a systematic formulation of the dynamic framework for the modi fied butterfly subdivision scheme is presented. Unlike the approximating schemes, the limit surface obtained via modified butterfly subdivision does not have any analytic expression even in a regular setting. Therefore, derivation of a local parameterization suitable for developing the dynamic framework for the limit surface is the key issue here which is discussed next. Figure 4.4. The smoothing effect of the subdivision process on the triangles of the initial mesh. 4.2.1 Local Parameterization The smooth limit surface defined by the modified butterfly subdivision technique is of arbitrary topology where a global parameterization is impossible. Nevertheless, the limit surface can be locally parameterized over the domain defined by the initial mesh following a similar approach described in Lounsbery et al. [53]. The idea is to track any arbitrary point on the initial mesh across the meshes obtained via the subdivision process (see Fig.4.4 and Fig.4.5), so that a correspondence can be established between the point being tracked in the initial mesh and its image on the limit surface. The modified butterfly subdivision scheme starts with an initial mesh consisting of a set of triangular faces. The recursive application of the subdivision rules smoothes out each triangular face, and in the limit, a smooth surface is obtained which can also be considered as a collection of smooth triangular patches. The subdivision process and the triangular decomposition of the limit surface is depicted in Fig.4.4. Note that, the limit surface can be represented by the same number of smooth triangular patches as that of the triangular faces in the initial mesh. Therefore, the limit surface  0 0 s can be expressed as s = Sk, (4.2) k=1 where n is the number of triangular faces in the initial mesh and sk is the smooth triangular patch in the limit surface corresponding to the kth triangular face in the initial mesh. The stage is now set for describing the parameterization of the limit surface over the initial mesh. The process is best explained through the following example. A simple planar mesh shown in Fig.4.5(a) is chosen as the initial mesh. An arbitrary point x inside the triangular face abc is tracked over the meshes obtained through subdivision. The vertices in the initial mesh are darkly shaded in Fig.4.5. After one step of subdivision, the initial mesh is refined by addition of new vertices which are lightly shaded. Another subdivision step on this refined mesh leads to a finer mesh with introduction of new vertices which are unshaded. It may be noted that any point inside the smooth triangular patch in the limit surface corresponding to the face abc in the initial mesh depends only on the vertices in the initial mesh which are within the 2neighborhood of the vertices a, b and c due to the local nature of the subdivision process. For example, the vertex d, introduced after first subdivision step, can be obtained using the 10 point stencil shown in Fig.4.3(a) on the edge ab. All the contributing vertices in the initial mesh are within the 1neighborhood of the vertices a and b. A 10 point stencil can be used again in the next subdivision step on the edge db to obtain the vertex g. Some of the contributing vertices at this level of subdivision, for example, the (lightly shaded) 1neighbors of the vertex b (except (a) (b) Figure 4.5. Tracking a point x through various levels of subdivision : (a) initial mesh, (b) the selected section (enclosed by dotted lines) of the mesh in (a), after one subdivision step, (c) the selected section of the mesh in (b), after another subdivision step. d and e) in Fig.4.5(b), depend on some vertices in the initial mesh which are within the 2neighborhood of the vertices a, b and c in the initial mesh. In the rest of the formulation, superscripts are used to indicate the subdivision level. For example, v7 denotes the collection of vertices at level j which control the smooth patch in the limit surface corresponding to the triangular face uvw at the jth level of subdivision. Let vob be the collection of vertices in the initial mesh which are within the 2neighborhood of the vertices a, b and c (marked black in Fig.4.5(a)). Let the number of such vertices be r. Then, the vector v0bc, which is the concatenation of the (x, y, z) positions for all the r vertices, is of dimension 3r. These r vertices control the smooth triangular patch in the limit surface corresponding to the triangular face abc in the initial mesh. Now, there exists four subdivision matrices (Aa,),, (Anac)t, (Aa.e), and (Ab), of dimension (3r, 3r) such that v = (abc),rVac Vbd = (Aoa)iVabc, vfe = (Aoate),ac, v, = (Aabc)mva (4.3) where the subscripts t, 1, r and m denote top, left, right and middle triangle positions respectively (indicating the relative position of the new triangle with respect to the original triangle), and Vid, V vi and vYef are the concatenation of the (x, y, z) positions for the vertices in the 2neighborhood of the corresponding triangle in the newly obtained subdivided mesh. The new vertices in this level of subdivision are lightly shaded in Fig.4.5(b). The 2neighborhood configuration of the vertices in the newly obtained triangles is exactly the same as that of the original triangle, hence local subdivision matrices are square and the vector dimensions on both sides of Eqn.4.3 are the same. This concept is further illustrated in Fig.4.6. Carrying out one more level of subdivision, a new set of vertices which are unshaded in Fig.4.5(c) are obtained along with the old vertices. Adopting a similar approach as in the derivation of Eqn.4.3, it can be shown that vd = (Abed)tVed Vbg = (Abed)ived 6vh = (Abed),vted vjh, = (Abed)mVed (4.4) The relative position of the triangular face dgi in Fig.4.5(c) with respect to the triangular face bed is topologically the same as of the triangular face adf in Fig.4.5(b) with respect to the triangular face abc. Therefore, one can write (Abed), = (Ab~)t,. Using similar reasoning, Eqn.4.4 can be rewritten as v~, = (AMd)tved = (Aabc)(tVed vbh = (Abed)Ved = (Aab)rVld v2 = (Abed)rVd = (Aabc)rVed Vethbe 67 .?'!' ,b , (b) 4 v Figure 4.6. Different set of control points at a subdivided level obtained by applying different subdivision matrices on a given set of control points in a coarser mesh. vg2 = (Ad)mVV = (Aak)mV d. (4.5) Combining Eqn.4.3 and Eqn.4.5, it can be shown that vd = (Ask)t(Aok),vo, = (Aa~)l(A.kc)iVO, vh = (Ao.),(Aai),v , vh = (A.)m(Aa ),v.bc (4.6) Let x be a point with barycentric coordinates (ac fac, Y7 ) inside the trian gular face abc. When the initial mesh is subdivided, x becomes a point inside the triangular face bed with barycentric coordinates (i,. r.,.,y*). Another level of subdivision causes x to be included in the triangular face dgi with barycentric coor dinates (admi, },"y/,i). Let sj denote the jth level approximation of the smooth triangular patch Sab in the limit surface corresponding to the triangular face abc in the initial mesh. Now v0 can be written as r r r c = C ., b, ,..., b,,c,,.. where the subscripts x, y and z indicate the x, y and z coordinates respectively of the corresponding vertex position. The expressions for v d and vi can also be written in a similar manner. Next, the matrix Bb can be constructed as follows: r r r "oa, b, 7, 0,... 0,0,...,0,0,..., 0 r r r Bob(x) = 0 .,, 0,1 ..., 0, bc r bc 7el ,0, .,0,0,..., 0 r r 0F. 0,0,0.. a5b 0 7l 0,.. , The matrices B1d and B' can also be constructed in a similar fashion. Now sobc(X), s9c(x), and sla(x) can be written as Sabc(X) B0 Sk(X) = Bn (x)vab, Sabc(x) = B' (x)v = Bd(x)(A ) = Bbed(X)Vbed = Bbed(X)(Aabc)ljVab, s(x) = B' (x)v, = B'd(X)(Aabc)tVd g B (X)(Aabc)t(Aasc)iVac. (4.7) Proceeding in a similar way, the expression for s8c(x), jth level approximation of sa(x), is given by sI(x) = B,,,(x)(Abc)... (Aak)t(Aabc) Ve c = Bu,(x)(Alb)vbc = Bjc(x)vObc, (4.8) where x is inside the triangular face uvw at level j (with an assumption that uvw is the triangular face in the middle with respect to its coarser level original triangular face in the previous level), (Ab) = (Ab,,)m... (A(Aa,,)(Aa,) and Bi (x) = Bj,(x)(A ). It may be noted that the sequence of applying (Aabs)t, (Aask), (Abc), and (Aabc)m depends on the triangle inside which the tracked point x falls after each subdivision step. Finally, the local parameterization process can be completed by writing Sab(x) = (lim BI.(x))v = Back(x)v (4.9) In the above equation, Bae is the collection of basis functions at the vertices of v b. It may also be noted that the modified butterfly subdivision scheme is a sta tionary subdivision process, and hence new vertex positions are obtained by affine combinations of nearby vertices. This guarantees that each row of the matrices (Anb),, (Aabe), (Aab), and (Aatc)m sums to one. The largest eigenvalue of such ma trices is 1 and therefore the limit in Eqn.4.9 exists. Now, assuming the triangular face abc is the kth face in the initial mesh, Eqn.4.9 can be rewritten as Sk(x) = Bk(x)v0 = Bk(x)Akp, (4.10) where p is the concatenation of the (x,y,z) positions of all the vertices in the initial mesh and the matrix Ak, when postmultiplied by p, selects the vertices vA controlling the kth smooth triangular patch in the limit surface. If there are t vertices in the initial mesh and r of them control the kth patch, then p is a vector of dimension 3t, Ak is a matrix of dimension (3r, 3t), and Bk(x) is a matrix of dimension (3,3r). Combining Eqn.4.2 and Eqn.4.10, it can be shown that s(x) =( Bk(x)Ak)p = J(x)p, (4.11) k=1 where J, a matrix of dimension (3,3t), is the collection of basis functions for the corresponding vertices in the initial mesh. The vector p is also known as the degrees of freedom vector of the smooth limit surface s. 4.2.2 Dynamics An expression of the limit surface obtained via modified butterfly subdivision process is derived in the last section. Now a dynamic framework for the limit surface can be derived using this expression. The derivation of the dynamic model from this point is exactly similar in spirit as that of the dynamic CatmullClark subdivision model presented in Section 3.2.4. The vertex positions in the initial mesh defining the smooth limit surface s are treated as a function of time in order to embed the modified butterfly subdivision model in a dynamic framework. The velocity of this surface model can be expressed as 9(x,p) = J(x)p, (4.12) where an overstruck dot denotes a time derivative and x e So, So being the domain defined by the initial mesh. As in the case of the dynamic CatmullClark subdivision surfaces, the Lagrangian equation of motion (Eqn.2.1) is used to derive the dynamic modified butterfly subdivision surface model. Let p be the mass density function of the surface. Then the kinetic energy of the surface is T= j ( p(x)s(x)x) = p Mp, (4.13) where (using Eqn.4.12) M = fx5so p(x)JT(x)J(x)dx is the mass matrix of dimen sion (3t, 3t). Similarly, let 7 be the damping density function of the surface. The dissipation energy is F = 1 y(x)s(x)s(x)dx = TDp, (4.14) 2 JeSO 2 where D= JxES y(x)JT(x)J(x)dx is the damping matrix of dimension (3t, 3t). The potential energy of the smooth limit surface can be expressed as U pTKp, (4.15) where K is the stiffness matrix of dimension (3t, 3t) obtained by assigning various internal energies to a discretized approximation of the limit surface and is detailed in Section 4.3. Using the expressions for the kinetic, dissipation and potential energy in Eqn.2.1, the same motion equation as in Eqn.2.5 can be obtained. The generalized force vector f, can be obtained is a similar fashion described in Section 3.2.4 and is given by f,= feso J'(x)f(x, t)dx. (4.16) 4.2.3 Multilevel Dynamics The initial mesh of the dynamic modified butterfly subdivision surface model can be subdivided to increase the degrees of freedom for model representation. The situation is exactly similar to that of the dynamic CatmullClark subdivision surface model as described in Section 3.2.5 and a similar example is used for illustration. After one step of modified butterfly subdivision, the initial degrees of freedom p (refer to Eqn.4.11 and Eqn.4.12) in the dynamic system is replaced by a larger number of degrees of freedom q, where q = Bp. B is a global subdivision matrix of size (3s, 3t) whose entries are uniquely determined by the weights used in the modified butterfly subdivision scheme (see Section 4.1 for the weights). Thus, p, expressed as a function of q, can be written as p = (BTB)IBTq = Btq, (4.17) where Bt = (BTB)'BT. Therefore, Eqn.4.11 and Eqn.4.12 is modified as s(x) = (J(x)Bt)q, (4.18) and s(x,q) = (J(x)Bt)il, (4.19) (b) Figure 4.7. (a) An initial mesh, and (b) the corresponding limit surface. The domains of the shaded elements in the limit surface are the corresponding triangular faces in the initial mesh. The encircled vertices in (a) are the degrees of freedom for the corresponding element. respectively. The motion equation, explicitly expressed as a function of q, can be written as Mq + Dq + Kq = fq, (4.20) where Mq = fxEsl p(x)(Bt)TJT(x)J(x)Btdx, S' being the domain defined by the newly obtained subdivided mesh. The derivations of Dq, K, and f, are similar. t  4.3 Finite Element Implementation In the previous section, the smooth limit surface obtained via modified butterfly subdivision process is expressed as a function of the control vertex positions in its initial mesh; mass and damping distribution, internal deformation energy and forces are assigned to the limit surface in order to develop the corresponding physical model. In this section, the implementation of this physical model is detailed. It may be noted that even though the formulation of the dynamic framework for modified butterfly subdivision scheme is similar to that of CatmullClark subdivision scheme once the local parameterization is derived, the finite element implementation is totally different as the limit surface does not have any analytic expression in case of the modified butterfly subdivision scheme. In Section 4.2 it was pointed out that the smooth limit surface obtained by the recursive application of the modified butterfly subdivision rules can be represented by a set of smooth triangular patches, each of which is represented by a finite element. The shape (basis) function of this finite element is obtained by smoothing a hat function through repeated application of the modified butterfly subdivision rules. The number of finite elements in the smooth limit surface is equal to the number of triangular faces in the initial mesh as mentioned earlier (refer Fig.4.4 and 4.7). A detailed discussion on how to derive the mass, damping and stiffness matrices for these elements is presented next. These elemental matrices can be assembled to obtain the global physical matrices M, D and K, and a numerical solution to the governing secondorder differential equation as given by Eqn.2.5 can be obtained using finite element analysis techniques [42]. The same example as in Section 4.2 (refer Fig.4.5) is used to develop the related concepts. The concept of elements along with the control vertices and their corresponding domains is further illustrated in Fig.4.7. Now it will be shown how to derive the mass, damping and stiffness matrices for the element corresponding to the triangular face abc in Fig.4.5. The derivations hold for any element in general. 4.3.1 Elemental Mass and Damping matrices The mass matrix for the element given by Sbc (Eqn.4.9) can be written as Mb = f p(x) {B.(x)}T{Ba(x)}dx. (4.21) However, from Eqn.4.9 it is known that the basis functions corresponding to the vertices in vk which are stored as entries in Ba, are obtained as a limiting process. These basis functions do not have any analytic form, hence computing the inner product of such basis functions as needed in Eqn.4.21 is a challenging problem. In Lounsbery et al. [53], an outline is provided on the computation of these inner products without performing any integration. In this dissertation, a different yet simpler approach is developed to solve this problem. The smooth triangular patch in the limit surface corresponding to the face abc in the initial mesh is approximated by a triangular mesh with 43 faces obtained after j levels of subdivision of the original triangular face abc (each subdivision step splits one triangular face into 4 triangular faces). Then the mass matrix can be expressed as Mbc= f L (x){BJi(x)}T{B,(x)}dx. (4.22) The jth level approximation of the corresponding basis functions can be explicitly evaluated (refer Eqn.4.8 for an expression of B1 ). An important point to note is that Eqn.4.8 involves several matrix multiplications and hence can be very expensive to evaluate. However, the matrix (AJ()(= (Aab) ... (Aoac)t(Aa.)t) in Eqn.4.8 encodes how vertices in the 2neighborhood of the triangular face uvw at level j are related to the vertices in the 2neighborhood of the triangular face abc in the initial mesh. In the implementation, how a new vertex is obtained from the contributing vertices in its immediate predecessor level is tracked. The information stored in (Abc) can be obtained by tracking the contributions from level j to level 0 recursively. Thus the entries of (Abc) can be determined without forming any local subdivision matrices and thereby avoiding subsequent matrix multiplications. By choosing a sufficiently high value of j, a reasonably good approximation of the elemental mass matrices is achieved. The computations involved in evaluating the integrals in Eqn.4.22 are eliminated by using discrete mass density function which has nonzero values only at the vertex positions of the jth subdivision level mesh. Therefore, the approximation of the mass matrix for the element can be written as Ma,. = k ),(V, T B v) B )}, (4.23) i=1 where k is the number of vertices in the triangular mesh with 43 faces. This approxi mation has been found to be very effective and efficient for implementation purposes. The elemental damping matrix can be obtained in a similar fashion. 4.3.2 Elemental Stiffness Matrix The internal energy is assigned to each element in the limit surface, thereby defining the internal energy of the smooth subdivision surface model. A similar approach as in the derivation of the elemental mass and damping matrix is taken and the internal energy is assigned to a jth level approximation of the element. Three types of internal energy are used tension, stiffness and spring energy. For the examples used throughout the paper, this energy at the jth level of approx imation can be written as Ekac Ej = (Ek), + (Ek),, + (Ef )p,, (4.24) where the subscripts t, st and sp denote tension, stiffness and spring respectively. The expression for the tension energy, which is essentially equivalent to the first order strain (membrane) energy [96], is (E )t = ktEvi vI2 2 0 = 1ktv .T (K )Jt(Vk}, (4.25) where kt is a constant, vf and vi, the 1th and mth vertex in the jth level mesh, are in the 1neighborhood of each other, Q is the domain defined by all such vertex pairs, and vb is the concatenation of the (x,y,z) positions of all the vertices in the jth subdivision level of the triangular face abc in the initial mesh. Similarly, the expression for stiffness energy, which is equivalent to the second order strain (thin plate) energy [96], can be written as (E ) k., Iv 2v + +vi2 n 1= k{v (L v }, (4.26) where vJ, vi and v are the three vertices of a triangular face. The summation in volves three terms corresponding to each triangular face, and is over all the triangular faces in the mesh at the jth level of subdivision. Finally, a spring energy term is added which is given by 1 (kim) (I vi v e1) (E bc)sp 2 i  vm (vj VIm) S2{v } (Ki ),,{v }, (4.27) where (km),, is the spring constant, Q is the domain as in Eqn.4.25 and ,im is the natural length of the spring connected between vi and vi. It may be noted that the entries in (Kb ), depend on the distance between the connected vertices and hence, (Kib),, unlike other elemental matrices, is a function of time which needs to be recomputed in each time step. Combining the expressions for tension, stiffness and spring energy, it can be shown that Ekb ab c} {kt(Ka + ket(K c),+(K ),j}{v } = {v }(K abc) tva T1 j T j j 2{(A f){v }} (K bc)(A bc){V } I T KT j )j (4.28) = {v } {(AIse) (K )(A ac)} vac}, (4.28) where (Abe) and vb are same as in Eqn.4.8. Therefore, the expression for the elemental stiffness matrix is given by Ka = (Ab)T(Ka)(Ab). (4.29) It may be noted that the matrix multiplications for constructing Kabc are avoided in the implementation by following the same technique described in Section 4.3.1. 4.3.3 Other Implementation Issues The techniques of applying synthesized forces on the smooth limit surface ob tained via modified butterfly subdivision process are similar to those described in Section 3.3.4 in the context of dynamic CatmullClark subdivision surface model. The techniques of model subdivision and discrete dynamic equation derivation are also the same. 82 4.4 Generalization of the Approach The proposed approach of deriving the dynamic framework for modified butter fly subdivision scheme is generalized for other interpolatory subdivision schemes in the next chapter, Chapter 5. CHAPTER 5 UNIFIED DYNAMIC FRAMEWORK FOR SUBDIVISION SURFACES Dynamic framework for specific subdivision schemes was presented in the last two chapters. In Chapter 3, CatmullClark subdivision scheme, a representative of the approximating subdivision schemes, is embedded in a physicsbased modeling en vironment. In Chapter 4, a dynamic framework is provided for the modified butterfly subdivision scheme, a popular interpolatory subdivision technique. In this chapter, a unified approach is presented for embedding any subdivision scheme in a dynamic framework. 5.1 Overview of the Unified Approach The key concept of the unified approach is the ability to express the smooth limit surface generated by a subdivision scheme as a collection of a single type of finite elements. The type of this finite element depends on the subdivision scheme involved. For example, the limit surface generated by the CatmullClark subdivision scheme can be expressed as a collection a quadrilateral finite element patches. This approach differs form the one taken in Chapter 3 where the limit surface was a collection of normal elements (corresponding to the quadrilateral bicubic Bspline patches away from the extraordinary points) and special elements (corresponding to the nsided patches near extraordinary vertices of degree n). This concept will be further elaborated in the next few sections. The unified approach adopts two different strategies depending on the subdivi sion scheme involved. One strategy is for the approximating schemes, whose limit surface is some type of Bspline surface away from the extraordinary points, and the other is for the interpolatory subdivision schemes, where the limit surface does not have any analytical expression. It may be recalled that recursive subdivision is a method for smoothing polyhedra, and hence corresponding to each nsided face in the control mesh, there would be a smooth nsided patch in the limit surface in general. For example, the control mesh (after at most one subdivision step) of the CatmullClark subdivision scheme consists of quadrilateral faces, which when smoothed through subdivision lead to quadrilateral patches in the limit surface (see Fig.5.1(b)). In the case of approximating subdivision schemes, the faces which do not contain any extraordinary vertex lead to spline patches with analytic expressions. On the other hand, faces which have one extraordinary vertex, lead to patches whose analytic expressions are difficult to determine. Continuing the example with the CatmullClark subdivision scheme, the faces without extraordinary vertices lead to quadrilateral bicubic Bspline patches in the limit, whereas the faces with one extraor dinary vertex lead to quadrilateral patches whose analytic expressions are difficult to obtain. In Chapter 3, n such adjacent patches corresponding to an extraordinary vertex of degree n are effectively combined together to obtain an expression for the re sulting nsided patch (see Fig.5.1(a)). Recently, Stam [90, 91] has introduced schemes for exactly evaluating the limit surface, even near the extraordinary points, for ap proximating subdivision schemes. Analytic expressions for the patches resulting from smoothing faces with an extraordinary vertex can be obtained using these schemes. Therefore, the limit surface can be expressed as a collection of a single patch type, and consequently a single type of finite element can be used to provide the dynamic framework for the approximating subdivision schemes. This unified framework for approximating subdivision schemes is fully worked out in this chapter for Catmull Clark and Loop subdivision schemes and a general outline is also presented on how to carry out the same strategy for other approximating subdivision schemes. The limit surface resulting from an interpolatory subdivision scheme does not have an analytic expression in general, and hence a strategy similar to the one adopted in Chapter 4 needs to be used. An outline of the strategy to be used for interpolatory subdivision schemes is also presented in this chapter. 5.2 Unified Dynamic Framework for CatmullClark Subdivision Surfaces A systematic formulation along with the implementation details of the unified dynamic framework for CatmullClark subdivision surfaces is presented in this sec tion. The key difference between the framework developed in Chapter 3 and the one presented here is the representation of the limit surface, which leads to different types of finite elements used for developing the dynamic framework. This is illustrated with a schematic diagram in Fig.5.1. 
Full Text 
xml version 1.0 encoding UTF8
REPORT xmlns http:www.fcla.edudlsmddaitss xmlns:xsi http:www.w3.org2001XMLSchemainstance xsi:schemaLocation http:www.fcla.edudlsmddaitssdaitssReport.xsd INGEST IEID EZKAYFT03_SDXRFY INGEST_TIME 20130123T14:39:46Z PACKAGE AA00012984_00001 AGREEMENT_INFO ACCOUNT UF PROJECT UFDC FILES 