ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/sm.c
Revision: 3.3
Committed: Mon Aug 24 12:38:57 1998 UTC (25 years, 8 months ago) by gwlarson
Content type: text/plain
Branch: MAIN
Changes since 3.2: +66 -18 lines
Log Message:
optimized triangle addition/deletion during edge swapping

File Contents

# Content
1 /* Copyright (c) 1998 Silicon Graphics, Inc. */
2
3 #ifndef lint
4 static char SCCSid[] = "$SunId$ SGI";
5 #endif
6
7 /*
8 * sm.c
9 */
10 #include "standard.h"
11
12 #include "object.h"
13
14 #include "sm_list.h"
15 #include "sm_geom.h"
16 #include "sm.h"
17
18
19 SM *smMesh = NULL;
20 double smDist_sum=0;
21 int smNew_tri_cnt=0;
22
23 #ifdef TEST_DRIVER
24 VIEW View = {0,{0,0,0},{0,0,-1},{0,1,0},60,60,0};
25 VIEW Current_View = {0,{0,0,0},{0,0,-1},{0,1,0},60,60,0};
26 int Pick_cnt =0;
27 int Pick_tri = -1,Picking = FALSE;
28 FVECT Pick_point[500],Pick_origin,Pick_dir;
29 FVECT Pick_v0[500], Pick_v1[500], Pick_v2[500];
30 FVECT P0,P1,P2;
31 FVECT FrustumNear[4],FrustumFar[4];
32 double dev_zmin=.01,dev_zmax=1000;
33 #endif
34
35 smDir(sm,ps,id)
36 SM *sm;
37 FVECT ps;
38 int id;
39 {
40 FVECT p;
41
42 VCOPY(p,SM_NTH_WV(sm,id));
43 point_on_sphere(ps,p,SM_VIEW_CENTER(sm));
44 }
45
46 smClear_mesh(sm)
47 SM *sm;
48 {
49 /* Reset the triangle counters */
50 SM_TRI_CNT(sm) = 0;
51 SM_NUM_TRIS(sm) = 0;
52 SM_FREE_TRIS(sm) = -1;
53 }
54
55 smClear_flags(sm,which)
56 SM *sm;
57 int which;
58 {
59 int i;
60
61 if(which== -1)
62 for(i=0; i < T_FLAGS;i++)
63 bzero(SM_NTH_FLAGS(sm,i),T_TOTAL_FLAG_BYTES(SM_MAX_TRIS(sm)));
64 else
65 bzero(SM_NTH_FLAGS(sm,which),T_TOTAL_FLAG_BYTES(SM_MAX_TRIS(sm)));
66 }
67
68 smClear_locator(sm)
69 SM *sm;
70 {
71 STREE *st;
72
73 st = SM_LOCATOR(sm);
74
75 stClear(st);
76 }
77
78 smInit_locator(sm,center,base)
79 SM *sm;
80 FVECT center,base[4];
81 {
82 STREE *st;
83
84 st = SM_LOCATOR(sm);
85
86 stInit(st,center,base);
87
88 }
89
90 smClear(sm)
91 SM *sm;
92 {
93 smClear_samples(sm);
94 smClear_mesh(sm);
95 smClear_locator(sm);
96 }
97
98 int
99 smAlloc_tris(sm,max_verts,max_tris)
100 SM *sm;
101 int max_verts,max_tris;
102 {
103 int i,nbytes,vbytes,fbytes;
104
105 vbytes = max_verts*sizeof(VERT);
106 fbytes = T_TOTAL_FLAG_BYTES(max_tris);
107 nbytes = vbytes + max_tris*sizeof(TRI) +T_FLAGS*fbytes + 8;
108 for(i = 1024; nbytes > i; i <<= 1)
109 ;
110 /* check if casting works correctly */
111 max_tris = (i-vbytes-8)/(sizeof(TRI) + T_FLAG_BYTES);
112 fbytes = T_TOTAL_FLAG_BYTES(max_tris);
113
114 SM_BASE(sm)=(char *)malloc(vbytes+max_tris*sizeof(TRI)+T_FLAGS*fbytes);
115
116 if (SM_BASE(sm) == NULL)
117 return(0);
118
119 SM_TRIS(sm) = (TRI *)SM_BASE(sm);
120 SM_VERTS(sm) = (VERT *)(SM_TRIS(sm) + max_tris);
121
122 SM_NTH_FLAGS(sm,0) = (int4 *)(SM_VERTS(sm) + max_verts);
123 for(i=1; i < T_FLAGS; i++)
124 SM_NTH_FLAGS(sm,i) = (int4 *)(SM_NTH_FLAGS(sm,i-1)+fbytes/sizeof(int4));
125
126 SM_MAX_VERTS(sm) = max_verts;
127 SM_MAX_TRIS(sm) = max_tris;
128
129 smClear_mesh(sm);
130
131 return(max_tris);
132 }
133
134
135
136 int
137 smAlloc_locator(sm)
138 SM *sm;
139 {
140 STREE *st;
141
142 st = SM_LOCATOR(sm);
143
144 st = stAlloc(st);
145
146 if(st)
147 return(TRUE);
148 else
149 return(FALSE);
150 }
151
152 /* Initialize/clear global smL sample list for at least n samples */
153 smAlloc(max_samples)
154 register int max_samples;
155 {
156 unsigned nbytes;
157 register unsigned i;
158 int total_points;
159 int max_tris;
160
161 /* If this is the first call, allocate sample,vertex and triangle lists */
162 if(!smMesh)
163 {
164 if(!(smMesh = (SM *)malloc(sizeof(SM))))
165 error(SYSTEM,"smAlloc():Unable to allocate memory");
166 bzero(smMesh,sizeof(SM));
167 }
168 else
169 { /* If existing structure: first deallocate */
170 if(SM_BASE(smMesh))
171 free(SM_BASE(smMesh));
172 if(SM_SAMP_BASE(smMesh))
173 free(SM_SAMP_BASE(smMesh));
174 }
175
176 /* First allocate at least n samples + extra points:at least enough
177 necessary to form the BASE MESH- Default = 4;
178 */
179 max_samples = smAlloc_samples(smMesh, max_samples,SM_EXTRA_POINTS);
180
181 total_points = max_samples + SM_EXTRA_POINTS;
182 max_tris = total_points*2;
183
184 /* Now allocate space for mesh vertices and triangles */
185 max_tris = smAlloc_tris(smMesh, total_points, max_tris);
186
187 /* Initialize the structure for point,triangle location.
188 */
189 smAlloc_locator(smMesh);
190
191 }
192
193
194
195 smInit_mesh(sm,vp)
196 SM *sm;
197 FVECT vp;
198 {
199
200 /* NOTE: Should be elsewhere?*/
201 smDist_sum = 0;
202 smNew_tri_cnt = 0;
203
204 VCOPY(SM_VIEW_CENTER(smMesh),vp);
205 smClear_locator(sm);
206 smInit_locator(sm,vp,0);
207 smClear_aux_samples(sm);
208 smClear_mesh(sm);
209 smCreate_base_mesh(sm,SM_DEFAULT);
210 }
211
212 /*
213 * int
214 * smInit(n) : Initialize/clear data structures for n entries
215 * int n;
216 *
217 * This routine allocates/initializes the sample, mesh, and point-location
218 * structures for at least n samples.
219 * If n is <= 0, then clear data structures. Returns number samples
220 * actually allocated.
221 */
222
223 int
224 smInit(n)
225 register int n;
226 {
227 int max_vertices;
228
229 sleep(10);
230
231 /* If n <=0, Just clear the existing structures */
232 if(n <= 0)
233 {
234 smClear(smMesh);
235 return(0);
236 }
237
238 /* Total mesh vertices includes the sample points and the extra vertices
239 to form the base mesh
240 */
241 max_vertices = n + SM_EXTRA_POINTS;
242
243 /* If the current mesh contains enough room, clear and return */
244 if(smMesh && max_vertices <= SM_MAX_VERTS(smMesh))
245 {
246 smClear(smMesh);
247 return(SM_MAX_SAMP(smMesh));
248 }
249 /* Otherwise- mesh must be allocated with the appropriate number of
250 samples
251 */
252 smAlloc(n);
253
254 return(SM_MAX_SAMP(smMesh));
255 }
256
257
258 int
259 smLocator_apply_func(sm,v0,v1,v2,func,arg)
260 SM *sm;
261 FVECT v0,v1,v2;
262 int (*func)();
263 char *arg;
264 {
265 STREE *st;
266 char found;
267 FVECT p0,p1,p2;
268
269 st = SM_LOCATOR(sm);
270
271 point_on_sphere(p0,v0,SM_VIEW_CENTER(sm));
272 point_on_sphere(p1,v1,SM_VIEW_CENTER(sm));
273 point_on_sphere(p2,v2,SM_VIEW_CENTER(sm));
274
275 found = stApply_to_tri_cells(st,p0,p1,p2,func,arg);
276
277 return(found);
278 }
279
280
281 int
282 smLocator_add_tri(sm,t_id,v0_id,v1_id,v2_id)
283 SM *sm;
284 int t_id;
285 int v0_id,v1_id,v2_id;
286 {
287 STREE *st;
288 FVECT p0,p1,p2;
289
290 st = SM_LOCATOR(sm);
291
292 smDir(sm,p0,v0_id);
293 smDir(sm,p1,v1_id);
294 smDir(sm,p2,v2_id);
295
296 stAdd_tri(st,t_id,p0,p1,p2);
297
298 return(1);
299 }
300
301 /* Add a triangle to the base array with vertices v1-v2-v3 */
302 int
303 smAdd_tri(sm, v0_id,v1_id,v2_id,tptr)
304 SM *sm;
305 int v0_id,v1_id,v2_id;
306 TRI **tptr;
307 {
308 int t_id;
309 TRI *t;
310
311
312 if(SM_TRI_CNT(sm)+1 > SM_MAX_TRIS(sm))
313 error(SYSTEM,"smAdd_tri():Too many triangles");
314
315 /* Get the id for the next available triangle */
316 SM_FREE_TRI_ID(sm,t_id);
317 if(t_id == -1)
318 return(t_id);
319
320 t = SM_NTH_TRI(sm,t_id);
321 T_CLEAR_NBRS(t);
322
323 if(SM_BASE_ID(sm,v0_id) || SM_BASE_ID(sm,v1_id) || SM_BASE_ID(sm,v2_id))
324 {
325 smClear_tri_flags(sm,t_id);
326 SM_SET_NTH_T_BASE(sm,t_id);
327 }
328 else
329 {
330 SM_CLEAR_NTH_T_BASE(sm,t_id);
331 SM_SET_NTH_T_ACTIVE(sm,t_id);
332 SM_SET_NTH_T_LRU(sm,t_id);
333 SM_SET_NTH_T_NEW(sm,t_id);
334 SM_NUM_TRIS(sm)++;
335 smNew_tri_cnt++;
336 }
337 /* set the triangle vertex ids */
338 T_NTH_V(t,0) = v0_id;
339 T_NTH_V(t,1) = v1_id;
340 T_NTH_V(t,2) = v2_id;
341
342 SM_NTH_VERT(sm,v0_id) = t_id;
343 SM_NTH_VERT(sm,v1_id) = t_id;
344 SM_NTH_VERT(sm,v2_id) = t_id;
345
346 if(t)
347 *tptr = t;
348 /* return initialized triangle */
349 return(t_id);
350 }
351
352 int
353 smClosest_vertex_in_tri(sm,v0_id,v1_id,v2_id,p,eps)
354 SM *sm;
355 int v0_id,v1_id,v2_id;
356 FVECT p;
357 double eps;
358 {
359 FVECT v;
360 double d,d0,d1,d2;
361 int closest = -1;
362
363 if(v0_id != -1)
364 {
365 smDir(sm,v,v0_id);
366 d0 = DIST(p,v);
367 if(eps < 0 || d0 < eps)
368 {
369 closest = v0_id;
370 d = d0;
371 }
372 }
373 if(v1_id != -1)
374 {
375 smDir(sm,v,v1_id);
376 d1 = DIST(p,v);
377 if(closest== -1)
378 {
379 if(eps < 0 || d1 < eps)
380 {
381 closest = v1_id;
382 d = d1;
383 }
384 }
385 else
386 if(d1 < d2)
387 {
388 if((eps < 0) || d1 < eps)
389 {
390 closest = v1_id;
391 d = d1;
392 }
393 }
394 else
395 if((eps < 0) || d2 < eps)
396 {
397 closest = v2_id;
398 d = d2;
399 }
400 }
401 if(v2_id != -1)
402 {
403 smDir(sm,v,v2_id);
404 d2 = DIST(p,v);
405 if((eps < 0) || d2 < eps)
406 if(closest== -1 ||(d2 < d))
407 return(v2_id);
408 }
409 return(closest);
410 }
411
412
413 int
414 smClosest_vertex_in_w_tri(sm,v0_id,v1_id,v2_id,p)
415 SM *sm;
416 int v0_id,v1_id,v2_id;
417 FVECT p;
418 {
419 FVECT v;
420 double d,d0,d1,d2;
421 int closest;
422
423 VCOPY(v,SM_NTH_WV(sm,v0_id));
424 d = d0 = DIST(p,v);
425 closest = v0_id;
426
427 VCOPY(v,SM_NTH_WV(sm,v1_id));
428 d1 = DIST(p,v);
429 if(d1 < d2)
430 {
431 closest = v1_id;
432 d = d1;
433 }
434 VCOPY(v,SM_NTH_WV(sm,v2_id));
435 d2 = DIST(p,v);
436 if(d2 < d)
437 return(v2_id);
438 else
439 return(closest);
440 }
441
442 void
443 smTris_swap_edge(sm,t_id,t1_id,e,e1,tn_id,tn1_id,add,del)
444 SM *sm;
445 int t_id,t1_id;
446 int e,e1;
447 int **tn_id,**tn1_id;
448 LIST **add,**del;
449
450 {
451 TRI *t,*t1;
452 TRI *ta,*tb;
453 int verts[3],enext,eprev,e1next,e1prev;
454 TRI *n;
455 FVECT p1,p2,p3;
456 int ta_id,tb_id;
457 /* swap diagonal (e relative to t, and e1 relative to t1)
458 defined by quadrilateral
459 formed by t,t1- swap for the opposite diagonal
460 */
461 t = SM_NTH_TRI(sm,t_id);
462 t1 = SM_NTH_TRI(sm,t1_id);
463 enext = (e+1)%3;
464 eprev = (e+2)%3;
465 e1next = (e1+1)%3;
466 e1prev = (e1+2)%3;
467 verts[e] = T_NTH_V(t,e);
468 verts[enext] = T_NTH_V(t1,e1prev);
469 verts[eprev] = T_NTH_V(t,eprev);
470 ta_id = smAdd_tri(sm,verts[0],verts[1],verts[2],&ta);
471 *add = push_data(*add,ta_id);
472
473 verts[e1] = T_NTH_V(t1,e1);
474 verts[e1next] = T_NTH_V(t,eprev);
475 verts[e1prev] = T_NTH_V(t1,e1prev);
476 tb_id = smAdd_tri(sm,verts[0],verts[1],verts[2],&tb);
477 *add = push_data(*add,tb_id);
478
479 /* set the neighbors */
480 T_NTH_NBR(ta,e) = T_NTH_NBR(t1,e1next);
481 T_NTH_NBR(tb,e1) = T_NTH_NBR(t,enext);
482 T_NTH_NBR(ta,enext) = tb_id;
483 T_NTH_NBR(tb,e1next) = ta_id;
484 T_NTH_NBR(ta,eprev) = T_NTH_NBR(t,eprev);
485 T_NTH_NBR(tb,e1prev) = T_NTH_NBR(t1,e1prev);
486
487 /* Reset neighbor pointers of original neighbors */
488 n = SM_NTH_TRI(sm,T_NTH_NBR(t,enext));
489 T_NTH_NBR(n,T_NTH_NBR_PTR(t_id,n)) = tb_id;
490 n = SM_NTH_TRI(sm,T_NTH_NBR(t,eprev));
491 T_NTH_NBR(n,T_NTH_NBR_PTR(t_id,n)) = ta_id;
492
493 n = SM_NTH_TRI(sm,T_NTH_NBR(t1,e1next));
494 T_NTH_NBR(n,T_NTH_NBR_PTR(t1_id,n)) = ta_id;
495 n = SM_NTH_TRI(sm,T_NTH_NBR(t1,e1prev));
496 T_NTH_NBR(n,T_NTH_NBR_PTR(t1_id,n)) = tb_id;
497
498 /* Delete two parent triangles */
499 *del = push_data(*del,t_id);
500 if(SM_IS_NTH_T_NEW(sm,t_id))
501 SM_CLEAR_NTH_T_NEW(sm,t_id);
502 else
503 SM_CLEAR_NTH_T_BASE(sm,t_id);
504 *del = push_data(*del,t1_id);
505 if(SM_IS_NTH_T_NEW(sm,t1_id))
506 SM_CLEAR_NTH_T_NEW(sm,t1_id);
507 else
508 SM_CLEAR_NTH_T_BASE(sm,t1_id);
509 *tn_id = ta_id;
510 *tn1_id = tb_id;
511 }
512
513 smUpdate_locator(sm,add_list,del_list)
514 SM *sm;
515 LIST *add_list,*del_list;
516 {
517 int t_id;
518 TRI *t;
519 while(add_list)
520 {
521 t_id = pop_list(&add_list);
522 if(!SM_IS_NTH_T_NEW(sm,t_id) && !SM_IS_NTH_T_BASE(sm,t_id))
523 {
524 SM_SET_NTH_T_NEW(sm,t_id);
525 continue;
526 }
527 t = SM_NTH_TRI(sm,t_id);
528 smLocator_add_tri(sm,t_id,T_NTH_V(t,0),T_NTH_V(t,1),T_NTH_V(t,2));
529 }
530
531 while(del_list)
532 {
533 t_id = pop_list(&del_list);
534 if(SM_IS_NTH_T_NEW(sm,t_id))
535 {
536 smDelete_tri(sm,t_id);
537 continue;
538 }
539 smLocator_remove_tri(sm,t_id);
540 smDelete_tri(sm,t_id);
541 }
542 }
543 /* MUST add check for constrained edges */
544 int
545 smFix_tris(sm,id,tlist)
546 SM *sm;
547 int id;
548 LIST *tlist;
549 {
550 TRI *t,*t_opp;
551 FVECT p,p1,p2,p3;
552 int e,e1,swapped = 0;
553 int t_id,t_opp_id;
554 LIST *add_list,*del_list;
555
556
557 add_list = del_list = NULL;
558 VSUB(p,SM_NTH_WV(sm,id),SM_VIEW_CENTER(sm));
559 while(tlist)
560 {
561 t_id = pop_list(&tlist);
562 t = SM_NTH_TRI(sm,t_id);
563 e = (T_WHICH_V(t,id)+1)%3;
564 t_opp_id = T_NTH_NBR(t,e);
565 t_opp = SM_NTH_TRI(sm,t_opp_id);
566
567 smDir(sm,p1,T_NTH_V(t_opp,0));
568 smDir(sm,p2,T_NTH_V(t_opp,1));
569 smDir(sm,p3,T_NTH_V(t_opp,2));
570 if(point_in_cone(p,p1,p2,p3))
571 {
572 swapped = 1;
573 e1 = T_NTH_NBR_PTR(t_id,t_opp);
574 /* check list for t_opp and Remove if there */
575 remove_from_list(t_opp_id,&tlist);
576 smTris_swap_edge(sm,t_id,t_opp_id,e,e1,&t_id,&t_opp_id,
577 &add_list,&del_list);
578 tlist = push_data(tlist,t_id);
579 tlist = push_data(tlist,t_opp_id);
580 }
581 }
582 smUpdate_locator(sm,add_list,del_list);
583
584 return(swapped);
585 }
586
587 /* Give the vertex "id" and a triangle "t" that it belongs to- return the
588 next nbr in a counter clockwise order about vertex "id"
589 */
590 int
591 smTri_next_ccw_nbr(sm,t,id)
592 SM *sm;
593 TRI *t;
594 int id;
595 {
596 int t_id;
597 int tri;
598
599 /* Want the edge for which "id" is the destination */
600 t_id = (T_WHICH_V(t,id)+ 2)% 3;
601 tri = T_NTH_NBR(t,t_id);
602 return(tri);
603 }
604
605 void
606 smReplace_point(sm,tri,id,nid)
607 SM *sm;
608 TRI *tri;
609 int id,nid;
610 {
611 TRI *t;
612 int t_id;
613
614 T_NTH_V(tri,T_WHICH_V(tri,id)) = nid;
615
616 t_id = smTri_next_ccw_nbr(sm,t,nid);
617 while((t= SM_NTH_TRI(sm,t_id)) != tri)
618 {
619 T_NTH_V(t,T_WHICH_V(t,id)) = nid;
620 t_id = smTri_next_ccw_nbr(sm,t,nid);
621 }
622 }
623
624
625 smClear_tri_flags(sm,id)
626 SM *sm;
627 int id;
628 {
629 int i;
630
631 for(i=0; i < T_FLAGS; i++)
632 SM_CLEAR_NTH_T_FLAG(sm,id,i);
633
634 }
635
636 /* Locate the point-id in the point location structure: */
637 int
638 smReplace_vertex(sm,c,dir,p,tri_id,snew_id,type,which)
639 SM *sm;
640 COLR c;
641 FVECT dir,p;
642 int tri_id,snew_id;
643 char type,which;
644 {
645 int n_id,s_id;
646 TRI *tri;
647
648 tri = SM_NTH_TRI(sm,tri_id);
649 /* Get the sample that corresponds to the "which" vertex of "tri" */
650 /* NEED: have re-init that sets clock flag */
651 /* If this is a base-sample: create a new sample and replace
652 all references to the base sample with references to the
653 new sample
654 */
655 s_id = T_NTH_V(tri,which);
656 if(SM_BASE_ID(sm,s_id))
657 {
658 if(snew_id != -1)
659 n_id = smAdd_sample_point(sm,c,dir,p);
660 else
661 n_id = snew_id;
662 smReplace_point(sm,tri,s_id,n_id);
663 s_id = n_id;
664 }
665 else /* If the sample exists, reset the values */
666 /* NOTE: This test was based on the SPHERICAL coordinates
667 of the point: If we are doing a multiple-sample-per
668 SPHERICAL pixel: we will want to test for equality-
669 and do other processing: for now: SINGLE SAMPLE PER
670 PIXEL
671 */
672 /* NOTE: snew_id needs to be marked as invalid?*/
673 if(snew_id == -1)
674 smInit_sample(sm,s_id,c,dir,p);
675 else
676 smReset_sample(sm,s_id,snew_id);
677 return(s_id);
678 }
679
680
681 /* Locate the point-id in the point location structure: */
682 int
683 smInsert_point_in_tri(sm,c,dir,p,s_id,tri_id)
684 SM *sm;
685 COLR c;
686 FVECT dir,p;
687 int s_id,tri_id;
688 {
689 TRI *tri,*t0,*t1,*t2,*nbr;
690 int v0_id,v1_id,v2_id,n_id;
691 int t0_id,t1_id,t2_id;
692 LIST *tlist;
693 FVECT npt;
694
695 if(s_id == SM_INVALID)
696 s_id = smAdd_sample_point(sm,c,dir,p);
697
698 tri = SM_NTH_TRI(sm,tri_id);
699 v0_id = T_NTH_V(tri,0);
700 v1_id = T_NTH_V(tri,1);
701 v2_id = T_NTH_V(tri,2);
702
703 n_id = -1;
704 if(SM_BASE_ID(sm,v0_id)||SM_BASE_ID(sm,v1_id)||SM_BASE_ID(sm,v2_id))
705 {
706 smDir(sm,npt,s_id);
707 /* Change to an add and delete */
708 t0_id = (SM_BASE_ID(sm,v0_id))?v0_id:-1;
709 t1_id = (SM_BASE_ID(sm,v1_id))?v1_id:-1;
710 t2_id = (SM_BASE_ID(sm,v2_id))?v2_id:-1;
711 n_id = smClosest_vertex_in_tri(sm,t0_id,t1_id,t2_id,npt,P_REPLACE_EPS);
712 }
713 t0_id = smAdd_tri(sm,s_id,v0_id,v1_id,&t0);
714 /* Add triangle to the locator */
715 smLocator_add_tri(sm,t0_id,s_id,v0_id,v1_id);
716
717 t1_id = smAdd_tri(sm,s_id,v1_id,v2_id,&t1);
718 smLocator_add_tri(sm,t1_id,s_id,v1_id,v2_id);
719 t2_id = smAdd_tri(sm,s_id,v2_id,v0_id,&t2);
720 smLocator_add_tri(sm,t2_id,s_id,v2_id,v0_id);
721
722 /* Set the neighbor pointers for the new tris */
723 T_NTH_NBR(t0,0) = t2_id;
724 T_NTH_NBR(t0,1) = T_NTH_NBR(tri,0);
725 T_NTH_NBR(t0,2) = t1_id;
726 T_NTH_NBR(t1,0) = t0_id;
727 T_NTH_NBR(t1,1) = T_NTH_NBR(tri,1);
728 T_NTH_NBR(t1,2) = t2_id;
729 T_NTH_NBR(t2,0) = t1_id;
730 T_NTH_NBR(t2,1) = T_NTH_NBR(tri,2);
731 T_NTH_NBR(t2,2) = t0_id;
732
733 /* Reset the neigbor pointers for the neighbors of the original */
734 nbr = SM_NTH_TRI(sm,T_NTH_NBR(tri,0));
735 T_NTH_NBR(nbr,T_NTH_NBR_PTR(tri_id,nbr)) = t0_id;
736 nbr = SM_NTH_TRI(sm,T_NTH_NBR(tri,1));
737 T_NTH_NBR(nbr,T_NTH_NBR_PTR(tri_id,nbr)) = t1_id;
738 nbr = SM_NTH_TRI(sm,T_NTH_NBR(tri,2));
739 T_NTH_NBR(nbr,T_NTH_NBR_PTR(tri_id,nbr)) = t2_id;
740
741 smLocator_remove_tri(sm,tri_id);
742 smDelete_tri(sm,tri_id);
743
744 /* Fix up the new triangles*/
745 tlist = push_data(NULL,t0_id);
746 tlist = push_data(tlist,t1_id);
747 tlist = push_data(tlist,t2_id);
748
749 smFix_tris(sm,s_id,tlist);
750
751 if(n_id != -1)
752 smDelete_point(sm,n_id);
753
754 return(s_id);
755 }
756
757
758 int
759 smPointLocate(sm,pt,type,which,norm)
760 SM *sm;
761 FVECT pt;
762 char *type,*which;
763 char norm;
764 {
765 STREE *st;
766 int tri;
767 FVECT npt;
768
769 st = SM_LOCATOR(sm);
770 if(norm)
771 {
772 point_on_sphere(npt,pt,SM_VIEW_CENTER(sm));
773 tri = stPoint_locate(st,npt,type,which);
774 }
775 else
776 tri = stPoint_locate(st,pt,type,which);
777 return(tri);
778 }
779
780 QUADTREE
781 smPointLocateCell(sm,pt,type,which,norm)
782 SM *sm;
783 FVECT pt;
784 char *type,*which;
785 char norm;
786 {
787 STREE *st;
788 QUADTREE qt;
789 FVECT npt;
790
791 st = SM_LOCATOR(sm);
792 if(norm)
793 {
794 point_on_sphere(npt,pt,SM_VIEW_CENTER(sm));
795
796 qt = stPoint_locate_cell(st,npt,type,which);
797 }
798 else
799 qt = stPoint_locate_cell(st,pt,type,which);
800
801 return(qt);
802 }
803
804 int
805 smAdd_sample_to_mesh(sm,c,dir,pt,s_id)
806 SM *sm;
807 COLR c;
808 FVECT dir,pt;
809 int s_id;
810 {
811 int t_id;
812 char type,which;
813 double d;
814 FVECT p;
815
816 /* If new, foreground pt */
817 if(pt)
818 {
819 /* NOTE: This should be elsewhere! */
820 d = DIST(pt,SM_VIEW_CENTER(smMesh));
821 smDist_sum += 1.0/d;
822 /************************************/
823 t_id = smPointLocate(smMesh,pt,&type,&which,TRUE);
824 if(type==GT_FACE)
825 s_id = smInsert_point_in_tri(smMesh,c,dir,pt,s_id,t_id);
826 else
827 if(type==GT_VERTEX)
828 s_id = smReplace_vertex(smMesh,c,dir,pt,t_id,s_id,type,which);
829 #ifdef DEBUG
830 else
831 eputs("smAdd_sample_to_mesh(): unrecognized type\n");
832 #endif
833 }
834 else if(s_id != -1)
835 {
836 VCOPY(p,SM_NTH_WV(sm,s_id));
837 if(SM_NTH_W_DIR(sm,s_id) != -1)
838 {
839 /* NOTE: This should be elsewhere! */
840 d = DIST(p,SM_VIEW_CENTER(smMesh));
841 smDist_sum += 1.0/d;
842 /************************************/
843 }
844 t_id = smPointLocate(smMesh,p,&type,&which,TRUE);
845 if(type==GT_FACE)
846 s_id = smInsert_point_in_tri(smMesh,c,dir,p,s_id,t_id);
847 else
848 if(type==GT_VERTEX)
849 s_id = smReplace_vertex(smMesh,c,dir,p,t_id,s_id,type,which);
850 #ifdef DEBUG
851 else
852 eputs("smAdd_sample_to_mesh(): unrecognized type\n");
853 #endif
854 }
855 /* Is a BG(sky point) */
856 else
857 {
858 t_id = smPointLocate(smMesh,dir,&type,&which,FALSE);
859 if(type==GT_FACE)
860 s_id = smInsert_point_in_tri(smMesh,c,dir,NULL,s_id,t_id);
861 else
862 if(type==GT_VERTEX)
863 s_id = smReplace_vertex(smMesh,c,dir,NULL,t_id,s_id,type,which);
864 #ifdef DEBUG
865 else
866 eputs("smAdd_sample_to_mesh(): unrecognized type\n");
867 #endif
868 }
869 return(s_id);
870 }
871
872 /*
873 * int
874 * smNewSamp(c, dir, p) : register new sample point and return index
875 * COLR c; : pixel color (RGBE)
876 * FVECT dir; : ray direction vector
877 * FVECT p; : world intersection point
878 *
879 * Add new sample point to data structures, removing old values as necessary.
880 * New sample representation will be output in next call to smUpdate().
881 * If the point is a sky point: then v=NULL
882 */
883 int
884 smNewSamp(c,dir,p)
885 COLR c;
886 FVECT dir;
887 FVECT p;
888
889 {
890 int s_id;
891
892 /* First check if this the first sample: if so initialize mesh */
893 if(SM_NUM_SAMP(smMesh) == 0)
894 #ifdef TEST_DRIVER
895 smInit_mesh(smMesh,View.vp);
896 #else
897 smInit_mesh(smMesh,odev.v.vp);
898 #endif
899 s_id = smAdd_sample_to_mesh(smMesh,c,dir,p,-1);
900 #if 0
901 {
902 char buff[100];
903 sprintf(buff,"Added sample %d\n",s_id);
904 eputs(buff);
905 }
906 #endif
907 return(s_id);
908
909 }
910 /*
911 * int
912 * smFindsamp(orig, dir): intersect ray with 3D rep. and find closest sample
913 * FVECT orig, dir;
914 *
915 * Find the closest sample to the given ray. Return -1 on failure.
916 */
917
918 /*
919 * smClean() : display has been wiped clean
920 *
921 * Called after display has been effectively cleared, meaning that all
922 * geometry must be resent down the pipeline in the next call to smUpdate().
923 */
924
925
926 /*
927 * smUpdate(vp, qua) : update OpenGL output geometry for view vp
928 * VIEW *vp; : desired view
929 * int qua; : quality level (percentage on linear time scale)
930 *
931 * Draw new geometric representation using OpenGL calls. Assume that the
932 * view has already been set up and the correct frame buffer has been
933 * selected for drawing. The quality level is on a linear scale, where 100%
934 * is full (final) quality. It is not necessary to redraw geometry that has
935 * been output since the last call to smClean().
936 */
937
938
939 int
940 smClear_vert(sm,id)
941 SM *sm;
942 int id;
943 {
944 if(SM_INVALID_POINT_ID(sm,id))
945 return(FALSE);
946
947 SM_NTH_VERT(sm,id) = SM_INVALID;
948
949 return(TRUE);
950 }
951
952 int
953 smAdd_base_vertex(sm,v,d)
954 SM *sm;
955 FVECT v,d;
956 {
957 int id;
958
959 /* First add coordinate to the sample array */
960 id = smAdd_aux_point(sm,v,d);
961 if(id == -1)
962 return(SM_INVALID);
963 /* Initialize triangle pointer to -1 */
964 smClear_vert(sm,id);
965 return(id);
966 }
967
968
969
970 /* Initialize a the point location DAG based on a 6 triangle tesselation
971 of the unit sphere centered on the view center. The DAG structure
972 contains 6 roots: one for each initial base triangle
973 */
974
975 int
976 smCreate_base_mesh(sm,type)
977 SM *sm;
978 int type;
979 {
980 int i,id;
981 int p[4],ids[4];
982 int v0_id,v1_id,v2_id;
983 TRI *tris[4];
984 FVECT d,pt,cntr;
985
986 /* First insert the base vertices into the sample point array */
987
988 for(i=0; i < 4; i++)
989 {
990 VADD(cntr,stDefault_base[i],SM_VIEW_CENTER(sm));
991 point_on_sphere(d,cntr,SM_VIEW_CENTER(sm));
992 id = smAdd_base_vertex(sm,cntr,d);
993 /* test to make sure vertex allocated */
994 if(id != -1)
995 p[i] = id;
996 else
997 return(0);
998 }
999 /* Create the base triangles */
1000 for(i=0;i < 4; i++)
1001 {
1002 v0_id = p[stTri_verts[i][0]];
1003 v1_id = p[stTri_verts[i][1]];
1004 v2_id = p[stTri_verts[i][2]];
1005 if((ids[i] = smAdd_tri(sm, v0_id,v1_id,v2_id,&(tris[i])))== -1)
1006 return(0);
1007 smLocator_add_tri(sm,ids[i],v0_id,v1_id,v2_id);
1008 }
1009 /* Set neighbors */
1010
1011 T_NTH_NBR(tris[0],0) = ids[3];
1012 T_NTH_NBR(tris[0],1) = ids[2];
1013 T_NTH_NBR(tris[0],2) = ids[1];
1014
1015 T_NTH_NBR(tris[1],0) = ids[3];
1016 T_NTH_NBR(tris[1],1) = ids[0];
1017 T_NTH_NBR(tris[1],2) = ids[2];
1018
1019 T_NTH_NBR(tris[2],0) = ids[3];
1020 T_NTH_NBR(tris[2],1) = ids[1];
1021 T_NTH_NBR(tris[2],2) = ids[0];
1022
1023 T_NTH_NBR(tris[3],0) = ids[1];
1024 T_NTH_NBR(tris[3],1) = ids[2];
1025 T_NTH_NBR(tris[3],2) = ids[0];
1026 return(1);
1027
1028 }
1029
1030
1031 int
1032 smNext_tri_flag_set(sm,i,which,b)
1033 SM *sm;
1034 int i,which;
1035 char b;
1036 {
1037
1038 for(; i < SM_TRI_CNT(sm);i++)
1039 {
1040 if(!SM_IS_NTH_T_FLAG(sm,i,which))
1041 continue;
1042
1043 if(!b)
1044 break;
1045 if((b==1) && !SM_BG_TRI(sm,i))
1046 break;
1047 if((b==2) && SM_BG_TRI(sm,i))
1048 break;
1049 }
1050
1051 return(i);
1052 }
1053
1054
1055 int
1056 smNext_valid_tri(sm,i)
1057 SM *sm;
1058 int i;
1059 {
1060
1061 while( i < SM_TRI_CNT(sm) && !T_IS_VALID(SM_NTH_TRI(sm,i)))
1062 i++;
1063
1064 return(i);
1065 }
1066
1067
1068
1069 qtTri_verts_from_id(t_id,v0,v1,v2)
1070 int t_id;
1071 FVECT v0,v1,v2;
1072 {
1073 TRI *t;
1074 int v0_id,v1_id,v2_id;
1075
1076 t = SM_NTH_TRI(smMesh,t_id);
1077 v0_id = T_NTH_V(t,0);
1078 v1_id = T_NTH_V(t,1);
1079 v2_id = T_NTH_V(t,2);
1080
1081 smDir(smMesh,v0,v0_id);
1082 smDir(smMesh,v1,v1_id);
1083 smDir(smMesh,v2,v2_id);
1084
1085 }
1086
1087
1088 int
1089 smIntersectTriSet(sm,t_set,orig,dir,pt)
1090 SM *sm;
1091 OBJECT *t_set;
1092 FVECT orig,dir,pt;
1093 {
1094 OBJECT *optr;
1095 int i,t_id,v_id;
1096 TRI *tri;
1097 FVECT p0,p1,p2;
1098 char type,which;
1099 int p0_id,p1_id,p2_id;
1100
1101 for(optr = QT_SET_PTR(t_set),i = QT_SET_CNT(t_set); i > 0; i--)
1102 {
1103 t_id = QT_SET_NEXT_ELEM(optr);
1104 tri = SM_NTH_TRI(sm,t_id);
1105 p0_id = T_NTH_V(tri,0);
1106 p1_id = T_NTH_V(tri,1);
1107 p2_id = T_NTH_V(tri,2);
1108 VCOPY(p0,SM_NTH_WV(sm,p0_id));
1109 VCOPY(p1,SM_NTH_WV(sm,p1_id));
1110 VCOPY(p2,SM_NTH_WV(sm,p2_id));
1111 if(type = ray_intersect_tri(orig,dir,p0,p1,p2,pt,&which))
1112 {
1113 if(type==GT_VERTEX)
1114 return(T_NTH_V(tri,which));
1115 v_id = smClosest_vertex_in_w_tri(sm,p0_id,p1_id,p2_id,pt);
1116 return(v_id);
1117 }
1118 }
1119 return(-1);
1120 }
1121
1122 /*
1123 * int
1124 * smTraceRay(SM *sm,FVECT orig, FVECT dir,FVECT v0,FVECT v1,FVECT v2,FVECT r)
1125 *
1126 * Intersect the ray with triangle v0v1v2, return intersection point in r
1127 *
1128 * Assumes orig,v0,v1,v2 are in spherical coordinates, and orig is
1129 * unit
1130 */
1131 char
1132 smTraceRay(sm,orig,dir,v0,v1,v2,r)
1133 SM *sm;
1134 FVECT orig,dir;
1135 FVECT v0,v1,v2;
1136 FVECT r;
1137 {
1138 FVECT n,p[3],d;
1139 double pt[3],r_eps;
1140 char i;
1141 int which;
1142
1143 /* Find the plane equation for the triangle defined by the edge v0v1 and
1144 the view center*/
1145 VCROSS(n,v0,v1);
1146 /* Intersect the ray with this plane */
1147 i = intersect_ray_plane(orig,dir,n,0.0,&(pt[0]),p[0]);
1148 if(i)
1149 which = 0;
1150 else
1151 which = -1;
1152
1153 VCROSS(n,v1,v2);
1154 i = intersect_ray_plane(orig,dir,n,0.0,&(pt[1]),p[1]);
1155 if(i)
1156 if((which==-1) || pt[1] < pt[0])
1157 which = 1;
1158
1159 VCROSS(n,v2,v0);
1160 i = intersect_ray_plane(orig,dir,n,0.0,&(pt[2]),p[2]);
1161 if(i)
1162 if((which==-1) || pt[2] < pt[which])
1163 which = 2;
1164
1165 if(which != -1)
1166 {
1167 /* Return point that is K*eps along projection of the ray on
1168 the sphere to push intersection point p[which] into next cell
1169 */
1170 normalize(p[which]);
1171 /* Calculate the ray perpendicular to the intersection point: approx
1172 the direction moved along the sphere a distance of K*epsilon*/
1173 r_eps = -DOT(p[which],dir);
1174 VSUM(n,dir,p[which],r_eps);
1175 /* Calculate the length along ray p[which]-dir needed to move to
1176 cause a move along the sphere of k*epsilon
1177 */
1178 r_eps = DOT(n,dir);
1179 VSUM(r,p[which],dir,(20*FTINY)/r_eps);
1180 normalize(r);
1181 return(TRUE);
1182 }
1183 else
1184 {
1185 /* Unable to find intersection: move along ray and try again */
1186 r_eps = -DOT(orig,dir);
1187 VSUM(n,dir,orig,r_eps);
1188 r_eps = DOT(n,dir);
1189 VSUM(r,orig,dir,(20*FTINY)/r_eps);
1190 normalize(r);
1191 #ifdef DEBUG
1192 eputs("smTraceRay:Ray does not intersect triangle");
1193 #endif
1194 return(FALSE);
1195 }
1196 }
1197
1198
1199 /*
1200 * int
1201 * smFindSamp(FVECT orig, FVECT dir)
1202 *
1203 * Find the closest sample to the given ray. Returns sample id, -1 on failure.
1204 * "dir" is assumed to be normalized
1205 */
1206 int
1207 smFindSamp(orig,dir)
1208 FVECT orig,dir;
1209 {
1210 FVECT r,v0,v1,v2,a,b,p;
1211 OBJECT os[MAXCSET+1],t_set[MAXSET+1];
1212 QUADTREE qt;
1213 int s_id;
1214 double d;
1215
1216 /* r is the normalized vector from the view center to the current
1217 * ray point ( starting with "orig"). Find the cell that r falls in,
1218 * and test the ray against all triangles stored in the cell. If
1219 * the test fails, trace the projection of the ray across to the
1220 * next cell it intersects: iterate until either an intersection
1221 * is found, or the projection ray is // to the direction. The sample
1222 * corresponding to the triangle vertex closest to the intersection
1223 * point is returned.
1224 */
1225
1226 /* First test if "orig" coincides with the View_center or if "dir" is
1227 parallel to r formed by projecting "orig" on the sphere. In
1228 either case, do a single test against the cell containing the
1229 intersection of "dir" and the sphere
1230 */
1231 point_on_sphere(b,orig,SM_VIEW_CENTER(smMesh));
1232 d = -DOT(b,dir);
1233 if(EQUAL_VEC3(orig,SM_VIEW_CENTER(smMesh)) || EQUAL(fabs(d),1.0))
1234 {
1235 qt = smPointLocateCell(smMesh,dir,NULL,NULL,FALSE);
1236 /* Test triangles in the set for intersection with Ray:returns
1237 first found
1238 */
1239 qtgetset(t_set,qt);
1240 s_id = smIntersectTriSet(smMesh,t_set,orig,dir,p);
1241 #ifdef TEST_DRIVER
1242 VCOPY(Pick_point[0],p);
1243 #endif
1244 return(s_id);
1245 }
1246 else
1247 {
1248 /* Starting with orig, Walk along projection of ray onto sphere */
1249 point_on_sphere(r,orig,SM_VIEW_CENTER(smMesh));
1250 qt = smPointLocateCell(smMesh,r,NULL,NULL,FALSE);
1251 qtgetset(t_set,qt);
1252 /* os will contain all triangles seen thus far */
1253 setcopy(os,t_set);
1254
1255 /* Calculate ray perpendicular to dir: when projection ray is // to dir,
1256 the dot product will become negative.
1257 */
1258 VSUM(a,b,dir,d);
1259 d = DOT(a,b);
1260 while(d > 0)
1261 {
1262 s_id = smIntersectTriSet(smMesh,t_set,orig,dir,p);
1263 #ifdef TEST_DRIVER
1264 VCOPY(Pick_point[0],p);
1265 #endif
1266 if(s_id != EMPTY)
1267 return(s_id);
1268 /* Find next cell that projection of ray intersects */
1269 smTraceRay(smMesh,r,dir,v0,v1,v2,r);
1270 qt = smPointLocateCell(smMesh,r,NULL,NULL,FALSE);
1271 qtgetset(t_set,qt);
1272 /* Check triangles in set against those seen so far(os):only
1273 check new triangles for intersection (t_set')
1274 */
1275 check_set(t_set,os);
1276 d = DOT(a,r);
1277 }
1278 }
1279 #ifdef DEBUG
1280 eputs("smFindSamp():Pick Ray did not intersect mesh");
1281 #endif
1282 return(EMPTY);
1283 }
1284
1285
1286 smRebuild_mesh(sm,vptr)
1287 SM *sm;
1288 VIEW *vptr;
1289 {
1290 int i;
1291 FVECT dir;
1292 COLR c;
1293 FVECT p,ov;
1294
1295 /* Clear the mesh- and rebuild using the current sample array */
1296 #ifdef TEST_DRIVER
1297 View = *vptr;
1298 #endif
1299
1300 VSUB(ov,vptr->vp,SM_VIEW_CENTER(sm));
1301 smInit_mesh(sm,vptr->vp);
1302
1303 SM_FOR_ALL_SAMPLES(sm,i)
1304 {
1305 if(SM_NTH_W_DIR(sm,i)==-1)
1306 VADD(SM_NTH_WV(sm,i),SM_NTH_WV(sm,i),ov);
1307 smAdd_sample_to_mesh(sm,NULL,NULL,NULL,i);
1308 }
1309 }
1310
1311