ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/sm_del.c
(Generate patch)

Comparing ray/src/hd/sm_del.c (file contents):
Revision 3.6 by gwlarson, Tue Oct 6 18:16:53 1998 UTC vs.
Revision 3.14 by greg, Wed Apr 23 00:52:34 2003 UTC

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

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines