ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/sm.c
Revision: 3.4
Committed: Tue Aug 25 11:03:27 1998 UTC (25 years, 8 months ago) by gwlarson
Content type: text/plain
Branch: MAIN
Changes since 3.3: +15 -87 lines
Log Message:
fixed problem with picking (ray tracking) of tetrahedron

File Contents

# User Rev Content
1 gwlarson 3.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 gwlarson 3.3 smTris_swap_edge(sm,t_id,t1_id,e,e1,tn_id,tn1_id,add,del)
444 gwlarson 3.1 SM *sm;
445     int t_id,t1_id;
446     int e,e1;
447     int **tn_id,**tn1_id;
448 gwlarson 3.3 LIST **add,**del;
449 gwlarson 3.1
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 gwlarson 3.3 *add = push_data(*add,ta_id);
472 gwlarson 3.1
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 gwlarson 3.3 *add = push_data(*add,tb_id);
478    
479 gwlarson 3.1 /* 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 gwlarson 3.3 *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 gwlarson 3.1 *tn_id = ta_id;
510     *tn1_id = tb_id;
511     }
512    
513 gwlarson 3.3 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 gwlarson 3.4 smNew_tri_cnt--;
526 gwlarson 3.3 continue;
527     }
528     t = SM_NTH_TRI(sm,t_id);
529     smLocator_add_tri(sm,t_id,T_NTH_V(t,0),T_NTH_V(t,1),T_NTH_V(t,2));
530     }
531    
532     while(del_list)
533     {
534     t_id = pop_list(&del_list);
535     if(SM_IS_NTH_T_NEW(sm,t_id))
536     {
537     smDelete_tri(sm,t_id);
538     continue;
539     }
540     smLocator_remove_tri(sm,t_id);
541     smDelete_tri(sm,t_id);
542     }
543     }
544 gwlarson 3.1 /* MUST add check for constrained edges */
545 gwlarson 3.3 int
546 gwlarson 3.1 smFix_tris(sm,id,tlist)
547     SM *sm;
548     int id;
549     LIST *tlist;
550     {
551     TRI *t,*t_opp;
552 gwlarson 3.3 FVECT p,p1,p2,p3;
553     int e,e1,swapped = 0;
554 gwlarson 3.1 int t_id,t_opp_id;
555 gwlarson 3.3 LIST *add_list,*del_list;
556    
557    
558     add_list = del_list = NULL;
559     VSUB(p,SM_NTH_WV(sm,id),SM_VIEW_CENTER(sm));
560 gwlarson 3.1 while(tlist)
561     {
562 gwlarson 3.3 t_id = pop_list(&tlist);
563 gwlarson 3.1 t = SM_NTH_TRI(sm,t_id);
564     e = (T_WHICH_V(t,id)+1)%3;
565     t_opp_id = T_NTH_NBR(t,e);
566     t_opp = SM_NTH_TRI(sm,t_opp_id);
567    
568     smDir(sm,p1,T_NTH_V(t_opp,0));
569     smDir(sm,p2,T_NTH_V(t_opp,1));
570     smDir(sm,p3,T_NTH_V(t_opp,2));
571     if(point_in_cone(p,p1,p2,p3))
572     {
573     swapped = 1;
574     e1 = T_NTH_NBR_PTR(t_id,t_opp);
575     /* check list for t_opp and Remove if there */
576     remove_from_list(t_opp_id,&tlist);
577 gwlarson 3.3 smTris_swap_edge(sm,t_id,t_opp_id,e,e1,&t_id,&t_opp_id,
578     &add_list,&del_list);
579 gwlarson 3.1 tlist = push_data(tlist,t_id);
580     tlist = push_data(tlist,t_opp_id);
581     }
582     }
583 gwlarson 3.3 smUpdate_locator(sm,add_list,del_list);
584    
585 gwlarson 3.1 return(swapped);
586     }
587    
588     /* Give the vertex "id" and a triangle "t" that it belongs to- return the
589     next nbr in a counter clockwise order about vertex "id"
590     */
591     int
592     smTri_next_ccw_nbr(sm,t,id)
593     SM *sm;
594     TRI *t;
595     int id;
596     {
597     int t_id;
598     int tri;
599    
600     /* Want the edge for which "id" is the destination */
601     t_id = (T_WHICH_V(t,id)+ 2)% 3;
602     tri = T_NTH_NBR(t,t_id);
603     return(tri);
604     }
605    
606     void
607     smReplace_point(sm,tri,id,nid)
608     SM *sm;
609     TRI *tri;
610     int id,nid;
611     {
612     TRI *t;
613     int t_id;
614    
615     T_NTH_V(tri,T_WHICH_V(tri,id)) = nid;
616    
617     t_id = smTri_next_ccw_nbr(sm,t,nid);
618     while((t= SM_NTH_TRI(sm,t_id)) != tri)
619     {
620     T_NTH_V(t,T_WHICH_V(t,id)) = nid;
621     t_id = smTri_next_ccw_nbr(sm,t,nid);
622     }
623     }
624    
625    
626     smClear_tri_flags(sm,id)
627     SM *sm;
628     int id;
629     {
630     int i;
631    
632     for(i=0; i < T_FLAGS; i++)
633     SM_CLEAR_NTH_T_FLAG(sm,id,i);
634    
635     }
636    
637     /* Locate the point-id in the point location structure: */
638     int
639     smReplace_vertex(sm,c,dir,p,tri_id,snew_id,type,which)
640     SM *sm;
641     COLR c;
642     FVECT dir,p;
643     int tri_id,snew_id;
644     char type,which;
645     {
646     int n_id,s_id;
647     TRI *tri;
648    
649     tri = SM_NTH_TRI(sm,tri_id);
650     /* Get the sample that corresponds to the "which" vertex of "tri" */
651     /* NEED: have re-init that sets clock flag */
652     /* If this is a base-sample: create a new sample and replace
653     all references to the base sample with references to the
654     new sample
655     */
656     s_id = T_NTH_V(tri,which);
657     if(SM_BASE_ID(sm,s_id))
658     {
659     if(snew_id != -1)
660     n_id = smAdd_sample_point(sm,c,dir,p);
661     else
662     n_id = snew_id;
663     smReplace_point(sm,tri,s_id,n_id);
664     s_id = n_id;
665     }
666     else /* If the sample exists, reset the values */
667     /* NOTE: This test was based on the SPHERICAL coordinates
668     of the point: If we are doing a multiple-sample-per
669     SPHERICAL pixel: we will want to test for equality-
670     and do other processing: for now: SINGLE SAMPLE PER
671     PIXEL
672     */
673     /* NOTE: snew_id needs to be marked as invalid?*/
674     if(snew_id == -1)
675     smInit_sample(sm,s_id,c,dir,p);
676     else
677     smReset_sample(sm,s_id,snew_id);
678     return(s_id);
679     }
680    
681    
682     /* Locate the point-id in the point location structure: */
683     int
684     smInsert_point_in_tri(sm,c,dir,p,s_id,tri_id)
685     SM *sm;
686     COLR c;
687     FVECT dir,p;
688     int s_id,tri_id;
689     {
690     TRI *tri,*t0,*t1,*t2,*nbr;
691     int v0_id,v1_id,v2_id,n_id;
692     int t0_id,t1_id,t2_id;
693     LIST *tlist;
694     FVECT npt;
695    
696     if(s_id == SM_INVALID)
697     s_id = smAdd_sample_point(sm,c,dir,p);
698    
699     tri = SM_NTH_TRI(sm,tri_id);
700     v0_id = T_NTH_V(tri,0);
701     v1_id = T_NTH_V(tri,1);
702     v2_id = T_NTH_V(tri,2);
703    
704     n_id = -1;
705     if(SM_BASE_ID(sm,v0_id)||SM_BASE_ID(sm,v1_id)||SM_BASE_ID(sm,v2_id))
706     {
707     smDir(sm,npt,s_id);
708     /* Change to an add and delete */
709     t0_id = (SM_BASE_ID(sm,v0_id))?v0_id:-1;
710     t1_id = (SM_BASE_ID(sm,v1_id))?v1_id:-1;
711     t2_id = (SM_BASE_ID(sm,v2_id))?v2_id:-1;
712     n_id = smClosest_vertex_in_tri(sm,t0_id,t1_id,t2_id,npt,P_REPLACE_EPS);
713     }
714     t0_id = smAdd_tri(sm,s_id,v0_id,v1_id,&t0);
715 gwlarson 3.3 /* Add triangle to the locator */
716     smLocator_add_tri(sm,t0_id,s_id,v0_id,v1_id);
717    
718 gwlarson 3.1 t1_id = smAdd_tri(sm,s_id,v1_id,v2_id,&t1);
719 gwlarson 3.3 smLocator_add_tri(sm,t1_id,s_id,v1_id,v2_id);
720 gwlarson 3.1 t2_id = smAdd_tri(sm,s_id,v2_id,v0_id,&t2);
721 gwlarson 3.3 smLocator_add_tri(sm,t2_id,s_id,v2_id,v0_id);
722    
723 gwlarson 3.1 /* Set the neighbor pointers for the new tris */
724     T_NTH_NBR(t0,0) = t2_id;
725     T_NTH_NBR(t0,1) = T_NTH_NBR(tri,0);
726     T_NTH_NBR(t0,2) = t1_id;
727     T_NTH_NBR(t1,0) = t0_id;
728     T_NTH_NBR(t1,1) = T_NTH_NBR(tri,1);
729     T_NTH_NBR(t1,2) = t2_id;
730     T_NTH_NBR(t2,0) = t1_id;
731     T_NTH_NBR(t2,1) = T_NTH_NBR(tri,2);
732     T_NTH_NBR(t2,2) = t0_id;
733    
734     /* Reset the neigbor pointers for the neighbors of the original */
735     nbr = SM_NTH_TRI(sm,T_NTH_NBR(tri,0));
736     T_NTH_NBR(nbr,T_NTH_NBR_PTR(tri_id,nbr)) = t0_id;
737     nbr = SM_NTH_TRI(sm,T_NTH_NBR(tri,1));
738     T_NTH_NBR(nbr,T_NTH_NBR_PTR(tri_id,nbr)) = t1_id;
739     nbr = SM_NTH_TRI(sm,T_NTH_NBR(tri,2));
740     T_NTH_NBR(nbr,T_NTH_NBR_PTR(tri_id,nbr)) = t2_id;
741    
742 gwlarson 3.3 smLocator_remove_tri(sm,tri_id);
743 gwlarson 3.1 smDelete_tri(sm,tri_id);
744    
745     /* Fix up the new triangles*/
746     tlist = push_data(NULL,t0_id);
747     tlist = push_data(tlist,t1_id);
748     tlist = push_data(tlist,t2_id);
749    
750     smFix_tris(sm,s_id,tlist);
751    
752     if(n_id != -1)
753     smDelete_point(sm,n_id);
754    
755     return(s_id);
756     }
757    
758    
759     int
760     smPointLocate(sm,pt,type,which,norm)
761     SM *sm;
762     FVECT pt;
763     char *type,*which;
764     char norm;
765     {
766     STREE *st;
767     int tri;
768     FVECT npt;
769    
770     st = SM_LOCATOR(sm);
771     if(norm)
772     {
773     point_on_sphere(npt,pt,SM_VIEW_CENTER(sm));
774     tri = stPoint_locate(st,npt,type,which);
775     }
776     else
777     tri = stPoint_locate(st,pt,type,which);
778     return(tri);
779     }
780    
781     QUADTREE
782 gwlarson 3.4 smPointLocateCell(sm,pt,norm,v0,v1,v2)
783 gwlarson 3.1 SM *sm;
784 gwlarson 3.2 FVECT pt;
785 gwlarson 3.1 char norm;
786 gwlarson 3.4 FVECT v0,v1,v2;
787 gwlarson 3.1 {
788     STREE *st;
789 gwlarson 3.4 QUADTREE *qtptr;
790 gwlarson 3.1 FVECT npt;
791    
792     st = SM_LOCATOR(sm);
793     if(norm)
794     {
795     point_on_sphere(npt,pt,SM_VIEW_CENTER(sm));
796    
797 gwlarson 3.4 qtptr = stPoint_locate_cell(st,npt,v0,v1,v2);
798 gwlarson 3.1 }
799     else
800 gwlarson 3.4 qtptr = stPoint_locate_cell(st,pt,v0,v1,v2);
801 gwlarson 3.1
802 gwlarson 3.4 if(qtptr)
803     return(*qtptr);
804     else
805     return(EMPTY);
806 gwlarson 3.1 }
807    
808     int
809     smAdd_sample_to_mesh(sm,c,dir,pt,s_id)
810     SM *sm;
811     COLR c;
812     FVECT dir,pt;
813     int s_id;
814     {
815     int t_id;
816     char type,which;
817     double d;
818     FVECT p;
819    
820     /* If new, foreground pt */
821     if(pt)
822     {
823     /* NOTE: This should be elsewhere! */
824     d = DIST(pt,SM_VIEW_CENTER(smMesh));
825     smDist_sum += 1.0/d;
826     /************************************/
827     t_id = smPointLocate(smMesh,pt,&type,&which,TRUE);
828     if(type==GT_FACE)
829     s_id = smInsert_point_in_tri(smMesh,c,dir,pt,s_id,t_id);
830     else
831     if(type==GT_VERTEX)
832     s_id = smReplace_vertex(smMesh,c,dir,pt,t_id,s_id,type,which);
833     #ifdef DEBUG
834     else
835     eputs("smAdd_sample_to_mesh(): unrecognized type\n");
836     #endif
837     }
838     else if(s_id != -1)
839     {
840     VCOPY(p,SM_NTH_WV(sm,s_id));
841     if(SM_NTH_W_DIR(sm,s_id) != -1)
842     {
843     /* NOTE: This should be elsewhere! */
844     d = DIST(p,SM_VIEW_CENTER(smMesh));
845     smDist_sum += 1.0/d;
846     /************************************/
847     }
848     t_id = smPointLocate(smMesh,p,&type,&which,TRUE);
849     if(type==GT_FACE)
850     s_id = smInsert_point_in_tri(smMesh,c,dir,p,s_id,t_id);
851     else
852     if(type==GT_VERTEX)
853     s_id = smReplace_vertex(smMesh,c,dir,p,t_id,s_id,type,which);
854     #ifdef DEBUG
855     else
856     eputs("smAdd_sample_to_mesh(): unrecognized type\n");
857     #endif
858     }
859     /* Is a BG(sky point) */
860     else
861     {
862     t_id = smPointLocate(smMesh,dir,&type,&which,FALSE);
863     if(type==GT_FACE)
864     s_id = smInsert_point_in_tri(smMesh,c,dir,NULL,s_id,t_id);
865     else
866     if(type==GT_VERTEX)
867     s_id = smReplace_vertex(smMesh,c,dir,NULL,t_id,s_id,type,which);
868     #ifdef DEBUG
869     else
870     eputs("smAdd_sample_to_mesh(): unrecognized type\n");
871     #endif
872     }
873     return(s_id);
874     }
875    
876     /*
877     * int
878     * smNewSamp(c, dir, p) : register new sample point and return index
879     * COLR c; : pixel color (RGBE)
880     * FVECT dir; : ray direction vector
881     * FVECT p; : world intersection point
882     *
883     * Add new sample point to data structures, removing old values as necessary.
884     * New sample representation will be output in next call to smUpdate().
885     * If the point is a sky point: then v=NULL
886     */
887     int
888     smNewSamp(c,dir,p)
889     COLR c;
890     FVECT dir;
891     FVECT p;
892    
893     {
894     int s_id;
895    
896     /* First check if this the first sample: if so initialize mesh */
897     if(SM_NUM_SAMP(smMesh) == 0)
898     #ifdef TEST_DRIVER
899     smInit_mesh(smMesh,View.vp);
900     #else
901     smInit_mesh(smMesh,odev.v.vp);
902     #endif
903     s_id = smAdd_sample_to_mesh(smMesh,c,dir,p,-1);
904     #if 0
905     {
906     char buff[100];
907     sprintf(buff,"Added sample %d\n",s_id);
908     eputs(buff);
909     }
910     #endif
911     return(s_id);
912    
913     }
914     /*
915     * int
916     * smFindsamp(orig, dir): intersect ray with 3D rep. and find closest sample
917     * FVECT orig, dir;
918     *
919     * Find the closest sample to the given ray. Return -1 on failure.
920     */
921    
922     /*
923     * smClean() : display has been wiped clean
924     *
925     * Called after display has been effectively cleared, meaning that all
926     * geometry must be resent down the pipeline in the next call to smUpdate().
927     */
928    
929    
930     /*
931     * smUpdate(vp, qua) : update OpenGL output geometry for view vp
932     * VIEW *vp; : desired view
933     * int qua; : quality level (percentage on linear time scale)
934     *
935     * Draw new geometric representation using OpenGL calls. Assume that the
936     * view has already been set up and the correct frame buffer has been
937     * selected for drawing. The quality level is on a linear scale, where 100%
938     * is full (final) quality. It is not necessary to redraw geometry that has
939     * been output since the last call to smClean().
940     */
941    
942    
943     int
944     smClear_vert(sm,id)
945     SM *sm;
946     int id;
947     {
948     if(SM_INVALID_POINT_ID(sm,id))
949     return(FALSE);
950    
951     SM_NTH_VERT(sm,id) = SM_INVALID;
952    
953     return(TRUE);
954     }
955    
956     int
957     smAdd_base_vertex(sm,v,d)
958     SM *sm;
959     FVECT v,d;
960     {
961     int id;
962    
963     /* First add coordinate to the sample array */
964     id = smAdd_aux_point(sm,v,d);
965     if(id == -1)
966     return(SM_INVALID);
967     /* Initialize triangle pointer to -1 */
968     smClear_vert(sm,id);
969     return(id);
970     }
971    
972    
973    
974     /* Initialize a the point location DAG based on a 6 triangle tesselation
975     of the unit sphere centered on the view center. The DAG structure
976     contains 6 roots: one for each initial base triangle
977     */
978    
979     int
980     smCreate_base_mesh(sm,type)
981     SM *sm;
982     int type;
983     {
984     int i,id;
985     int p[4],ids[4];
986     int v0_id,v1_id,v2_id;
987     TRI *tris[4];
988     FVECT d,pt,cntr;
989    
990     /* First insert the base vertices into the sample point array */
991    
992     for(i=0; i < 4; i++)
993     {
994     VADD(cntr,stDefault_base[i],SM_VIEW_CENTER(sm));
995     point_on_sphere(d,cntr,SM_VIEW_CENTER(sm));
996     id = smAdd_base_vertex(sm,cntr,d);
997     /* test to make sure vertex allocated */
998     if(id != -1)
999     p[i] = id;
1000     else
1001     return(0);
1002     }
1003     /* Create the base triangles */
1004     for(i=0;i < 4; i++)
1005     {
1006     v0_id = p[stTri_verts[i][0]];
1007     v1_id = p[stTri_verts[i][1]];
1008     v2_id = p[stTri_verts[i][2]];
1009     if((ids[i] = smAdd_tri(sm, v0_id,v1_id,v2_id,&(tris[i])))== -1)
1010     return(0);
1011 gwlarson 3.3 smLocator_add_tri(sm,ids[i],v0_id,v1_id,v2_id);
1012 gwlarson 3.1 }
1013     /* Set neighbors */
1014    
1015     T_NTH_NBR(tris[0],0) = ids[3];
1016     T_NTH_NBR(tris[0],1) = ids[2];
1017     T_NTH_NBR(tris[0],2) = ids[1];
1018    
1019     T_NTH_NBR(tris[1],0) = ids[3];
1020     T_NTH_NBR(tris[1],1) = ids[0];
1021     T_NTH_NBR(tris[1],2) = ids[2];
1022    
1023     T_NTH_NBR(tris[2],0) = ids[3];
1024     T_NTH_NBR(tris[2],1) = ids[1];
1025     T_NTH_NBR(tris[2],2) = ids[0];
1026    
1027     T_NTH_NBR(tris[3],0) = ids[1];
1028     T_NTH_NBR(tris[3],1) = ids[2];
1029     T_NTH_NBR(tris[3],2) = ids[0];
1030     return(1);
1031    
1032     }
1033    
1034    
1035     int
1036     smNext_tri_flag_set(sm,i,which,b)
1037     SM *sm;
1038     int i,which;
1039     char b;
1040     {
1041    
1042     for(; i < SM_TRI_CNT(sm);i++)
1043     {
1044     if(!SM_IS_NTH_T_FLAG(sm,i,which))
1045     continue;
1046    
1047     if(!b)
1048     break;
1049     if((b==1) && !SM_BG_TRI(sm,i))
1050     break;
1051     if((b==2) && SM_BG_TRI(sm,i))
1052     break;
1053     }
1054    
1055     return(i);
1056     }
1057    
1058    
1059     int
1060     smNext_valid_tri(sm,i)
1061     SM *sm;
1062     int i;
1063     {
1064    
1065     while( i < SM_TRI_CNT(sm) && !T_IS_VALID(SM_NTH_TRI(sm,i)))
1066     i++;
1067    
1068     return(i);
1069     }
1070    
1071    
1072    
1073     qtTri_verts_from_id(t_id,v0,v1,v2)
1074     int t_id;
1075     FVECT v0,v1,v2;
1076     {
1077     TRI *t;
1078     int v0_id,v1_id,v2_id;
1079    
1080     t = SM_NTH_TRI(smMesh,t_id);
1081     v0_id = T_NTH_V(t,0);
1082     v1_id = T_NTH_V(t,1);
1083     v2_id = T_NTH_V(t,2);
1084    
1085     smDir(smMesh,v0,v0_id);
1086     smDir(smMesh,v1,v1_id);
1087     smDir(smMesh,v2,v2_id);
1088    
1089     }
1090    
1091    
1092     int
1093     smIntersectTriSet(sm,t_set,orig,dir,pt)
1094     SM *sm;
1095     OBJECT *t_set;
1096     FVECT orig,dir,pt;
1097     {
1098     OBJECT *optr;
1099     int i,t_id,v_id;
1100     TRI *tri;
1101     FVECT p0,p1,p2;
1102     char type,which;
1103     int p0_id,p1_id,p2_id;
1104    
1105     for(optr = QT_SET_PTR(t_set),i = QT_SET_CNT(t_set); i > 0; i--)
1106     {
1107     t_id = QT_SET_NEXT_ELEM(optr);
1108     tri = SM_NTH_TRI(sm,t_id);
1109     p0_id = T_NTH_V(tri,0);
1110     p1_id = T_NTH_V(tri,1);
1111     p2_id = T_NTH_V(tri,2);
1112     VCOPY(p0,SM_NTH_WV(sm,p0_id));
1113     VCOPY(p1,SM_NTH_WV(sm,p1_id));
1114     VCOPY(p2,SM_NTH_WV(sm,p2_id));
1115     if(type = ray_intersect_tri(orig,dir,p0,p1,p2,pt,&which))
1116     {
1117     if(type==GT_VERTEX)
1118     return(T_NTH_V(tri,which));
1119     v_id = smClosest_vertex_in_w_tri(sm,p0_id,p1_id,p2_id,pt);
1120     return(v_id);
1121     }
1122     }
1123     return(-1);
1124     }
1125    
1126    
1127     /*
1128     * int
1129     * smFindSamp(FVECT orig, FVECT dir)
1130     *
1131     * Find the closest sample to the given ray. Returns sample id, -1 on failure.
1132     * "dir" is assumed to be normalized
1133     */
1134     int
1135     smFindSamp(orig,dir)
1136     FVECT orig,dir;
1137     {
1138     FVECT r,v0,v1,v2,a,b,p;
1139     OBJECT os[MAXCSET+1],t_set[MAXSET+1];
1140     QUADTREE qt;
1141     int s_id;
1142     double d;
1143    
1144     /* r is the normalized vector from the view center to the current
1145     * ray point ( starting with "orig"). Find the cell that r falls in,
1146     * and test the ray against all triangles stored in the cell. If
1147     * the test fails, trace the projection of the ray across to the
1148     * next cell it intersects: iterate until either an intersection
1149     * is found, or the projection ray is // to the direction. The sample
1150     * corresponding to the triangle vertex closest to the intersection
1151     * point is returned.
1152     */
1153    
1154     /* First test if "orig" coincides with the View_center or if "dir" is
1155     parallel to r formed by projecting "orig" on the sphere. In
1156     either case, do a single test against the cell containing the
1157     intersection of "dir" and the sphere
1158     */
1159     point_on_sphere(b,orig,SM_VIEW_CENTER(smMesh));
1160     d = -DOT(b,dir);
1161     if(EQUAL_VEC3(orig,SM_VIEW_CENTER(smMesh)) || EQUAL(fabs(d),1.0))
1162     {
1163 gwlarson 3.4 qt = smPointLocateCell(smMesh,dir,FALSE,NULL,NULL,NULL);
1164 gwlarson 3.1 /* Test triangles in the set for intersection with Ray:returns
1165     first found
1166     */
1167     qtgetset(t_set,qt);
1168     s_id = smIntersectTriSet(smMesh,t_set,orig,dir,p);
1169     #ifdef TEST_DRIVER
1170     VCOPY(Pick_point[0],p);
1171     #endif
1172     return(s_id);
1173     }
1174     else
1175     {
1176     /* Starting with orig, Walk along projection of ray onto sphere */
1177     point_on_sphere(r,orig,SM_VIEW_CENTER(smMesh));
1178 gwlarson 3.4 qt = smPointLocateCell(smMesh,r,FALSE,v0,v1,v2);
1179 gwlarson 3.1 qtgetset(t_set,qt);
1180     /* os will contain all triangles seen thus far */
1181     setcopy(os,t_set);
1182    
1183     /* Calculate ray perpendicular to dir: when projection ray is // to dir,
1184     the dot product will become negative.
1185     */
1186     VSUM(a,b,dir,d);
1187     d = DOT(a,b);
1188     while(d > 0)
1189     {
1190     s_id = smIntersectTriSet(smMesh,t_set,orig,dir,p);
1191     #ifdef TEST_DRIVER
1192     VCOPY(Pick_point[0],p);
1193     #endif
1194     if(s_id != EMPTY)
1195     return(s_id);
1196     /* Find next cell that projection of ray intersects */
1197 gwlarson 3.4 traceRay(r,dir,v0,v1,v2,r);
1198     qt = smPointLocateCell(smMesh,r,FALSE,v0,v1,v2);
1199 gwlarson 3.1 qtgetset(t_set,qt);
1200     /* Check triangles in set against those seen so far(os):only
1201     check new triangles for intersection (t_set')
1202     */
1203     check_set(t_set,os);
1204 gwlarson 3.4 d = DOT(a,r);
1205 gwlarson 3.1 }
1206     }
1207     #ifdef DEBUG
1208     eputs("smFindSamp():Pick Ray did not intersect mesh");
1209     #endif
1210     return(EMPTY);
1211     }
1212    
1213    
1214     smRebuild_mesh(sm,vptr)
1215     SM *sm;
1216     VIEW *vptr;
1217     {
1218     int i;
1219     FVECT dir;
1220     COLR c;
1221     FVECT p,ov;
1222    
1223     /* Clear the mesh- and rebuild using the current sample array */
1224     #ifdef TEST_DRIVER
1225     View = *vptr;
1226     #endif
1227    
1228     VSUB(ov,vptr->vp,SM_VIEW_CENTER(sm));
1229     smInit_mesh(sm,vptr->vp);
1230    
1231     SM_FOR_ALL_SAMPLES(sm,i)
1232     {
1233     if(SM_NTH_W_DIR(sm,i)==-1)
1234     VADD(SM_NTH_WV(sm,i),SM_NTH_WV(sm,i),ov);
1235     smAdd_sample_to_mesh(sm,NULL,NULL,NULL,i);
1236     }
1237     }
1238    
1239