FreeFOAM The Cross-Platform CFD Toolkit
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. '^' defined in vector
182  float dotprod = u->face[i]->normal ^ sides[j]->normal;
183  mincurv = min(mincurv,(1-dotprod)/2.0f);
184  }
185  curvature = max(curvature,mincurv);
186  }
187  // the more coplanar the lower the curvature term
188  return edgelength * curvature;
189 }
190 
191 void ComputeEdgeCostAtVertex(Vertex *v) {
192  // compute the edge collapse cost for all edges that start
193  // from vertex v. Since we are only interested in reducing
194  // the object by selecting the min cost edge at each step, we
195  // only cache the cost of the least cost edge at this vertex
196  // (in member variable collapse) as well as the value of the
197  // cost (in member variable objdist).
198  if(v->neighbor.num==0) {
199  // v doesn't have neighbors so it costs nothing to collapse
200  v->collapse=NULL;
201  v->objdist=-0.01f;
202  return;
203  }
204  v->objdist = 1000000;
205  v->collapse=NULL;
206  // search all neighboring edges for "least cost" edge
207  for(int i=0;i<v->neighbor.num;i++) {
208  float dist;
209  dist = ComputeEdgeCollapseCost(v,v->neighbor[i]);
210  if(dist<v->objdist) {
211  v->collapse=v->neighbor[i]; // candidate for edge collapse
212  v->objdist=dist; // cost of the collapse
213  }
214  }
215 }
216 void ComputeAllEdgeCollapseCosts() {
217  // For all the edges, compute the difference it would make
218  // to the model if it was collapsed. The least of these
219  // per vertex is cached in each vertex object.
220  for(int i=0;i<vertices.num;i++) {
221  ComputeEdgeCostAtVertex(vertices[i]);
222  }
223 }
224 
225 void Collapse(Vertex *u,Vertex *v){
226  // Collapse the edge uv by moving vertex u onto v
227  // Actually remove tris on uv, then update tris that
228  // have u to have v, and then remove u.
229  if(!v) {
230  // u is a vertex all by itself so just delete it
231  delete u;
232  return;
233  }
234  int i;
235  List<Vertex *>tmp;
236  // make tmp a list of all the neighbors of u
237  for(i=0;i<u->neighbor.num;i++) {
238  tmp.Add(u->neighbor[i]);
239  }
240  // delete triangles on edge uv:
241  for(i=u->face.num-1;i>=0;i--) {
242  if(u->face[i]->HasVertex(v)) {
243  delete(u->face[i]);
244  }
245  }
246  // update remaining triangles to have v instead of u
247  for(i=u->face.num-1;i>=0;i--) {
248  u->face[i]->ReplaceVertex(u,v);
249  }
250  delete u;
251  // recompute the edge collapse costs for neighboring vertices
252  for(i=0;i<tmp.num;i++) {
253  ComputeEdgeCostAtVertex(tmp[i]);
254  }
255 }
256 
257 void AddVertex(List<Vector> &vert){
258  for(int i=0;i<vert.num;i++) {
259  new Vertex(vert[i],i);
260  }
261 }
262 void AddFaces(List<tridata> &tri){
263  for(int i=0;i<tri.num;i++) {
264  new Triangle(
265  vertices[tri[i].v[0]],
266  vertices[tri[i].v[1]],
267  vertices[tri[i].v[2]] );
268  }
269 }
270 
271 Vertex *MinimumCostEdge(){
272  // Find the edge that when collapsed will affect model the least.
273  // This funtion actually returns a Vertex, the second vertex
274  // of the edge (collapse candidate) is stored in the vertex data.
275  // Serious optimization opportunity here: this function currently
276  // does a sequential search through an unsorted list :-(
277  // Our algorithm could be O(n*lg(n)) instead of O(n*n)
278  Vertex *mn=vertices[0];
279  for(int i=0;i<vertices.num;i++) {
280  if(vertices[i]->objdist < mn->objdist) {
281  mn = vertices[i];
282  }
283  }
284  return mn;
285 }
286 
287 void ProgressiveMesh(List<Vector> &vert, List<tridata> &tri,
288  List<int> &map, List<int> &permutation)
289 {
290  AddVertex(vert); // put input data into our data structures
291  AddFaces(tri);
292  ComputeAllEdgeCollapseCosts(); // cache all edge collapse costs
293  permutation.SetSize(vertices.num); // allocate space
294  map.SetSize(vertices.num); // allocate space
295  // reduce the object down to nothing:
296  while(vertices.num > 0) {
297  // get the next vertex to collapse
298  Vertex *mn = MinimumCostEdge();
299  // keep track of this vertex, i.e. the collapse ordering
300  permutation[mn->id]=vertices.num-1;
301  // keep track of vertex to which we collapse to
302  map[vertices.num-1] = (mn->collapse)?mn->collapse->id:-1;
303  // Collapse this edge
304  Collapse(mn,mn->collapse);
305  }
306  // reorder the map list based on the collapse ordering
307  for(int i=0;i<map.num;i++) {
308  map[i] = (map[i]==-1)?0:permutation[map[i]];
309  }
310  // The caller of this function should reorder their vertices
311  // according to the returned "permutation".
312 }
313 
314 
315 // ************************ vim: set sw=4 sts=4 et: ************************ //