ViewVC Help

View File | Revision Log | Show Annotations | Download File | Root Listing

root/Development/ray/src/hd/sm_stree.c

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.10 by gwlarson, Sun Jan 10 10:27:44 1999 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"
17	<
16	>	#include "sm_list.h"
17	>	#include "sm_flag.h"
18		#include "sm_geom.h"
19	+	#include "sm_qtree.h"
20		#include "sm_stree.h"
21
22	+	#ifdef TEST_DRIVER
23	+	extern FVECT Pick_point[500],Pick_v0[500],Pick_v1[500],Pick_v2[500];
24	+	extern int Pick_cnt;
25	+	#endif
26	+	/* octahedron coordinates */
27	+	FVECT stDefault_base[6] = {  {1.,0.,0.},{0.,1.,0.}, {0.,0.,1.},
28	+	{-1.,0.,0.},{0.,-1.,0.},{0.,0.,-1.}};
29	+	/* octahedron triangle vertices */
30	+	int stBase_verts[8][3] = { {0,1,2},{3,1,2},{0,4,2},{3,4,2},
31	+	{0,1,5},{3,1,5},{0,4,5},{3,4,5}};
32	+	/* octahedron triangle nbrs ; nbr i is the face opposite vertex i*/
33	+	int stBase_nbrs[8][3] =  { {1,2,4},{0,3,5},{3,0,6},{2,1,7},
34	+	{5,6,0},{4,7,1},{7,4,2},{6,5,3}};
35	+	int stRoot_indices[8][3] = {{1,1,1},{-1,1,1},{1,-1,1},{-1,-1,1},
36	+	{1,1,-1},{-1,1,-1},{1,-1,-1},{-1,-1,-1}};
37	+	/*
38	+	+z   y                -z   y
39	+	|                     |
40	+	1    |   0             5   |   4
41	+	______|______ x      _______|______ x
42	+	3    |   2             7   |   6
43	+	|                     |
44
45	<	/* Define 4 vertices on the sphere to create a tetrahedralization on
46	<	the sphere: triangles are as follows:
47	<	(0,1,2),(0,2,3), (0,3,1), (1,3,2)
45	>	Nbrs
46	>	+z   y                -z   y
47	>	/0|1\                 /1|0\
48	>	5  /  |  \ 4             /  |  \
49	>	/(1)|(0)\           1 /(5)|(4)\ 0
50	>	/    |    \           /    |    \
51	>	/2   1|0   2\         /2   0|1   2\
52	>	/------|------\x      /------|------\x
53	>	\0    1|2    0/       \0    2|2    1/
54	>	\     |     /         \     |     /
55	>	7\ (3)|(2) / 6       3 \ (7)|(6) / 2
56	>	\   |   /             \   |   /
57	>	\ 2|1 /               \ 1|0 /
58		*/
59
22	–	FVECT stDefault_base[4] = { {SQRT3_INV, SQRT3_INV, SQRT3_INV},
23	–	{-SQRT3_INV, -SQRT3_INV, SQRT3_INV},
24	–	{-SQRT3_INV, SQRT3_INV, -SQRT3_INV},
25	–	{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}};
60
61	<	stNth_base_verts(st,i,v1,v2,v3)
61	>	stInit(st)
62		STREE *st;
33	–	int i;
34	–	FVECT v1,v2,v3;
63		{
64	<	VCOPY(v1,ST_NTH_BASE(st,stTri_verts[i][0]));
65	<	VCOPY(v2,ST_NTH_BASE(st,stTri_verts[i][1]));
66	<	VCOPY(v3,ST_NTH_BASE(st,stTri_verts[i][2]));
64	>	int i,j;
65	>
66	>	ST_TOP_QT(st) = qtAlloc();
67	>	ST_BOTTOM_QT(st) = qtAlloc();
68	>	/* Clear the children */
69	>
70	>	QT_CLEAR_CHILDREN(ST_TOP_QT(st));
71	>	QT_CLEAR_CHILDREN(ST_BOTTOM_QT(st));
72		}
73
74	<	/* Frees the 4 quadtrees rooted at st */
74	>	/* Frees the children of the 2 quadtrees rooted at st,
75	>	Does not free root nodes: just clears
76	>	*/
77		stClear(st)
78	+	STREE *st;
79	+	{
80	+	qtDone();
81	+	stInit(st);
82	+	}
83	+
84	+	/* Allocates a stree structure  and creates octahedron base */
85	+	STREE
86	+	*stAlloc(st)
87		STREE *st;
88		{
89	<	int i;
89	>	int i,m;
90	>	FVECT v0,v1,v2;
91	>	FVECT n;
92	>
93	>	if(!st)
94	>	if(!(st = (STREE *)malloc(sizeof(STREE))))
95	>	error(SYSTEM,"stAlloc(): Unable to allocate memory\n");
96
97	<	/* stree always has 4 children corresponding to the base tris
98	<	*/
99	<	for (i = 0; i < 4; i++)
100	<	qtFree(ST_NTH_ROOT(st, i));
97	>	/* Allocate the top and bottom quadtree root nodes */
98	>	stInit(st);
99	>
100	>
101	>	/* will go ********************************************/
102	>	/* Set the octahedron base */
103	>	ST_SET_BASE(st,stDefault_base);
104
105	<	QT_CLEAR_CHILDREN(ST_ROOT(st));
105	>	/* Calculate octahedron face and edge normals */
106	>	for(i=0; i < ST_NUM_ROOT_NODES; i++)
107	>	{
108	>	VCOPY(v0,ST_NTH_V(st,i,0));
109	>	VCOPY(v1,ST_NTH_V(st,i,1));
110	>	VCOPY(v2,ST_NTH_V(st,i,2));
111	>	tri_plane_equation(v0,v1,v2, &ST_NTH_PLANE(st,i),FALSE);
112	>	m = max_index(FP_N(ST_NTH_PLANE(st,i)),NULL);
113	>	FP_X(ST_NTH_PLANE(st,i)) = (m+1)%3;
114	>	FP_Y(ST_NTH_PLANE(st,i)) = (m+2)%3;
115	>	FP_Z(ST_NTH_PLANE(st,i)) = m;
116	>	VCROSS(ST_EDGE_NORM(st,i,0),v0,v1);
117	>	VCROSS(ST_EDGE_NORM(st,i,1),v1,v2);
118	>	VCROSS(ST_EDGE_NORM(st,i,2),v2,v0);
119	>	}
120
121	+	/*****************************************************************/
122	+	return(st);
123		}
124
125	+	#define BARY_INT(v,b)  if((v)>2.0) (b) = MAXBCOORD;else \
126	+	if((v)<-2.0) (b)=-MAXBCOORD;else (b)=(BCOORD)((v)*MAXBCOORD2);
127
128	<	STREE
129	<	*stInit(st,center,base)
130	<	STREE *st;
131	<	FVECT  center,base[4];
128	>	vert_to_qt_frame(root,v,b)
129	>	int root;
130	>	FVECT v;
131	>	BCOORD b[3];
132		{
133	+	int i;
134	+	double scale;
135	+	double d0,d1,d2;
136
137	<	if(base)
138	<	ST_SET_BASE(st,base);
137	>	if(STR_NTH_INDEX(root,0)==-1)
138	>	d0 = -v[0];
139		else
140	<	ST_SET_BASE(st,stDefault_base);
140	>	d0 = v[0];
141	>	if(STR_NTH_INDEX(root,1)==-1)
142	>	d1 = -v[1];
143	>	else
144	>	d1 = v[1];
145	>	if(STR_NTH_INDEX(root,2)==-1)
146	>	d2 = -v[2];
147	>	else
148	>	d2 = v[2];
149
150	<	ST_SET_CENTER(st,center);
151	<	stClear(st);
150	>	/* Plane is now x+y+z = 1 - intersection of pt ray is qtv/den */
151	>	scale = 1.0/ (d0 + d1 + d2);
152	>	d0 *= scale;
153	>	d1 *= scale;
154	>	d2 *= scale;
155
156	<	return(st);
156	>	BARY_INT(d0,b[0])
157	>	BARY_INT(d1,b[1])
158	>	BARY_INT(d2,b[2])
159		}
160
161
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.
162
163	<	*/
164	<	STREE
165	<	*stAlloc(st)
166	<	STREE *st;
163	>
164	>	ray_to_qt_frame(root,v,dir,b,db)
165	>	int root;
166	>	FVECT v,dir;
167	>	BCOORD b[3],db[3];
168		{
169		int i;
170	+	double scale;
171	+	double d0,d1,d2;
172	+	double dir0,dir1,dir2;
173
174	<	if(!st)
175	<	st = (STREE *)malloc(sizeof(STREE));
174	>	if(STR_NTH_INDEX(root,0)==-1)
175	>	{
176	>	d0 = -v[0];
177	>	dir0 = -dir[0];
178	>	}
179	>	else
180	>	{
181	>	d0 = v[0];
182	>	dir0 = dir[0];
183	>	}
184	>	if(STR_NTH_INDEX(root,1)==-1)
185	>	{
186	>	d1 = -v[1];
187	>	dir1 = -dir[1];
188	>	}
189	>	else
190	>	{
191	>	d1 = v[1];
192	>	dir1 = dir[1];
193	>	}
194	>	if(STR_NTH_INDEX(root,2)==-1)
195	>	{
196	>	d2 = -v[2];
197	>	dir2 = -dir[2];
198	>	}
199	>	else
200	>	{
201	>	d2 = v[2];
202	>	dir2 = dir[2];
203	>	}
204	>	/* Plane is now x+y+z = 1 - intersection of pt ray is qtv/den */
205	>	scale = 1.0/ (d0 + d1 + d2);
206	>	d0 *= scale;
207	>	d1 *= scale;
208	>	d2 *= scale;
209
210	<	ST_ROOT(st) = qtAlloc();
211	<
212	<	QT_CLEAR_CHILDREN(ST_ROOT(st));
210	>	/* Calculate intersection point of orig+dir: This calculation is done
211	>	after the origin is projected into the plane in order to constrain
212	>	the projection( i.e. the size of the projection of the unit direction
213	>	vector translated to the origin depends on how close
214	>	the origin is to the view center
215	>	*/
216	>	/* Must divide by at least root2 to insure that projection will fit
217	>	int [-2,2] bounds: assumed length is 1: therefore greatest projection
218	>	from endpoint of triangle is at 45 degrees or projected length of root2
219	>	*/
220	>	dir0 = d0 + dir0*0.5;
221	>	dir1 = d1 + dir1*0.5;
222	>	dir2 = d2 + dir2*0.5;
223
224	<	return(st);
225	<	}
224	>	scale = 1.0/ (dir0 + dir1 + dir2);
225	>	dir0 *= scale;
226	>	dir1 *= scale;
227	>	dir2 *= scale;
228
229	+	BARY_INT(d0,b[0])
230	+	BARY_INT(d1,b[1])
231	+	BARY_INT(d2,b[2])
232	+	BARY_INT(dir0,db[0])
233	+	BARY_INT(dir1,db[1])
234	+	BARY_INT(dir2,db[2])
235
236	<	/* Find location of sample point in the DAG and return lowest level
237	<	containing triangle. "type" indicates whether the point was found
238	<	to be in interior to the triangle: GT_FACE, on one of its
239	<	edges GT_EDGE or coinciding with one of its vertices GT_VERTEX.
240	<	"which" specifies which vertex (0,1,2) or edge (0=v0v1, 1 = v1v2, 2 = v21)
241	<	*/
242	<	int
243	<	stPoint_locate(st,npt,type,which)
244	<	STREE *st;
106	<	FVECT npt;
107	<	char *type,*which;
236	>	db[0]  -= b[0];
237	>	db[1]  -= b[1];
238	>	db[2]  -= b[2];
239	>	}
240	>
241	>	qt_frame_to_vert(root,b,v)
242	>	int root;
243	>	BCOORD b[3];
244	>	FVECT v;
245		{
246	<	int i,d,j,id;
247	<	QUADTREE *rootptr,*qtptr;
111	<	FVECT v1,v2,v3;
112	<	OBJECT os[MAXSET+1],*optr;
113	<	char w;
114	<	FVECT p0,p1,p2;
246	>	int i;
247	>	double d0,d1,d2;
248
249	<	/* Test each of the root triangles against point id */
250	<	for(i=0; i < 4; i++)
251	<	{
252	<	rootptr = ST_NTH_ROOT_PTR(st,i);
253	<	stNth_base_verts(st,i,v1,v2,v3);
254	<	/* Return tri that p falls in */
255	<	qtptr = qtRoot_point_locate(rootptr,v1,v2,v3,npt,NULL,NULL,NULL);
256	<	if(!qtptr)
257	<	continue;
258	<	/* Get the set */
259	<	qtgetset(os,*qtptr);
260	<	for (j = QT_SET_CNT(os),optr = QT_SET_PTR(os); j > 0; j--)
261	<	{
262	<	/* Find the first triangle that pt falls */
263	<	id = QT_SET_NEXT_ELEM(optr);
264	<	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);
249	>	d0 = b[0]/(double)MAXBCOORD2;
250	>	d1 = b[1]/(double)MAXBCOORD2;
251	>	d2 = b[2]/(double)MAXBCOORD2;
252	>
253	>	if(STR_NTH_INDEX(root,0)==-1)
254	>	v[0] = -d0;
255	>	else
256	>	v[0] = d0;
257	>	if(STR_NTH_INDEX(root,1)==-1)
258	>	v[1] = -d1;
259	>	else
260	>	v[1] = d1;
261	>	if(STR_NTH_INDEX(root,2)==-1)
262	>	v[2] = -d2;
263	>	else
264	>	v[2] = d2;
265		}
266
267	+
268	+	/* Return quadtree leaf node containing point 'p'*/
269		QUADTREE
270	<	*stPoint_locate_cell(st,p,t0,t1,t2)
270	>	stPoint_locate(st,p)
271		STREE *st;
272		FVECT p;
154	–	FVECT t0,t1,t2;
273		{
274	<	int i,d;
275	<	QUADTREE *rootptr,*qtptr;
276	<	FVECT v0,v1,v2;
274	>	QUADTREE qt;
275	>	BCOORD bcoordi[3];
276	>	int i;
277
278	+	/* Find root quadtree that contains p */
279	+	i = stLocate_root(p);
280	+	qt = ST_ROOT_QT(st,i);
281
282	<	/* Test each of the root triangles against point id */
283	<	for(i=0; i < 4; i++)
284	<	{
285	<	rootptr = ST_NTH_ROOT_PTR(st,i);
286	<	stNth_base_verts(st,i,v0,v1,v2);
287	<	/* Return tri that p falls in */
288	<	qtptr = qtRoot_point_locate(rootptr,v0,v1,v2,p,t0,t1,t2);
289	<	/* NOTE: For now return only one triangle */
290	<	if(qtptr)
291	<	return(qtptr);
292	<	}    /* Point not found */
293	<	return(NULL);
282	>	/* Will return lowest level triangle containing point: It the
283	>	point is on an edge or vertex: will return first associated
284	>	triangle encountered in the child traversal- the others can
285	>	be derived using triangle adjacency information
286	>	*/
287	>	if(QT_IS_TREE(qt))
288	>	{
289	>	vert_to_qt_frame(i,p,bcoordi);
290	>	i = bary_child(bcoordi);
291	>	return(qtLocate(QT_NTH_CHILD(qt,i),bcoordi));
292	>	}
293	>	else
294	>	return(qt);
295		}
296
297	+	static unsigned int nbr_b[8][3] ={{2,4,16},{1,8,32},{8,1,64},{4,2,128},
298	+	{32,64,1},{16,128,2},{128,16,4},{64,32,8}};
299	+	unsigned int
300	+	stTri_cells(st,v)
301	+	STREE *st;
302	+	FVECT v[3];
303	+	{
304	+	unsigned int cells,cross;
305	+	unsigned int vcell[3];
306	+	double t0,t1;
307	+	int i,inext;
308
309	<	QUADTREE
310	<	*stAdd_tri_from_pt(st,p,t_id)
311	<	STREE *st;
312	<	FVECT p;
313	<	int t_id;
309	>	/* First find base cells that tri vertices are in (0-7)*/
310	>	vcell[0] = stLocate_root(v[0]);
311	>	vcell[1] = stLocate_root(v[1]);
312	>	vcell[2] = stLocate_root(v[2]);
313	>
314	>	/* If all in same cell- return that bit only */
315	>	if(vcell[0] == vcell[1] && vcell[1] == vcell[2])
316	>	return( 1 << vcell[0]);
317	>
318	>	cells = 0;
319	>	for(i=0;i<3; i++)
320	>	{
321	>	if(i==2)
322	>	inext = 0;
323	>	else
324	>	inext = i+1;
325	>	/* Mark cell containing initial vertex */
326	>	cells |= 1 << vcell[i];
327	>
328	>	/* Take the exclusive or: will have bits set where edge crosses axis=0*/
329	>	cross = vcell[i] ^ vcell[inext];
330	>	/* If crosses 2 planes: then have 2 options for edge crossing-pick closest
331	>	otherwise just hits two*/
332	>	/* Neighbors are zyx */
333	>	switch(cross){
334	>	case 3: /* crosses x=0 and y=0 */
335	>	t0 = -v[i][0]/(v[inext][0]-v[i][0]);
336	>	t1 = -v[i][1]/(v[inext][1]-v[i][1]);
337	>	if(t0==t1)
338	>	break;
339	>	else if(t0 < t1)
340	>	cells |= nbr_b[vcell[i]][0];
341	>	else
342	>	cells |= nbr_b[vcell[i]][1];
343	>	break;
344	>	case 5: /* crosses x=0 and z=0 */
345	>	t0 = -v[i][0]/(v[inext][0]-v[i][0]);
346	>	t1 = -v[i][2]/(v[inext][2]-v[i][2]);
347	>	if(t0==t1)
348	>	break;
349	>	else if(t0 < t1)
350	>	cells |= nbr_b[vcell[i]][0];
351	>	else
352	>	cells |=nbr_b[vcell[i]][2];
353	>
354	>	break;
355	>	case 6:/* crosses  z=0 and y=0 */
356	>	t0 = -v[i][2]/(v[inext][2]-v[i][2]);
357	>	t1 = -v[i][1]/(v[inext][1]-v[i][1]);
358	>	if(t0==t1)
359	>	break;
360	>	else if(t0 < t1)
361	>	{
362	>	cells |= nbr_b[vcell[i]][2];
363	>	}
364	>	else
365	>	{
366	>	cells |=nbr_b[vcell[i]][1];
367	>	}
368	>	break;
369	>	case 7:
370	>	error(CONSISTENCY," Insert:Edge shouldnt be able to span 3 cells");
371	>	break;
372	>	}
373	>	}
374	>	return(cells);
375	>	}
376	>
377	>
378	>	stRoot_trace_ray(qt,root,orig,dir,nextptr,func,f)
379	>	QUADTREE qt;
380	>	int root;
381	>	FVECT orig,dir;
382	>	int *nextptr;
383	>	FUNC func;
384	>	int *f;
385		{
386	<	int i,d;
387	<	QUADTREE *rootptr,*qtptr;
388	<	FVECT v0,v1,v2;
386	>	double br[3];
387	>	BCOORD bi[3],dbi[3];
388	>
389	>	/* Project the origin onto the root node plane */
390	>	/* Find the intersection point of the origin */
391	>	ray_to_qt_frame(root,orig,dir,bi,dbi);
392
393	<
394	<	/* Test each of the root triangles against point id */
395	<	for(i=0; i < 4; i++)
396	<	{
397	<	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);
393	>	/* trace the ray starting with this node */
394	>	qtTrace_ray(qt,bi,dbi[0],dbi[1],dbi[2],nextptr,0,0,func,f);
395	>	if(!QT_FLAG_IS_DONE(*f))
396	>	qt_frame_to_vert(root,bi,orig);
397	>
398		}
399
400	+	/* Trace ray 'orig-dir' through stree and apply 'func(qtptr,f,argptr)' at each
401	+	node that it intersects
402	+	*/
403		int
404	<	stAdd_tri(st,id,v0,v1,v2)
404	>	stTrace_ray(st,orig,dir,func)
405	>	STREE *st;
406	>	FVECT orig,dir;
407	>	FUNC func;
408	>	{
409	>	int next,last,i,f=0;
410	>	QUADTREE qt;
411	>	FVECT o,n,v;
412	>	double pd,t,d;
413	>
414	>	VCOPY(o,orig);
415	>	#ifdef TEST_DRIVER
416	>	Pick_cnt=0;
417	>	#endif;
418	>	/* Find the qt node that o falls in */
419	>	i = stLocate_root(o);
420	>	qt = ST_ROOT_QT(st,i);
421	>
422	>	stRoot_trace_ray(qt,i,o,dir,&next,func,&f);
423	>
424	>	if(QT_FLAG_IS_DONE(f))
425	>	return(TRUE);
426	>	/*
427	>	d = DOT(orig,dir)/sqrt(DOT(orig,orig));
428	>	VSUM(v,orig,dir,-d);
429	>	*/
430	>	/* Crossed over to next cell: id = nbr */
431	>	while(1)
432	>	{
433	>	/* test if ray crosses plane between this quadtree triangle and
434	>	its neighbor- if it does then find intersection point with
435	>	ray and plane- this is the new origin
436	>	*/
437	>	if(next == INVALID)
438	>	return(FALSE);
439	>	/*
440	>	if(DOT(o,v) < 0.0)
441	>	return(FALSE);
442	>	*/
443	>	i = stBase_nbrs[i][next];
444	>	qt = ST_ROOT_QT(st,i);
445	>	stRoot_trace_ray(qt,i,o,dir,&next,func,&f);
446	>	if(QT_FLAG_IS_DONE(f))
447	>	return(TRUE);
448	>	}
449	>	}
450	>
451	>
452	>	stVisit_poly(st,verts,l,root,func)
453		STREE *st;
454	<	int id;
455	<	FVECT v0,v1,v2;
454	>	FVECT *verts;
455	>	LIST *l;
456	>	unsigned int root;
457	>	FUNC func;
458		{
459	<
460	<	int i,found;
209	<	QUADTREE *rootptr;
210	<	FVECT t0,t1,t2;
459	>	int id0,id1,id2;
460	>	FVECT tri[3];
461
462	<	found = 0;
463	<	for(i=0; i < 4; i++)
462	>	id0 = pop_list(&l);
463	>	id1 = pop_list(&l);
464	>	while(l)
465		{
466	<	rootptr = ST_NTH_ROOT_PTR(st,i);
467	<	stNth_base_verts(st,i,t0,t1,t2);
468	<	found |= qtRoot_add_tri(rootptr,id,v0,v1,v2,t0,t1,t2,0);
466	>	id2 = pop_list(&l);
467	>	VCOPY(tri[0],verts[id0]);
468	>	VCOPY(tri[1],verts[id1]);
469	>	VCOPY(tri[2],verts[id2]);
470	>	stRoot_visit_tri(st,root,tri,func);
471	>	id1 = id2;
472		}
219	–	return(found);
473		}
474
475	<	int
476	<	stApply_to_tri_cells(st,v0,v1,v2,func,arg)
477	<	STREE *st;
478	<	FVECT v0,v1,v2;
479	<	int (*func)();
480	<	char *arg;
475	>	stVisit_clip(st,i,verts,vcnt,l,cell,func)
476	>	STREE *st;
477	>	int i;
478	>	FVECT *verts;
479	>	int *vcnt;
480	>	LIST *l;
481	>	unsigned int cell;
482	>	FUNC func;
483		{
229	–	int i,found;
230	–	QUADTREE *rootptr;
231	–	FVECT t0,t1,t2;
484
485	<	found = 0;
486	<	for(i=0; i < 4; i++)
485	>	LIST *labove,*lbelow,*endb,*enda;
486	>	int last = -1;
487	>	int id,last_id;
488	>	int first,first_id;
489	>	unsigned int cellb;
490	>
491	>	labove = lbelow = NULL;
492	>	enda = endb = NULL;
493	>	while(l)
494		{
495	<	rootptr = ST_NTH_ROOT_PTR(st,i);
496	<	stNth_base_verts(st,i,t0,t1,t2);
497	<	found |= qtApply_to_tri_cells(rootptr,v0,v1,v2,t0,t1,t2,func,arg);
495	>	id = pop_list(&l);
496	>	if(ZERO(verts[id][i]))
497	>	{
498	>	if(last==-1)
499	>	{/* add below and above */
500	>	first = 2;
501	>	first_id= id;
502	>	}
503	>	lbelow=add_data(lbelow,id,&endb);
504	>	labove=add_data(labove,id,&enda);
505	>	last_id = id;
506	>	last = 2;
507	>	continue;
508	>	}
509	>	if(verts[id][i] < 0)
510	>	{
511	>	if(last != 1)
512	>	{
513	>	lbelow=add_data(lbelow,id,&endb);
514	>	if(last==-1)
515	>	{
516	>	first = 0;
517	>	first_id = id;
518	>	}
519	>	last_id = id;
520	>	last = 0;
521	>	continue;
522	>	}
523	>	/* intersect_edges */
524	>	intersect_edge_coord_plane(verts[last_id],verts[id],i,verts[*vcnt]);
525	>	/*newpoint goes to above and below*/
526	>	lbelow=add_data(lbelow,*vcnt,&endb);
527	>	lbelow=add_data(lbelow,id,&endb);
528	>	labove=add_data(labove,*vcnt,&enda);
529	>	last = 0;
530	>	last_id = id;
531	>	(*vcnt)++;
532	>	}
533	>	else
534	>	{
535	>	if(last != 0)
536	>	{
537	>	labove=add_data(labove,id,&enda);
538	>	if(last==-1)
539	>	{
540	>	first = 1;
541	>	first_id = id;
542	>	}
543	>	last_id = id;
544	>	last = 1;
545	>	continue;
546	>	}
547	>	/* intersect_edges */
548	>	/*newpoint goes to above and below*/
549	>	intersect_edge_coord_plane(verts[last_id],verts[id],i,verts[*vcnt]);
550	>	lbelow=add_data(lbelow,*vcnt,&endb);
551	>	labove=add_data(labove,*vcnt,&enda);
552	>	labove=add_data(labove,id,&enda);
553	>	last_id = id;
554	>	(*vcnt)++;
555	>	last = 1;
556	>	}
557		}
558	<	return(found);
558	>	if(first != 2 && first != last)
559	>	{
560	>	intersect_edge_coord_plane(verts[id],verts[first_id],i,verts[*vcnt]);
561	>	/*newpoint goes to above and below*/
562	>	lbelow=add_data(lbelow,*vcnt,&endb);
563	>	labove=add_data(labove,*vcnt,&enda);
564	>	(*vcnt)++;
565	>
566	>	}
567	>	if(i==2)
568	>	{
569	>	if(lbelow)
570	>	{
571	>	if(LIST_NEXT(lbelow) && LIST_NEXT(LIST_NEXT(lbelow)))
572	>	{
573	>	cellb = cell | (1 << i);
574	>	stVisit_poly(st,verts,lbelow,cellb,func);
575	>	}
576	>	else
577	>	free_list(lbelow);
578	>	}
579	>	if(labove)
580	>	{
581	>	if(LIST_NEXT(labove) && LIST_NEXT(LIST_NEXT(labove)))
582	>	stVisit_poly(st,verts,labove,cell,func);
583	>	else
584	>	free_list(labove);
585	>	}
586	>	}
587	>	else
588	>	{
589	>	if(lbelow)
590	>	{
591	>	if(LIST_NEXT(lbelow) && LIST_NEXT(LIST_NEXT(lbelow)))
592	>	{
593	>	cellb = cell | (1 << i);
594	>	stVisit_clip(st,i+1,verts,vcnt,lbelow,cellb,func);
595	>	}
596	>	else
597	>	free_list(lbelow);
598	>	}
599	>	if(labove)
600	>	{
601	>	if(LIST_NEXT(labove) && LIST_NEXT(LIST_NEXT(labove)))
602	>	stVisit_clip(st,i+1,verts,vcnt,labove,cell,func);
603	>	else
604	>	free_list(labove);
605	>	}
606	>	}
607	>
608		}
609
610	+	stVisit(st,tri,func)
611	+	STREE *st;
612	+	FVECT tri[3];
613	+	FUNC func;
614	+	{
615	+	int r0,r1,r2;
616	+	LIST *l;
617
618	+	r0 = stLocate_root(tri[0]);
619	+	r1 = stLocate_root(tri[1]);
620	+	r2 = stLocate_root(tri[2]);
621	+	if(r0 == r1 && r1==r2)
622	+	stRoot_visit_tri(st,r0,tri,func);
623	+	else
624	+	{
625	+	FVECT verts[ST_CLIP_VERTS];
626	+	int cnt;
627
628	+	VCOPY(verts[0],tri[0]);
629	+	VCOPY(verts[1],tri[1]);
630	+	VCOPY(verts[2],tri[2]);
631	+
632	+	l = add_data(NULL,0,NULL);
633	+	l = add_data(l,1,NULL);
634	+	l = add_data(l,2,NULL);
635	+	cnt = 3;
636	+	stVisit_clip(st,0,verts,&cnt,l,0,func);
637	+	}
638	+	}
639
640	<	int
641	<	stRemove_tri(st,id,v0,v1,v2)
642	<	STREE *st;
643	<	int id;
644	<	FVECT v0,v1,v2;
640	>
641	>	/* New Insertion code!!! */
642	>
643	>
644	>	BCOORD qtRoot[3][3] = { {MAXBCOORD2,0,0},{0,MAXBCOORD2,0},{0,0,MAXBCOORD2}};
645	>
646	>
647	>	convert_tri_to_frame(root,tri,b0,b1,b2,db10,db21,db02)
648	>	int root;
649	>	FVECT tri[3];
650	>	BCOORD b0[3],b1[3],b2[3];
651	>	BCOORD db10[3],db21[3],db02[3];
652		{
653	<
654	<	int i,found;
655	<	QUADTREE *rootptr;
656	<	FVECT t0,t1,t2;
653	>	/* Project the vertex into the qtree plane */
654	>	vert_to_qt_frame(root,tri[0],b0);
655	>	vert_to_qt_frame(root,tri[1],b1);
656	>	vert_to_qt_frame(root,tri[2],b2);
657
658	<	found = 0;
659	<	for(i=0; i < 4; i++)
660	<	{
661	<	rootptr = ST_NTH_ROOT_PTR(st,i);
261	<	stNth_base_verts(st,i,t0,t1,t2);
262	<	found |= qtRemove_tri(rootptr,id,v0,v1,v2,t0,t1,t2);
263	<	}
264	<	return(found);
658	>	/* calculate triangle edge differences in new frame */
659	>	db10[0] = b1[0] - b0[0]; db10[1] = b1[1] - b0[1]; db10[2] = b1[2] - b0[2];
660	>	db21[0] = b2[0] - b1[0]; db21[1] = b2[1] - b1[1]; db21[2] = b2[2] - b1[2];
661	>	db02[0] = b0[0] - b2[0]; db02[1] = b0[1] - b2[1]; db02[2] = b0[2] - b2[2];
662		}
663
664
665	+	QUADTREE
666	+	stRoot_insert_tri(st,root,tri,f)
667	+	STREE *st;
668	+	int root;
669	+	FVECT tri[3];
670	+	FUNC f;
671	+	{
672	+	BCOORD b0[3],b1[3],b2[3];
673	+	BCOORD db10[3],db21[3],db02[3];
674	+	unsigned int s0,s1,s2,sq0,sq1,sq2;
675	+	QUADTREE qt;
676
677	+	/* Map the triangle vertices into the canonical barycentric frame */
678	+	convert_tri_to_frame(root,tri,b0,b1,b2,db10,db21,db02);
679
680	+	/* Calculate initial sidedness info */
681	+	SIDES_GTR(b0,b1,b2,s0,s1,s2,qtRoot[1][0],qtRoot[0][1],qtRoot[0][2]);
682	+	SIDES_GTR(b0,b1,b2,sq0,sq1,sq2,qtRoot[0][0],qtRoot[1][1],qtRoot[2][2]);
683
684	<	#if 0
685	<	int
686	<	stAdd_tri_opt(st,id,v0,v1,v2)
687	<	STREE *st;
688	<	int id;
689	<	FVECT v0,v1,v2;
684	>	qt = ST_ROOT_QT(st,root);
685	>	/* Visit cells that triangle intersects */
686	>	qt = qtInsert_tri(root,qt,qtRoot[0],qtRoot[1],qtRoot[2],
687	>	b0,b1,b2,db10,db21,db02,MAXBCOORD2 >> 1,s0,s1,s2, sq0,sq1,sq2,f,0);
688	>
689	>	return(qt);
690	>	}
691	>
692	>	stRoot_visit_tri(st,root,tri,f)
693	>	STREE *st;
694	>	int root;
695	>	FVECT tri[3];
696	>	FUNC f;
697		{
698	+	BCOORD b0[3],b1[3],b2[3];
699	+	BCOORD db10[3],db21[3],db02[3];
700	+	unsigned int s0,s1,s2,sq0,sq1,sq2;
701	+	QUADTREE qt;
702	+
703	+	/* Map the triangle vertices into the canonical barycentric frame */
704	+	convert_tri_to_frame(root,tri,b0,b1,b2,db10,db21,db02);
705	+
706	+	/* Calculate initial sidedness info */
707	+	SIDES_GTR(b0,b1,b2,s0,s1,s2,qtRoot[1][0],qtRoot[0][1],qtRoot[0][2]);
708	+	SIDES_GTR(b0,b1,b2,sq0,sq1,sq2,qtRoot[0][0],qtRoot[1][1],qtRoot[2][2]);
709	+
710	+	qt = ST_ROOT_QT(st,root);
711	+	QT_SET_FLAG(ST_QT(st,root));
712	+	/* Visit cells that triangle intersects */
713	+	qtVisit_tri(root,qt,qtRoot[0],qtRoot[1],qtRoot[2],
714	+	b0,b1,b2,db10,db21,db02,MAXBCOORD2 >> 1,s0,s1,s2, sq0,sq1,sq2,f);
715	+
716	+	}
717	+
718	+	stInsert_tri(st,tri,f)
719	+	STREE *st;
720	+	FVECT tri[3];
721	+	FUNC f;
722	+	{
723	+	unsigned int cells,which;
724	+	int root;
725
279	–	int i,found;
280	–	QUADTREE *qtptr;
281	–	FVECT pt,t0,t1,t2;
726
727	<	/* First add all of the leaf cells lying on the triangle perimeter:
728	<	mark all cells seen on the way
285	<	*/
286	<	/* clear all of the flags */
287	<	qtClearAllFlags();            /* clear all quadtree branch flags */
727	>	/* calculate entry/exit points of edges through the cells */
728	>	cells = stTri_cells(st,tri);
729
730	<	/* Now trace each triangle edge-marking cells visited, and adding tri to
731	<	the leafs
732	<	*/
733	<	stAdd_tri_from_pt(st,v0,id,t0,t1,t2);
734	<	/* Find next cell that projection of ray intersects */
735	<	VCOPY(pt,v0);
295	<	/* NOTE: Check if in same cell */
296	<	while(traceEdge(pt,v1,t0,t1,t2,pt))
730	>	/* For each cell that quadtree intersects: Map the triangle vertices into
731	>	the canonical barycentric frame of (1,0,0), (0,1,0),(0,0,1). Insert
732	>	by first doing a trivial reject on the interior nodes, and then a
733	>	tri/tri intersection at the leaf nodes.
734	>	*/
735	>	for(root=0,which=1; root < ST_NUM_ROOT_NODES; root++,which <<= 1)
736		{
737	<	stAdd_tri_from_pt(st,pt,id,t0,t1,t2);
738	<	traceEdge(pt,v1,t0,t1,t2,pt);
737	>	/* For each of the quadtree roots: check if marked as intersecting tri*/
738	>	if(cells & which)
739	>	/* Visit tri cells */
740	>	ST_ROOT_QT(st,root) = stRoot_insert_tri(st,root,tri,f);
741		}
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	–
742		}
743
744	<	#endif
744	>
745	>
746	>
747	>
748	>
749	>
750	>
751	>
752	>
753	>

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines (old)
>	Changed lines (new)