--- ray/src/hd/sm_del.c 1998/09/16 18:16:28 3.5 +++ ray/src/hd/sm_del.c 1998/10/06 18:16:53 3.6 @@ -8,70 +8,57 @@ static char SCCSid[] = "$SunId$ SGI"; * sm_del.c */ #include "standard.h" - +#include "sm_flag.h" #include "sm_list.h" #include "sm_geom.h" +#include "sm_qtree.h" +#include "sm_stree.h" #include "sm.h" -static EDGE Edges[MAX_EDGES]; +static int Max_edges=200; +static EDGE *Edges=NULL; static int Ecnt=0; -int -remove_tri(qtptr,fptr,t_id) +#define remove_tri_compress remove_tri +remove_tri(qtptr,fptr,tptr) QUADTREE *qtptr; int *fptr; - int t_id; + int *tptr; { int n; - if(QT_IS_EMPTY(*qtptr)) - return(FALSE); - /* remove id from set */ - else - { - if(!qtinset(*qtptr,t_id)) - return(FALSE); - n = QT_SET_CNT(qtqueryset(*qtptr))-1; - *qtptr = qtdelelem(*qtptr,t_id); - if(n == 0) - (*fptr) |= QT_COMPRESS; - if(!QT_FLAG_FILL_TRI(*fptr)) - (*fptr)++; - } - return(TRUE); -} + if(QT_IS_EMPTY(*qtptr)) + return; + if(QT_LEAF_IS_FLAG(*qtptr)) + return; -int -remove_tri_compress(qtptr,q0,q1,q2,t0,t1,t2,n,arg,t_id) -QUADTREE *qtptr; -FVECT q0,q1,q2; -FVECT t0,t1,t2; -int n; -int *arg; -int t_id; -{ - int f = 0; - /* NOTE compress */ - return(remove_tri(qtptr,&f,t_id)); + n = QT_SET_CNT(qtqueryset(*qtptr))-1; + *qtptr = qtdelelem(*qtptr,*tptr); + if(n == 0) + (*fptr) |= QT_COMPRESS; + if(!QT_FLAG_FILL_TRI(*fptr)) + (*fptr)++; } -smLocator_remove_tri(sm,t_id) +smLocator_remove_tri(sm,t_id,v0_id,v1_id,v2_id) SM *sm; int t_id; +int v0_id,v1_id,v2_id; { STREE *st; - TRI *t; FVECT v0,v1,v2; - + st = SM_LOCATOR(sm); - t = SM_NTH_TRI(sm,t_id); + VSUB(v0,SM_NTH_WV(sm,v0_id),SM_VIEW_CENTER(sm)); + VSUB(v1,SM_NTH_WV(sm,v1_id),SM_VIEW_CENTER(sm)); + VSUB(v2,SM_NTH_WV(sm,v2_id),SM_VIEW_CENTER(sm)); - VSUB(v0,SM_T_NTH_WV(sm,t,0),SM_VIEW_CENTER(sm)); - VSUB(v1,SM_T_NTH_WV(sm,t,1),SM_VIEW_CENTER(sm)); - VSUB(v2,SM_T_NTH_WV(sm,t,2),SM_VIEW_CENTER(sm)); - stApply_to_tri(st,v0,v1,v2,remove_tri,remove_tri_compress,t_id,NULL); + qtClearAllFlags(); + + stApply_to_tri(st,v0,v1,v2,remove_tri,remove_tri_compress,&t_id); + } smFree_tri(sm,id) @@ -95,7 +82,6 @@ SM *sm; int t_id; { - /* NOTE: Assumes that a new triangle adjacent to each vertex has been added- before the deletion: replacing the old tri- and therefore dont need to dereference any pointers @@ -104,73 +90,97 @@ int t_id; */ if(!SM_IS_NTH_T_BASE(sm,t_id)) { - SM_NUM_TRIS(sm)--; + SM_SAMPLE_TRIS(sm)--; if(SM_IS_NTH_T_NEW(sm,t_id)) smNew_tri_cnt--; } smClear_tri_flags(sm,t_id); smFree_tri(sm,t_id); - } +int +eNew_edge() +{ + if(!Edges) + if(!(Edges = (EDGE *)realloc(NULL,(Max_edges+1)*sizeof(EDGE)))) + goto memerr; + + if(Ecnt >= Max_edges) + { + if(Max_edges > 10000) + error(CONSISTENCY,"Too many edges in vertex loop\n"); + Max_edges += 100; + if(!(Edges = (EDGE *)realloc(Edges,(Max_edges+1)*sizeof(EDGE)))) + goto memerr; + } + return(++Ecnt); + +memerr: + error(SYSTEM,"eNew_edge(): Unable to allocate memory"); +} + LIST -*smVertex_star_polygon(sm,id,del_set) +*smVertex_star_polygon(sm,id,delptr) SM *sm; int id; -OBJECT *del_set; +QUADTREE *delptr; { TRI *tri,*t_next; LIST *elist,*end; int t_id,v_next,t_next_id; int e; + OBJECT del_set[2]; elist = end = NULL; /* Get the first triangle adjacent to vertex id */ t_id = SM_NTH_VERT(sm,id); tri = SM_NTH_TRI(sm,t_id); - - if((e = eNew_edge()) == SM_INVALID) - { -#ifdef DEBUG - eputs("smVertex_star_polygon():Too many edges\n"); -#endif - return(NULL); - } - elist = add_data_to_circular_list(elist,&end,e); + if((e = eNew_edge()) == INVALID) + return(NULL); + v_next = (T_WHICH_V(tri,id)+1)%3; SET_E_NTH_VERT(e,0,T_NTH_V(tri,v_next)); - SET_E_NTH_TRI(e,0,SM_INVALID); + SET_E_NTH_TRI(e,0,INVALID); SET_E_NTH_TRI(e,1,T_NTH_NBR(tri,v_next)); v_next = (T_WHICH_V(tri,id)+2)%3; SET_E_NTH_VERT(e,1,T_NTH_V(tri,v_next)); + elist = add_data_to_circular_list(elist,&end,e); - t_next_id = t_id; t_next = tri; - insertelem(del_set,t_id); + del_set[0] =1; del_set[1] = t_id; + *delptr = qtnewleaf(del_set); - while((t_next_id = smTri_next_ccw_nbr(sm,t_next,id)) != t_id) + while((t_next_id = T_NTH_NBR(t_next,v_next)) != t_id) { - if((e = eNew_edge()) == SM_INVALID) + if((e = eNew_edge()) == INVALID) + return(NULL); + + t_next = SM_NTH_TRI(sm,t_next_id); + v_next = (T_WHICH_V(t_next,id)+1)%3; + + SET_E_NTH_VERT(e,0,T_NTH_V(t_next,v_next)); + SET_E_NTH_TRI(e,0,INVALID); + SET_E_NTH_TRI(e,1,T_NTH_NBR(t_next,v_next)); + v_next = (T_WHICH_V(t_next,id)+2)%3; + SET_E_NTH_VERT(e,1,T_NTH_V(t_next,v_next)); + elist = add_data_to_circular_list(elist,&end,e); + + + if(qtinset(*delptr,t_next_id)) { #ifdef DEBUG - eputs("smVertex_star_polygon():Too many edges\n"); -#endif - return(NULL); + eputs("smVertex_star_polygon(): id already in set\n"); +#endif + free_list(elist); + return(NULL); } - elist = add_data_to_circular_list(elist,&end,e); - t_next = SM_NTH_TRI(sm,t_next_id); - v_next = (T_WHICH_V(t_next,id)+1)%3; - SET_E_NTH_VERT(e,0,T_NTH_V(t_next,v_next)); - SET_E_NTH_TRI(e,0,SM_INVALID); - SET_E_NTH_TRI(e,1,T_NTH_NBR(t_next,v_next)); - v_next = (T_WHICH_V(t_next,id)+2)%3; - SET_E_NTH_VERT(e,1,T_NTH_V(t_next,v_next)); - insertelem(del_set,t_next_id); + else + qtaddelem(*delptr,t_next_id); } return(elist); } @@ -224,7 +234,7 @@ smFind_next_convex_vertex(sm,id0,id1,v0,v1,l) vertex v such that v0v1v is a convex angle, and the edge v1v does not intersect any other edges */ - id = SM_INVALID; + id = INVALID; el = l; while(id != id0) { @@ -240,7 +250,7 @@ smFind_next_convex_vertex(sm,id0,id1,v0,v1,l) if(el == l) break; } - return(SM_INVALID); + return(INVALID); } int @@ -251,10 +261,10 @@ LIST **l,**lnew; LIST *list,*lptr,*end; int e,e1,e2,new_e; - e2 = SM_INVALID; + e2 = INVALID; list = lptr = *l; - if((new_e = eNew_edge())==SM_INVALID) + if((new_e = eNew_edge())==INVALID) { #ifdef DEBUG eputs("split_edge_list():Too many edges\n"); @@ -263,8 +273,8 @@ LIST **l,**lnew; } SET_E_NTH_VERT(new_e,0,id0); SET_E_NTH_VERT(new_e,1,id_new); - SET_E_NTH_TRI(new_e,0,SM_INVALID); - SET_E_NTH_TRI(new_e,1,SM_INVALID); + SET_E_NTH_TRI(new_e,0,INVALID); + SET_E_NTH_TRI(new_e,1,INVALID); while(e2 != id_new) { @@ -289,7 +299,7 @@ LIST **l,**lnew; /* now follow other cycle */ list = lptr; - e2 = SM_INVALID; + e2 = INVALID; while(e2 != id0) { lptr = LIST_NEXT(lptr); @@ -317,7 +327,6 @@ smTriangulate_convex(sm,plist,add_ptr) SM *sm; LIST *plist,**add_ptr; { - TRI *tri; int t_id,e_id0,e_id1,e_id2; int v_id0,v_id1,v_id2; LIST *lptr; @@ -333,7 +342,7 @@ LIST *plist,**add_ptr; v_id1 = E_NTH_VERT(e_id1,0); v_id2 = E_NTH_VERT(e_id1,1); /* form a triangle for each triple of with v0 as base of star */ - t_id = smAdd_tri(sm,v_id0,v_id1,v_id2,&tri); + t_id = smAdd_tri(sm,v_id0,v_id1,v_id2); *add_ptr = push_data(*add_ptr,t_id); /* add which pointer?*/ @@ -395,7 +404,7 @@ LIST *plist,**add_ptr; } /* if concave: add edge and recurse on two sub polygons */ id_next = smFind_next_convex_vertex(sm,id0,id1,v0,v1,LIST_NEXT(l)); - if(id_next == SM_INVALID) + if(id_next == INVALID) { #ifdef DEBUG eputs("smTriangulate_elist():Unable to find convex vertex\n"); @@ -419,37 +428,163 @@ LIST *plist,**add_ptr; } int -smTriangulate(sm,plist,add_ptr) +smTriangulate_add_tri(sm,id0,id1,id2,e0,e1,e2ptr) SM *sm; +int id0,id1,id2,e0,e1,*e2ptr; +{ + int t_id; + int e2; + + t_id = smAdd_tri(sm,id0,id1,id2); + if(*e2ptr == 0) + { + e2 = eNew_edge(); + SET_E_NTH_VERT(e2,0,id2); + SET_E_NTH_VERT(e2,1,id0); + } + else + e2 = *e2ptr; + /* set appropriate tri for each edge*/ + SET_E_NTH_TRI(e0,0,t_id); + SET_E_NTH_TRI(e1,0,t_id); + SET_E_NTH_TRI(e2,0,t_id); + + *e2ptr = e2; + return(t_id); +} +int +smTriangulate_elist_new(sm,id,plist,add_ptr) +SM *sm; +int id; LIST *plist,**add_ptr; { + LIST *l,*prev,*t; + FVECT v0,v1,v2,n,p; + int is_tri,loop,t_id,id0,id1,id2,e2,eprev,enext; + double dp; + + smDir(sm,p,id); + enext=0; + is_tri= loop = FALSE; + l = prev = plist; + /* get v0,v1,v2 */ + eprev = (int)LIST_DATA(l); + id0 = E_NTH_VERT(eprev,0); + id1 = E_NTH_VERT(eprev,1); + smDir(sm,v0,id0); + smDir(sm,v1,id1); + while(l) + { + l = LIST_NEXT(l); + e2 = (int)LIST_DATA(l); + id2 = E_NTH_VERT(e2,1); + /* Check if have a triangle */ + if(LIST_NEXT(LIST_NEXT(l)) == prev) + { + is_tri = TRUE; + break; + } + if(LIST_NEXT(l) == plist) + { + if(!loop) + loop = 1; + else + loop++; + if(loop > 3) + break; + } + smDir(sm,v2,id2); + /* determine if convex (left turn), or concave(right turn) angle */ + if(!convex_angle(v0,v1,v2)) + { + VCOPY(v0,v1); + VCOPY(v1,v2); + id0 = id1; + id1 = id2; + prev = l; + eprev = e2; + continue; + } + VCROSS(n,v0,v2); + dp = DOT(n,p); + if(loop <=1 && (!ZERO(dp) && dp < 0.0)) + { + VCOPY(v0,v1); + VCOPY(v1,v2); + id0 = id1; + id1 = id2; + eprev = e2; + prev = l; + continue; + } + loop = FALSE; + + enext = 0; + t_id = smTriangulate_add_tri(sm,id0,id1,id2,eprev,e2,&enext); + *add_ptr = push_data(*add_ptr,t_id); + + LIST_NEXT(prev) = LIST_NEXT(l); + LIST_DATA(prev) = eprev = -enext; + LIST_NEXT(l)=NULL; + if(l== plist) + plist = prev; + free_list(l); + l = prev; + VCOPY(v1,v2); + id1 = id2; + } + if(is_tri) + { + l = LIST_NEXT(l); + enext = (int)LIST_DATA(l); + t_id = smTriangulate_add_tri(sm,id0,id1,id2,eprev,e2,&enext); + *add_ptr = push_data(*add_ptr,t_id); + free_list(l); + } + else + { +#ifdef DEBUG + eputs("smTriangulate_elist()Unable to triangulate\n"); +#endif + return(FALSE); + } + return(TRUE); +} + +int +smTriangulate(sm,p_id,plist,add_ptr) +SM *sm; +int p_id; +LIST *plist,**add_ptr; +{ int e,id_t0,id_t1,e0,e1; - TRI *t0,*t1; int test; - + + test = smTriangulate_elist_new(sm,p_id,plist,add_ptr); +#if 0 test = smTriangulate_elist(sm,plist,add_ptr); +#endif if(!test) return(test); + FOR_ALL_EDGES(e) { id_t0 = E_NTH_TRI(e,0); id_t1 = E_NTH_TRI(e,1); - if((id_t0==SM_INVALID) || (id_t1==SM_INVALID)) + if((id_t0==INVALID) || (id_t1==INVALID)) { #ifdef DEBUG eputs("smTriangulate(): Unassigned edge neighbor\n"); #endif continue; } - t0 = SM_NTH_TRI(sm,id_t0); - t1 = SM_NTH_TRI(sm,id_t1); - e0 = T_WHICH_V(t0,E_NTH_VERT(e,0)); - T_NTH_NBR(t0,e0) = id_t1; + e0 = T_WHICH_V(SM_NTH_TRI(sm,id_t0),E_NTH_VERT(e,0)); + T_NTH_NBR(SM_NTH_TRI(sm,id_t0),e0) = id_t1; - e1 = T_WHICH_V(t1,E_NTH_VERT(e,1)); - T_NTH_NBR(t1,e1) = id_t0; + e1 = T_WHICH_V(SM_NTH_TRI(sm,id_t1),E_NTH_VERT(e,1)); + T_NTH_NBR(SM_NTH_TRI(sm,id_t1),e1) = id_t0; } return(test); } @@ -468,80 +603,76 @@ TRI *t; return(T_NTH_V(t,0)==E_NTH_VERT(e,1)||T_NTH_V(t,1)==E_NTH_VERT(e,1)); return(FALSE); } -smFix_edges(sm,add_list,del_set) + +smFix_edges(sm,add_list,delptr) SM *sm; LIST *add_list; - OBJECT *del_set; + QUADTREE *delptr; + { - int e,id_t0,id_t1,e_new,e0,e1,e0_next,e1_next; - TRI *t0,*t1,*nt0,*nt1; - int i,id_v0,id_v1,id_v2,id_p,nid_t0,nid_t1; + int e,t0_id,t1_id,e_new,e0,e1,e0_next,e1_next; + int i,v0_id,v1_id,v2_id,p_id,t0_nid,t1_nid; FVECT v0,v1,v2,p,np,v; FOR_ALL_EDGES(e) { - id_t0 = E_NTH_TRI(e,0); - id_t1 = E_NTH_TRI(e,1); - if((id_t0==SM_INVALID) || (id_t1==SM_INVALID)) + t0_id = E_NTH_TRI(e,0); + t1_id = E_NTH_TRI(e,1); + if((t0_id==INVALID) || (t1_id==INVALID)) { #ifdef DEBUG eputs("smFix_edges: Unassigned edge nbr\n"); #endif continue; } - t0 = SM_NTH_TRI(sm,id_t0); - t1 = SM_NTH_TRI(sm,id_t1); - - e0 = T_WHICH_V(t0,E_NTH_VERT(e,0)); - e1 = T_WHICH_V(t1,E_NTH_VERT(-e,0)); + e0 = T_WHICH_V(SM_NTH_TRI(sm,t0_id),E_NTH_VERT(e,0)); + e1 = T_WHICH_V(SM_NTH_TRI(sm,t1_id),E_NTH_VERT(-e,0)); e0_next = (e0+2)%3; e1_next = (e1+2)%3; - id_v0 = E_NTH_VERT(e,0); - id_v1 = E_NTH_VERT(e,1); - id_v2 = T_NTH_V(t0,e0_next); - id_p = T_NTH_V(t1,e1_next); + v0_id = E_NTH_VERT(e,0); + v1_id = E_NTH_VERT(e,1); + v2_id = T_NTH_V(SM_NTH_TRI(sm,t0_id),e0_next); + p_id = T_NTH_V(SM_NTH_TRI(sm,t1_id),e1_next); - smDir_in_cone(sm,v0,id_v0); - smDir_in_cone(sm,v1,id_v1); - smDir_in_cone(sm,v2,id_v2); + smDir_in_cone(sm,v0,v0_id); + smDir_in_cone(sm,v1,v1_id); + smDir_in_cone(sm,v2,v2_id); - VCOPY(p,SM_NTH_WV(sm,id_p)); + VCOPY(p,SM_NTH_WV(sm,p_id)); VSUB(p,p,SM_VIEW_CENTER(sm)); if(point_in_cone(p,v0,v1,v2)) { - smTris_swap_edge(sm,id_t0,id_t1,e0,e1,&nid_t0,&nid_t1,&add_list, - del_set); + smTris_swap_edge(sm,t0_id,t1_id,e0,e1,&t0_nid,&t1_nid,&add_list, + delptr); - nt0 = SM_NTH_TRI(sm,nid_t0); - nt1 = SM_NTH_TRI(sm,nid_t1); FOR_ALL_EDGES_FROM(e,i) { - if(E_NTH_TRI(i,0)==id_t0 || E_NTH_TRI(i,0)==id_t1) + if(E_NTH_TRI(i,0)==t0_id || E_NTH_TRI(i,0)==t1_id) { - if(eIn_tri(i,nt0)) - SET_E_NTH_TRI(i,0,nid_t0); + if(eIn_tri(i,SM_NTH_TRI(sm,t0_nid))) + SET_E_NTH_TRI(i,0,t0_nid); else - SET_E_NTH_TRI(i,0,nid_t1); + SET_E_NTH_TRI(i,0,t1_nid); } - if(E_NTH_TRI(i,1)==id_t0 || E_NTH_TRI(i,1)==id_t1) + if(E_NTH_TRI(i,1)==t0_id || E_NTH_TRI(i,1)==t1_id) { - if(eIn_tri(i,nt0)) - SET_E_NTH_TRI(i,1,nid_t0); + if(eIn_tri(i,SM_NTH_TRI(sm,t0_nid))) + SET_E_NTH_TRI(i,1,t0_nid); else - SET_E_NTH_TRI(i,1,nid_t1); + SET_E_NTH_TRI(i,1,t1_nid); } } - id_t0 = nid_t0; - id_t1 = nid_t1; + t0_id = t0_nid; + t1_id = t1_nid; e_new = eNew_edge(); - SET_E_NTH_VERT(e_new,0,id_p); - SET_E_NTH_VERT(e_new,1,id_v2); - SET_E_NTH_TRI(e_new,0,id_t0); - SET_E_NTH_TRI(e_new,1,id_t1); + SET_E_NTH_VERT(e_new,0,p_id); + SET_E_NTH_VERT(e_new,1,v2_id); + SET_E_NTH_TRI(e_new,0,t0_id); + SET_E_NTH_TRI(e_new,1,t1_id); } } - smUpdate_locator(sm,add_list,del_set); + smUpdate_locator(sm,add_list,qtqueryset(*delptr)); } int @@ -549,48 +680,44 @@ smMesh_remove_vertex(sm,id) SM *sm; int id; { - int tri; + int t_id; LIST *elist,*add_list; int cnt,debug; - OBJECT del_set[QT_MAXSET +1]; - + QUADTREE delnode; /* generate list of vertices that form the boundary of the star polygon formed by vertex id and all of its adjacent triangles */ eClear_edges(); - QT_CLEAR_SET(del_set); - elist = smVertex_star_polygon(sm,id,del_set); - if(!elist) - return(FALSE); + elist = smVertex_star_polygon(sm,id,&delnode); + if(!elist) + { +#ifdef DEBUG + eputs("smMesh_remove_vertex(): Unable to remove vertex"); +#endif + qtfreeleaf(delnode); + return(FALSE); + } add_list = NULL; /* Triangulate spherical polygon */ - smTriangulate(sm,elist,&add_list); - - + if(!smTriangulate(sm,id,elist,&add_list)) + { + while(add_list) + { + t_id = pop_list(&add_list); + smDelete_tri(sm,t_id); + } + qtfreeleaf(delnode); + return(FALSE); + } /* Fix up new triangles to be Delaunay */ - smFix_edges(sm,add_list,del_set); + smFix_edges(sm,add_list,&delnode); + qtfreeleaf(delnode); return(TRUE); } -/* Remove point from samples, and from mesh. Delete any triangles - adjacent to the point and re-triangulate the hole - Return TRUE is point found , FALSE otherwise -*/ -int -smDelete_point(sm,id) -SM *sm; -int id; -{ - - /* Remove the corresponding vertex from the mesh */ - smMesh_remove_vertex(sm,id); - /* Free the sample point */ - smDelete_sample(sm,id); - return(TRUE); -}