Citation |

- Permanent Link:
- http://ufdc.ufl.edu/UF00101371/00001
## Material Information- Title:
- Dynamic framework for subdivision surfaces
- Creator:
- Mandal, Chhandomay (
*Dissertant*) Vemuri, Baba C. (*Thesis advisor*) Qin, Hong (*Thesis advisor*) Sahni, Sartaj (*Reviewer*) Fishwick, Paul A. (*Reviewer*) Cenzer, Douglas A. (*Reviewer*) - Place of Publication:
- Gainesville, Fla.
- Publisher:
- Department of Computer and Information Science and Engineering, University of Florida
- Publication Date:
- 1998
- Copyright Date:
- 1998
- Language:
- English
- Physical Description:
- xv, 168 leaves : ill. ; 29 cm.
## Subjects- Subjects / Keywords:
- Butterflies ( jstor )
Damping ( jstor ) Data ranges ( jstor ) Datasets ( jstor ) Degrees of freedom ( jstor ) Dynamic modeling ( jstor ) Geometric shapes ( jstor ) Matrices ( jstor ) Topology ( jstor ) Vertices ( jstor ) Computer and Information Science and Engineering thesis, Ph.D Computer graphics -- Mathematical models ( lcsh ) Computer vision -- Mathematical models ( lcsh ) Dissertations, Academic -- UF -- Computer and Information Science and Engineering Image processing -- Digital techniques ( lcsh ) - Genre:
- bibliography ( marcgt )
theses ( marcgt )
## Notes- Abstract:
- Subdivision surfaces are extensively used to model smooth shapes of arbitrary topology. Recursive subdivision on a user-defined initial control mesh generates a visually pleasing smooth surface in the limit. However, users have to carefully select the intiail 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 fromt eh lack of direct manipulation tools for the limit surface. In this dissertation, techniques from physics-based 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 physics based â€œ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 surface 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 damping distribution functions, internal deformation energy, forces, and other physical quantities required to develop a physics-based subdivision surface model. In this dissertation, local parameterization techniques suitable for embedding the geometric subdivision surface model in a physics-based modeling framework have been developed. Specific local parameterization techniques have been fully developed for the Catmull-Clark, 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 physics-based deformable models are incorporated 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 non-rigid 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.
- General Note:
- Typescript.
- General Note:
- Vita.
- Thesis:
- Thesis (Ph.D.)--University of Florida, 1998.
- Bibliography:
- Includes bibliographical references (leaves 159-167).
## Record Information- Source Institution:
- University of Florida
- Holding Location:
- University of Florida
- Rights Management:
- Copyright [name of dissertation author]. Permission granted to the University of Florida to digitize, archive and distribute this item for non-profit research and educational purposes. Any reuse of this item in excess of fair use or other copyright exemptions requires permission of the copyright holder.
- Resource Identifier:
- 41869700 ( OCLC )
002427664 ( ALEPH )
## UFDC Membership |

Downloads |

## This item has the following downloads: |

Full Text |

A DYNA i',C FRA -.i1l WORK FOR SUBDIVISION SURFACES CHHANDO.iAY I. iANDAL A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFIl I'-.i.l:NT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 1998 Copyright 1998 by C !,!,,i ,,l, I_.,-, ', I,,i ,l I To iy Parents, and .liy Brother ACKNOWLEDGE; 1: NTS 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 non-academic 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 co-advisor, 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. Fi-hl, i, 1 and Dr. Douglas A. Center for serving in my supervisory committee and advising me on various aspects of this dissertation. .iy dissertation research was supported in part by the NSF grant ECS-9210648 and the NIH grant ROl-i:'. 'i944 of Dr. Vemuri, and by the NSF CAREER award CCR-9702103 and DI 1i-9700129 of Dr. Qin. I wish to acknowledge Dr. Tim :.1. il,- erney, Dr. Gregoire M.il.,iiil.iii Dr. Hughes Hoppe, Dr. Kari Pulli, Dr. Dimitry Goldgof, Dr. C!, i-i ii., Leonard and Dr. S!,.li.-TII, ig 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 C11 ii, Li Wang, 'l illi.iying Huang, Clii.iiir Kapoor, Fengting C11 i Arun Srini- vasan, Jundong Liu and Jun Ye. Special thanks goes to my long time roommates, Raja Cli.ii i rjee and Kingshuk .ijii inl.ii. for helping me in various ways during my stay at Gainesville. Finilly, 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 ACKNOW LEDG ( il NTS ............................. 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 Physics-based Deformable Surface M'.idels ....... ...... 17 2.3 S~11 ie '.ibdeling Using Physics-based Subdivision-surface '.idel . 20 2.4 11 1iip and '. ition Estimation Using Physics-based Subdivision-surface '.idel ....... ... .... ....... ...... 22 2.5 Unified Framework for !,liip Recovery and S!liipi dealing g . . 24 3 DYNA' iiC CAT'.IULL-CLARK SUBDIVISION SURFACES ...... ..27 3.1 Overview of the Catmull-Clark 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 '.i ili!; vel Dynamics ................ .... .. 46 3.3 Finiii 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 .i, del Subdivision .................. .. .. .. 55 3.4 Generalization of the Approach ............... .. .. 56 4 DYNA-.IIC BUTTER FLY SUBDIVISION SURFACES . . . .... 57 4.1 Overview of the ('.1>dified) Butterfly Subdivision Scheme ...... . 57 4.2 Formulation .................. ........... .. .. 61 4.2.1 Local Parameterization ................. . . .. 62 4.2.2 Dynamics .................. .......... .. 72 4.2.3 `.ililk!vel Dynamics ................ .... .. 74 4.3 Filiii Element Implementation ............... .. ... 76 4.3.1 Elemental '.I1--, and Damping matrices . . . . .... 77 4.3.2 Elemental Stiffness '.,i.i.\ .................. . 79 4.3.3 Other Implementation Issues ..... . . . . . ..... 81 4.4 Generalization of the Approach ............... .. .. 82 5 UNIFIED DYNAI .iC FRAI ii;WORK FOR SUBDIVISION SURFACES 83 5.1 Overview of the Unified Approach . . . . . . . .. .. 83 5.2 Unified Dynamic Framework for Catmull-Clark Subdivision Surfaces 85 5.2.1 Local Parameterization ................. . . .. 88 5.2.2 Dynamics .................. .......... .. 92 5.2.3 Fi;iii 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 Fi;iii 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 'iULTIRESOLUTION DYNAIiCS ............. . . .. .. 103 6.1 Overview of ', ltl i resolution Analysis and Wavelets . . . ... 109 6.2 ':. li ii resolution Analysis for Surfaces of Arbitrary Topology . . 114 6.2.1 .-1 ii Spaces through Subdivision . . . . . . 116 6.2.2 Inner Product ........ . . . . . ........ 118 6.2.3 Biorthogonal Surface Wavelets on Arbitrary .' li l 1 . . .119 6.3 '. li it resolution Representation ...... . . . . . ..125 6.4 Dynamics .................. . . . . ... 128 6.5 Implementation Details .................. .. ...... 129 7 SYSTEM ARCHITECTURE ............ 7.1 Topological Information Processing ':. dule ............. 135 7.2 Subdivision '.1, dule ........... ................ 135 7.3 Fiiii Element Analysis '.idule ...... ....... . . . .. 135 7.4 Force Synthesis '.ibdule .................. .. ..... . 136 7.5 Update Engine .................. ......... .. .. 136 7.6 Display '.iodule ............... . . . . . .. . 137 8 APPLICATIONS ................. . . . . ... 138 8.1 Geometric ':. Ideling ......... . . . . . . ... 138 8.2 SliipeI Recovery from Range and Volume Data. . . . . . . 1 1 8.3 Non-rigid M.i )tion Estimation. ............. .. . .. .. 148 8.4 'liili resolution Editing .................. .. ..... . 150 8.5 .'; 1 i 1.! Terrain ':.MIdeling ................... . . 151 9 CONCLUSIONS AND FUTURE WORK .. . . . . . ..... 155 9.1 Conclusions ................. . . . . . ... 155 9.2 Future Directions ............. . . . . . ...... 156 9.2.1 Automatic Cn111, of Topology . . . . . ..... 156 9.2.2 Local Refinement ........ . . . . . . ... 157 9.2.3 Automatic i. del Parameter Selection . . . . .... 157 9.2.4 Constraint Imposition ....... . . . . . ... 157 9.2.5 Recovery of S'1,irp Features ..... . . . . . ..... 158 9.2.6 Automatic o. del 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 Catmull-Clark subdivision : (a) initial mesh, (b) mesh obtained after one step of Catmull-Clark 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 T!i, data structure used for dynamic subdivision surface implementation. 50 3.7 '1ldel 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) Ti, control polygon with triangular faces; (b) mesh obtained after one subdivision step using butterfly subdivision rules. . . . ... 58 4.2 (a) T!i, contributing vertices in the j-th level for the vertex in the j+1- th level corresponding to the edge between vi and vi; (b) the weighing factors for different vertices. .................. .... 58 4.3 (a) Tli, 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 Ti, smoothing effect of the subdivision process on the triangles of the initial m esh. . . . . . . . . . .. ... . . . 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. Tli, do- mains of the shaded elements in the limit surface are the corresponding triangular faces in the initial mesh. Tli, 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 Cliipi r 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 Catmull-Clark subdivision, each non-boundary 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) Tli, 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) Tli, 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. T!i, domains of the shaded triangular patches in the limit surface are the corresponding triangular faces in the initial mesh. Tli, 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 non-boundary 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 '.i llil, 1vel dynamics vs. multiresolution dynamics . . . . . 106 6.3 Representation of the degrees of freedom in multilevel dynamics and multiresolution dynamics approach. ..... . . . . . ...... 108 6.4 Tli. 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 Catmull-Clark 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 '.11:1 (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 Catmull-Clark subdivision surface model fitted to the cerebellum. . . . . ... 146 8.7 (a) Data points identifying the boundary of the caudate nucleus on a I ii: 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) Ti, initialized model along with the data set; (b) the fitted model with two subdivisions on the initial mesh, along with attached springs for editing. Tli~ 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 DYNAI ilC FRA. 1il:WORK FOR SUBDIVISION SURFACES By n1 1 i i ii Id I l l December 1998 C11.!l ,, ini Dr. Baba C. Vemuri Cochairman : Dr. Hong Qin ., ij, r Department: Computer and Information Science and Engineering Subdivision surfaces are extensively used to model smooth shapes of arbitrary topology. Recursive subdivision on a user-defined 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 physics-based 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 physics-based "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. '.i, I- 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 physics-based subdivision surface model. In this dissertation, local parameterization techniques suitable for embedding the geometric subdivision surface model in a physics-based modeling framework have been devel- oped. Specific local parameterization techniques have been fully developed for the Catmull-Clark, 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 physics-based deformable models are in- corporated into the conventional subdivision schemes, and a dynamic hierarchical control of this model is introduced. Fi illy, a multiresolution representation of the control mesh is developed and a unified approach for deriving subdivision surface- based finite elements is presented. Ti, 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 non-rigid motion tracking from a temporal sequence of data sets. ':. ill resolution 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, computer-aided geometric design and scientific visualization researchers. Tli, 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. Ti, 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 user-defined initial control mesh. Tli, 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. Ti, current state-of-the-art 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 surface-based 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. t' ., i cheless, 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 Physics-based 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 physics-based modeling paradigm, a deformable model is derived by assigning mass, damping, stiffness and other physical properties to a surface model. Ti, 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 physics-based modeling paradigm, then the smooth limit surface can be deformed directly by using physics-based "force" tools. However, this procedure-based surface model obtained through the subdivision pro- cess does not have a closed-form 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 physics-based 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. Tli1 motion of this physics- based deformable subdivision surface model is governed by a second-order differential equation, which is solved numerically using the finite element method. N"- types of finite elements for the chosen subdivision scheme are also presented for representing the smooth limit surface. Tli, 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. Ti, 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. Tli, 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 physics-based 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 .iiit!" clay/play-dough modeling environment. Ti, proposed model can be edited directly in a hierarchical fashion using synthesized forces. Also, this physics-based 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 surface-based and physics-based modeling techniques to solve important theoretical and practical problems. Although the principles of physics-based modeling are well understood by the graphics experts and modeling researchers, this dissertation will greatly advance the state of the art in physics-based 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 physics-based smooth limit surface. Tli, smooth dynamic subdivision surface in the limit is treated as a i !" physical object represented by a set of novel finite elements. Ti, basis (shape) functions of these new variety of finite elements are derived using the subdivision schemes. Ti, 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. Ti, subdivision techniques are used to create a surface model that incorporates mass and damping distribution functions, internal deformation energy, forces, and other physical quantities. Tli, motion equations are also systematically derived for this dynamic subdivision surface model whose degrees of freedom are the collection of initial user-specified control vertices. Tli, i, fi, e, the advantages of both the physics-based modeling philosophy and the geometric subdivision schemes are incorporated within a single unified framework. Users will be able to manipulate this physics-based model in an arbitrary region, and the model will respond naturally (just like a real-world object would) to this force application. This shape deformation is quantitatively characterized by the time-varying displacement of control points that uniquely define the geometry of the limit surface. Tlii 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 physics-based 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 force-based 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 hierarchically-structured approximation can satisfy any user-specified error tol- erance. Ti, 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. Fini lly, 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 Tli, rest of the dissertation is organized as follows : Cli~,ipi r 2 contains an overview of the subdivision surfaces along with a review of the related literature. Ti, motivation for embedding the subdivision surface models in a physics-based modeling framework is discussed, and the proposed model is com- pared with the existing physics-based models. Tli, advantages of a unified framework for shape modeling and shape recovery are also pointed out in this chapter. Cllip, ir 3 provides a dynamic framework for the Catmull-Clark 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 Catmull-Clark subdivision scheme is derived and the p!1i, !" quantities required to develop the dynamic model are introduced [57, 60, 76]. Til, governing dynamic differential equation is derived using Lagrangian mechanics and is imple- mented using a finite element technique. In Cl11,lpi r 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. Tlii physics-based 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]. Cli,,1pi r 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 Cliipi, r 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 Clilpi, r 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. Ti, proposed modeling framework is used in various application areas and the results are presented in Clilipi r 8. Ti, proposed dynamic model has been success- fully used in geometric modeling, shape recovery and non-rigid 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. Finilly, conclusions are drawn in Cliipi, r 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 physics-based models along with a review of the re- lated literature. Tli, motivating factors for embedding geometric subdivision surface models in a physics-based modeling paradigm are presented in the next section. TI, advantages of the proposed dynamic models over the existing physics-based 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 Tli, input to any subdivision scheme is an initial mesh (also known as control mesh), \1" = (V, Fo), which is a collection of vertices V0 and a collection of faces F. Ti, subdivision surface S(.\l") associated with the initial mesh \'1" is defined as the limit of the recursive application of the refinement R as shown in Fig.2.1. TIi, refinement R, when applied on a mesh I.k = (V, Fk), creates a refined mesh \/k+1 = (Vk+l, Fk+l) where the vertices in Vk+1 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 sub-faces. (c) Mesh M' = R(R(M)) (d) Limit surface S(M) = R"(M) Figure 2.1. Refinement of an initial control mesh to obtain the limit surface. (Cour- tesy : H. Hoppe) 'b) Mesh M-' = R(1M) 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). Ti, vertices which do not have this degree are called extraordinary vertices. Ti, number of extraordinary vertices are very few as these vertices are not introduced via subdivision refinement in general. For example, in the Catmull-Clark 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. Ti, limit surface defined by the subdivision process is most often C1(first derivative continuous) or C2 (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 physics-based modeling paradigm to enhance the applicabil- ity of these surface models. In the rest of this section, these i iliid .1 !" techniques of subdivision surfaces are reviewed. Cli.il.i i [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 Cli;i .i i'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. Tih 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 sub-vertices, 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 B-splines, 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. Tih Catmull-Clark subdivision surface, defined by an arbitrary initial mesh, can be reduced to a set of standard bicubic B-spline patches except at a finite number of degenerate points. This scheme is an example of face subdivision scheme, where a k-sided face is split into k sub- faces. Tl, details of this subdivision scheme are discussed in Section 3.1. Halstead et al. [38] proposed an algorithm to construct a Catmull-Clark 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 B-splines 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. T11 -- 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. Tli, ii refinement rules yield a C1 surface that has a piecewise quadratic pa- rameterization except at isolated extraordinary points. ':.1- recently, non-uniform Doo-Sabin and Catmull-Clark surfaces that generalize non-uniform tensor product B-spline 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]. Ti, most well-known interpolation-based subdivision scheme is the ill i, ill ' 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 C1 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. T111 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. Ti, i, f; e, these schemes are global, i.e., every new vertex position depends on all the vertex positions of the coarser level mesh. Ti, local refinement property which makes the subdivision schemes attractive for implementation in the computer graphics applications is not retained in the variational approach. Tli, 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 eigen-analysis 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]. '.1i -, recently, Stam [90] presented a method for exact evaluation of Catmull-Clark 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 Physics-based Deformable Surface .odels SIilipe 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 physics-based framework to allow only a pli -i, .llI mean- i.hlil model behavior. Various types of energies are assigned to the model using the degrees of freedom, and the model "d(, I .i ii-" to find an equilibrium state with min- imal energy. Ti, motion equation is derived using Lagrangian dynamics [36], where various energies associated with model gives rise to internal and external forces. Ti, equilibrium state of the model is a model position where the internal deformation force becomes equal to the externally applied force. Tli, physics-based 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). Ti, 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. Tir, i, the Lagrangian equation of motion for the model can be expressed as d OT OT OF OU + --+ f=. (2.1) dt 0 /; 0 1 (0/ ( 0/, Let i(u, v) be the mass density function of the surface. Til 11, the kinetic energy of the surface is T = f _s'sdudv = pMp, (2.2) 2 i2 where M is called the mass matrix. Similarly, let j(u, v) be the damping density function of the surface. Tli, I the dissipation energy is F = 1 'sdudv = pDp, (2.3) 2JJ 2 where D is called the damping matrix. Tli, potential energy of the model can be defined using the thin-plate-under-tension energy model [96], and is given by S s' Os Os' Os 02s' 02s U (a + O22 P11 U S( u u ov ov o2u o2u 02T 02s 02s' 2s 1 + 312 v + 322 )dudv = p Kp, (2.4) uLiv uLiv d2- d2v 2 ' where aii(u, v) and 3ij (u, v) are elasticity functions which control tension and rigidity respectively of the deformable surface model, and K is known as the stiffness matrix. Ti, 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. Tlii free-form 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], i i-.11.1,x and Terzopou- los [65] and Vemuri and Radisavljevic [104]. Celniker and Gossard [11] developed a system for interactive free-form 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 B-spline 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 (D-NURBS) which are very sophisticated models suit- able for representing a wide variety of free-form as well as standard analytic shapes. T1i, D-NURBS 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 ':.1 II1i Iny and Terzopoulos [64] and the references therein. Some specific existing deformable models used for shape recovery and non-rigid motion estimation will be reviewed in Section 2.4. 2.3 1,ip,1 .ideling Using Physics-based Subdivision-surface :iodel 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 surface-based modeling schemes. Despite the advent of advanced 3D graphics interaction tools, these indirect geometric operations remain non-intuitive 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, physics-based modeling provides a superior approach to shape modeling that can overcome most of the limitations associated with traditional geometric modeling approaches. Free-form deformable models governed by the laws of continuum mechanics are of particular interest in this context. Ti, -, dynamic models respond to externally applied forces in a very intuitive manner. Tlii 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. Tli, deformation energy functionals can be formulated to satisfy local and global modeling criteria, and geometric constraints relevant to shape design can also be imposed. Tli, dynamic approach subsumes all of the aforementioned modeling capabilities in a formulation which grounds everything in real-world physical behavior. A severe limitation of the existing deformable models, including D-NURBS, 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. Ti, 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 real-world !Li A dynamic framework for subdivision surface-based 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. TI formulation and implementation details of these dynamic frameworks are discussed in subsequent chapters. 2.4 11i,11i and otion Estimation Using Physics-based Subdivision-surface : .odel Tlii 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 non-rigid 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. Tli, 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 physics-based 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 physics-based modeling paradigm. T111 models involve the use of either fixed size [21, 66, 71, 98, 104] or adaptive size [15, 41, 47, 62, 95, 101] grids. Tli, 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 .1 i, i i iy 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. i'.1 .reover, almost all of these schemes require a globally parameterized mesh as their input which may be infeasible at times. Tli, 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. Ti, 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. Ti, process of subdivision after achieving an approximate fit is continued till a pre- scribed error criteria for fitting the data points is achieved. Tli, 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 Sliple Recovery and ,l,1. 1i 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 fine-tuned 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 state-of-the-art methods yield large polygonal meshes which are very difficult to manipulate, especially for global changes in shape. .- often the resulting meshes from a shape recovery application are not directly amenable to multiresolution analysis. Computationally expensive re-meshing 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 lidivision-( .iI II,, 1 i.;l 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. Tli, 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 physics-based force tools. T!i11- the modeler is re- lieved of the burden of both building an initial model and editing a cumbersome huge polygonal mesh. Tli, 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 pre-defined error in fit criteria. Ti, initial mesh of the final recovered smooth shape has subdivi- sion connectivity in-built, as it is obtained by subdivision refinement, and therefore no re-meshing is needed for multiresolution analysis. Tli, dynamic framework for subdivision surface-based wavelets makes multiresolution editing using physics-based force tools easy to perform. Ti, -, advantages of shape modeling and shape recovery in a unified framework will be further illustrated in later chapters. CHAPTER 3 DYNA-i, C CAT'.IULL-CLARK SUBDIVISION SURFACES Subdivision surfaces, as mentioned earlier, can be broadly classified into two categories -approximating and interpolatory subdivision schemes. Ti, 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 Catmull-Clark subdivision scheme, a popular approx- imating subdivision technique. T!i, Catmull-Clark subdivision scheme generalizes the idea of obtaining a bicubic B-spline surface by recursive subdivision of a rectangu- lar initial mesh. Before discussing the dynamic framework, first the Catmull-Clark subdivision scheme is briefly reviewed. 3.1 Overview of the Catmull-Clark Subdivision Scheme Catmull-Clark subdivision scheme, like any other subdivision scheme, starts with a user-defined 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. Ti, subdivision rules are as follows: (a) (b) Figure 3.1. Catmull-Clark subdivision : (a) initial mesh, (b) mesh obtained after one step of Catmull-Clark 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 (non-boundary) 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 (non-boundary) vertex V, a new vertex is introduced whose position is + 2E- (n 3) where F is the average of the new face vertices of all faces n n n 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. -, "i edges are formed by connecting each newx face vertex to the newx edge vertices of the edges defining the old face and by connecting each newx vertex corresponding to an old vertex to the new edge vertices of all old edges incident on the old vertex. * i faces are defined as faces enclosed by new edges. An example of Catmull-Clark subdivision on an initial mesh is shown in Fig.3.1. Tli, most important property of the Catmull-Clark subdivision surfaces is that a smooth surface can be generated from any control mesh of arbitrary topology. Tl ,i 1 - fore, this subdivision scheme is extremely valuable for modeling various complicated geometric objects of arbitrary topology. Catmull-Clark subdivision surfaces include standard bicubic B-spline surfaces as their special case (i.e., the limit surface is a bicubic B-spline surface for a rectangular mesh with all non-boundary vertices of degree 4). In addition, the aforementioned subdivision rules generalize the recur- sive bicubic B-spline patch subdivision algorithm. For non-rectangular meshes, the limit surface converges to a bicubic B-spline surface except at a finite number of ex- traordinary points. T111 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 newx vertices created subsequently will have four incident edges. Tli, 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 Catmull-Clark subdivision on the initial mesh. Tli, limit surface is curvature-continuous everywhere except at extraordinary vertices, where only tangent plane continuity is achieved. In spite of the popularity of Catmull-Clark subdivision surfaces for representing complex geometric shapes of arbitrary topology, these subdivision surfaces may not be easily parameterizable and deriving a closed-form analytic formulation of the limit surface can be very difficult. T111 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 Catmull-Clark subdivision surfaces which offers a closed-form 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 Catmull-Clark subdivision surfaces at arbitrary parameter values, but the work presented here on dynamic Catmull-Clark subdivision surfaces [57, 60, 76] was published prior to the publication of the above-mentioned work. 3.2 Formulation In this section, a systematic formulation of the new dynamic model based on Catmull-Clark 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. TI 1, i .,_re, 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 non-rectangular 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 B-spline 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. Ti, 1, F; re, 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. Tli, 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. Tli! i, F, re, face F1 can be thought of as the portion of the parametric domain over which the patch S1 is defined, i.e., has non-zero values. t' ,. i. heless, each rectangular face (e.g. Fi) can be parametrically I I -------- ,------- -- ----- r II I I -- -" ----- S S 4 --_--- --------------------------------- _ ---- ---- - S2 S3 r----- - - - 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 B-spline 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. T11, faces F0, F1, ..., F8 are assigned to bicubic patches So, S1,..., Sg respectively (as they all have vertices of degree 4) following the aforementioned scheme. However, the central smooth surface enclosed by the patches So, Si,..., Sg 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 Catmull-Clark subdivision surface model. Tli1 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 F0, F1, .., F8, Po, P1, 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 F II II 1 I I S2 S4 30Z) 8 S 6 \ \ S S Si F - - - - - - - - - -- i o Fu i Figure 3.3. A mesh xith an extraordinary point of degree 3 and its limit surface. selected mesh for -- local subdivision t 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 B-spline patch is obtained as the limiting process of the Catmull-Clark subdivision algorithm applied to an initial 4 by 4 rectangular control mesh. Let sp(u, v), where (u, v) c [0, 1]2, denote this bicubic B-spline patch which can be expressed analytically as 3 3 Sp(u, v) (x(u, v), y(u, v), z(u, v))T= E di,jBi,4(u)Bj,4(v), (3.1) i=0 j=0 where dj represents a 3-dimensional position vector at the (i,j)th control point location and Bi,4(u) and Bj,4(v) are the cubic B-spline basis functions. Ti, subscript p on s denotes the patch under consideration. Expressing Eqn.3.1 in a generalized coordinate system we have Sp = Jpqp, (3.2) where Jp is a transformation matrix storing the basis functions of a bicubic B-spline patch, and is of size (3, 48). Vector q, is the concatenation of all control points defining a B-spline 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) = J(1, k + 1) = Jp(2,k + 2) = Bi,4 ()Bj,4(v) and Jp(0, k + 1) = Jp(0, k + 2) = Jp((1, k) = Jp(, 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 (u, v) = s, (2u,2v)+s2(2(u- -),2v)+s, (2(u- ),2(v )) +s, (2u, 2(v )), 2 2 2 2 (3.3) where sm,(2u, 2v) = sm(u,v) for 0 < u, v < and 0 otherwise. Similarly, sm,,Sm3 and s,4 are also equal to sm(u, v) for an appropriate range of values of u, v and 0 outside. It may be noted that sm, sm s,, m ,4 correspond to patches S1, S2, S3, S4 respectively in Fig.3.2. Eqn.3.3 in generalized coordinates becomes 4 Sm= Jlql + J2q2 + J3q3 + J4q4 = iqi, (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 = JiAlq + J2A2qm + J3A3q[m + J4A4qm = JiAiqm = Jmqm. (3.5) i=1 Where, q, is the 75-component vector of 3D positions of the 25 vertex control mesh defining the limit surface s,. ':., 111 i'es Ai, 1 < i < 4, are of size (48, 75) each row consisting of a single nonzero entry (= 1) and the (3, 75)-sized matrix J = Zi=1 JiAi. Limit surface with many bicubic patches from an arbitrary initial mesh Tli, stage is now set to define the limit surface s using the vertices of initial mesh AM 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 AM be a, and let I of these be the extraordinary vertices. Tli, 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 111 I i1i, 1. I ) and each of the remaining (b- k) faces have one of the I extraordinary vertices (henceforth termed a -"1' i! 1, 1. "). 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 I s = ni + sj, (3.6) i=1 j=l where ni 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 ni 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 n, = (Ji)("pi) = (("Ji)(nAi))p = ("J)p, (3.7) i=1 i=1 i=1 Where "Ji," pi and "Ai are the equivalent of Ji, pi in Eqn.3.4 and Ai in Eqn.3.5 respec- tively. Ti, pre-superscript n is used to indicate that these mathematical quantities describe bicubic patch in the limit surface corresponding to normal faces. Ti, 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. Tli, pre-superscript s is used to represent a collection of bicubic patches around an extraordinary vertex, the subscript j is used to indicate the j-th extraordinary point, the post-superscript 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. Ti, 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 sj = {SJj}{Spj} + {sj}l, (3.8) where the first term of Eqn.3.8 is the generalized coordinate representation of the bicubic B-spline 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}1 represents the rest of the infinite bicubic B-spline patches surrounding the extraordinary point (similar to the central patch enclosed by patches marked n in Fig.3.4). Tli, vertices in the newly obtained lo- cal subdivision mesh {'"pj} 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 }. Furthermore, there exists a matrix {'Bj}I of size (3c, 3d), such that {'Bj}{pjJ}0o = {'pj}1 where {"pJ}, and {(PJ}0 are vectors of dimension 3c and 3d respectively. Applying the idea of recursive local subdivision again on {sj}I, sj can be further expanded as sj {'Jj}l{'Bj}l{'pj}o + fJJ2{ Bj}2 + {s}2. (3.9) In the above derivation, {'Pj}I is a vector of dimension 3d, comprising of a subset of the vertices defining the 3c dimensional vector {'pj}I. Note that, {'Pj}I has the same structure as {Spj}o, therefore, there exists a (3d,3d) matrix {'Cj}1 such that {SCj} {spj } = {sPj}l. Each subdivision of a local mesh with d vertices creates a new local mesh with c vertices which contributes a fixed number of bicubic B-spline patches. Proceeding one step further, it can be shown that Sj = {sJj}l{Bj}lp J+pJ}o JJ2{z"BJC}2 Cj SPj}o +{"J}{SBy C}{s}o + {s,}3. (3.10) Because of the intrinsic property of the local recursive subdivision around the ex- traordinary point, we have {"J} = {sJj}2 ... = }{Jj = ... = .Jj. In addition, the subdivision rules remain the same throughout the refinement process, we also have {(Bj} = {"Bj}2 = ... = {(Bj}, ... = { Bj},. So, the above equations can be further simplified leading to Sj = {sJj} {Bj}l {spj}o {SJjlSBjl {sCj} {spj}o +{sJj}{SBj "sC} 2 sP }o ... = sJjl{sBj}l (E {Cj} ){"Spj}o. (3.11) i=0 Vector Sj can be written as sj = (SJj)(Spj), (3.12) where "Jj = {(JjJ}I"BN}I(E o {"CN}i) and "pj = {"pJ}0. Tl, idea of local 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 {(CjN} sums to 1. Tlii 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]. Ti, 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 Ssj = (sJj)(spj) = ( s(JjA)(sAj))p = (sJ)p. (3.13) j=1 j=1 j=1 From Eqn.3.6, Eqn.3.7 and Eqn.3.13, s = (nJ)p + (SJ)p. (3.14) Let J = ("J) + ("J), and hence s = Jp. (3.15) fIspi ------ fsB. 1sj 2s so Figure 3.5. Local subdivision around the extraordinary point and the corresponding 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 Catmull-Clark subdivision model in a dynamic framework. Tlii velocity of this surface model can be expressed as s(u,v,p) = Jp, (3.16) where an overstruck dot denotes a time derivative. Ti, physics of the dynamic subdivision surface model is based on the work-energy version of Lagrangian dynamics [36] and is formulated in an analogous way to the dynamic framework presented in Section 2.2. Let it(u, v) be the mass density function of the surface. T 11, the kinetic energy of the surface is T 2= I ,IsTsdudv 2= TMp, (3.17) where (using Eqn.3.16) M = if Jdudv (3.18) is the mass matrix of size (N, N). Similarly, let '(u, v) be the damping density function of the surface. Ti 1, the dissipation energy is F I I 'sdudv = Dp, (3.19) where D= J JJdudv (3.20) is the damping matrix of size (N, N). Ti, potential energy of the smooth limit surface can be expressed as 1 U = p'Kp, (3.21) 2 where K is the stiffness matrix of size (N, N). A thin-plate-under-tension energy model [96] is used to compute the elastic potential energy of the dynamic subdivision surface. Ti, corresponding expression for the stiffness matrix K is K = (/ (a(IJo f + OJJ2 +v 11 + iT 12 J + 322 VJ,)dudv. (3.22) where the subscripts on J denote the parametric partial derivatives. Ti11 aii(u,v) and 3ij (u, v)s are elasticity functions controlling local tension and rigidity in the two parametric coordinate directions. Note that the thin-plate energy expression might diverge at extraordinary points on the limit surface for Catmull-Clark subdivision scheme as shown in Halstead et al. [38]. Several methods have been -il-.-. -i ,1 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. Tli11 i fi .e, the thin-plate energy at extraordinary points is zero. Tli, effect is negligible to the overall thin-plate 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. Ti, generalized force vector fp, which can be obtained through the principle of virtual work [36], is expressed as fp = Jff j (, 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 ',lii i!lvel Dynamics At this point, a dynamic framework is developed for the Catmull-Clark subdi- vision scheme where the smooth limit surface evolves over time in response to the applied forces. Ti, 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. Ti, deformation of the surface in response to the applied forces is governed by the motion equation (Eqn.2.5). Within the physics-based 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. Tli, 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. T!, ,, f, re, the dynamic Catmull-Clark 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 (.~, N) whose entries are uniquely determined by Catmull-Clark subdivision rules (see Section 3.1 for the details about the rules). Tlil-, p, expressed as a function of q, can be written as p = (ATA) -ATq = Atq, (3.24) where At = (A A)- AT. Tli, i, fI .re, 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. Ti, mass, damping and stiffness matrices need to be recomputed for this "li "I level. T!i, 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 = ff J/(At)'J' JAtdudv and the derivation of Dq, Kq and fq follow suit. It may be noted that further subdivision, if necessary, can be carried out in a similar fashion. T, i i, ,re, 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 application-specific criteria. 3.3 Finii Element Implementation Ti, evolution of the generalized coordinates for our new dynamic surface model can be determined by the second-order 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]. Tli, limit surface of the dynamic Catmull-Clark 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). Tli, 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. Tli, detailed implementation is explained in the following sections. 3.3.1 Data Structures Tli, limit surface of the dynamic Catmull-Clark 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 NLL list of faces List of edges Subdiv object -Suv o t list of vertices list of normal elements 1 c 1 subdivision around extraordinary vertices List of faces list of faces list of edges -u v oj list of edges -ubdiv object 3ubdiv object list of vertices list of vertices List of normal elements -list of normal elements Figure 3.6. T!i 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. Ti, 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 user-defined 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 8-connected 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. Ti, matrix J in Eqn.3.18, 3.20 and 3.22 needs to be replaced by Jp (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. Ti, mass matrix M, of this patch can be written as M, = Mp, (3.28) where Mp is the normal element mass matrix (scaled by a factor of to take into account of the area shrinkage in bicubic patches at higher level of subdivision) and Q2, is the transformation matrix of the control points of that arbitrary patch from the corresponding control points in the initial mesh. Tli, damping and stiffness matrices for the given bicubic patch can be derived in a similar fashion. Tli1 -i 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. Tli1 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 Ti!, force f(u, v, t) in Eqn.3.23 represents the net effect of all externally applied forces. Ti, current implementation supports spring, inflation as well as image-based 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 j-th level approximation mesh), the net applied spring force being f(x, t) = (do s(x,t))(x- xo)dx, (3.29) where 8 is the unit impulse function implying f(xo, t) = / (,, 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. Tli, 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 image-based forces. A 3D edge detection is performed on a volume data set using the 3D '-.1 n;.i-Deriche(. Ii)) operator [69] to produce a 3D potential field P(x, y, z), which we use as an external potential for the model. Tli, force distribution is then computed as VP(x" y, z) f(x, y,z)= A P( (3.30) || VPr y, ) || where A controls the strength of the force. Tli, applied force on each vertex at the j- th approximation level is computed by trilinear interpolation for evaluating Eqn.3.23 in Cartesian coordinates. '.i ,re sophisticated image-based forces which incorporate region-based 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 image-based forces by placing points on the boundary of the region of interest in each slices of the 3D volume (' i : CT etc.) image data. 3.3.5 Discrete Dynamic Equation Ti, differential equation given by Eqn.2.5 is integrated through time by dis- cretizing the time derivative of p over time steps At. Tli, 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 (t + At) = pp(t + At) 2p(t) + p(t At) (3.31) p (t +A t) = 2 (3.31) At2 and p(t + At) = p(t + At) p(t At) (3.32) p(t t)- =t> (3.32) 2At Ti, 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) = 2At2f,(t + At) + (DAt 2M)p(t At) + 4Mp(t). (3.33) evolLtion evolution subdiv sion same limit surface at equili rium with more patches evol tion evol tion (a) (b) Figure 3.7. .idel 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 '.iodel Subdivision Tli, initialized model grows dynamically according to the equation of motion (Eqn.2.5). Tli, 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. '- idel 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. Tl1I 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 image-based 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 above-mentioned replacement scheme until the model energy is sufficiently small and the change in model energy between successive iterations becomes less than a pre-specified tolerance. 3.4 Generalization of the Approach Tli, proposed approach can be generalized for other approximating subdivision schemes. However, a more general approach is presented in Clplii, r, 5 for deriving the dynamic framework. Dynamic Catmull-Clark 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 Cl.,pli r 5. CHAPTER 4 DYNAI, iC BUTTERFLY SUBDIVISION SURFACES In the previous chapter, a dynamic framework has been presented for an approx- imating subdivision scheme namely, Catmull-Clark 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. ':i \, 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. Tli, imple- mentation details are also discussed. 4.1 Overview of the ( modifiedd) Butterfly Subdivision Scheme Tli, 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. Ti, 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) (b) Figure 4.1. (a) Ti, control polygon with triangular faces; (b) mesh obtained after one subdivision step using butterfly subdivision rules. -W -W -W -W Figure 4.2. (a) Tli, contributing vertices in the j-th level for the vertex in the j+1- th level corresponding to the edge between vi and v ; (b) the weighing factors for different vertices. 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 eight-point 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). Tli, name of the scheme originated from the ill i ii Il'-like configuration of the contributing vertices. T11, weighing factors for different contributing vertex positions are shown in Fig.4.2(b). Ti, vertex el in the j 1-th level of subdivision, corresponding to the edge connecting vertices v] and vi at level j, is obtained by e1 = 0.5(v1 + vj) + 2w(v- + vj) w(vj + v- + v' + vJ), (4.1) where 0 < w < 1, and vJ denotes the position of the i-th vertex at the j-th level. Tli, butterfly subdivision scheme produces a smooth C1 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 -0.06 0625-w n-2 w (a) (b) Figure 4.3. (a) Tli, 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 = 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 $ 6, and si = (0.25 + cos(27ri/n) + 0.5cos(47i/n))/n, i = 0, 1,... 1, are the weights associated with the vertices of degree 6), and (iii) edges connecting two vertices of degree n $ 6. Ti, 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. Ti, iI fi 'e, derivation of a local parameterization suitable for developing the dynamic framework for the limit surface is the key issue here which is discussed next. * 0- Figure 4.4. Tli, smoothing effect of the subdivision process on the triangles of the initial mesh. 4.2.1 Local Parameterization Ti, smooth limit surface defined by the modified butterfly subdivision technique is of arbitrary topology where a global parameterization is impossible. t'. '. i cheless, 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]. Tli, 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. Tli, modified butterfly subdivision scheme starts with an initial mesh consisting of a set of triangular faces. T!i, 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. Ti, 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. Tli, i fI .e, the limit surface s can be expressed as s = Sk, (4.2) k=l 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 k-th triangular face in the initial mesh. Tih stage is now set for describing the parameterization of the limit surface over the initial mesh. Tih 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. Th, 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 newx vertices which are lightly shaded. Another subdivision step on this refined mesh leads to a finer mesh with introduction of newx 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 2-neighborhood 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 1-neighborhood 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) 1-neighbors 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 2-neighborhood 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, vuvw, 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 j-th level of subdivision. Let vobe be the collection of vertices in the initial mesh which are within the 2-neighborhood of the vertices a, b and c (marked black in Fig.4.5(a)). Let the number of such vertices be r. Ti,11 the vector Vab which is the concatenation of the (x, y, z) positions for all the r vertices, is of dimension 3r. Til, -, 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 (Aabc)t, (Aabc)l, (Aabc)r and (Aabc)m of dimension (3r, 3r) such that Vadf = (Abc\Vc Vbed = (Aabc)lVabc f = (eAabc)rVbc Vde = (Aabc),mVb0, (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 V Vd, ved, Vfe and vef are the concatenation of the (x, y, z) positions for the vertices in the 2-neighborhood of the corresponding triangle in the newly obtained subdivided mesh. Ti. new vertices in this level of subdivision are lightly shaded in Fig.4.5(b). Ti, 2-neighborhood 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 V gi = (Abed)tV d v 2g = (Abed)V Ied Veih = (Abed)Ved Vghi = (Abed)mVbed (4.4) Tli, 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. Ti, I r ;.e, one can write (Abed)t = (Aabc)t. Using similar reasoning, Eqn.4.4 can be rewritten as V gi = (Abed)tV d = (Aabc)tV d V hg = (Abed)V Ied = (Aabc)V Ied Vi = (Abed)rVed = (Aabc)rV1ed Veih __bd 67 IV,-bc abc l I 1 eb ------------------------------ (a) (Aabc) 1 / (b) Figure 4.6. Different set of control points at different subdivision matrices on a given set a subdivided level obtained by applying of control points in a coarser mesh. vi = (Abed),Vd = (Aabc),Ved. (4.5) Combining Eqn.4.3 and Eqn.4.5, it can be shown that dgi = (Aabc)t(Aabc)tVbc, Vhg = (Aabc)t(Aabc)tVbc, 2 0 Veih = (Aabc)r(Aabc)lVabc, Vhi = (Aabc)m(Aabc)tVbc. (4.6) Let x be a point with barycentric coordinates (Oabc( 0bc 7abc) inside the trian- gular face abc. When the initial mesh is subdivided, x becomes a point inside the triangular face bed with barycentric coordinates ( bd, led, Tbed). Another level of subdivision causes x to be included in the triangular face dgi with barycentric coor- dinates (Ogi,, q2g 7ygi). Let sbc denote the j-th level approximation of the smooth triangular patch Sabc in the limit surface corresponding to the triangular face abc in the initial mesh. Now vbc can be written as r r r oabc [=ax, bx, C, ..., ay, by, cy, . ., az, bz, c, .... where the subscripts x, y and z indicate the x, y and z coordinates respectively of the corresponding vertex position. Tli, expressions for vd and v 2 can also be written ~h~1~311V1 1 II bed dg in a similar manner. 'i \ the matrix BOb can be constructed as follows: r r r bc ybc, 0,..., 0, 0,... 0,..., 0 r r r B xabc, ( 0,.,00 1 abc,, c bc, O, O, 0...,0 0 0, "", "' r r r 0,...,0,0,...,0, ab bc bc, ... (2.)ab, V ^c. .... abc, b Tli, matrices Bed and B' can also be constructed in a similar fashion. Now so~ (x), sbc((x), and sb c(x) can be written as sabc(x) 0 Bbc(X)Vbc, bc(x) = Bbed(x)vbed = B (x)(Aabc)V bc, sbc(x) = B ( (x)v = B 0 j(x)(Aabc,)tved = B j (x)(Abc)(Aabc)V b (4.7) Proceeding in a similar way, the expression for sbc(x), j-th level approximation of Sabc (x), is given by sbc(x) = BL,,(x) (Aabc)m ... (Aabc)t(Aabc); Vabc = BI(x) (Abc)Vabc S bc (X) V abc, (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), (A.bc) = (A.abc (Aabc),(Aabc) and Bbc (x) = BW,(x)(AbJ). It may be noted that the sequence of applying (Aabc)t, (Aabc)l, (Aabc), and (Aabc), depends on the triangle inside which the tracked point x falls after each subdivision step. Fin.illy, the local parameterization process can be completed by writing Sabc(x) = (lim Bjbc(X))V bc = Babc(X)V.be. (4.9) In the above equation, Babc is the collection of basis functions at the vertices of vabc. 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 (Aabc)t, (Aabc)l, (Aabc)r and (Aab(c), sums to one. Tli, 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 k-th face in the initial mesh, Eqn.4.9 can be rewritten as Sk(x) = Bk(x)v = 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 post-multiplied by p, selects the vertices vk controlling the k-th smooth triangular patch in the limit surface. If there are t vertices in the initial mesh and r of them control the k-th 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) = (E Bk(x)Ak)p = J(x)p, (4.11) k=l where J, a matrix of dimension (3, 3t), is the collection of basis functions for the corresponding vertices in the initial mesh. Ti, 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. Tli, derivation of the dynamic model from this point is exactly similar in spirit as that of the dynamic Catmull-Clark subdivision model presented in Section 3.2.4. Ti, 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. Ti, velocity of this surface model can be expressed as s(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 Catmull-Clark 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. Ti 1, the kinetic energy of the surface is r 1. T (x)s(x)s(x)x Mp, (4.13) 2 xESO 2 where (using Eqn.4.12) M = fXo t(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. T!i, dissipation energy is 1r 1. F I (x)s (x)s(x)dx = Dp, (4.14) 2 Jxso 2 where D = fXEso y(x)JT(x)J(x)dx is the damping matrix of dimension (3t, 3t). Tli, potential energy of the smooth limit surface can be expressed as 1 U = p'Kp, (4.15) 2 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. Tli, generalized force vector fp can be obtained is a similar fashion described in Section 3.2.4 and is given by fp= I J (x)f(x, t)dx. (4.16) 4.2.3 Ai li ilvel Dynamics Ti, initial mesh of the dynamic modified butterfly subdivision surface model can be subdivided to increase the degrees of freedom for model representation. Ti11 situation is exactly similar to that of the dynamic Catmull-Clark 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). Tlit-, p, expressed as a function of q, can be written as p = (BTB)-IBTq = B'q, (4.17) where Bt = (B 'B)-1B'. Tli, i, ._re, Eqn.4.11 and Eqn.4.12 is modified as s(x) = (J(x)Bt)q, (4.18) and (4.19) s(x, q) = (J(x)Bt)q, Figure 4.7. (a) An initial mesh, and (b) the corresponding limit surface. Ti, domains of the shaded elements in the limit surface are the corresponding triangular faces in the initial mesh. Ti, encircled vertices in (a) are the degrees of freedom for the corresponding element. respectively. Tli, motion equation, explicitly expressed as a function of q, can be written as M~iq + Dq + Kqq = fL, where Mq (4.20) fxst 1(x)(Bt)TJT(x)J(x)Btdx, S1 being the domain defined by the newly obtained subdivided mesh. Ti, derivations of D,, K, and f, are similar. 4.3 Fini, 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 Catmull-Clark 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. Tli, shape (basis) function of this finite element is obtained by smoothing a hat function through repeated application of the modified butterfly subdivision rules. Tli, 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. Tli, elemental matrices can be assembled to obtain the global physical matrices M, D and K, and a numerical solution to the governing second-order differential equation as given by Eqn.2.5 can be obtained using finite element analysis techniques [42]. Tli, same example as in Section 4.2 (refer Fig.4.5) is used to develop the related concepts. Tli, 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. Tli, derivations hold for any element in general. 4.3.1 Elemental .I,--, and Damping matrices Ti, mass matrix for the element given by Sab, (Eqn.4.9) can be written as Mabc = i (x) s Babc(x) (B(x)T Bab (x) }dx. (4.21) However, from Eqn.4.9 it is known that the basis functions corresponding to the vertices in vabe which are stored as entries in Babc are obtained as a limiting process. Ti, -, 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. Tli, smooth triangular patch in the limit surface corresponding to the face abc in the initial mesh is approximated by a triangular mesh with 4J faces obtained after j levels of subdivision of the original triangular face abc (each subdivision step splits one triangular face into 4 triangular faces). Ti, ii the mass matrix can be expressed as 4J Mabc = : ( t(x){BNbc(x) } { Bbc(X) }dx. (4.22) i=1 xEAi Ti, j-th level approximation of the corresponding basis functions can be explicitly evaluated (refer Eqn.4.8 for an expression of B bc). 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 (Ab,)(= (Aabc)m .. (Aabc),(Aabc)1) in Eqn.4.8 encodes how vertices in the 2-neighborhood of the triangular face uvw at level j are related to the vertices in the 2-neighborhood 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. TI,, information stored in (Aabc) can be obtained by tracking the contributions from level j to level 0 recursively. TI [ -. 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. Tli, computations involved in evaluating the integrals in Eqn.4.22 are eliminated by using discrete mass density function which has non-zero values only at the vertex positions of the j-th subdivision level mesh. TI! 1, fi .e, the approximation of the mass matrix for the element can be written as k Mabc= E (v)L{Bbc(v )} b { c(vi)}, (4.23) i=1 where k is the number of vertices in the triangular mesh with 4- faces. This approxi- mation has been found to be very effective and efficient for implementation purposes. Tli, elemental damping matrix can be obtained in a similar fashion. 4.3.2 Elemental Stiffness '.i i7,\ Ti, 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 j-th 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 j-th level of approx- imation can be written as Eab, Eb, = (E b,)t + (Eabt + (Eabc)sp, (4.24) where the subscripts t, st and sp denote tension, stiffness and spring respectively. Tli, expression for the tension energy, which is essentially equivalent to the first order strain (membrane) energy [96], is 1 2 (abc) 2 v -vI = 2I b} (Kbc)IVbc}, (4.25) where kt is a constant, v, and v the 1-th and m-th vertex in the j-th level mesh, are in the 1-neighborhood 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 j-th 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 1 j Eab)St 2kst IV] 2v + v | 2 n = 1vk,,u} (Ku )s{v,}, (4.26) where v{, vi and vi are the three vertices of a triangular face. Ti, summation in- volves three terms corresponding to each triangular face, and is over all the triangular faces in the mesh at the j-th level of subdivision. Fiiilly, a spring energy term is added which is given by 1 (km),,p(I|v v~| 1m,) (EIp = (vE v )- v (abcp 2 vm = 1{v,} (K ){vb}, (4.27) where (kim), is the spring constant, Q is the domain as in Eqn.4.25 and tm is the natural length of the spring connected between v1 and v-. It may be noted that the entries in (K b,) depend on the distance between the connected vertices and hence, (K b,) 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 EJb =c {vb( {(K b)t + k. I (K~).(t+ (K bc) sp Vbc} = {vbc } (Kabca){v bta 1 2 I = {(Ab)vb.}} (Kb)(Ab bc}J = {vb}T {(Ajb)Tc(K 'b)(Abc)}{ }, (fI where (Aabc) and vabe are same as in Eqn.4.8. T, i, f; re, the expression for the elemental stiffness matrix is given by Kabc = (AbT (Kjb)(Ajb). (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 Tli, 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 Catmull-Clark subdivision surface model. Ti, techniques of model subdivision and discrete dynamic equation derivation are also the same. 4.4 Generalization of the Approach Ti, proposed approach of deriving the dynamic framework for modified butter- fly subdivision scheme is generalized for other interpolatory subdivision schemes in the next chapter, Cli.,11i r 5. CHAPTER 5 UNIFIED DYNA I C FRA 1ii;WORK FOR SUBDIVISION SURFACES Dynamic framework for specific subdivision schemes was presented in the last two chapters. In Clilipi r 3, Catmull-Clark subdivision scheme, a representative of the approximating subdivision schemes, is embedded in a physics-based modeling en- vironment. In C(ri _iir 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 Tli, 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. Ti, type of this finite element depends on the subdivision scheme involved. For example, the limit surface generated by the Catmull-Clark subdivision scheme can be expressed as a collection a quadrilateral finite element patches. This approach differs form the one taken in Cliip i, 3 where the limit surface was a collection of normal elements (corresponding to the quadrilateral bicubic B-spline patches away from the extraordinary points) and special elements (corresponding to the n-sided patches near extraordinary vertices of degree n). This concept will be further elaborated in the next few sections. Tli, 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 B-spline 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 n-sided face in the control mesh, there would be a smooth n-sided patch in the limit surface in general. For example, the control mesh (after at most one subdivision step) of the Catmull-Clark 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 Catmull-Clark subdivision scheme, the faces without extraordinary vertices lead to quadrilateral bicubic B-spline patches in the limit, whereas the faces with one extraor- dinary vertex lead to quadrilateral patches whose analytic expressions are difficult to obtain. In Cliipl r 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 n-sided 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. T!i, 1, fi .e, 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. Tli1 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 Cliipi, r 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 Catmull-Clark Subdivision Surfaces A systematic formulation along with the implementation details of the unified dynamic framework for Catmull-Clark subdivision surfaces is presented in this sec- tion. Tli1 key difference between the framework developed in Clilipi r 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. |