ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/sm_del.c
Revision: 3.12
Committed: Thu Jun 10 15:22:22 1999 UTC (24 years, 10 months ago) by gwlarson
Content type: text/plain
Branch: MAIN
Changes since 3.11: +460 -184 lines
Log Message: 
Implemented sample quadtree in place of triangle quadtree
Made geometric predicates more robust
Added #define LORES which utilizes a single precision floating point
  sample array, the default is a double sample array
Added topology DEBUG commands (for DEBUG > 1)
Made code optimizations

File Contents

# User Rev Content
1 gwlarson 3.1 /* Copyright (c) 1998 Silicon Graphics, Inc. */
2    
3     #ifndef lint
4     static char SCCSid[] = "$SunId$ SGI";
5     #endif
6    
7     /*
8     * sm_del.c
9     */
10     #include "standard.h"
11 gwlarson 3.6 #include "sm_flag.h"
12 gwlarson 3.1 #include "sm_list.h"
13     #include "sm_geom.h"
14     #include "sm.h"
15    
16 gwlarson 3.12 #define MAX_EDGE_INIT 100
17     static int Max_edges= MAX_EDGE_INIT;
18 gwlarson 3.6 static EDGE *Edges=NULL;
19 gwlarson 3.1 static int Ecnt=0;
20    
21 gwlarson 3.12 smFree_tri(sm,id,t)
22 gwlarson 3.1 SM *sm;
23     int id;
24 gwlarson 3.12 TRI *t;
25 gwlarson 3.1 {
26     /* Add to the free_list */
27 gwlarson 3.12 T_NEXT_FREE(t) = SM_FREE_TRIS(sm);
28 gwlarson 3.1 SM_FREE_TRIS(sm) = id;
29 gwlarson 3.12 #ifdef DEBUG
30     T_VALID_FLAG(t) = INVALID;
31     #endif
32 gwlarson 3.1 }
33    
34     /* Assumes mesh pointers have been cleaned up appropriately: just deletes from
35     Point location and triangle data structure
36     */
37 gwlarson 3.12 smDelete_tri(sm,t_id,t)
38 gwlarson 3.1 SM *sm;
39     int t_id;
40 gwlarson 3.12 TRI *t;
41 gwlarson 3.1 {
42    
43     /* NOTE: Assumes that a new triangle adjacent to each vertex
44     has been added- before the deletion: replacing
45     the old tri- and therefore dont need to dereference any pointers
46     to id because the vertices can no longer
47     point to tri id as being the first triangle pointer
48     */
49 gwlarson 3.12 SM_CLR_NTH_T_ACTIVE(sm,t_id);
50     smFree_tri(sm,t_id,t);
51 gwlarson 3.9
52 gwlarson 3.6 }
53 gwlarson 3.1
54 gwlarson 3.6 int
55     eNew_edge()
56     {
57     if(!Edges)
58     if(!(Edges = (EDGE *)realloc(NULL,(Max_edges+1)*sizeof(EDGE))))
59     goto memerr;
60    
61     if(Ecnt >= Max_edges)
62     {
63 gwlarson 3.12 #ifdef DEBUG
64 gwlarson 3.6 if(Max_edges > 10000)
65 gwlarson 3.12 error(CONSISTENCY,"Too many edges in vertex loop\n");
66     #endif
67 gwlarson 3.6 Max_edges += 100;
68     if(!(Edges = (EDGE *)realloc(Edges,(Max_edges+1)*sizeof(EDGE))))
69     goto memerr;
70     }
71 gwlarson 3.12 #ifdef DEBUG
72     SET_E_NTH_TRI(Ecnt+1,0,INVALID);
73     SET_E_NTH_TRI(Ecnt+1,1,INVALID);
74     #endif
75 gwlarson 3.6 return(++Ecnt);
76    
77     memerr:
78 gwlarson 3.12 error(SYSTEM,"eNew_edge): Unable to allocate memory");
79 gwlarson 3.1 }
80    
81 gwlarson 3.7 /* Return list of edges defining polygon formed by boundary of triangles
82     adjacent to id. Return set of triangles adjacent to id to delete in delptr
83     */
84 gwlarson 3.5 LIST
85 gwlarson 3.12 *smVertexStar(sm,id)
86 gwlarson 3.1 SM *sm;
87     int id;
88     {
89     TRI *tri,*t_next;
90 gwlarson 3.5 LIST *elist,*end;
91 gwlarson 3.12 int e,t_id,v_next,t_next_id,b_id,v_id,t_last_id;
92 gwlarson 3.1
93     elist = end = NULL;
94 gwlarson 3.7
95 gwlarson 3.1 /* Get the first triangle adjacent to vertex id */
96     t_id = SM_NTH_VERT(sm,id);
97     tri = SM_NTH_TRI(sm,t_id);
98    
99 gwlarson 3.7 e = eNew_edge();
100     /* Get the next vertex on the polygon boundary */
101     v_id = T_WHICH_V(tri,id);
102     b_id = (v_id + 1)%3;
103     /* Create an edge */
104     SET_E_NTH_VERT(e,0,T_NTH_V(tri,b_id));
105 gwlarson 3.6 SET_E_NTH_TRI(e,0,INVALID);
106 gwlarson 3.7 SET_E_NTH_TRI(e,1,T_NTH_NBR(tri,v_id));
107 gwlarson 3.12 v_next = T_NTH_V(tri,(b_id+1)%3);
108 gwlarson 3.7 SET_E_NTH_VERT(e,1,v_next);
109 gwlarson 3.12
110 gwlarson 3.6 elist = add_data_to_circular_list(elist,&end,e);
111 gwlarson 3.1 t_next_id = t_id;
112     t_next = tri;
113 gwlarson 3.12 t_last_id = -1;
114    
115 gwlarson 3.7 /* Create a set to hold all of the triangles for deletion later */
116 gwlarson 3.1
117 gwlarson 3.7 while((t_next_id = T_NTH_NBR(t_next,b_id)) != t_id)
118     {
119     e = eNew_edge();
120 gwlarson 3.12 if(t_last_id != -1)
121     {
122     smDelete_tri(sm,t_last_id,t_next);
123     }
124 gwlarson 3.6 t_next = SM_NTH_TRI(sm,t_next_id);
125 gwlarson 3.12 t_last_id = t_next_id;
126 gwlarson 3.7 SET_E_NTH_VERT(e,0,v_next);
127 gwlarson 3.6 SET_E_NTH_TRI(e,0,INVALID);
128 gwlarson 3.7 v_id = T_WHICH_V(t_next,id);
129     b_id = (v_id + 1)%3;
130     SET_E_NTH_TRI(e,1,T_NTH_NBR(t_next,v_id));
131     v_next = T_NTH_V(t_next,(b_id+1)%3);
132     SET_E_NTH_VERT(e,1,v_next);
133 gwlarson 3.6 elist = add_data_to_circular_list(elist,&end,e);
134 gwlarson 3.12
135 gwlarson 3.1 }
136 gwlarson 3.12 smDelete_tri(sm,t_last_id,t_next);
137     smDelete_tri(sm,t_id,tri);
138     return(elist);
139 gwlarson 3.1 }
140    
141     int
142 gwlarson 3.7 smTriangulate_add_tri(sm,id0,id1,id2,e0,e1,e2ptr)
143 gwlarson 3.1 SM *sm;
144 gwlarson 3.7 int id0,id1,id2,e0,e1,*e2ptr;
145 gwlarson 3.1 {
146 gwlarson 3.12 int t_id,e2;
147     TRI *t;
148 gwlarson 3.1
149 gwlarson 3.12 t_id = smAdd_tri(sm,id0,id1,id2,&t);
150 gwlarson 3.7 if(*e2ptr == 0)
151     {
152     e2 = eNew_edge();
153     SET_E_NTH_VERT(e2,0,id2);
154     SET_E_NTH_VERT(e2,1,id0);
155     }
156     else
157     e2 = *e2ptr;
158     /* set appropriate tri for each edge*/
159     SET_E_NTH_TRI(e0,0,t_id);
160     SET_E_NTH_TRI(e1,0,t_id);
161     SET_E_NTH_TRI(e2,0,t_id);
162 gwlarson 3.12 #ifdef DEBUG
163     #if DEBUG > 1
164     if(E_NTH_TRI(e0,1) == E_NTH_TRI(e0,0) ||
165     E_NTH_TRI(e1,1) == E_NTH_TRI(e1,0) ||
166     E_NTH_TRI(e2,1) == E_NTH_TRI(e2,0))
167     error(CONSISTENCY,"Non manifold\n");
168     #endif
169     #endif
170 gwlarson 3.7 *e2ptr = e2;
171     return(t_id);
172 gwlarson 3.12
173 gwlarson 3.1 }
174    
175 gwlarson 3.12 int
176     smTriangulate_quad(sm,l)
177     SM *sm;
178     LIST *l;
179 gwlarson 3.1 {
180 gwlarson 3.12 int e1,e2,e3,e4,e,t_id,id0,id1,id2,id3;
181     FVECT v0,v1,v2,v3,n;
182     double dp;
183     TRI *tri;
184     LIST *lptr;
185 gwlarson 3.1
186 gwlarson 3.12 #ifdef DEBUG
187     if(LIST_NEXT(LIST_NEXT(LIST_NEXT(LIST_NEXT(l)))) != l)
188     {
189     eputs("smTriangulate_quad(): not a quadrilateral\n");
190     return(FALSE);
191     }
192     eputs("smTriangulate_quad():\n");
193     #endif
194     lptr=l;
195     e1 = (int)LIST_DATA(l);
196     id0 = E_NTH_VERT(e1,0);
197     id1 = E_NTH_VERT(e1,1);
198     VSUB(v0,SM_NTH_WV(sm,id0),SM_VIEW_CENTER(sm));
199     VSUB(v1,SM_NTH_WV(sm,id1),SM_VIEW_CENTER(sm));
200     /* Get v2 */
201     l = LIST_NEXT(l);
202     e2 = (int)LIST_DATA(l);
203     id2 = E_NTH_VERT(e2,1);
204     VSUB(v2,SM_NTH_WV(sm,id2),SM_VIEW_CENTER(sm));
205     /* Get v3 */
206     l = LIST_NEXT(l);
207     e3 = (int)LIST_DATA(l);
208     id3 = E_NTH_VERT(e3,1);
209     VSUB(v3,SM_NTH_WV(sm,id3),SM_VIEW_CENTER(sm));
210     l = LIST_NEXT(l);
211     e4 = (int)LIST_DATA(l);
212     free_list(lptr);
213    
214     VCROSS(n,v0,v2);
215     normalize(n);
216     normalize(v1);
217     dp = DOT(n,v1);
218     e = 0;
219     if(dp > 0)
220     {
221     if(dp >= EV_EPS)
222 gwlarson 3.1 {
223 gwlarson 3.12 normalize(v3);
224     if(DOT(n,v3) < 0)
225     {
226     t_id = smTriangulate_add_tri(sm,id0,id1,id2,e1,e2,&e);
227     e *= -1;
228     t_id = smTriangulate_add_tri(sm,id2,id3,id0,e3,e4,&e);
229     return(TRUE);
230     }
231     }
232 gwlarson 3.1
233 gwlarson 3.12 }
234     #ifdef DEBUG
235     #if DEBUG > 1
236     VCROSS(n,v1,v3);
237     normalize(n);
238     normalize(v0);
239     normalize(v2);
240     dp = DOT(n,v2);
241     if((dp < 0) || (dp < EV_EPS) || (DOT(n,v0) >= 0.0))
242     error(CONSISTENCY,"smTriangulate_quad: cannot triangulate\n");
243     #endif
244     #endif
245     t_id = smTriangulate_add_tri(sm,id1,id2,id3,e2,e3,&e);
246     e *= -1;
247     t_id = smTriangulate_add_tri(sm,id3,id0,id1,e4,e1,&e);
248     return(TRUE);
249 gwlarson 3.1 }
250    
251 gwlarson 3.7 /* Triangulate the polygon defined by plist, and generating vertex p_id.
252     Return list of added triangles in list add_ptr. Returns TRUE if
253     successful, FALSE otherwise. This is NOT a general triangulation routine,
254     assumes polygon star relative to id
255     */
256 gwlarson 3.6
257     int
258 gwlarson 3.12 smTriangulate(sm,id,plist)
259 gwlarson 3.6 SM *sm;
260     int id;
261 gwlarson 3.12 LIST *plist;
262 gwlarson 3.1 {
263 gwlarson 3.6 LIST *l,*prev,*t;
264     FVECT v0,v1,v2,n,p;
265 gwlarson 3.12 int is_tri,cut,t_id,id0,id1,id2,e2,e1,enew;
266     double dp1,dp2;
267 gwlarson 3.6
268 gwlarson 3.7 VSUB(p,SM_NTH_WV(sm,id),SM_VIEW_CENTER(sm));
269 gwlarson 3.12 normalize(p);
270     l = plist;
271    
272 gwlarson 3.7 enew = 0;
273 gwlarson 3.12 cut = is_tri= FALSE;
274     prev = l;
275 gwlarson 3.7 /* get v0,v1 */
276     e1 = (int)LIST_DATA(l);
277     id0 = E_NTH_VERT(e1,0);
278     id1 = E_NTH_VERT(e1,1);
279     VSUB(v0,SM_NTH_WV(sm,id0),SM_VIEW_CENTER(sm));
280 gwlarson 3.12 normalize(v0);
281 gwlarson 3.7 VSUB(v1,SM_NTH_WV(sm,id1),SM_VIEW_CENTER(sm));
282 gwlarson 3.12 normalize(v1);
283 gwlarson 3.6 while(l)
284     {
285     l = LIST_NEXT(l);
286 gwlarson 3.7 /* Get v2 */
287 gwlarson 3.6 e2 = (int)LIST_DATA(l);
288     id2 = E_NTH_VERT(e2,1);
289 gwlarson 3.7 VSUB(v2,SM_NTH_WV(sm,id2),SM_VIEW_CENTER(sm));
290 gwlarson 3.12 normalize(v2);
291 gwlarson 3.6 if(LIST_NEXT(LIST_NEXT(l)) == prev)
292 gwlarson 3.7 {/* Check if have a triangle */
293     is_tri = TRUE;
294     break;
295     }
296     /* determine if v0-v1-v2 is convex:defined clockwise on the sphere-
297     so switch orientation
298     */
299 gwlarson 3.12 VCROSS(n,v0,v2);
300     normalize(n);
301     dp1 = DOT(n,p);
302     if(dp1 <= 0.0)
303 gwlarson 3.6 {
304 gwlarson 3.7 /* test if safe to cut off v0-v1-v2 by testing if p lies outside of
305     triangle v0-v1-v2: if so, because plist is the star polygon around p,
306     the new edge v2-v0 cannot intersect any existing edges
307     */
308 gwlarson 3.12 dp1 = DOT(n,v1);
309     /* Distance from point s to plane is d = |N.p0s|/||N|| */
310     /* sp1 = s-p0 for point on plane p0, a.(b+c) = a.b + a.c */
311     if(dp1 >= EV_EPS)
312     {
313 gwlarson 3.7 /* remove edges e1,e2 and add triangle id0,id1,id2 */
314 gwlarson 3.12 enew = 0;
315     t_id = smTriangulate_add_tri(sm,id0,id1,id2,e1,e2,&enew);
316     cut = TRUE;
317     /* Insert edge enew into the list, reuse prev list element */
318     LIST_NEXT(prev) = LIST_NEXT(l);
319     LIST_DATA(prev) = e1 = -enew;
320     /* If removing head of list- reset plist pointer */
321     if(l== plist)
322     plist = prev;
323     /* free list element for e2 */
324     LIST_NEXT(l)=NULL;
325     free_list(l);
326     l = prev;
327     VCOPY(v1,v2);
328     id1 = id2;
329     continue;
330 gwlarson 3.7 }
331 gwlarson 3.6 }
332 gwlarson 3.7 VCOPY(v0,v1);
333     VCOPY(v1,v2);
334     id0 = id1;
335     id1 = id2;
336     e1 = e2;
337     /* check if gone around circular list without adding any
338     triangles: prevent infinite loop */
339     if(l == plist)
340 gwlarson 3.6 {
341 gwlarson 3.7 if(LIST_NEXT(LIST_NEXT(l)) == prev)
342     {/* Check if have a triangle */
343     is_tri = TRUE;
344     break;
345     }
346 gwlarson 3.12 if(!cut)
347 gwlarson 3.6 break;
348 gwlarson 3.7 cut = FALSE;
349 gwlarson 3.6 }
350 gwlarson 3.7 prev = l;
351 gwlarson 3.6 }
352     if(is_tri)
353     {
354     l = LIST_NEXT(l);
355 gwlarson 3.7 enew = (int)LIST_DATA(l);
356     t_id = smTriangulate_add_tri(sm,id0,id1,id2,e1,e2,&enew);
357 gwlarson 3.6 free_list(l);
358 gwlarson 3.7 }
359 gwlarson 3.6 else
360 gwlarson 3.12 if(!smTriangulate_quad(sm,l))
361     goto Terror;
362     #if 0
363     {int i;
364     eputs("\n\n");
365     for(i=1;i<=Ecnt;i++)
366     fprintf(stderr,"%d verts %d %d tris %d %d\n",
367     i,Edges[i].verts[0],Edges[i].verts[1],
368     Edges[i].tris[0],Edges[i].tris[1]);
369     }
370     #endif
371 gwlarson 3.6
372 gwlarson 3.7 /* Set triangle adjacencies based on edge adjacencies */
373     FOR_ALL_EDGES(enew)
374 gwlarson 3.1 {
375 gwlarson 3.7 id0 = E_NTH_TRI(enew,0);
376     id1 = E_NTH_TRI(enew,1);
377 gwlarson 3.1
378 gwlarson 3.7 e1 = (T_WHICH_V(SM_NTH_TRI(sm,id0),E_NTH_VERT(enew,0))+2)%3;
379     T_NTH_NBR(SM_NTH_TRI(sm,id0),e1) = id1;
380    
381     e2 = (T_WHICH_V(SM_NTH_TRI(sm,id1),E_NTH_VERT(enew,1))+2)%3;
382     T_NTH_NBR(SM_NTH_TRI(sm,id1),e2) = id0;
383 gwlarson 3.1 }
384 gwlarson 3.12 #if 0
385     {int i;
386     eputs("\n\n");
387     for(i=1;i<=Ecnt;i++)
388     fprintf(stderr,"%d verts %d %d tris %d %d\n",
389     i,Edges[i].verts[0],Edges[i].verts[1],
390     Edges[i].tris[0],Edges[i].tris[1]);
391     }
392     #endif
393    
394     #ifdef DEBUG
395     #if DEBUG > 1
396     {
397     TRI *t0,*t1,*n;
398    
399     FOR_ALL_EDGES(enew)
400     {
401     id0 = E_NTH_TRI(enew,0);
402     id1 = E_NTH_TRI(enew,1);
403     t0 = SM_NTH_TRI(sm,id0);
404     t1 = SM_NTH_TRI(sm,id1);
405     if(T_NTH_NBR(t0,0) == T_NTH_NBR(t0,1) ||
406     T_NTH_NBR(t0,1) == T_NTH_NBR(t0,2) ||
407     T_NTH_NBR(t0,0) == T_NTH_NBR(t0,2))
408     error(CONSISTENCY,"Invalid topology\n");
409     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(t0,0))) ||
410     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(t0,1))) ||
411     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(t0,2))))
412     error(CONSISTENCY,"Invalid topology\n");
413     n = SM_NTH_TRI(sm,T_NTH_NBR(t0,0));
414     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
415     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
416     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
417     error(CONSISTENCY,"Invalid topology\n");
418     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
419     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
420     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
421     error(CONSISTENCY,"Invalid topology\n");
422     n = SM_NTH_TRI(sm,T_NTH_NBR(t0,1));
423     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
424     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
425     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
426     error(CONSISTENCY,"Invalid topology\n");
427     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
428     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
429     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
430     error(CONSISTENCY,"Invalid topology\n");
431     n = SM_NTH_TRI(sm,T_NTH_NBR(t0,2));
432     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
433     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
434     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
435     error(CONSISTENCY,"Invalid topology\n");
436     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
437     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
438     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
439     error(CONSISTENCY,"Invalid topology\n");
440    
441     if(T_NTH_NBR(t1,0) == T_NTH_NBR(t1,1) ||
442     T_NTH_NBR(t1,1) == T_NTH_NBR(t1,2) ||
443     T_NTH_NBR(t1,0) == T_NTH_NBR(t1,2))
444     error(CONSISTENCY,"Invalid topology\n");
445     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(t1,0))) ||
446     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(t1,1))) ||
447     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(t1,2))))
448     error(CONSISTENCY,"Invalid topology\n");
449     n = SM_NTH_TRI(sm,T_NTH_NBR(t1,0));
450     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
451     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
452     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
453     error(CONSISTENCY,"Invalid topology\n");
454     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
455     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
456     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
457     error(CONSISTENCY,"Invalid topology\n");
458     n = SM_NTH_TRI(sm,T_NTH_NBR(t1,1));
459     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
460     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
461     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
462     error(CONSISTENCY,"Invalid topology\n");
463     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
464     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
465     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
466     error(CONSISTENCY,"Invalid topology\n");
467     n = SM_NTH_TRI(sm,T_NTH_NBR(t1,2));
468     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
469     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
470     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
471     error(CONSISTENCY,"Invalid topology\n");
472     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
473     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
474     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
475     error(CONSISTENCY,"Invalid topology\n");
476     }
477     }
478     #endif
479     #endif
480 gwlarson 3.7 return(TRUE);
481 gwlarson 3.12
482     Terror:
483     #ifdef DEBUG
484     error(CONSISTENCY,"smTriangulate():Unable to triangulate\n");
485     #endif
486 gwlarson 3.1 }
487    
488     eIn_tri(e,t)
489     int e;
490     TRI *t;
491     {
492    
493     if(T_NTH_V(t,0)==E_NTH_VERT(e,0))
494     return(T_NTH_V(t,1)==E_NTH_VERT(e,1)||T_NTH_V(t,2)==E_NTH_VERT(e,1));
495     else
496     if(T_NTH_V(t,1)==E_NTH_VERT(e,0))
497     return(T_NTH_V(t,0)==E_NTH_VERT(e,1)||T_NTH_V(t,2)==E_NTH_VERT(e,1));
498     else if(T_NTH_V(t,2)==E_NTH_VERT(e,0))
499     return(T_NTH_V(t,0)==E_NTH_VERT(e,1)||T_NTH_V(t,1)==E_NTH_VERT(e,1));
500 gwlarson 3.7
501 gwlarson 3.1 return(FALSE);
502     }
503 gwlarson 3.6
504 gwlarson 3.12
505     void
506     smTris_swap_edge(sm,t_id,t1_id,e,e1,tn_id,tn1_id)
507     SM *sm;
508     int t_id,t1_id;
509     int e,e1;
510     int *tn_id,*tn1_id;
511     {
512     int verts[3],enext,eprev,e1next,e1prev;
513     TRI *n,*ta,*tb,*t,*t1;
514     FVECT p1,p2,p3;
515     int ta_id,tb_id;
516     /* form new diagonal (e relative to t, and e1 relative to t1)
517     defined by quadrilateral formed by t,t1- swap for the opposite diagonal
518     */
519     t = SM_NTH_TRI(sm,t_id);
520     t1 = SM_NTH_TRI(sm,t1_id);
521     enext = (e+1)%3;
522     eprev = (e+2)%3;
523     e1next = (e1+1)%3;
524     e1prev = (e1+2)%3;
525     verts[e] = T_NTH_V(t,e);
526     verts[enext] = T_NTH_V(t,enext);
527     verts[eprev] = T_NTH_V(t1,e1);
528     ta_id = smAdd_tri(sm,verts[0],verts[1],verts[2],&ta);
529     #if 0
530     fprintf(stderr,"Added tri %d %d %d %d\n",ta_id,T_NTH_V(ta,0),
531     T_NTH_V(ta,1), T_NTH_V(ta,2));
532     #endif
533     verts[e1] = T_NTH_V(t1,e1);
534     verts[e1next] = T_NTH_V(t1,e1next);
535     verts[e1prev] = T_NTH_V(t,e);
536     tb_id = smAdd_tri(sm,verts[0],verts[1],verts[2],&tb);
537     #if 0
538     fprintf(stderr,"Added tri %d %d %d %d\n",tb_id,T_NTH_V(tb,0),
539     T_NTH_V(tb,1), T_NTH_V(tb,2));
540     #endif
541     /* set the neighbors */
542     T_NTH_NBR(ta,e) = T_NTH_NBR(t1,e1next);
543     T_NTH_NBR(tb,e1) = T_NTH_NBR(t,enext);
544     T_NTH_NBR(ta,enext)= tb_id;
545     T_NTH_NBR(tb,e1next)= ta_id;
546     T_NTH_NBR(ta,eprev)=T_NTH_NBR(t,eprev);
547     T_NTH_NBR(tb,e1prev)=T_NTH_NBR(t1,e1prev);
548    
549     /* Reset neighbor pointers of original neighbors */
550     n = SM_NTH_TRI(sm,T_NTH_NBR(t,enext));
551     T_NTH_NBR(n,T_NTH_NBR_PTR(t_id,n)) = tb_id;
552     n = SM_NTH_TRI(sm,T_NTH_NBR(t,eprev));
553     T_NTH_NBR(n,T_NTH_NBR_PTR(t_id,n)) = ta_id;
554     n = SM_NTH_TRI(sm,T_NTH_NBR(t1,e1next));
555     T_NTH_NBR(n,T_NTH_NBR_PTR(t1_id,n)) = ta_id;
556     n = SM_NTH_TRI(sm,T_NTH_NBR(t1,e1prev));
557     T_NTH_NBR(n,T_NTH_NBR_PTR(t1_id,n)) = tb_id;
558    
559     smDelete_tri(sm,t_id,t);
560     smDelete_tri(sm,t1_id,t1);
561    
562     #ifdef DEBUG
563     #if DEBUG > 1
564     if(T_NTH_NBR(ta,0) == T_NTH_NBR(ta,1) ||
565     T_NTH_NBR(ta,1) == T_NTH_NBR(ta,2) ||
566     T_NTH_NBR(ta,0) == T_NTH_NBR(ta,2))
567     error(CONSISTENCY,"Invalid topology\n");
568     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(ta,0))) ||
569     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(ta,1))) ||
570     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(ta,2))))
571     error(CONSISTENCY,"Invalid topology\n");
572     n = SM_NTH_TRI(sm,T_NTH_NBR(ta,0));
573     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
574     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
575     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
576     error(CONSISTENCY,"Invalid topology\n");
577     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
578     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
579     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
580     error(CONSISTENCY,"Invalid topology\n");
581     n = SM_NTH_TRI(sm,T_NTH_NBR(ta,1));
582     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
583     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
584     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
585     error(CONSISTENCY,"Invalid topology\n");
586     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
587     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
588     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
589     error(CONSISTENCY,"Invalid topology\n");
590     n = SM_NTH_TRI(sm,T_NTH_NBR(ta,2));
591     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
592     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
593     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
594     error(CONSISTENCY,"Invalid topology\n");
595     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
596     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
597     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
598     error(CONSISTENCY,"Invalid topology\n");
599     if(T_NTH_NBR(ta,0) == T_NTH_NBR(ta,1) ||
600     T_NTH_NBR(ta,1) == T_NTH_NBR(ta,2) ||
601     T_NTH_NBR(ta,0) == T_NTH_NBR(ta,2))
602     error(CONSISTENCY,"Invalid topology\n");
603    
604     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(tb,0))) ||
605     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(tb,1))) ||
606     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(tb,2))))
607     error(CONSISTENCY,"Invalid topology\n");
608     n = SM_NTH_TRI(sm,T_NTH_NBR(tb,0));
609     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
610     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
611     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
612     error(CONSISTENCY,"Invalid topology\n");
613     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
614     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
615     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
616     error(CONSISTENCY,"Invalid topology\n");
617     n = SM_NTH_TRI(sm,T_NTH_NBR(tb,1));
618     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
619     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
620     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
621     error(CONSISTENCY,"Invalid topology\n");
622     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
623     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
624     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
625     error(CONSISTENCY,"Invalid topology\n");
626     n = SM_NTH_TRI(sm,T_NTH_NBR(tb,2));
627     if(T_NTH_NBR(n,0) == T_NTH_NBR(n,1) ||
628     T_NTH_NBR(n,1) == T_NTH_NBR(n,2) ||
629     T_NTH_NBR(n,0) == T_NTH_NBR(n,2))
630     error(CONSISTENCY,"Invalid topology\n");
631     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,0))) ||
632     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,1))) ||
633     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(n,2))))
634     error(CONSISTENCY,"Invalid topology\n");
635     #endif
636     #endif
637     *tn_id = ta_id;
638     *tn1_id = tb_id;
639    
640     return;
641     }
642    
643 gwlarson 3.7 /* Test the new set of triangles for Delaunay condition. 'Edges' contains
644     all of the new edges added. The CCW triangle assoc with each edge is
645     tested against the opposite vertex of the CW triangle. If the vertex
646     lies inside the circle defined by the CCW triangle- the edge is swapped
647     for the opposite diagonal
648     */
649 gwlarson 3.12 smFixEdges(sm)
650 gwlarson 3.1 SM *sm;
651     {
652 gwlarson 3.6 int e,t0_id,t1_id,e_new,e0,e1,e0_next,e1_next;
653     int i,v0_id,v1_id,v2_id,p_id,t0_nid,t1_nid;
654 gwlarson 3.1 FVECT v0,v1,v2,p,np,v;
655 gwlarson 3.7 TRI *t0,*t1;
656 gwlarson 3.5
657 gwlarson 3.1 FOR_ALL_EDGES(e)
658     {
659 gwlarson 3.6 t0_id = E_NTH_TRI(e,0);
660     t1_id = E_NTH_TRI(e,1);
661 gwlarson 3.12 #ifdef DEBUG
662 gwlarson 3.6 if((t0_id==INVALID) || (t1_id==INVALID))
663 gwlarson 3.12 error(CONSISTENCY,"smFix_edges: Unassigned edge nbr\n");
664 gwlarson 3.1 #endif
665 gwlarson 3.7 t0 = SM_NTH_TRI(sm,t0_id);
666     t1 = SM_NTH_TRI(sm,t1_id);
667     e0 = T_NTH_NBR_PTR(t1_id,t0);
668     e1 = T_NTH_NBR_PTR(t0_id,t1);
669    
670 gwlarson 3.6 v0_id = E_NTH_VERT(e,0);
671     v1_id = E_NTH_VERT(e,1);
672 gwlarson 3.7 v2_id = T_NTH_V(t0,e0);
673     p_id = T_NTH_V(t1,e1);
674 gwlarson 3.1
675 gwlarson 3.11 smDir(sm,v0,v0_id);
676     smDir(sm,v1,v1_id);
677     smDir(sm,v2,v2_id);
678 gwlarson 3.1
679 gwlarson 3.12 VSUB(p,SM_NTH_WV(sm,p_id),SM_VIEW_CENTER(sm));
680     normalize(p);
681     if(pt_in_cone(p,v2,v1,v0))
682 gwlarson 3.1 {
683 gwlarson 3.12 smTris_swap_edge(sm,t0_id,t1_id,e0,e1,&t0_nid,&t1_nid);
684 gwlarson 3.7 /* Adjust the triangle pointers of the remaining edges to be
685     processed
686     */
687 gwlarson 3.1 FOR_ALL_EDGES_FROM(e,i)
688     {
689 gwlarson 3.6 if(E_NTH_TRI(i,0)==t0_id || E_NTH_TRI(i,0)==t1_id)
690 gwlarson 3.1 {
691 gwlarson 3.6 if(eIn_tri(i,SM_NTH_TRI(sm,t0_nid)))
692     SET_E_NTH_TRI(i,0,t0_nid);
693 gwlarson 3.1 else
694 gwlarson 3.6 SET_E_NTH_TRI(i,0,t1_nid);
695 gwlarson 3.1 }
696    
697 gwlarson 3.6 if(E_NTH_TRI(i,1)==t0_id || E_NTH_TRI(i,1)==t1_id)
698 gwlarson 3.1 {
699 gwlarson 3.6 if(eIn_tri(i,SM_NTH_TRI(sm,t0_nid)))
700     SET_E_NTH_TRI(i,1,t0_nid);
701 gwlarson 3.1 else
702 gwlarson 3.6 SET_E_NTH_TRI(i,1,t1_nid);
703 gwlarson 3.1 }
704 gwlarson 3.12 #ifdef DEBUG
705     if(E_NTH_TRI(i,1) == E_NTH_TRI(i,0) )
706     error(CONSISTENCY,"invalid edge\n");
707     #endif
708 gwlarson 3.1 }
709 gwlarson 3.6 t0_id = t0_nid;
710     t1_id = t1_nid;
711 gwlarson 3.1 e_new = eNew_edge();
712 gwlarson 3.6 SET_E_NTH_VERT(e_new,0,p_id);
713     SET_E_NTH_VERT(e_new,1,v2_id);
714     SET_E_NTH_TRI(e_new,0,t0_id);
715     SET_E_NTH_TRI(e_new,1,t1_id);
716 gwlarson 3.12 #ifdef DEBUG
717     if(E_NTH_TRI(i,1) == E_NTH_TRI(i,0) )
718     error(CONSISTENCY,"invalid edge\n");
719     #endif
720 gwlarson 3.1 }
721     }
722     }
723    
724 gwlarson 3.12
725     smDelete_samp(sm,s_id)
726     SM *sm;
727     int s_id;
728     {
729     QUADTREE qt;
730     OBJECT *os;
731    
732     /* Mark as free */
733     smUnalloc_samp(sm,s_id);
734    
735     #ifdef DEBUG
736     SM_NTH_VERT(sm,s_id) = INVALID;
737     /* fprintf(stderr,"deleting sample %d\n",s_id); */
738     #endif
739     /* remove from its set */
740     qt = SM_S_NTH_QT(sm,s_id);
741     os = qtqueryset(qt);
742     deletuelem(os, s_id); /* delete obj from unsorted os, no questions */
743     }
744 gwlarson 3.7 /* Remove vertex "id" from the mesh- and retriangulate the resulting
745     hole: Returns TRUE if successful, FALSE otherwise.
746     */
747 gwlarson 3.1 int
748 gwlarson 3.7 smRemoveVertex(sm,id)
749 gwlarson 3.1 SM *sm;
750     int id;
751     {
752 gwlarson 3.12 LIST *b_list;
753 gwlarson 3.7 /* generate list of edges that form the boundary of the
754 gwlarson 3.12 polygon formed by the triangles adjacent to vertex 'id'*/
755     b_list = smVertexStar(sm,id);
756     #if 0
757     {int i;
758     eputs("\n\n");
759     for(i=1;i<=Ecnt;i++)
760     fprintf(stderr,"%d verts %d %d tris %d %d\n",
761     i,Edges[i].verts[0],Edges[i].verts[1],
762     Edges[i].tris[0],Edges[i].tris[1]);
763     }
764     #endif
765    
766 gwlarson 3.7 /* Triangulate polygonal hole */
767 gwlarson 3.12 smTriangulate(sm,id,b_list);
768    
769     /* Fix up new triangles to be Delaunay*/
770 gwlarson 3.7
771 gwlarson 3.12 smFixEdges(sm);
772     smDelete_samp(sm,id);
773     eClear_edges();
774 gwlarson 3.1 return(TRUE);
775     }
776    
777    
778    
779 gwlarson 3.8
780 gwlarson 3.1
781    
782    
783    
784    
785    
786    
787    
788    
789    
790