progmesh.C
Go to the documentation of this file.
1 /*
2  * Progressive Mesh type Polygon Reduction Algorithm
3  * by Stan Melax (c) 1998
4  * Permission to use any of this code wherever you want is granted..
5  * Although, please do acknowledge authorship if appropriate.
6  *
7  * See the header file progmesh.h for a description of this module
8  */
9 
10 #include <stdio.h>
11 #include <math.h>
12 #include <stdlib.h>
13 #include <assert.h>
14 //#include <windows.h>
15 
16 #include "vector.h"
17 #include "list.h"
18 #include "progmesh.h"
19 
20 #define min(x,y) (((x) <= (y)) ? (x) : (y))
21 #define max(x,y) (((x) >= (y)) ? (x) : (y))
22 
23 
24 /*
25  * For the polygon reduction algorithm we use data structures
26  * that contain a little bit more information than the usual
27  * indexed face set type of data structure.
28  * From a vertex we wish to be able to quickly get the
29  * neighboring faces and vertices.
30  */
31 class Triangle;
32 class Vertex;
33 
34 class Triangle {
35  public:
36  Vertex * vertex[3]; // the 3 points that make this tri
37  Vector normal; // unit vector othogonal to this face
38  Triangle(Vertex *v0,Vertex *v1,Vertex *v2);
39  ~Triangle();
40  void ComputeNormal();
41  void ReplaceVertex(Vertex *vold,Vertex *vnew);
42  int HasVertex(Vertex *v);
43 };
44 class Vertex {
45  public:
46  Vector position; // location of point in euclidean space
47  int id; // place of vertex in original list
48  List<Vertex *> neighbor; // adjacent vertices
49  List<Triangle *> face; // adjacent triangles
50  float objdist; // cached cost of collapsing edge
51  Vertex * collapse; // candidate vertex for collapse
52  Vertex(Vector v,int _id);
53  ~Vertex();
54  void RemoveIfNonNeighbor(Vertex *n);
55 };
56 List<Vertex *> vertices;
57 List<Triangle *> triangles;
58 
59 
60 Triangle::Triangle(Vertex *v0,Vertex *v1,Vertex *v2){
61  assert(v0!=v1 && v1!=v2 && v2!=v0);
62  vertex[0]=v0;
63  vertex[1]=v1;
64  vertex[2]=v2;
65  ComputeNormal();
66  triangles.Add(this);
67  for(int i=0;i<3;i++) {
68  vertex[i]->face.Add(this);
69  for(int j=0;j<3;j++) if(i!=j) {
70  vertex[i]->neighbor.AddUnique(vertex[j]);
71  }
72  }
73 }
74 Triangle::~Triangle(){
75  int i;
76  triangles.Remove(this);
77  for(i=0;i<3;i++) {
78  if(vertex[i]) vertex[i]->face.Remove(this);
79  }
80  for(i=0;i<3;i++) {
81  int i2 = (i+1)%3;
82  if(!vertex[i] || !vertex[i2]) continue;
83  vertex[i ]->RemoveIfNonNeighbor(vertex[i2]);
84  vertex[i2]->RemoveIfNonNeighbor(vertex[i ]);
85  }
86 }
87 int Triangle::HasVertex(Vertex *v) {
88  return (v==vertex[0] ||v==vertex[1] || v==vertex[2]);
89 }
90 void Triangle::ComputeNormal(){
91  Vector v0=vertex[0]->position;
92  Vector v1=vertex[1]->position;
93  Vector v2=vertex[2]->position;
94  normal = (v1-v0)*(v2-v1);
95  if(magnitude(normal)==0)return;
96  normal = normalize(normal);
97 }
98 void Triangle::ReplaceVertex(Vertex *vold,Vertex *vnew) {
99  assert(vold && vnew);
100  assert(vold==vertex[0] || vold==vertex[1] || vold==vertex[2]);
101  assert(vnew!=vertex[0] && vnew!=vertex[1] && vnew!=vertex[2]);
102  if(vold==vertex[0]){
103  vertex[0]=vnew;
104  }
105  else if(vold==vertex[1]){
106  vertex[1]=vnew;
107  }
108  else {
109  assert(vold==vertex[2]);
110  vertex[2]=vnew;
111  }
112  int i;
113  vold->face.Remove(this);
114  assert(!vnew->face.Contains(this));
115  vnew->face.Add(this);
116  for(i=0;i<3;i++) {
117  vold->RemoveIfNonNeighbor(vertex[i]);
118  vertex[i]->RemoveIfNonNeighbor(vold);
119  }
120  for(i=0;i<3;i++) {
121  assert(vertex[i]->face.Contains(this)==1);
122  for(int j=0;j<3;j++) if(i!=j) {
123  vertex[i]->neighbor.AddUnique(vertex[j]);
124  }
125  }
126  ComputeNormal();
127 }
128 
129 Vertex::Vertex(Vector v,int _id) {
130  position =v;
131  id=_id;
132  vertices.Add(this);
133 }
134 
135 Vertex::~Vertex(){
136  assert(face.num==0);
137  while(neighbor.num) {
138  neighbor[0]->neighbor.Remove(this);
139  neighbor.Remove(neighbor[0]);
140  }
141  vertices.Remove(this);
142 }
143 void Vertex::RemoveIfNonNeighbor(Vertex *n) {
144  // removes n from neighbor list if n isn't a neighbor.
145  if(!neighbor.Contains(n)) return;
146  for(int i=0;i<face.num;i++) {
147  if(face[i]->HasVertex(n)) return;
148  }
149  neighbor.Remove(n);
150 }
151 
152 
153 float ComputeEdgeCollapseCost(Vertex *u,Vertex *v) {
154  // if we collapse edge uv by moving u to v then how
155  // much different will the model change, i.e. how much "error".
156  // Texture, vertex normal, and border vertex code was removed
157  // to keep this demo as simple as possible.
158  // The method of determining cost was designed in order
159  // to exploit small and coplanar regions for
160  // effective polygon reduction.
161  // Is is possible to add some checks here to see if "folds"
162  // would be generated. i.e. normal of a remaining face gets
163  // flipped. I never seemed to run into this problem and
164  // therefore never added code to detect this case.
165  int i;
166  float edgelength = magnitude(v->position - u->position);
167  float curvature=0;
168 
169  // find the "sides" triangles that are on the edge uv
170  List<Triangle *> sides;
171  for(i=0;i<u->face.num;i++) {
172  if(u->face[i]->HasVertex(v)){
173  sides.Add(u->face[i]);
174  }
175  }
176  // use the triangle facing most away from the sides
177  // to determine our curvature term
178  for(i=0;i<u->face.num;i++) {
179  float mincurv=1; // curve for face i and closer side to it
180  for(int j=0;j<sides.num;j++) {
181  // use dot product of face normals. '^'
182  // defined in vector
183  float dotprod = u->face[i]->normal ^ sides[j]->normal;
184  mincurv = min(mincurv,(1-dotprod)/2.0f);
185  }
186  curvature = max(curvature,mincurv);
187  }
188  // the more coplanar the lower the curvature term
189  return edgelength * curvature;
190 }
191 
192 void ComputeEdgeCostAtVertex(Vertex *v) {
193  // compute the edge collapse cost for all edges that start
194  // from vertex v. Since we are only interested in reducing
195  // the object by selecting the min cost edge at each step, we
196  // only cache the cost of the least cost edge at this vertex
197  // (in member variable collapse) as well as the value of the
198  // cost (in member variable objdist).
199  if(v->neighbor.num==0) {
200  // v doesn't have neighbors so it costs nothing to collapse
201  v->collapse=NULL;
202  v->objdist=-0.01f;
203  return;
204  }
205  v->objdist = 1000000;
206  v->collapse=NULL;
207  // search all neighboring edges for "least cost" edge
208  for(int i=0;i<v->neighbor.num;i++) {
209  float dist;
210  dist = ComputeEdgeCollapseCost(v,v->neighbor[i]);
211  if(dist<v->objdist) {
212  // candidate for edge collapse
213  v->collapse=v->neighbor[i];
214  // cost of the collapse
215  v->objdist=dist;
216  }
217  }
218 }
219 void ComputeAllEdgeCollapseCosts() {
220  // For all the edges, compute the difference it would make
221  // to the model if it was collapsed. The least of these
222  // per vertex is cached in each vertex object.
223  for(int i=0;i<vertices.num;i++) {
224  ComputeEdgeCostAtVertex(vertices[i]);
225  }
226 }
227 
228 void Collapse(Vertex *u,Vertex *v){
229  // Collapse the edge uv by moving vertex u onto v
230  // Actually remove tris on uv, then update tris that
231  // have u to have v, and then remove u.
232  if(!v) {
233  // u is a vertex all by itself so just delete it
234  delete u;
235  return;
236  }
237  int i;
238  List<Vertex *>tmp;
239  // make tmp a list of all the neighbors of u
240  for(i=0;i<u->neighbor.num;i++) {
241  tmp.Add(u->neighbor[i]);
242  }
243  // delete triangles on edge uv:
244  for(i=u->face.num-1;i>=0;i--) {
245  if(u->face[i]->HasVertex(v)) {
246  delete(u->face[i]);
247  }
248  }
249  // update remaining triangles to have v instead of u
250  for(i=u->face.num-1;i>=0;i--) {
251  u->face[i]->ReplaceVertex(u,v);
252  }
253  delete u;
254  // recompute the edge collapse costs for neighboring vertices
255  for(i=0;i<tmp.num;i++) {
256  ComputeEdgeCostAtVertex(tmp[i]);
257  }
258 }
259 
260 void AddVertex(List<Vector> &vert){
261  for(int i=0;i<vert.num;i++) {
262  new Vertex(vert[i],i);
263  }
264 }
265 void AddFaces(List<tridata> &tri){
266  for(int i=0;i<tri.num;i++) {
267  new Triangle(
268  vertices[tri[i].v[0]],
269  vertices[tri[i].v[1]],
270  vertices[tri[i].v[2]] );
271  }
272 }
273 
274 Vertex *MinimumCostEdge(){
275  // Find the edge that when collapsed will affect model the least.
276  // This funtion actually returns a Vertex, the second vertex
277  // of the edge (collapse candidate) is stored in the vertex data.
278  // Serious optimization opportunity here: this function currently
279  // does a sequential search through an unsorted list :-(
280  // Our algorithm could be O(n*lg(n)) instead of O(n*n)
281  Vertex *mn=vertices[0];
282  for(int i=0;i<vertices.num;i++) {
283  if(vertices[i]->objdist < mn->objdist) {
284  mn = vertices[i];
285  }
286  }
287  return mn;
288 }
289 
290 void ProgressiveMesh(List<Vector> &vert, List<tridata> &tri,
291  List<int> &map, List<int> &permutation)
292 {
293  AddVertex(vert); // put input data into our data structures
294  AddFaces(tri);
295  ComputeAllEdgeCollapseCosts(); // cache all edge collapse costs
296  permutation.SetSize(vertices.num); // allocate space
297  map.SetSize(vertices.num); // allocate space
298  // reduce the object down to nothing:
299  while(vertices.num > 0) {
300  // get the next vertex to collapse
301  Vertex *mn = MinimumCostEdge();
302  // keep track of this vertex, i.e. the collapse ordering
303  permutation[mn->id]=vertices.num-1;
304  // keep track of vertex to which we collapse to
305  map[vertices.num-1] = (mn->collapse)?mn->collapse->id:-1;
306  // Collapse this edge
307  Collapse(mn,mn->collapse);
308  }
309  // reorder the map list based on the collapse ordering
310  for(int i=0;i<map.num;i++) {
311  map[i] = (map[i]==-1)?0:permutation[map[i]];
312  }
313  // The caller of this function should reorder their vertices
314  // according to the returned "permutation".
315 }
graph_traits< Graph >::vertex_descriptor Vertex
Definition: SloanRenumber.C:72
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
quaternion normalize(const quaternion &q)
Return the normalized (unit) quaternion of the given quaternion.
Definition: quaternionI.H:609
labelList f(nPoints)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
label n