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

Comparing ray/src/hd/sm_stree.c (file contents):
Revision 3.3 by gwlarson, Tue Aug 25 11:03:27 1998 UTC vs.
Revision 3.7 by gwlarson, Wed Nov 11 12:05:41 1998 UTC

# Line 6 | Line 6 | static char SCCSid[] = "$SunId$ SGI";
6  
7   /*
8   * sm_stree.c
9 + *  An stree (spherical quadtree) is defined by an octahedron in
10 + *  canonical form,and a world center point. Each face of the
11 + *  octahedron is adaptively subdivided as a planar triangular quadtree.
12 + *  World space geometry is projected onto the quadtree faces from the
13 + *  sphere center.
14   */
15   #include "standard.h"
16 < #include "object.h"
12 <
16 > #include "sm_flag.h"
17   #include "sm_geom.h"
18 + #include "sm_qtree.h"
19   #include "sm_stree.h"
20  
21 + #ifdef TEST_DRIVER
22 + extern FVECT Pick_point[500],Pick_v0[500],Pick_v1[500],Pick_v2[500];
23 + extern int Pick_cnt;
24 + #endif
25 + /* octahedron coordinates */
26 + FVECT stDefault_base[6] = {  {1.,0.,0.},{0.,1.,0.}, {0.,0.,1.},
27 +                            {-1.,0.,0.},{0.,-1.,0.},{0.,0.,-1.}};
28 + /* octahedron triangle vertices */
29 + int stBase_verts[8][3] = { {2,1,0},{1,5,0},{5,1,3},{1,2,3},
30 +                          {4,2,0},{4,0,5},{3,4,5},{4,3,2}};
31 + /* octahedron triangle nbrs ; nbr i is the face opposite vertex i*/
32 + int stBase_nbrs[8][3] =  { {1,4,3},{5,0,2},{3,6,1},{7,2,0},
33 +                          {0,5,7},{1,6,4},{5,2,7},{3,4,6}};
34 + /* look up table for octahedron point location */
35 + int stlocatetbl[8] = {6,7,2,3,5,4,1,0};
36  
17 /* Define 4 vertices on the sphere to create a tetrahedralization on
18   the sphere: triangles are as follows:
19        (0,1,2),(0,2,3), (0,3,1), (1,3,2)
20 */
37  
38 < FVECT stDefault_base[4] = { {SQRT3_INV, SQRT3_INV, SQRT3_INV},
39 <                          {-SQRT3_INV, -SQRT3_INV, SQRT3_INV},
40 <                          {-SQRT3_INV, SQRT3_INV, -SQRT3_INV},
41 <                          {SQRT3_INV, -SQRT3_INV, -SQRT3_INV}};
26 < int stTri_verts[4][3] = { {2,1,0},
27 <                          {3,2,0},
28 <                          {1,3,0},
29 <                          {2,3,1}};
30 <
31 < stNth_base_verts(st,i,v1,v2,v3)
38 > /* Initializes an stree structure with origin 'center':
39 >   Frees existing quadtrees hanging off of the roots
40 > */
41 > stInit(st)
42   STREE *st;
33 int i;
34 FVECT v1,v2,v3;
43   {
44 <  VCOPY(v1,ST_NTH_BASE(st,stTri_verts[i][0]));
45 <  VCOPY(v2,ST_NTH_BASE(st,stTri_verts[i][1]));
46 <  VCOPY(v3,ST_NTH_BASE(st,stTri_verts[i][2]));
44 >    ST_TOP_ROOT(st) = qtAlloc();
45 >    ST_BOTTOM_ROOT(st) = qtAlloc();
46 >    ST_INIT_ROOT(st);
47   }
48  
49 < /* Frees the 4 quadtrees rooted at st */
49 > /* Frees the children of the 2 quadtrees rooted at st,
50 >   Does not free root nodes: just clears
51 > */
52   stClear(st)
53 < STREE *st;
53 >   STREE *st;
54   {
55 <  int i;
56 <
47 <  /* stree always has 4 children corresponding to the base tris
48 <  */
49 <  for (i = 0; i < 4; i++)
50 <    qtFree(ST_NTH_ROOT(st, i));
51 <
52 <  QT_CLEAR_CHILDREN(ST_ROOT(st));
53 <
55 >    qtDone();
56 >    stInit(st);
57   }
58  
59 <
59 > /* Allocates a stree structure  and creates octahedron base */
60   STREE
58 *stInit(st,center,base)
59 STREE *st;
60 FVECT  center,base[4];
61 {
62
63  if(base)
64    ST_SET_BASE(st,base);
65  else
66    ST_SET_BASE(st,stDefault_base);
67
68  ST_SET_CENTER(st,center);
69  stClear(st);
70
71  return(st);
72 }
73
74
75 /* "base" defines 4 vertices on the sphere to create a tetrahedralization on
76     the sphere: triangles are as follows:(0,1,2),(0,2,3), (0,3,1), (1,3,2)
77     if base is null: does default.
78
79 */
80 STREE
61   *stAlloc(st)
62   STREE *st;
63   {
64 <  int i;
65 <
64 >  int i,m;
65 >  FVECT v0,v1,v2;
66 >  FVECT n;
67 >  
68    if(!st)
69 <    st = (STREE *)malloc(sizeof(STREE));
69 >    if(!(st = (STREE *)malloc(sizeof(STREE))))
70 >       error(SYSTEM,"stAlloc(): Unable to allocate memory\n");
71  
72 <  ST_ROOT(st) = qtAlloc();
73 <    
74 <  QT_CLEAR_CHILDREN(ST_ROOT(st));
72 >  /* Allocate the top and bottom quadtree root nodes */
73 >  stInit(st);
74 >  
75 >  /* Set the octahedron base */
76 >  ST_SET_BASE(st,stDefault_base);
77  
78 +  /* Calculate octahedron face and edge normals */
79 +  for(i=0; i < ST_NUM_ROOT_NODES; i++)
80 +  {
81 +      VCOPY(v0,ST_NTH_V(st,i,0));
82 +      VCOPY(v1,ST_NTH_V(st,i,1));
83 +      VCOPY(v2,ST_NTH_V(st,i,2));
84 +      tri_plane_equation(v0,v1,v2, &ST_NTH_PLANE(st,i),FALSE);
85 +      m = max_index(FP_N(ST_NTH_PLANE(st,i)),NULL);
86 +      FP_X(ST_NTH_PLANE(st,i)) = (m+1)%3;
87 +      FP_Y(ST_NTH_PLANE(st,i)) = (m+2)%3;
88 +      FP_Z(ST_NTH_PLANE(st,i)) = m;
89 +      VCROSS(ST_EDGE_NORM(st,i,0),v0,v1);
90 +      VCROSS(ST_EDGE_NORM(st,i,1),v1,v2);
91 +      VCROSS(ST_EDGE_NORM(st,i,2),v2,v0);
92 +  }
93    return(st);
94   }
95  
96  
97 < /* Find location of sample point in the DAG and return lowest level
98 <   containing triangle. "type" indicates whether the point was found
99 <   to be in interior to the triangle: GT_FACE, on one of its
100 <   edges GT_EDGE or coinciding with one of its vertices GT_VERTEX.
101 <   "which" specifies which vertex (0,1,2) or edge (0=v0v1, 1 = v1v2, 2 = v21)
102 < */
103 < int
104 < stPoint_locate(st,npt,type,which)
105 <    STREE *st;
106 <    FVECT npt;
107 <    char *type,*which;
108 < {
109 <    int i,d,j,id;
110 <    QUADTREE *rootptr,*qtptr;
111 <    FVECT v1,v2,v3;
112 <    OBJECT os[MAXSET+1],*optr;
113 <    char w;
114 <    FVECT p0,p1,p2;
115 <
116 <    /* Test each of the root triangles against point id */
117 <    for(i=0; i < 4; i++)
118 <     {
119 <         rootptr = ST_NTH_ROOT_PTR(st,i);
120 <         stNth_base_verts(st,i,v1,v2,v3);
121 <         /* Return tri that p falls in */
122 <         qtptr = qtRoot_point_locate(rootptr,v1,v2,v3,npt,NULL,NULL,NULL);
123 <         if(!qtptr)
124 <            continue;
125 <         /* Get the set */
126 <         qtgetset(os,*qtptr);
127 <         for (j = QT_SET_CNT(os),optr = QT_SET_PTR(os); j > 0; j--)
128 <         {
129 <             /* Find the first triangle that pt falls */
130 <             id = QT_SET_NEXT_ELEM(optr);
131 <             qtTri_verts_from_id(id,p0,p1,p2);
132 <             d = test_single_point_against_spherical_tri(p0,p1,p2,npt,&w);  
133 <             if(d)
134 <                {
135 <                    if(type)
136 <                       *type = d;
137 <                    if(which)
138 <                       *which = w;
139 <                    return(id);
140 <                }
141 <         }
142 <     }
143 <    if(which)
144 <      *which = 0;
145 <    if(type)
146 <      *type = 0;
147 <    return(EMPTY);
148 < }
149 <
97 > /* Return quadtree leaf node containing point 'p'*/
98   QUADTREE
99 < *stPoint_locate_cell(st,p,t0,t1,t2)
99 > stPoint_locate(st,p)
100      STREE *st;
101      FVECT p;
154    FVECT t0,t1,t2;
102   {
103 <    int i,d;
104 <    QUADTREE *rootptr,*qtptr;
158 <    FVECT v0,v1,v2;
103 >    int i;
104 >    QUADTREE root,qt;
105  
106 +    /* Find root quadtree that contains p */
107 +    i = stPoint_in_root(p);
108 +    root = ST_NTH_ROOT(st,i);
109      
110 <    /* Test each of the root triangles against point id */
111 <    for(i=0; i < 4; i++)
112 <     {
113 <         rootptr = ST_NTH_ROOT_PTR(st,i);
165 <         stNth_base_verts(st,i,v0,v1,v2);
166 <         /* Return tri that p falls in */
167 <         qtptr = qtRoot_point_locate(rootptr,v0,v1,v2,p,t0,t1,t2);
168 <         /* NOTE: For now return only one triangle */
169 <         if(qtptr)
170 <            return(qtptr);
171 <     }    /* Point not found */
172 <    return(NULL);
110 >    /* Traverse quadtree to leaf level */
111 >    qt = qtRoot_point_locate(root,ST_NTH_V(st,i,0),ST_NTH_V(st,i,1),
112 >                        ST_NTH_V(st,i,2),ST_NTH_PLANE(st,i),p);
113 >    return(qt);
114   }
115  
116 + /* Add triangle 'id' with coordinates 't0,t1,t2' to the stree: returns
117 +   FALSE on error, TRUE otherwise
118 + */
119  
120 < QUADTREE
177 < *stAdd_tri_from_pt(st,p,t_id)
178 <    STREE *st;
179 <    FVECT p;
180 <    int t_id;
181 < {
182 <    int i,d;
183 <    QUADTREE *rootptr,*qtptr;
184 <    FVECT v0,v1,v2;
185 <
186 <    
187 <    /* Test each of the root triangles against point id */
188 <    for(i=0; i < 4; i++)
189 <     {
190 <         rootptr = ST_NTH_ROOT_PTR(st,i);
191 <         stNth_base_verts(st,i,v0,v1,v2);
192 <         /* Return tri that p falls in */
193 <         qtptr = qtRoot_add_tri_from_point(rootptr,v0,v1,v2,p,t_id);
194 <         /* NOTE: For now return only one triangle */
195 <         if(qtptr)
196 <            return(qtptr);
197 <     }    /* Point not found */
198 <    return(NULL);
199 < }
200 <
201 < int
202 < stAdd_tri(st,id,v0,v1,v2)
120 > stAdd_tri(st,id,t0,t1,t2)
121   STREE *st;
122   int id;
123 < FVECT v0,v1,v2;
123 > FVECT t0,t1,t2;
124   {
125 <  
126 <  int i,found;
209 <  QUADTREE *rootptr;
210 <  FVECT t0,t1,t2;
125 >  int i;
126 >  QUADTREE root;
127  
128 <  found = 0;
213 <  for(i=0; i < 4; i++)
128 >  for(i=0; i < ST_NUM_ROOT_NODES; i++)
129    {
130 <    rootptr = ST_NTH_ROOT_PTR(st,i);
131 <    stNth_base_verts(st,i,t0,t1,t2);
132 <    found |= qtRoot_add_tri(rootptr,id,v0,v1,v2,t0,t1,t2,0);
130 >    root = ST_NTH_ROOT(st,i);
131 >    ST_NTH_ROOT(st,i) = qtRoot_add_tri(root,ST_NTH_V(st,i,0),ST_NTH_V(st,i,1),
132 >                            ST_NTH_V(st,i,2),t0,t1,t2,id,0);
133    }
219  return(found);
134   }
135  
136 < int
137 < stApply_to_tri_cells(st,v0,v1,v2,func,arg)
136 > /* Remove triangle 'id' with coordinates 't0,t1,t2' to the stree: returns
137 >   FALSE on error, TRUE otherwise
138 > */
139 >
140 > stRemove_tri(st,id,t0,t1,t2)
141   STREE *st;
142 < FVECT v0,v1,v2;
143 < int (*func)();
227 < char *arg;
142 > int id;
143 > FVECT t0,t1,t2;
144   {
145 <  int i,found;
146 <  QUADTREE *rootptr;
231 <  FVECT t0,t1,t2;
145 >  int i;
146 >  QUADTREE root;
147  
148 <  found = 0;
234 <  for(i=0; i < 4; i++)
148 >  for(i=0; i < ST_NUM_ROOT_NODES; i++)
149    {
150 <    rootptr = ST_NTH_ROOT_PTR(st,i);
151 <    stNth_base_verts(st,i,t0,t1,t2);
152 <    found |= qtApply_to_tri_cells(rootptr,v0,v1,v2,t0,t1,t2,func,arg);
150 >    root = ST_NTH_ROOT(st,i);
151 >    ST_NTH_ROOT(st,i)=qtRoot_remove_tri(root,id,ST_NTH_V(st,i,0),ST_NTH_V(st,i,1),
152 >                          ST_NTH_V(st,i,2),t0,t1,t2);
153    }
240  return(found);
154   }
155  
156 + /* Visit all nodes that are intersected by the edges of triangle 't0,t1,t2'
157 +   and apply 'func'
158 + */
159  
160 + stVisit_tri_edges(st,t0,t1,t2,func,fptr,argptr)
161 +   STREE *st;
162 +   FVECT t0,t1,t2;
163 +   int (*func)(),*fptr;
164 +   int *argptr;
165 + {
166 +    int id,i,w,next;
167 +    QUADTREE root;
168 +    FVECT v[3],i_pt;
169  
170 +    VCOPY(v[0],t0); VCOPY(v[1],t1); VCOPY(v[2],t2);
171 +    w = -1;
172  
173 +    /* Locate the root containing triangle vertex v0 */
174 +    i = stPoint_in_root(v[0]);
175 +    /* Mark the root node as visited */
176 +    QT_SET_FLAG(ST_ROOT(st,i));
177 +    root = ST_NTH_ROOT(st,i);
178 +    
179 +    ST_NTH_ROOT(st,i) = qtRoot_visit_tri_edges(root,ST_NTH_V(st,i,0),
180 +       ST_NTH_V(st,i,1),ST_NTH_V(st,i,2),ST_NTH_PLANE(st,i),v,i_pt,&w,
181 +                                               &next,func,fptr,argptr);
182 +    if(QT_FLAG_IS_DONE(*fptr))
183 +      return;
184 +        
185 +    /* Crossed over to next node: id = nbr */
186 +    while(1)
187 +    {
188 +        /* test if ray crosses plane between this quadtree triangle and
189 +           its neighbor- if it does then find intersection point with
190 +           ray and plane- this is the new start point
191 +           */
192 +        i = stBase_nbrs[i][next];
193 +        root = ST_NTH_ROOT(st,i);
194 +        ST_NTH_ROOT(st,i) =
195 +          qtRoot_visit_tri_edges(root,ST_NTH_V(st,i,0),ST_NTH_V(st,i,1),
196 +         ST_NTH_V(st,i,2),ST_NTH_PLANE(st,i),v,i_pt,&w,&next,func,fptr,argptr);
197 +        if(QT_FLAG_IS_DONE(*fptr))
198 +          return;
199 +    }
200 + }
201 +
202 + /* Trace ray 'orig-dir' through stree and apply 'func(qtptr,f,argptr)' at each
203 +   node that it intersects
204 + */
205   int
206 < stRemove_tri(st,id,v0,v1,v2)
207 < STREE *st;
208 < int id;
209 < FVECT v0,v1,v2;
206 > stTrace_ray(st,orig,dir,func,argptr)
207 >   STREE *st;
208 >   FVECT orig,dir;
209 >   int (*func)();
210 >   int *argptr;
211   {
212 <  
213 <  int i,found;
214 <  QUADTREE *rootptr;
215 <  FVECT t0,t1,t2;
212 >    int next,last,i,f=0;
213 >    QUADTREE root;
214 >    FVECT o,n,v;
215 >    double pd,t,d;
216  
217 <  found = 0;
218 <  for(i=0; i < 4; i++)
219 <  {
220 <    rootptr = ST_NTH_ROOT_PTR(st,i);
221 <    stNth_base_verts(st,i,t0,t1,t2);
222 <   found |= qtRemove_tri(rootptr,id,v0,v1,v2,t0,t1,t2);
223 <  }
224 <  return(found);
217 >    VCOPY(o,orig);
218 > #ifdef TEST_DRIVER
219 >       Pick_cnt=0;
220 > #endif;
221 >    /* Find the root node that o falls in */
222 >    i = stPoint_in_root(o);
223 >    root = ST_NTH_ROOT(st,i);
224 >    
225 >    ST_NTH_ROOT(st,i) =
226 >      qtRoot_trace_ray(root,ST_NTH_V(st,i,0),ST_NTH_V(st,i,1),
227 >          ST_NTH_V(st,i,2),ST_NTH_PLANE(st,i),o,dir,&next,func,&f,argptr);
228 >
229 >    if(QT_FLAG_IS_DONE(f))
230 >      return(TRUE);
231 >    
232 >    d = DOT(orig,dir)/sqrt(DOT(orig,orig));
233 >    VSUM(v,orig,dir,-d);
234 >    /* Crossed over to next cell: id = nbr */
235 >    while(1)
236 >      {
237 >        /* test if ray crosses plane between this quadtree triangle and
238 >           its neighbor- if it does then find intersection point with
239 >           ray and plane- this is the new origin
240 >           */
241 >        if(next == INVALID)
242 >          return(FALSE);
243 > #if 0
244 >        if(!intersect_ray_oplane(o,dir,ST_EDGE_NORM(st,i,(next+1)%3),NULL,o))
245 > #endif
246 >           if(DOT(o,v) < 0.0)
247 >           /* Ray does not cross into next cell: done and tri not found*/
248 >           return(FALSE);
249 >
250 >        VSUM(o,o,dir,10*FTINY);
251 >        i = stBase_nbrs[i][next];
252 >        root = ST_NTH_ROOT(st,i);
253 >        
254 >        ST_NTH_ROOT(st,i) =
255 >          qtRoot_trace_ray(root, ST_NTH_V(st,i,0),ST_NTH_V(st,i,1),
256 >               ST_NTH_V(st,i,2),ST_NTH_PLANE(st,i),o,dir,&next,func,&f,argptr);
257 >        if(QT_FLAG_IS_DONE(f))
258 >          return(TRUE);
259 >      }
260   }
261  
262  
263 + /* Visit nodes intersected by tri 't0,t1,t2' and apply 'func(arg1,arg2,arg3):
264 +   assumes that stVisit_tri_edges has already been called such that all nodes
265 +   intersected by tri edges are already marked as visited
266 + */
267 + stVisit_tri(st,t0,t1,t2,func,f,argptr)
268 +   STREE *st;
269 +   FVECT t0,t1,t2;
270 +   int (*func)(),*f;
271 +   int *argptr;
272 + {
273 +  int i;
274 +  QUADTREE root;
275 +  FVECT n0,n1,n2;
276  
277 +  /* Calcuate the edge normals for tri */
278 +  VCROSS(n0,t0,t1);
279 +  VCROSS(n1,t1,t2);
280 +  VCROSS(n2,t2,t0);
281  
282 +  for(i=0; i < ST_NUM_ROOT_NODES; i++)
283 +  {
284 +    root = ST_NTH_ROOT(st,i);
285 +    ST_NTH_ROOT(st,i) = qtVisit_tri_interior(root,ST_NTH_V(st,i,0),
286 +        ST_NTH_V(st,i,1),ST_NTH_V(st,i,2),t0,t1,t2,n0,n1,n2,0,func,f,argptr);
287  
288 < #if 0
288 >  }
289 > }
290 >
291 > /* Visit nodes intersected by tri 't0,t1,t2'.Apply 'edge_func(arg1,arg2,arg3)',
292 >   to those nodes intersected by edges, and interior_func to ALL nodes:
293 >   ie some Nodes  will be visited more than once
294 > */
295   int
296 < stAdd_tri_opt(st,id,v0,v1,v2)
297 < STREE *st;
298 < int id;
299 < FVECT v0,v1,v2;
296 > stApply_to_tri(st,t0,t1,t2,edge_func,tri_func,argptr)
297 >   STREE *st;
298 >   FVECT t0,t1,t2;
299 >   int (*edge_func)(),(*tri_func)();
300 >   int *argptr;
301   {
302 <  
303 <  int i,found;
280 <  QUADTREE *qtptr;
281 <  FVECT pt,t0,t1,t2;
302 >    int f;
303 >    FVECT dir;
304  
305 + #ifdef TEST_DRIVER
306 +    Pick_cnt=0;
307 + #endif;
308    /* First add all of the leaf cells lying on the triangle perimeter:
309       mark all cells seen on the way
310     */
311 <  /* clear all of the flags */
312 <  qtClearAllFlags();            /* clear all quadtree branch flags */
311 >    f = 0;
312 >    /* Visit cells along edges of the tri */
313 >    stVisit_tri_edges(st,t0,t1,t2,edge_func,&f,argptr);
314  
315 <  /* Now trace each triangle edge-marking cells visited, and adding tri to
316 <     the leafs
317 <   */
292 <  stAdd_tri_from_pt(st,v0,id,t0,t1,t2);
293 <  /* Find next cell that projection of ray intersects */
294 <  VCOPY(pt,v0);
295 <  /* NOTE: Check if in same cell */
296 <  while(traceEdge(pt,v1,t0,t1,t2,pt))
297 <  {
298 <      stAdd_tri_from_pt(st,pt,id,t0,t1,t2);
299 <      traceEdge(pt,v1,t0,t1,t2,pt);
300 <  }
301 <  while(traceEdge(pt,v2,t0,t1,t2,pt))
302 <  {
303 <      stAdd_tri_from_pt(st,pt,id,t0,t1,t2);
304 <      traceEdge(pt,v2,t0,t1,t2,pt);
305 <  }
306 <  while(traceEdge(pt,v0,t0,t1,t2,pt))
307 <  {
308 <      stAdd_tri_from_pt(st,pt,id,t0,t1,t2);
309 <      traceEdge(pt,v2,t0,t1,t2,pt);
310 <  }
311 <
312 <  /* NOTE: Optimization: if <= 2 cells added: dont need to fill */
313 <  /* Traverse: follow nodes with flag set or one vertex in triangle */
314 <  
315 >    /* Now visit All cells interior */
316 >    if(QT_FLAG_FILL_TRI(f) || QT_FLAG_UPDATE(f))
317 >       stVisit_tri(st,t0,t1,t2,tri_func,&f,argptr);
318   }
319  
320 < #endif
320 >
321 >
322 >
323 >
324 >
325 >
326 >
327 >
328 >
329 >
330 >
331 >

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines