ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/sm.c
Revision: 3.10
Committed: Mon Dec 28 19:30:27 1998 UTC (25 years, 3 months ago) by gwlarson
Content type: text/plain
Branch: MAIN
Changes since 3.9: +1148 -425 lines
Log Message:
New insertion routine
New Culling routine based on insertion algorithm
Adapted old insertion code: now used by picking
Point location code returns on-vertex,on-edge, or in-triangle
Added on_edge case for subdivision
Implemented unordered sets
Removed deletion from quadtree- added set compression to replace functionality

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 gwlarson 3.8 #include "sm_flag.h"
12 gwlarson 3.1 #include "sm_list.h"
13     #include "sm_geom.h"
14 gwlarson 3.8 #include "sm_qtree.h"
15     #include "sm_stree.h"
16 gwlarson 3.1 #include "sm.h"
17    
18    
19     SM *smMesh = NULL;
20     double smDist_sum=0;
21     int smNew_tri_cnt=0;
22 gwlarson 3.10 double smMinSampDiff = 4.25e-4; /* min edge length in radians */
23 gwlarson 3.1
24 gwlarson 3.10 /* Each edge extends .5536 radians - 31.72 degrees */
25     static FVECT icosa_verts[162] =
26     {{-0.018096,0.495400,0.868477},{0.018614,-0.554780,0.831789},
27     {0.850436,0.011329,0.525956},{0.850116,0.048016,-0.524402},
28     {-0.850116,-0.048016,0.524402},{-0.850436,-0.011329,-0.525956},
29     {0.495955,0.867825,0.030147},{-0.555329,0.831118,0.029186},
30     {0.555329,-0.831118,-0.029186},{-0.495955,-0.867825,-0.030147},
31     {-0.018614,0.554780,-0.831789},{0.018096,-0.495400,-0.868477},
32     {0.000305,-0.034906,0.999391},{0.489330,0.297837,0.819664},
33     {0.510900,-0.319418,0.798094},{-0.488835,-0.354308,0.797186},
34     {-0.510409,0.262947,0.818744},{0.999391,0.034875,0.000914},
35     {0.791335,0.516672,0.326862},{0.791142,0.538234,-0.290513},
36     {0.826031,-0.460219,-0.325378},{0.826216,-0.481780,0.291985},
37     {-0.999391,-0.034875,-0.000914},{-0.826216,0.481780,-0.291985},
38     {-0.826031,0.460219,0.325378},{-0.791142,-0.538234,0.290513},
39     {-0.791335,-0.516672,-0.326862},{-0.034911,0.998782,0.034881},
40     {0.280899,0.801316,0.528194},{-0.337075,0.779739,0.527624},
41     {-0.337377,0.814637,-0.471746},{0.280592,0.836213,-0.471186},
42     {0.034911,-0.998782,-0.034881},{-0.280592,-0.836213,0.471186},
43     {0.337377,-0.814637,0.471746},{0.337075,-0.779739,-0.527624},
44     {-0.280899,-0.801316,-0.528194},{-0.000305,0.034906,-0.999391},
45     {0.510409,-0.262947,-0.818744},{0.488835,0.354308,-0.797186},
46     {-0.510900,0.319418,-0.798094},{-0.489330,-0.297837,-0.819664},
47     {0.009834,-0.306509,0.951817},{0.268764,-0.186284,0.945021},
48     {0.275239,-0.454405,0.847207},{-0.009247,0.239357,0.970888},
49     {0.244947,0.412323,0.877491},{0.257423,0.138235,0.956360},
50     {0.525850,-0.011345,0.850502},{0.696394,0.160701,0.699436},
51     {0.707604,-0.160141,0.688223},{-0.268186,0.119892,0.955878},
52     {-0.274715,0.394186,0.877011},{-0.244420,-0.472542,0.846737},
53     {-0.256843,-0.204627,0.944542},{-0.525332,-0.048031,0.849541},
54     {-0.695970,-0.209123,0.686945},{-0.707182,0.111718,0.698149},
55     {0.961307,0.043084,-0.272090},{0.941300,0.301289,-0.152245},
56     {0.853076,0.304715,-0.423568},{0.961473,0.024015,0.273848},
57     {0.853344,0.274439,0.443269},{0.941400,0.289953,0.172315},
58     {0.831918,0.554570,0.019109},{0.669103,0.719629,0.185565},
59     {0.669003,0.730836,-0.135332},{0.959736,-0.234941,0.153979},
60     {0.871472,-0.244526,0.425140},{0.871210,-0.214250,-0.441690},
61     {0.959638,-0.223606,-0.170573},{0.868594,-0.495213,-0.017555},
62     {0.717997,-0.671205,-0.184294},{0.718092,-0.682411,0.136596},
63     {-0.961307,-0.043084,0.272090},{-0.959638,0.223606,0.170573},
64     {-0.871210,0.214250,0.441690},{-0.961473,-0.024015,-0.273848},
65     {-0.871472,0.244526,-0.425140},{-0.959736,0.234941,-0.153979},
66     {-0.868594,0.495213,0.017555},{-0.718092,0.682411,-0.136596},
67     {-0.717997,0.671205,0.184294},{-0.941400,-0.289953,-0.172315},
68     {-0.853344,-0.274439,-0.443269},{-0.853076,-0.304715,0.423568},
69     {-0.941300,-0.301289,0.152245},{-0.831918,-0.554570,-0.019109},
70     {-0.669003,-0.730836,0.135332},{-0.669103,-0.719629,-0.185565},
71     {-0.306809,0.951188,0.033302},{-0.195568,0.935038,0.295731},
72     {-0.463858,0.837300,0.289421},{0.239653,0.970270,0.033802},
73     {0.403798,0.867595,0.290218},{0.129326,0.946383,0.296031},
74     {-0.029535,0.831255,0.555106},{0.136603,0.674019,0.725974},
75     {-0.184614,0.662804,0.725678},{0.129164,0.964728,-0.229382},
76     {0.403637,0.885733,-0.229245},{-0.464015,0.855438,-0.230036},
77     {-0.195726,0.953383,-0.229677},{-0.029855,0.867949,-0.495755},
78     {-0.185040,0.711807,-0.677562},{0.136173,0.723022,-0.677271},
79     {0.306809,-0.951188,-0.033302},{0.195726,-0.953383,0.229677},
80     {0.464015,-0.855438,0.230036},{-0.239653,-0.970270,-0.033802},
81     {-0.403637,-0.885733,0.229245},{-0.129164,-0.964728,0.229382},
82     {0.029855,-0.867949,0.495755},{-0.136173,-0.723022,0.677271},
83     {0.185040,-0.711807,0.677562},{-0.129326,-0.946383,-0.296031},
84     {-0.403798,-0.867595,-0.290218},{0.463858,-0.837300,-0.289421}
85     ,{0.195568,-0.935038,-0.295731},{0.029535,-0.831255,-0.555106},
86     {0.184614,-0.662804,-0.725678},{-0.136603,-0.674019,-0.725974},
87     {-0.009834,0.306509,-0.951817},{0.256843,0.204627,-0.944542},
88     {0.244420,0.472542,-0.846737},{0.009247,-0.239357,-0.970888},
89     {0.274715,-0.394186,-0.877011},{0.268186,-0.119892,-0.955878},
90     {0.525332,0.048031,-0.849541},{0.707182,-0.111718,-0.698149},
91     {0.695970,0.209123,-0.686945},{-0.257423,-0.138235,-0.956360},
92     {-0.244947,-0.412323,-0.877491},{-0.275239,0.454405,-0.847207},
93     {-0.268764,0.186284,-0.945021},{-0.525850,0.011345,-0.850502},
94     {-0.707604,0.160141,-0.688223},{-0.696394,-0.160701,-0.699436},
95     {0.563717,0.692920,0.449538},{0.673284,0.428212,0.602763},
96     {0.404929,0.577853,0.708603},{0.445959,-0.596198,0.667584},
97     {0.702960,-0.421213,0.573086},{0.611746,-0.681577,0.401523},
98     {-0.445543,0.548166,0.707818},{-0.702606,0.380190,0.601498},
99     {-0.611490,0.651895,0.448456},{-0.563454,-0.722602,0.400456},
100     {-0.672921,-0.469235,0.571835},{-0.404507,-0.625886,0.666814},
101     {0.404507,0.625886,-0.666814},{0.672921,0.469235,-0.571835},
102     {0.563454,0.722602,-0.400456},{0.611490,-0.651895,-0.448456},
103     {0.702606,-0.380190,-0.601498},{0.445543,-0.548166,-0.707818},
104     {-0.611746,0.681577,-0.401523},{-0.702960,0.421213,-0.573086},
105     {-0.445959,0.596198,-0.667584},{-0.404929,-0.577853,-0.708603},
106     {-0.673284,-0.428212,-0.602763},{-0.563717,-0.692920,-0.449538}};
107     static int icosa_tris[320][3] =
108     { {1,42,44},{42,12,43},{44,43,14},{42,43,44},{12,45,47},{45,0,46},{47,46,13},
109     {45,46,47},{14,48,50},{48,13,49},{50,49,2},{48,49,50},{12,47,43},{47,13,48},
110     {43,48,14},{47,48,43},{0,45,52},{45,12,51},{52,51,16},{45,51,52},{12,42,54},
111     {42,1,53},{54,53,15},{42,53,54},{16,55,57},{55,15,56},{57,56,4},{55,56,57},
112     {12,54,51},{54,15,55},{51,55,16},{54,55,51},{3,58,60},{58,17,59},{60,59,19},
113     {58,59,60},{17,61,63},{61,2,62},{63,62,18},{61,62,63},{19,64,66},{64,18,65},
114     {66,65,6},{64,65,66},{17,63,59},{63,18,64},{59,64,19},{63,64,59},{2,61,68},
115     {61,17,67},{68,67,21},{61,67,68},{17,58,70},{58,3,69},{70,69,20},{58,69,70},
116     {21,71,73},{71,20,72},{73,72,8},{71,72,73},{17,70,67},{70,20,71},{67,71,21},
117     {70,71,67},{4,74,76},{74,22,75},{76,75,24},{74,75,76},{22,77,79},{77,5,78},
118     {79,78,23},{77,78,79},{24,80,82},{80,23,81},{82,81,7},{80,81,82},{22,79,75},
119     {79,23,80},{75,80,24},{79,80,75},{5,77,84},{77,22,83},{84,83,26},{77,83,84},
120     {22,74,86},{74,4,85},{86,85,25},{74,85,86},{26,87,89},{87,25,88},{89,88,9},
121     {87,88,89},{22,86,83},{86,25,87},{83,87,26},{86,87,83},{7,90,92},{90,27,91},
122     {92,91,29},{90,91,92},{27,93,95},{93,6,94},{95,94,28},{93,94,95},{29,96,98},
123     {96,28,97},{98,97,0},{96,97,98},{27,95,91},{95,28,96},{91,96,29},{95,96,91},
124     {6,93,100},{93,27,99},{100,99,31},{93,99,100},{27,90,102},{90,7,101},
125     {102,101,30},{90,101,102},{31,103,105},{103,30,104},{105,104,10},{103,104,105},
126     {27,102,99},{102,30,103},{99,103,31},{102,103,99},{8,106,108},{106,32,107},
127     {108,107,34},{106,107,108},{32,109,111},{109,9,110},{111,110,33},{109,110,111},
128     {34,112,114},{112,33,113},{114,113,1},{112,113,114},{32,111,107},{111,33,112},
129     {107,112,34},{111,112,107},{9,109,116},{109,32,115},{116,115,36},{109,115,116},
130     {32,106,118},{106,8,117},{118,117,35},{106,117,118},{36,119,121},{119,35,120},
131     {121,120,11},{119,120,121},{32,118,115},{118,35,119},{115,119,36},{118,119,115},
132     {10,122,124},{122,37,123},{124,123,39},{122,123,124},{37,125,127},{125,11,126},
133     {127,126,38},{125,126,127},{39,128,130},{128,38,129},{130,129,3},{128,129,130},
134     {37,127,123},{127,38,128},{123,128,39},{127,128,123},{11,125,132},{125,37,131},
135     {132,131,41},{125,131,132},{37,122,134},{122,10,133},{134,133,40},{122,133,134},
136     {41,135,137},{135,40,136},{137,136,5},{135,136,137},{37,134,131},{134,40,135},
137     {131,135,41},{134,135,131},{6,65,94},{65,18,138},{94,138,28},{65,138,94},
138     {18,62,139},{62,2,49},{139,49,13},{62,49,139},{28,140,97},{140,13,46},{97,46,0},
139     {140,46,97},{18,139,138},{139,13,140},{138,140,28},{139,140,138},{1,44,114},
140     {44,14,141},{114,141,34},{44,141,114},{14,50,142},{50,2,68},{142,68,21},
141     {50,68,142},{34,143,108},{143,21,73},{108,73,8},{143,73,108},{14,142,141},
142     {142,21,143},{141,143,34},{142,143,141},{0,52,98},{52,16,144},{98,144,29},
143     {52,144,98},{16,57,145},{57,4,76},{145,76,24},{57,76,145},{29,146,92},
144     {146,24,82},{92,82,7},{146,82,92},{16,145,144},{145,24,146},{144,146,29},
145     {145,146,144},{9,88,110},{88,25,147},{110,147,33},{88,147,110},{25,85,148},
146     {85,4,56},{148,56,15},{85,56,148},{33,149,113},{149,15,53},{113,53,1},
147     {149,53,113},{25,148,147},{148,15,149},{147,149,33},{148,149,147},{10,124,105},
148     {124,39,150},{105,150,31},{124,150,105},{39,130,151},{130,3,60},{151,60,19},
149     {130,60,151},{31,152,100},{152,19,66},{100,66,6},{152,66,100},{39,151,150},
150     {151,19,152},{150,152,31},{151,152,150},{8,72,117},{72,20,153},{117,153,35},
151     {72,153,117},{20,69,154},{69,3,129},{154,129,38},{69,129,154},{35,155,120},
152     {155,38,126},{120,126,11},{155,126,120},{20,154,153},{154,38,155},{153,155,35},
153     {154,155,153},{7,81,101},{81,23,156},{101,156,30},{81,156,101},{23,78,157},
154     {78,5,136},{157,136,40},{78,136,157},{30,158,104},{158,40,133},{104,133,10},
155     {158,133,104},{23,157,156},{157,40,158},{156,158,30},{157,158,156},{11,132,121},
156     {132,41,159},{121,159,36},{132,159,121},{41,137,160},{137,5,84},{160,84,26},
157     {137,84,160},{36,161,116},{161,26,89},{116,89,9},{161,89,116},{41,160,159},
158     {160,26,161},{159,161,36},{160,161,159}};
159     static int icosa_nbrs[320][3] =
160     {{3,208,21},{12,3,20},{14,209,3},{2,0,1},{7,12,17},{202,7,16},{201,13,7},
161     {6,4,5},{11,212,14},{198,11,13},{197,213,11},{10,8,9},{15,1,4},{9,15,6},
162     {8,2,15},{14,12,13},{19,224,5},{28,19,4},{30,225,19},{18,16,17},{23,28,1},
163     {250,23,0},{249,29,23},{22,20,21},{27,228,30},{246,27,29},{245,229,27},
164     {26,24,25},{31,17,20},{25,31,22},{24,18,31},{30,28,29},{35,261,53},{44,35,52},
165     {46,262,35},{34,32,33},{39,44,49},{197,39,48},{196,45,39},{38,36,37},
166     {43,265,46},{193,43,45},{192,266,43},{42,40,41},{47,33,36},{41,47,38},
167     {40,34,47},{46,44,45},{51,213,37},{60,51,36},{62,214,51},{50,48,49},{55,60,33},
168     {277,55,32},{276,61,55},{54,52,53},{59,217,62},{273,59,61},{272,218,59},
169     {58,56,57},{63,49,52},{57,63,54},{56,50,63},{62,60,61},{67,229,85},{76,67,84},
170     {78,230,67},{66,64,65},{71,76,81},{293,71,80},{292,77,71},{70,68,69},
171     {75,233,78},{289,75,77},{288,234,75},{74,72,73},{79,65,68},{73,79,70},
172     {72,66,79},{78,76,77},{83,309,69},{92,83,68},{94,310,83},{82,80,81},{87,92,65},
173     {245,87,64},{244,93,87},{86,84,85},{91,313,94},{241,91,93},{240,314,91},
174     {90,88,89},{95,81,84},{89,95,86},{88,82,95},{94,92,93},{99,234,117},
175     {108,99,116},{110,232,99},{98,96,97},{103,108,113},{192,103,112},{194,109,103},
176     {102,100,101},{107,226,110},{200,107,109},{202,224,107},{106,104,105},
177     {111,97,100},{105,111,102},{104,98,111},{110,108,109},{115,266,101},
178     {124,115,100},{126,264,115},{114,112,113},{119,124,97},{288,119,96},
179     {290,125,119},{118,116,117},{123,258,126},{296,123,125},{298,256,123},
180     {122,120,121},{127,113,116},{121,127,118},{120,114,127},{126,124,125},
181     {131,218,149},{140,131,148},{142,216,131},{130,128,129},{135,140,145},
182     {240,135,144},{242,141,135},{134,132,133},{139,210,142},{248,139,141},
183     {250,208,139},{138,136,137},{143,129,132},{137,143,134},{136,130,143},
184     {142,140,141},{147,314,133},{156,147,132},{158,312,147},{146,144,145},
185     {151,156,129},{272,151,128},{274,157,151},{150,148,149},{155,306,158},
186     {280,155,157},{282,304,155},{154,152,153},{159,145,148},{153,159,150},
187     {152,146,159},{158,156,157},{163,256,181},{172,163,180},{174,257,163},
188     {162,160,161},{167,172,177},{282,167,176},{281,173,167},{166,164,165},
189     {171,260,174},{278,171,173},{277,261,171},{170,168,169},{175,161,164},
190     {169,175,166},{168,162,175},{174,172,173},{179,304,165},{188,179,164},
191     {190,305,179},{178,176,177},{183,188,161},{298,183,160},{297,189,183},
192     {182,180,181},{187,308,190},{294,187,189},{293,309,187},{186,184,185},
193     {191,177,180},{185,191,182},{184,178,191},{190,188,189},{195,101,42},
194     {204,195,41},{206,102,195},{194,192,193},{199,204,38},{10,199,37},{9,205,199},
195     {198,196,197},{203,105,206},{6,203,205},{5,106,203},{202,200,201},{207,193,196},
196     {201,207,198},{200,194,207},{206,204,205},{211,138,0},{220,211,2},{222,136,211},
197     {210,208,209},{215,220,8},{48,215,10},{50,221,215},{214,212,213},{219,130,222},
198     {56,219,221},{58,128,219},{218,216,217},{223,209,212},{217,223,214},
199     {216,210,223},{222,220,221},{227,106,16},{236,227,18},{238,104,227},
200     {226,224,225},{231,236,24},{64,231,26},{66,237,231},{230,228,229},{235,98,238},
201     {72,235,237},{74,96,235},{234,232,233},{239,225,228},{233,239,230},
202     {232,226,239},{238,236,237},{243,133,90},{252,243,89},{254,134,243},
203     {242,240,241},{247,252,86},{26,247,85},{25,253,247},{246,244,245},{251,137,254},
204     {22,251,253},{21,138,251},{250,248,249},{255,241,244},{249,255,246},
205     {248,242,255},{254,252,253},{259,122,160},{268,259,162},{270,120,259},
206     {258,256,257},{263,268,168},{32,263,170},{34,269,263},{262,260,261},
207     {267,114,270},{40,267,269},{42,112,267},{266,264,265},{271,257,260},
208     {265,271,262},{264,258,271},{270,268,269},{275,149,58},{284,275,57},
209     {286,150,275},{274,272,273},{279,284,54},{170,279,53},{169,285,279},
210     {278,276,277},{283,153,286},{166,283,285},{165,154,283},{282,280,281},
211     {287,273,276},{281,287,278},{280,274,287},{286,284,285},{291,117,74},
212     {300,291,73},{302,118,291},{290,288,289},{295,300,70},{186,295,69},
213     {185,301,295},{294,292,293},{299,121,302},{182,299,301},{181,122,299},
214     {298,296,297},{303,289,292},{297,303,294},{296,290,303},{302,300,301},
215     {307,154,176},{316,307,178},{318,152,307},{306,304,305},{311,316,184},
216     {80,311,186},{82,317,311},{310,308,309},{315,146,318},{88,315,317},
217     {90,144,315},{314,312,313},{319,305,308},{313,319,310},{312,306,319},
218     {318,316,317}};
219 gwlarson 3.8
220 gwlarson 3.5
221 gwlarson 3.1 #ifdef TEST_DRIVER
222     VIEW Current_View = {0,{0,0,0},{0,0,-1},{0,1,0},60,60,0};
223     int Pick_cnt =0;
224 gwlarson 3.5 int Pick_tri = -1,Picking = FALSE,Pick_samp=-1;
225 gwlarson 3.1 FVECT Pick_point[500],Pick_origin,Pick_dir;
226     FVECT Pick_v0[500], Pick_v1[500], Pick_v2[500];
227 gwlarson 3.9 int Pick_q[500];
228 gwlarson 3.1 FVECT P0,P1,P2;
229     FVECT FrustumNear[4],FrustumFar[4];
230     double dev_zmin=.01,dev_zmax=1000;
231     #endif
232    
233     smDir(sm,ps,id)
234     SM *sm;
235     FVECT ps;
236     int id;
237     {
238 gwlarson 3.8 VSUB(ps,SM_NTH_WV(sm,id),SM_VIEW_CENTER(sm));
239     normalize(ps);
240 gwlarson 3.1 }
241 gwlarson 3.8
242 gwlarson 3.7 smDir_in_cone(sm,ps,id)
243     SM *sm;
244     FVECT ps;
245     int id;
246     {
247    
248 gwlarson 3.8 VSUB(ps,SM_NTH_WV(sm,id),SM_VIEW_CENTER(sm));
249     normalize(ps);
250 gwlarson 3.7 }
251 gwlarson 3.1
252     smClear_flags(sm,which)
253     SM *sm;
254     int which;
255     {
256     int i;
257    
258     if(which== -1)
259     for(i=0; i < T_FLAGS;i++)
260 gwlarson 3.8 bzero(SM_NTH_FLAGS(sm,i),FLAG_BYTES(SM_MAX_TRIS(sm)));
261 gwlarson 3.1 else
262 gwlarson 3.8 bzero(SM_NTH_FLAGS(sm,which),FLAG_BYTES(SM_MAX_TRIS(sm)));
263 gwlarson 3.1 }
264    
265 gwlarson 3.8 /* Given an allocated mesh- initialize */
266     smInit_mesh(sm)
267     SM *sm;
268 gwlarson 3.1 {
269 gwlarson 3.8 /* Reset the triangle counters */
270     SM_NUM_TRI(sm) = 0;
271     SM_SAMPLE_TRIS(sm) = 0;
272     SM_FREE_TRIS(sm) = -1;
273 gwlarson 3.10 SM_AVAILABLE_TRIS(sm) = -1;
274 gwlarson 3.8 smClear_flags(sm,-1);
275 gwlarson 3.1 }
276    
277    
278 gwlarson 3.8 /* Clear the SM for reuse: free any extra allocated memory:does initialization
279     as well
280     */
281 gwlarson 3.1 smClear(sm)
282     SM *sm;
283     {
284     smClear_samples(sm);
285     smClear_locator(sm);
286 gwlarson 3.8 smInit_mesh(sm);
287 gwlarson 3.1 }
288    
289 gwlarson 3.8 static int realloc_cnt=0;
290    
291 gwlarson 3.1 int
292 gwlarson 3.8 smRealloc_mesh(sm)
293 gwlarson 3.1 SM *sm;
294 gwlarson 3.8 {
295     int fbytes,i,max_tris,max_verts;
296    
297     if(realloc_cnt <=0)
298     realloc_cnt = SM_MAX_TRIS(sm);
299    
300     max_tris = SM_MAX_TRIS(sm) + realloc_cnt;
301     fbytes = FLAG_BYTES(max_tris);
302    
303     if(!(SM_TRIS(sm) = (TRI *)realloc(SM_TRIS(sm),max_tris*sizeof(TRI))))
304     goto memerr;
305    
306     for(i=0; i< T_FLAGS; i++)
307     if(!(SM_NTH_FLAGS(sm,i)=(int4 *)realloc((char *)SM_NTH_FLAGS(sm,i),fbytes)))
308     goto memerr;
309    
310     SM_MAX_TRIS(sm) = max_tris;
311     return(max_tris);
312    
313     memerr:
314     error(SYSTEM, "out of memory in smRealloc_mesh()");
315     }
316    
317     /* Allocate and initialize a new mesh with max_verts and max_tris */
318     int
319     smAlloc_mesh(sm,max_verts,max_tris)
320     SM *sm;
321 gwlarson 3.1 int max_verts,max_tris;
322     {
323 gwlarson 3.8 int fbytes,i;
324 gwlarson 3.1
325 gwlarson 3.8 fbytes = FLAG_BYTES(max_tris);
326 gwlarson 3.1
327 gwlarson 3.8 if(!(SM_TRIS(sm) = (TRI *)realloc(NULL,max_tris*sizeof(TRI))))
328     goto memerr;
329 gwlarson 3.1
330 gwlarson 3.8 if(!(SM_VERTS(sm) = (VERT *)realloc(NULL,max_verts*sizeof(VERT))))
331     goto memerr;
332 gwlarson 3.1
333 gwlarson 3.8 for(i=0; i< T_FLAGS; i++)
334     if(!(SM_NTH_FLAGS(sm,i)=(int4 *)realloc(NULL,fbytes)))
335     goto memerr;
336    
337 gwlarson 3.1 SM_MAX_VERTS(sm) = max_verts;
338     SM_MAX_TRIS(sm) = max_tris;
339    
340 gwlarson 3.8 realloc_cnt = max_verts >> 1;
341 gwlarson 3.1
342 gwlarson 3.8 smInit_mesh(sm);
343    
344 gwlarson 3.1 return(max_tris);
345 gwlarson 3.8 memerr:
346     error(SYSTEM, "out of memory in smAlloc_mesh()");
347 gwlarson 3.1 }
348    
349    
350 gwlarson 3.10 int
351     compress_set(os)
352     OBJECT *os;
353     {
354     int i,j;
355    
356     for (i = 1, j = os[0]; i <= j; i++)
357     {
358     if(SM_T_ID_VALID(smMesh,os[i]))
359     continue;
360     if(i==j)
361     break;
362     while(!SM_T_ID_VALID(smMesh,os[j]) && (i < j))
363     j--;
364     if(i==j)
365     break;
366     os[i] = os[j--];
367     }
368     i--;
369    
370     os[0] = i;
371     return(i);
372    
373     }
374    
375 gwlarson 3.1 int
376 gwlarson 3.8 smAlloc_tri(sm)
377 gwlarson 3.1 SM *sm;
378     {
379 gwlarson 3.8 int id;
380    
381     /* First check if there are any tris on the free list */
382     /* Otherwise, have we reached the end of the list yet?*/
383     if(SM_NUM_TRI(sm) < SM_MAX_TRIS(sm))
384     return(SM_NUM_TRI(sm)++);
385 gwlarson 3.10 if((id = SM_AVAILABLE_TRIS(sm)) != -1)
386     {
387     SM_AVAILABLE_TRIS(sm) = T_NEXT_AVAILABLE(SM_NTH_TRI(sm,id));
388     return(id);
389     }
390 gwlarson 3.8 if((id = SM_FREE_TRIS(sm)) != -1)
391     {
392 gwlarson 3.10 dosets(compress_set);
393     SM_AVAILABLE_TRIS(sm) = T_NEXT_FREE(SM_NTH_TRI(sm,id));
394     SM_FREE_TRIS(sm) = -1;
395 gwlarson 3.8 return(id);
396     }
397     /*Else: realloc */
398     smRealloc_mesh(sm);
399     return(SM_NUM_TRI(sm)++);
400 gwlarson 3.1 }
401    
402 gwlarson 3.8 smFree_mesh(sm)
403     SM *sm;
404     {
405     int i;
406    
407     if(SM_TRIS(sm))
408     free(SM_TRIS(sm));
409     if(SM_VERTS(sm))
410     free(SM_VERTS(sm));
411     for(i=0; i< T_FLAGS; i++)
412     if(SM_NTH_FLAGS(sm,i))
413     free(SM_NTH_FLAGS(sm,i));
414     }
415    
416    
417 gwlarson 3.1 /* Initialize/clear global smL sample list for at least n samples */
418     smAlloc(max_samples)
419     register int max_samples;
420     {
421     unsigned nbytes;
422     register unsigned i;
423     int total_points;
424 gwlarson 3.8 int max_tris,n;
425 gwlarson 3.1
426 gwlarson 3.8 n = max_samples;
427 gwlarson 3.1 /* If this is the first call, allocate sample,vertex and triangle lists */
428     if(!smMesh)
429     {
430     if(!(smMesh = (SM *)malloc(sizeof(SM))))
431 gwlarson 3.8 error(SYSTEM,"smAlloc():Unable to allocate memory\n");
432 gwlarson 3.1 bzero(smMesh,sizeof(SM));
433     }
434     else
435     { /* If existing structure: first deallocate */
436 gwlarson 3.8 smFree_mesh(smMesh);
437     smFree_samples(smMesh);
438     smFree_locator(smMesh);
439     }
440 gwlarson 3.1
441     /* First allocate at least n samples + extra points:at least enough
442     necessary to form the BASE MESH- Default = 4;
443     */
444 gwlarson 3.8 SM_SAMP(smMesh) = sAlloc(&n,SM_EXTRA_POINTS);
445 gwlarson 3.1
446 gwlarson 3.8 total_points = n + SM_EXTRA_POINTS;
447 gwlarson 3.1
448 gwlarson 3.8 max_tris = total_points*4;
449 gwlarson 3.1 /* Now allocate space for mesh vertices and triangles */
450 gwlarson 3.8 max_tris = smAlloc_mesh(smMesh, total_points, max_tris);
451 gwlarson 3.1
452     /* Initialize the structure for point,triangle location.
453     */
454     smAlloc_locator(smMesh);
455     }
456    
457    
458    
459 gwlarson 3.8 smInit_sm(sm,vp)
460 gwlarson 3.1 SM *sm;
461     FVECT vp;
462     {
463    
464     smDist_sum = 0;
465     smNew_tri_cnt = 0;
466    
467 gwlarson 3.8 VCOPY(SM_VIEW_CENTER(sm),vp);
468     smInit_locator(sm,vp);
469     smInit_samples(sm);
470     smInit_mesh(sm);
471 gwlarson 3.1 smCreate_base_mesh(sm,SM_DEFAULT);
472     }
473    
474     /*
475     * int
476     * smInit(n) : Initialize/clear data structures for n entries
477     * int n;
478     *
479     * This routine allocates/initializes the sample, mesh, and point-location
480     * structures for at least n samples.
481     * If n is <= 0, then clear data structures. Returns number samples
482     * actually allocated.
483     */
484    
485     int
486     smInit(n)
487     register int n;
488     {
489     int max_vertices;
490    
491 gwlarson 3.8
492 gwlarson 3.1 /* If n <=0, Just clear the existing structures */
493     if(n <= 0)
494     {
495     smClear(smMesh);
496     return(0);
497     }
498    
499     /* Total mesh vertices includes the sample points and the extra vertices
500     to form the base mesh
501     */
502     max_vertices = n + SM_EXTRA_POINTS;
503    
504     /* If the current mesh contains enough room, clear and return */
505 gwlarson 3.8 if(smMesh && (max_vertices <= SM_MAX_VERTS(smMesh)) && SM_MAX_SAMP(smMesh) <=
506     n && SM_MAX_POINTS(smMesh) <= max_vertices)
507 gwlarson 3.1 {
508     smClear(smMesh);
509     return(SM_MAX_SAMP(smMesh));
510     }
511     /* Otherwise- mesh must be allocated with the appropriate number of
512     samples
513     */
514     smAlloc(n);
515    
516     return(SM_MAX_SAMP(smMesh));
517     }
518    
519    
520 gwlarson 3.10 smLocator_apply(sm,v0,v1,v2,func)
521 gwlarson 3.7 SM *sm;
522     FVECT v0,v1,v2;
523 gwlarson 3.10 FUNC func;
524 gwlarson 3.1 {
525     STREE *st;
526 gwlarson 3.10 FVECT tri[3];
527 gwlarson 3.1
528     st = SM_LOCATOR(sm);
529    
530 gwlarson 3.10 VSUB(tri[0],v0,SM_VIEW_CENTER(sm));
531     VSUB(tri[1],v1,SM_VIEW_CENTER(sm));
532     VSUB(tri[2],v2,SM_VIEW_CENTER(sm));
533     stVisit(st,tri,func);
534 gwlarson 3.1
535     }
536    
537    
538 gwlarson 3.10 /* NEW INSERTION!!! ******************************************************/
539 gwlarson 3.8 QUADTREE
540 gwlarson 3.10 insert_tri(argptr,root,qt,q0,q1,q2,t0,t1,t2,dt10,dt21,dt02,scale,
541     s0,s1,s2,sq0,sq1,sq2,rev,f,n)
542     int *argptr;
543     int root;
544     QUADTREE qt;
545     BCOORD q0[3],q1[3],q2[3];
546     BCOORD t0[3],t1[3],t2[3];
547     BCOORD dt10[3],dt21[3],dt02[3];
548     unsigned int scale,s0,s1,s2,sq0,sq1,sq2,rev;
549     FUNC f;
550     int n;
551 gwlarson 3.8 {
552 gwlarson 3.10 OBJECT *optr,*tptr,t_set[QT_MAXSET+1];
553     int i,t_id;
554 gwlarson 3.8 TRI *t;
555 gwlarson 3.10 FVECT tri[3];
556     BCOORD b0[3],b1[3],b2[3];
557     BCOORD tb0[3],tb1[3],tb2[3];
558     BCOORD db10[3],db21[3],db02[3];
559     BCOORD a,b,c,qa,qb,qc;
560     FUNC tfunc;
561 gwlarson 3.8
562 gwlarson 3.10 t_id = *argptr;
563 gwlarson 3.8
564 gwlarson 3.10 /* If the set is empty - just add */
565     if(QT_IS_EMPTY(qt))
566     return(qtadduelem(qt,t_id));
567    
568     /* If the set size is below expansion threshold,or at maximum
569     number of quadtree levels : just add */
570     optr = qtqueryset(qt);
571     if(QT_SET_CNT(optr) < QT_SET_THRESHOLD || (n >= QT_MAX_LEVELS))
572     return(qtadduelem(qt,t_id));
573     /* otherwise: expand node- and insert tri, and reinsert existing tris
574     in set to children of new node
575     */
576     /* First tri compressing set */
577     #ifdef NEWSETS0
578     if(qtcompressuelem(qt,compress_set) < QT_SET_THRESHOLD)
579     /* Get set: If greater than standard size: allocate memory to hold */
580     return(qtadduelem(qt,t_id));
581     #endif
582     if(QT_SET_CNT(optr) > QT_MAXSET)
583     {
584     if(!(tptr = (OBJECT *)malloc((QT_SET_CNT(optr)+1)*sizeof(OBJECT))))
585     goto memerr;
586     }
587     else
588     tptr = t_set;
589     qtgetset(tptr,qt);
590    
591     /* subdivide node */
592 gwlarson 3.8 qtfreeleaf(qt);
593     qtSubdivide(qt);
594    
595 gwlarson 3.10 a = q1[0]; b = q0[1]; c = q0[2];
596     qa = q0[0]; qb = q1[1]; qc = q2[2];
597     /* First add current tri: have all of the information */
598     if(rev)
599     qt = qtInsert_tri_rev(root,qt,q0,q1,q2,t0,t1,t2,dt10,dt21,dt02,scale,
600     s0,s1,s2,sq0,sq1,sq2,f,n);
601     else
602     qt = qtInsert_tri(root,qt,q0,q1,q2,t0,t1,t2,dt10,dt21,dt02,scale,
603     s0,s1,s2,sq0,sq1,sq2,f,n);
604    
605     /* Now add in all of the rest; */
606     F_FUNC(tfunc) = F_FUNC(f);
607 gwlarson 3.8 for(optr = QT_SET_PTR(tptr),i=QT_SET_CNT(tptr); i > 0; i--)
608     {
609     t_id = QT_SET_NEXT_ELEM(optr);
610     t = SM_NTH_TRI(smMesh,t_id);
611     if(!T_IS_VALID(t))
612     continue;
613 gwlarson 3.10 F_ARGS(tfunc) = &t_id;
614     VSUB(tri[0],SM_T_NTH_WV(smMesh,t,0),SM_VIEW_CENTER(smMesh));
615     VSUB(tri[1],SM_T_NTH_WV(smMesh,t,1),SM_VIEW_CENTER(smMesh));
616     VSUB(tri[2],SM_T_NTH_WV(smMesh,t,2),SM_VIEW_CENTER(smMesh));
617     convert_tri_to_frame(root,tri,b0,b1,b2,db10,db21,db02);
618    
619     /* Calculate appropriate sidedness values */
620     SIDES_GTR(b0,b1,b2,s0,s1,s2,a,b,c);
621     SIDES_GTR(b0,b1,b2,sq0,sq1,sq2,qa,qb,qc);
622     if(rev)
623     qt = qtInsert_tri_rev(root,qt,q0,q1,q2,b0,b1,b2,db10,db21,db02,scale,
624     s0,s1,s2,sq0,sq1,sq2,tfunc,n);
625     else
626     qt = qtInsert_tri(root,qt,q0,q1,q2,b0,b1,b2,db10,db21,db02,scale,
627     s0,s1,s2,sq0,sq1,sq2,tfunc,n);
628 gwlarson 3.8 }
629     /* If we allocated memory: free it */
630     if( tptr != t_set)
631     free(tptr);
632    
633     return(qt);
634     memerr:
635     error(SYSTEM,"expand_node():Unable to allocate memory");
636    
637 gwlarson 3.5 }
638    
639    
640 gwlarson 3.8
641 gwlarson 3.10 smLocator_add_tri(sm,t_id,v0_id,v1_id,v2_id)
642 gwlarson 3.1 SM *sm;
643     int t_id;
644 gwlarson 3.7 int v0_id,v1_id,v2_id;
645 gwlarson 3.1 {
646 gwlarson 3.10 FVECT tri[3];
647     FUNC f;
648 gwlarson 3.8
649 gwlarson 3.10 #ifdef DEBUG
650     if(T_NTH_NBR(SM_NTH_TRI(sm,t_id),0)== -1 ||
651     T_NTH_NBR(SM_NTH_TRI(sm,t_id),1)== -1 ||
652     T_NTH_NBR(SM_NTH_TRI(sm,t_id),2)== -1)
653     error("Invalid tri: smLocator_add_tri()\n");
654 gwlarson 3.1
655 gwlarson 3.10 if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,t_id),0))) ||
656     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,t_id),1))) ||
657     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,t_id),2))))
658     error("Invalid tri: smLocator_add_tri()\n");
659 gwlarson 3.9 #endif
660 gwlarson 3.8
661 gwlarson 3.10 VSUB(tri[0],SM_NTH_WV(sm,v0_id),SM_VIEW_CENTER(sm));
662     VSUB(tri[1],SM_NTH_WV(sm,v1_id),SM_VIEW_CENTER(sm));
663     VSUB(tri[2],SM_NTH_WV(sm,v2_id),SM_VIEW_CENTER(sm));
664 gwlarson 3.1
665 gwlarson 3.10 F_FUNC(f) = insert_tri;
666     F_ARGS(f) = &t_id;
667     stInsert_tri(SM_LOCATOR(sm),tri,f);
668 gwlarson 3.1 }
669    
670     /* Add a triangle to the base array with vertices v1-v2-v3 */
671     int
672 gwlarson 3.8 smAdd_tri(sm, v0_id,v1_id,v2_id)
673 gwlarson 3.1 SM *sm;
674     int v0_id,v1_id,v2_id;
675     {
676     int t_id;
677     TRI *t;
678 gwlarson 3.9 #ifdef DEBUG
679     if(v0_id==v1_id || v0_id==v2_id || v1_id==v2_id)
680 gwlarson 3.10 error(CONSISTENCY,"smAdd_tri: invalid vertex ids\n");
681 gwlarson 3.9 #endif
682 gwlarson 3.8 t_id = smAlloc_tri(sm);
683 gwlarson 3.1
684     if(t_id == -1)
685     return(t_id);
686    
687     t = SM_NTH_TRI(sm,t_id);
688 gwlarson 3.8
689 gwlarson 3.1 T_CLEAR_NBRS(t);
690 gwlarson 3.8 /* set the triangle vertex ids */
691     T_NTH_V(t,0) = v0_id;
692     T_NTH_V(t,1) = v1_id;
693     T_NTH_V(t,2) = v2_id;
694 gwlarson 3.1
695 gwlarson 3.8 SM_NTH_VERT(sm,v0_id) = t_id;
696     SM_NTH_VERT(sm,v1_id) = t_id;
697     SM_NTH_VERT(sm,v2_id) = t_id;
698    
699 gwlarson 3.1 if(SM_BASE_ID(sm,v0_id) || SM_BASE_ID(sm,v1_id) || SM_BASE_ID(sm,v2_id))
700     {
701     smClear_tri_flags(sm,t_id);
702     SM_SET_NTH_T_BASE(sm,t_id);
703     }
704     else
705     {
706 gwlarson 3.8 SM_CLR_NTH_T_BASE(sm,t_id);
707 gwlarson 3.1 SM_SET_NTH_T_ACTIVE(sm,t_id);
708     SM_SET_NTH_T_NEW(sm,t_id);
709 gwlarson 3.8 S_SET_FLAG(T_NTH_V(t,0));
710     S_SET_FLAG(T_NTH_V(t,1));
711     S_SET_FLAG(T_NTH_V(t,2));
712     SM_SAMPLE_TRIS(sm)++;
713 gwlarson 3.1 smNew_tri_cnt++;
714     }
715    
716     /* return initialized triangle */
717     return(t_id);
718     }
719    
720    
721     void
722 gwlarson 3.10 smTris_swap_edge(sm,t_id,t1_id,e,e1,tn_id,tn1_id,add_ptr)
723 gwlarson 3.1 SM *sm;
724     int t_id,t1_id;
725     int e,e1;
726 gwlarson 3.5 int *tn_id,*tn1_id;
727 gwlarson 3.7 LIST **add_ptr;
728 gwlarson 3.1 {
729     int verts[3],enext,eprev,e1next,e1prev;
730 gwlarson 3.9 TRI *n,*ta,*tb,*t,*t1;
731 gwlarson 3.1 FVECT p1,p2,p3;
732     int ta_id,tb_id;
733 gwlarson 3.9 /* form new diagonal (e relative to t, and e1 relative to t1)
734     defined by quadrilateral formed by t,t1- swap for the opposite diagonal
735 gwlarson 3.1 */
736     enext = (e+1)%3;
737     eprev = (e+2)%3;
738     e1next = (e1+1)%3;
739     e1prev = (e1+2)%3;
740 gwlarson 3.8 verts[e] = T_NTH_V(SM_NTH_TRI(sm,t_id),e);
741 gwlarson 3.9 verts[enext] = T_NTH_V(SM_NTH_TRI(sm,t_id),enext);
742     verts[eprev] = T_NTH_V(SM_NTH_TRI(sm,t1_id),e1);
743 gwlarson 3.8 ta_id = smAdd_tri(sm,verts[0],verts[1],verts[2]);
744 gwlarson 3.7 *add_ptr = push_data(*add_ptr,ta_id);
745 gwlarson 3.8 verts[e1] = T_NTH_V(SM_NTH_TRI(sm,t1_id),e1);
746 gwlarson 3.9 verts[e1next] = T_NTH_V(SM_NTH_TRI(sm,t1_id),e1next);
747     verts[e1prev] = T_NTH_V(SM_NTH_TRI(sm,t_id),e);
748 gwlarson 3.8 tb_id = smAdd_tri(sm,verts[0],verts[1],verts[2]);
749 gwlarson 3.10
750 gwlarson 3.7 *add_ptr = push_data(*add_ptr,tb_id);
751 gwlarson 3.9 ta = SM_NTH_TRI(sm,ta_id);
752     tb = SM_NTH_TRI(sm,tb_id);
753     t = SM_NTH_TRI(sm,t_id);
754     t1 = SM_NTH_TRI(sm,t1_id);
755 gwlarson 3.1 /* set the neighbors */
756 gwlarson 3.9 T_NTH_NBR(ta,e) = T_NTH_NBR(t1,e1next);
757     T_NTH_NBR(tb,e1) = T_NTH_NBR(t,enext);
758     T_NTH_NBR(ta,enext)= tb_id;
759     T_NTH_NBR(tb,e1next)= ta_id;
760     T_NTH_NBR(ta,eprev)=T_NTH_NBR(t,eprev);
761     T_NTH_NBR(tb,e1prev)=T_NTH_NBR(t1,e1prev);
762 gwlarson 3.1
763 gwlarson 3.10
764     #ifdef DEBUG
765     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(ta,0))) ||
766     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(ta,1))) ||
767     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(ta,2))) ||
768     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(tb,0))) ||
769     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(tb,1))) ||
770     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(tb,2))))
771     goto Ltri_error;
772     if(SM_NTH_TRI(sm,T_NTH_NBR(ta,0))==SM_NTH_TRI(sm,T_NTH_NBR(ta,1)) ||
773     SM_NTH_TRI(sm,T_NTH_NBR(ta,0))==SM_NTH_TRI(sm,T_NTH_NBR(ta,2)) ||
774     SM_NTH_TRI(sm,T_NTH_NBR(ta,1))==SM_NTH_TRI(sm,T_NTH_NBR(ta,2)) ||
775     SM_NTH_TRI(sm,T_NTH_NBR(tb,0))==SM_NTH_TRI(sm,T_NTH_NBR(tb,1)) ||
776     SM_NTH_TRI(sm,T_NTH_NBR(tb,0))==SM_NTH_TRI(sm,T_NTH_NBR(tb,2)) ||
777     SM_NTH_TRI(sm,T_NTH_NBR(tb,1))==SM_NTH_TRI(sm,T_NTH_NBR(tb,2)) )
778     goto Ltri_error;
779     #endif
780 gwlarson 3.1 /* Reset neighbor pointers of original neighbors */
781 gwlarson 3.9 n = SM_NTH_TRI(sm,T_NTH_NBR(t,enext));
782 gwlarson 3.10 #ifdef DEBUG
783     if(!T_IS_VALID(n))
784     goto Ltri_error;
785     #endif
786 gwlarson 3.1 T_NTH_NBR(n,T_NTH_NBR_PTR(t_id,n)) = tb_id;
787 gwlarson 3.9 n = SM_NTH_TRI(sm,T_NTH_NBR(t,eprev));
788 gwlarson 3.10 #ifdef DEBUG
789     if(!T_IS_VALID(n))
790     goto Ltri_error;
791     #endif
792 gwlarson 3.1 T_NTH_NBR(n,T_NTH_NBR_PTR(t_id,n)) = ta_id;
793    
794 gwlarson 3.9 n = SM_NTH_TRI(sm,T_NTH_NBR(t1,e1next));
795 gwlarson 3.10 #ifdef DEBUG
796     if(!T_IS_VALID(n))
797     goto Ltri_error;
798     #endif
799 gwlarson 3.1 T_NTH_NBR(n,T_NTH_NBR_PTR(t1_id,n)) = ta_id;
800 gwlarson 3.9 n = SM_NTH_TRI(sm,T_NTH_NBR(t1,e1prev));
801 gwlarson 3.10 #ifdef DEBUG
802     if(!T_IS_VALID(n))
803     goto Ltri_error;
804     #endif
805 gwlarson 3.1 T_NTH_NBR(n,T_NTH_NBR_PTR(t1_id,n)) = tb_id;
806    
807 gwlarson 3.10 #ifdef DEBUG
808     T_VALID_FLAG(SM_NTH_TRI(sm,t_id))=-1;
809     T_VALID_FLAG(SM_NTH_TRI(sm,t1_id))=-1;
810     #endif
811 gwlarson 3.5
812 gwlarson 3.10 remove_from_list(t_id,add_ptr);
813     remove_from_list(t1_id,add_ptr);
814     #if 1
815     smDelete_tri(sm,t_id);
816     smDelete_tri(sm,t1_id);
817     #else
818     *del_ptr = push_data(*del_ptr,t_id);
819     *del_ptr = push_data(*del_ptr,t1_id);
820     #endif
821 gwlarson 3.7
822 gwlarson 3.10
823    
824 gwlarson 3.1 *tn_id = ta_id;
825     *tn1_id = tb_id;
826 gwlarson 3.10
827     #ifdef DEBUG
828     if(T_NTH_NBR(SM_NTH_TRI(sm,ta_id),0)== -1 ||
829     T_NTH_NBR(SM_NTH_TRI(sm,ta_id),1)== -1 ||
830     T_NTH_NBR(SM_NTH_TRI(sm,ta_id),2)== -1 ||
831     T_NTH_NBR(SM_NTH_TRI(sm,tb_id),0)== -1 ||
832     T_NTH_NBR(SM_NTH_TRI(sm,tb_id),1)== -1 ||
833     T_NTH_NBR(SM_NTH_TRI(sm,tb_id),2)== -1)
834     goto Ltri_error;
835     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,ta_id),0))) ||
836     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,ta_id),1))) ||
837     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,ta_id),2))))
838     goto Ltri_error;
839     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,tb_id),0))) ||
840     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,tb_id),1))) ||
841     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,tb_id),2))))
842     goto Ltri_error;
843     #endif
844     return;
845     Ltri_error:
846     error(CONSISTENCY,"Invalid tri: smTris_swap_edge()\n");
847 gwlarson 3.1 }
848    
849 gwlarson 3.10 smUpdate_locator(sm,add_list)
850 gwlarson 3.3 SM *sm;
851 gwlarson 3.7 LIST *add_list;
852 gwlarson 3.3 {
853 gwlarson 3.7 int t_id,i;
854 gwlarson 3.3 TRI *t;
855 gwlarson 3.7 OBJECT *optr;
856    
857 gwlarson 3.3 while(add_list)
858     {
859     t_id = pop_list(&add_list);
860 gwlarson 3.6 t = SM_NTH_TRI(sm,t_id);
861 gwlarson 3.10 smLocator_add_tri(sm,t_id,T_NTH_V(t,0),T_NTH_V(t,1),T_NTH_V(t,2));
862 gwlarson 3.3 }
863 gwlarson 3.10 }
864 gwlarson 3.7
865 gwlarson 3.3 int
866 gwlarson 3.10 smFix_tris(sm,id,tlist,add_list)
867 gwlarson 3.1 SM *sm;
868     int id;
869 gwlarson 3.8 LIST *tlist,*add_list;
870 gwlarson 3.1 {
871     TRI *t,*t_opp;
872 gwlarson 3.9 FVECT p,p0,p1,p2;
873 gwlarson 3.3 int e,e1,swapped = 0;
874 gwlarson 3.1 int t_id,t_opp_id;
875 gwlarson 3.10 LIST *del_list=NULL;
876 gwlarson 3.3
877     VSUB(p,SM_NTH_WV(sm,id),SM_VIEW_CENTER(sm));
878 gwlarson 3.1 while(tlist)
879     {
880 gwlarson 3.3 t_id = pop_list(&tlist);
881 gwlarson 3.9 #ifdef DEBUG
882     if(t_id==INVALID || t_id > smMesh->num_tri)
883     error(CONSISTENCY,"Invalid tri id smFix_tris()\n");
884     #endif
885 gwlarson 3.1 t = SM_NTH_TRI(sm,t_id);
886 gwlarson 3.9 e = T_WHICH_V(t,id);
887 gwlarson 3.1 t_opp_id = T_NTH_NBR(t,e);
888 gwlarson 3.9 #ifdef DEBUG
889     if(t_opp_id==INVALID || t_opp_id > smMesh->num_tri)
890     error(CONSISTENCY,"Invalid tri id smFix_tris()\n");
891     #endif
892     t_opp = SM_NTH_TRI(sm,t_opp_id);
893 gwlarson 3.10 #ifdef DEBUG
894     if(!T_IS_VALID(t_opp))
895     error(CONSISTENCY,"Invalid trismFix_tris()\n");
896     #endif
897 gwlarson 3.9 smDir_in_cone(sm,p0,T_NTH_V(t_opp,0));
898     smDir_in_cone(sm,p1,T_NTH_V(t_opp,1));
899     smDir_in_cone(sm,p2,T_NTH_V(t_opp,2));
900     if(point_in_cone(p,p0,p1,p2))
901 gwlarson 3.1 {
902     swapped = 1;
903     e1 = T_NTH_NBR_PTR(t_id,t_opp);
904     /* check list for t_opp and Remove if there */
905     remove_from_list(t_opp_id,&tlist);
906 gwlarson 3.3 smTris_swap_edge(sm,t_id,t_opp_id,e,e1,&t_id,&t_opp_id,
907 gwlarson 3.10 &add_list);
908 gwlarson 3.1 tlist = push_data(tlist,t_id);
909     tlist = push_data(tlist,t_opp_id);
910     }
911     }
912 gwlarson 3.10 #if 0
913     while(del_list)
914     {
915     t_id = pop_list(&del_list);
916     smDelete_tri(sm,t_id);
917     }
918     #endif
919     smUpdate_locator(sm,add_list);
920    
921 gwlarson 3.1 return(swapped);
922     }
923    
924     /* Give the vertex "id" and a triangle "t" that it belongs to- return the
925     next nbr in a counter clockwise order about vertex "id"
926     */
927     int
928     smTri_next_ccw_nbr(sm,t,id)
929     SM *sm;
930     TRI *t;
931     int id;
932     {
933 gwlarson 3.10 int which;
934 gwlarson 3.8 int nbr_id;
935 gwlarson 3.10
936     /* Want the edge for which "id" is the destination */
937     which = T_WHICH_V(t,id);
938     if(which== INVALID)
939     return(which);
940 gwlarson 3.1
941 gwlarson 3.10 nbr_id = T_NTH_NBR(t,(which + 1)% 3);
942 gwlarson 3.8 return(nbr_id);
943 gwlarson 3.1 }
944    
945     smClear_tri_flags(sm,id)
946     SM *sm;
947     int id;
948     {
949     int i;
950    
951     for(i=0; i < T_FLAGS; i++)
952 gwlarson 3.8 SM_CLR_NTH_T_FLAG(sm,id,i);
953 gwlarson 3.1
954     }
955    
956     /* Locate the point-id in the point location structure: */
957     int
958 gwlarson 3.10 smInsert_samp(sm,s_id,tri_id,on,which)
959 gwlarson 3.1 SM *sm;
960     int s_id,tri_id;
961 gwlarson 3.10 int on,which;
962 gwlarson 3.1 {
963 gwlarson 3.10 int v_id[3],topp_id,i;
964     int t0_id,t1_id,t2_id,t3_id,replaced;
965     int e0,e1,e2,opp0,opp1,opp2,opp_id;
966 gwlarson 3.7 LIST *tlist,*add_list;
967 gwlarson 3.1 FVECT npt;
968 gwlarson 3.10 TRI *tri,*nbr,*topp;
969 gwlarson 3.1
970 gwlarson 3.10
971 gwlarson 3.7 add_list = NULL;
972 gwlarson 3.10 for(i=0; i<3;i++)
973     v_id[i] = T_NTH_V(SM_NTH_TRI(sm,tri_id),i);
974 gwlarson 3.1
975 gwlarson 3.10 /* If falls on existing edge */
976     if(on == ON_E)
977     {
978     tri = SM_NTH_TRI(sm,tri_id);
979     topp_id = T_NTH_NBR(tri,which);
980     e0 = which;
981     e1 = (which+1)%3;
982     e2 = (which+2)%3;
983     topp = SM_NTH_TRI(sm,topp_id);
984     opp0 = T_NTH_NBR_PTR(tri_id,topp);
985     opp_id = T_NTH_V(topp,opp0);
986     opp1 = (opp0+1)%3;
987     opp2 = (opp0+2)%3;
988    
989     /* Create 4 triangles */
990     /* /|\ /|\
991     / | \ / | \
992     / * \ ----> /__|__\
993     \ | / \ | /
994     \ | / \ | /
995     \|/ \|/
996     */
997     t0_id = smAdd_tri(sm,s_id,v_id[e0],v_id[e1]);
998     add_list = push_data(add_list,t0_id);
999     t1_id = smAdd_tri(sm,s_id,v_id[e2],v_id[e0]);
1000     add_list = push_data(add_list,t1_id);
1001     t2_id = smAdd_tri(sm,s_id,v_id[e1],opp_id);
1002     add_list = push_data(add_list,t2_id);
1003     t3_id = smAdd_tri(sm,s_id,opp_id,v_id[e2]);
1004     add_list = push_data(add_list,t3_id);
1005     /* CAUTION: tri must be assigned after Add_tri: pointers may change*/
1006     tri = SM_NTH_TRI(sm,tri_id);
1007     topp =SM_NTH_TRI(sm,topp_id);
1008 gwlarson 3.8
1009 gwlarson 3.10 /* Set the neighbor pointers for the new tris */
1010     T_NTH_NBR(SM_NTH_TRI(sm,t0_id),0) = T_NTH_NBR(tri,e2);
1011     T_NTH_NBR(SM_NTH_TRI(sm,t0_id),1) = t2_id;
1012     T_NTH_NBR(SM_NTH_TRI(sm,t0_id),2) = t1_id;
1013     T_NTH_NBR(SM_NTH_TRI(sm,t1_id),0) = T_NTH_NBR(tri,e1);
1014     T_NTH_NBR(SM_NTH_TRI(sm,t1_id),1) = t0_id;
1015     T_NTH_NBR(SM_NTH_TRI(sm,t1_id),2) = t3_id;
1016     T_NTH_NBR(SM_NTH_TRI(sm,t2_id),0) = T_NTH_NBR(topp,opp1);
1017     T_NTH_NBR(SM_NTH_TRI(sm,t2_id),1) = t3_id;
1018     T_NTH_NBR(SM_NTH_TRI(sm,t2_id),2) = t0_id;
1019     T_NTH_NBR(SM_NTH_TRI(sm,t3_id),0) = T_NTH_NBR(topp,opp2);
1020     T_NTH_NBR(SM_NTH_TRI(sm,t3_id),1) = t1_id;
1021     T_NTH_NBR(SM_NTH_TRI(sm,t3_id),2) = t2_id;
1022 gwlarson 3.3
1023 gwlarson 3.10 /* Reset the neigbor pointers for the neighbors of the original */
1024     nbr = SM_NTH_TRI(sm,T_NTH_NBR(tri,e1));
1025     T_NTH_NBR(nbr,T_NTH_NBR_PTR(tri_id,nbr)) = t1_id;
1026     nbr = SM_NTH_TRI(sm,T_NTH_NBR(tri,e2));
1027     T_NTH_NBR(nbr,T_NTH_NBR_PTR(tri_id,nbr)) = t0_id;
1028     nbr = SM_NTH_TRI(sm,T_NTH_NBR(topp,opp1));
1029     T_NTH_NBR(nbr,T_NTH_NBR_PTR(topp_id,nbr)) = t2_id;
1030     nbr = SM_NTH_TRI(sm,T_NTH_NBR(topp,opp2));
1031     T_NTH_NBR(nbr,T_NTH_NBR_PTR(topp_id,nbr)) = t3_id;
1032     #ifdef DEBUG
1033     T_VALID_FLAG(SM_NTH_TRI(sm,tri_id))=-1;
1034     T_VALID_FLAG(SM_NTH_TRI(sm,topp_id))=-1;
1035     #endif
1036     tlist = push_data(NULL,t0_id);
1037     tlist = push_data(tlist,t1_id);
1038     tlist = push_data(tlist,t2_id);
1039     tlist = push_data(tlist,t3_id);
1040     }
1041     else
1042     {
1043     /* Create 3 triangles */
1044     /* / \ /|\
1045     / \ / | \
1046     / * \ ----> / | \
1047     / \ / / \ \
1048     / \ / / \ \
1049     /___________\ //_________\\
1050     */
1051 gwlarson 3.6
1052 gwlarson 3.10 t0_id = smAdd_tri(sm,s_id,v_id[0],v_id[1]);
1053     add_list = push_data(add_list,t0_id);
1054     t1_id = smAdd_tri(sm,s_id,v_id[1],v_id[2]);
1055     add_list = push_data(add_list,t1_id);
1056     t2_id = smAdd_tri(sm,s_id,v_id[2],v_id[0]);
1057     add_list = push_data(add_list,t2_id);
1058 gwlarson 3.3
1059 gwlarson 3.10 /* CAUTION: tri must be assigned after Add_tri: pointers may change*/
1060     tri = SM_NTH_TRI(sm,tri_id);
1061     /* Set the neighbor pointers for the new tris */
1062     T_NTH_NBR(SM_NTH_TRI(sm,t0_id),0) = T_NTH_NBR(tri,2);
1063     T_NTH_NBR(SM_NTH_TRI(sm,t0_id),1) = t1_id;
1064     T_NTH_NBR(SM_NTH_TRI(sm,t0_id),2) = t2_id;
1065     T_NTH_NBR(SM_NTH_TRI(sm,t1_id),0) = T_NTH_NBR(tri,0);
1066     T_NTH_NBR(SM_NTH_TRI(sm,t1_id),1) = t2_id;
1067     T_NTH_NBR(SM_NTH_TRI(sm,t1_id),2) = t0_id;
1068     T_NTH_NBR(SM_NTH_TRI(sm,t2_id),0) = T_NTH_NBR(tri,1);
1069     T_NTH_NBR(SM_NTH_TRI(sm,t2_id),1) = t0_id;
1070     T_NTH_NBR(SM_NTH_TRI(sm,t2_id),2) = t1_id;
1071 gwlarson 3.1
1072 gwlarson 3.10 /* Reset the neigbor pointers for the neighbors of the original */
1073     nbr = SM_NTH_TRI(sm,T_NTH_NBR(tri,0));
1074     T_NTH_NBR(nbr,T_NTH_NBR_PTR(tri_id,nbr)) = t1_id;
1075     nbr = SM_NTH_TRI(sm,T_NTH_NBR(tri,1));
1076     T_NTH_NBR(nbr,T_NTH_NBR_PTR(tri_id,nbr)) = t2_id;
1077     nbr = SM_NTH_TRI(sm,T_NTH_NBR(tri,2));
1078     T_NTH_NBR(nbr,T_NTH_NBR_PTR(tri_id,nbr)) = t0_id;
1079 gwlarson 3.6
1080 gwlarson 3.10 tlist = push_data(NULL,t0_id);
1081     tlist = push_data(tlist,t1_id);
1082     tlist = push_data(tlist,t2_id);
1083     }
1084 gwlarson 3.1 /* Fix up the new triangles*/
1085 gwlarson 3.10 smFix_tris(sm,s_id,tlist,add_list);
1086 gwlarson 3.1
1087 gwlarson 3.10 smDelete_tri(sm,tri_id);
1088     if(on==ON_E)
1089     smDelete_tri(sm,topp_id);
1090 gwlarson 3.1
1091 gwlarson 3.8 return(TRUE);
1092 gwlarson 3.1
1093 gwlarson 3.10 #ifdef DEBUG
1094     Ladderror:
1095     error(CONSISTENCY,"Invalid tri: smInsert_samp()\n");
1096     #endif
1097     }
1098    
1099 gwlarson 3.1 int
1100 gwlarson 3.10 smTri_in_set(sm,p,optr,onptr,whichptr)
1101 gwlarson 3.1 SM *sm;
1102 gwlarson 3.8 FVECT p;
1103     OBJECT *optr;
1104 gwlarson 3.10 int *onptr,*whichptr;
1105 gwlarson 3.1 {
1106 gwlarson 3.10 int i,t_id,on0,on1,on2;
1107 gwlarson 3.8 FVECT n,v0,v1,v2;
1108     TRI *t;
1109 gwlarson 3.10 double side;
1110    
1111 gwlarson 3.8 for (i = QT_SET_CNT(optr),optr = QT_SET_PTR(optr);i > 0; i--)
1112 gwlarson 3.1 {
1113 gwlarson 3.8 /* Find the first triangle that pt falls */
1114     t_id = QT_SET_NEXT_ELEM(optr);
1115     t = SM_NTH_TRI(sm,t_id);
1116     if(!T_IS_VALID(t))
1117     continue;
1118     VSUB(v0,SM_T_NTH_WV(sm,t,0),SM_VIEW_CENTER(sm));
1119     VSUB(v1,SM_T_NTH_WV(sm,t,1),SM_VIEW_CENTER(sm));
1120     VSUB(v2,SM_T_NTH_WV(sm,t,2),SM_VIEW_CENTER(sm));
1121    
1122 gwlarson 3.10 if(EQUAL_VEC3(v0,p))
1123     {
1124     *onptr = ON_V;
1125     *whichptr = 0;
1126     return(t_id);
1127     }
1128     if(EQUAL_VEC3(v1,p))
1129     {
1130     *onptr = ON_V;
1131     *whichptr = 1;
1132     return(t_id);
1133     }
1134     if(EQUAL_VEC3(v2,p))
1135     {
1136     *onptr = ON_V;
1137     *whichptr = 2;
1138     return(t_id);
1139     }
1140 gwlarson 3.8
1141 gwlarson 3.10 on0 = on1 =on2 = 0;
1142 gwlarson 3.9 VCROSS(n,v0,v1);
1143 gwlarson 3.10 side = DOT(n,p);
1144     if(ZERO(side))
1145     on2 = TRUE;
1146     else
1147     if(side >0.0)
1148     continue;
1149    
1150 gwlarson 3.9 VCROSS(n,v1,v2);
1151 gwlarson 3.10 side = DOT(n,p);
1152     if(ZERO(side))
1153     on0 = TRUE;
1154     else
1155     if(side >0.0)
1156     continue;
1157 gwlarson 3.8
1158 gwlarson 3.9 VCROSS(n,v2,v0);
1159 gwlarson 3.10 side = DOT(n,p);
1160     if(ZERO(side))
1161     on1 = TRUE;
1162     else
1163     if(side >0.0)
1164     continue;
1165     if(on0)
1166     {
1167     if(on1)
1168     {
1169     *onptr = ON_P;
1170     *whichptr = 2;
1171     }
1172     else
1173     if(on2)
1174     {
1175     *onptr = ON_P;
1176     *whichptr = 1;
1177     }
1178     else
1179     {
1180     *onptr = ON_E;
1181     *whichptr = 0;
1182     }
1183     return(t_id);
1184     }
1185     if(on1)
1186     {
1187     if(on2)
1188     {
1189     *onptr = ON_P;
1190     *whichptr = 0;
1191     }
1192     else
1193     {
1194     *onptr = ON_E;
1195     *whichptr = 1;
1196     }
1197     return(t_id);
1198     }
1199     if(on2)
1200     {
1201     *onptr = ON_E;
1202     *whichptr = 2;
1203     return(t_id);
1204     }
1205     *onptr = IN_T;
1206 gwlarson 3.8 return(t_id);
1207 gwlarson 3.1 }
1208 gwlarson 3.8 return(INVALID);
1209 gwlarson 3.1 }
1210    
1211 gwlarson 3.8 int
1212 gwlarson 3.10 smPointLocateTri(sm,p,onptr,whichptr)
1213 gwlarson 3.1 SM *sm;
1214 gwlarson 3.9 FVECT p;
1215 gwlarson 3.10 int *onptr,*whichptr;
1216 gwlarson 3.1 {
1217 gwlarson 3.9 QUADTREE qt,*optr;
1218 gwlarson 3.8 FVECT tpt;
1219 gwlarson 3.10 int tri_id;
1220 gwlarson 3.9
1221     VSUB(tpt,p,SM_VIEW_CENTER(sm));
1222 gwlarson 3.8 qt = smPointLocateCell(sm,tpt);
1223     if(QT_IS_EMPTY(qt))
1224     return(INVALID);
1225    
1226     optr = qtqueryset(qt);
1227 gwlarson 3.10 tri_id = smTri_in_set(sm,tpt,optr,onptr,whichptr);
1228    
1229     #ifdef DEBUG
1230     if(!T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,tri_id),0))) ||
1231     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,tri_id),1))) ||
1232     !T_IS_VALID(SM_NTH_TRI(sm,T_NTH_NBR(SM_NTH_TRI(sm,tri_id),2))))
1233     eputs("Invalid tri nbrs smPointLocateTri()\n");
1234     #endif
1235     return(tri_id);
1236 gwlarson 3.8 }
1237 gwlarson 3.1
1238 gwlarson 3.8
1239     /*
1240     Determine whether this sample should:
1241     a) be added to the mesh by subdividing the triangle
1242     b) ignored
1243     c) replace values of an existing mesh vertex
1244    
1245     In case a, the routine will return TRUE, for b,c will return FALSE
1246     In case a, also determine if this sample should cause the deletion of
1247     an existing vertex. If so *rptr will contain the id of the sample to delete
1248     after the new sample has been added.
1249    
1250     ASSUMPTION: this will not be called with a new sample that is
1251     a BASE point.
1252    
1253     The following tests are performed (in order) to determine the fate
1254     of the sample:
1255    
1256 gwlarson 3.10 1) If the new sample is close in ws, and close in the spherical projection
1257 gwlarson 3.8 to one of the triangle vertex samples:
1258     pick the point with dir closest to that of the canonical view
1259     If it is the new point: mark existing point for deletion,and return
1260     TRUE,else return FALSE
1261 gwlarson 3.10 2) If the spherical projection of new is < REPLACE_EPS from a base point:
1262 gwlarson 3.8 add new and mark the base for deletion: return TRUE
1263 gwlarson 3.10 3) If the addition of the new sample point would introduce a "puncture"
1264 gwlarson 3.8 or cause new triangles with large depth differences:return FALSE
1265     This attempts to throw out points that should be occluded
1266     */
1267     int
1268 gwlarson 3.10 smTest_sample(sm,tri_id,dir,p,rptr)
1269 gwlarson 3.8 SM *sm;
1270 gwlarson 3.9 int tri_id;
1271     FVECT dir,p;
1272 gwlarson 3.8 int *rptr;
1273     {
1274     TRI *tri;
1275     double d,d2,dnear,dspt,d01,d12,d20,s0,s1,s2,ds,dv;
1276     int vid[3],i,nearid,norm,dcnt,bcnt;
1277     FVECT diff[3],spt,npt;
1278    
1279     *rptr = INVALID;
1280     bcnt = dcnt = 0;
1281     tri = SM_NTH_TRI(sm,tri_id);
1282     vid[0] = T_NTH_V(tri,0);
1283     vid[1] = T_NTH_V(tri,1);
1284     vid[2] = T_NTH_V(tri,2);
1285    
1286     for(i=0; i<3; i++)
1287     {
1288     if(SM_BASE_ID(sm,vid[i]))
1289     {
1290     bcnt++;
1291     continue;
1292     }
1293     if(SM_DIR_ID(sm,vid[i]))
1294     dcnt++;
1295     VSUB(diff[i],SM_NTH_WV(sm,vid[i]),p);
1296     /* If same world point: replace */
1297     }
1298 gwlarson 3.10 /* TEST 1: If the new sample is close in ws, and close in the spherical
1299 gwlarson 3.8 projection to one of the triangle vertex samples
1300     */
1301     norm = FALSE;
1302     if(bcnt + dcnt != 3)
1303     {
1304     VSUB(spt,p,SM_VIEW_CENTER(sm));
1305     ds = DOT(spt,spt);
1306     dnear = FHUGE;
1307     for(i=0; i<3; i++)
1308     {
1309     if(SM_BASE_ID(sm,vid[i]) || SM_DIR_ID(sm,vid[i]))
1310     continue;
1311     d = DOT(diff[i],diff[i])/ds;
1312     if(d < dnear)
1313     {
1314     dnear = d;
1315     nearid=vid[i];
1316     }
1317     }
1318    
1319     if(dnear <= smMinSampDiff*smMinSampDiff)
1320     {
1321     /* Pick the point with dir closest to that of the canonical view
1322     if it is the new sample: mark existing point for deletion
1323     */
1324     normalize(spt);
1325     norm = TRUE;
1326     VSUB(npt,SM_NTH_WV(sm,nearid),SM_VIEW_CENTER(sm));
1327     normalize(npt);
1328     d = fdir2diff(SM_NTH_W_DIR(sm,nearid), npt);
1329 gwlarson 3.9 d2 = 2. - 2.*DOT(dir,spt);
1330 gwlarson 3.8 /* The existing sample is a better sample:punt */
1331     if(d2 > d)
1332     return(FALSE);
1333     else
1334     {
1335     /* The new sample is better: mark the existing one
1336     for deletion after the new one is added*/
1337     *rptr = nearid;
1338     return(TRUE);
1339     }
1340     }
1341     }
1342     /* TEST 3: If the spherical projection of new is < S_REPLACE_EPS
1343 gwlarson 3.10 from a base point: Edge length is constrained to subtend <45 degrees:
1344     original base mesh edges are approx 32 degrees- so have about 13 degrees
1345     to work in: S_REPLACE_EPS is the square of the radian value
1346 gwlarson 3.8 */
1347     if(bcnt)
1348     {
1349     dnear = FHUGE;
1350     if(bcnt + dcnt ==3)
1351     VSUB(spt,p,SM_VIEW_CENTER(sm));
1352     if(!norm)
1353     normalize(spt);
1354    
1355     for(i=0; i<3; i++)
1356     {
1357     if(!SM_BASE_ID(sm,vid[i]))
1358     continue;
1359     VSUB(npt,SM_NTH_WV(sm,vid[i]),SM_VIEW_CENTER(sm));
1360     d = DIST_SQ(npt,spt);
1361     if(d < S_REPLACE_EPS && d < dnear)
1362     {
1363     dnear = d;
1364     nearid = vid[i];
1365     }
1366     }
1367     if(dnear != FHUGE)
1368     {
1369     /* add new and mark the base for deletion */
1370     *rptr = nearid;
1371     return(TRUE);
1372     }
1373     }
1374    
1375     /* TEST 4:
1376     If the addition of the new sample point would introduce a "puncture"
1377     or cause new triangles with large depth differences:dont add:
1378     */
1379     if(bcnt || dcnt)
1380     return(TRUE);
1381     /* If the new point is in front of the existing points- add */
1382     dv = DIST_SQ(SM_NTH_WV(sm,vid[0]),SM_VIEW_CENTER(sm));
1383     if(ds < dv)
1384     return(TRUE);
1385    
1386     d01 = DIST_SQ(SM_NTH_WV(sm,vid[1]),SM_NTH_WV(sm,vid[0]));
1387     s0 = DOT(diff[0],diff[0]);
1388     if(s0 < S_REPLACE_SCALE*d01)
1389     return(TRUE);
1390     d12 = DIST_SQ(SM_NTH_WV(sm,vid[2]),SM_NTH_WV(sm,vid[1]));
1391     if(s0 < S_REPLACE_SCALE*d12)
1392     return(TRUE);
1393     d20 = DIST_SQ(SM_NTH_WV(sm,vid[0]),SM_NTH_WV(sm,vid[2]));
1394     if(s0 < S_REPLACE_SCALE*d20)
1395     return(TRUE);
1396     d = MIN3(d01,d12,d20);
1397     s1 = DOT(diff[1],diff[1]);
1398     if(s1 < S_REPLACE_SCALE*d)
1399     return(TRUE);
1400     s2 = DOT(diff[2],diff[2]);
1401     if(s2 < S_REPLACE_SCALE*d)
1402     return(TRUE);
1403    
1404 gwlarson 3.9 /* TEST 5:
1405     Check to see if triangle is relatively large and should therefore
1406     be subdivided anyway.
1407     */
1408     dv *= DIST_SQ(SM_NTH_WV(sm,vid[1]),SM_VIEW_CENTER(sm));
1409     dv *= DIST_SQ(SM_NTH_WV(sm,vid[2]),SM_VIEW_CENTER(sm));
1410     if (d01*d12*d20/dv > S_REPLACE_TRI)
1411     return(TRUE);
1412 gwlarson 3.10
1413 gwlarson 3.8 return(FALSE);
1414     }
1415    
1416    
1417     int
1418 gwlarson 3.10 smAlloc_samp(sm)
1419 gwlarson 3.8 SM *sm;
1420     {
1421     int s_id,replaced,cnt;
1422     SAMP *s;
1423     FVECT p;
1424    
1425     s = SM_SAMP(sm);
1426     s_id = sAlloc_samp(s,&replaced);
1427    
1428     cnt=0;
1429     while(replaced)
1430 gwlarson 3.1 {
1431 gwlarson 3.9 if(smRemoveVertex(sm,s_id))
1432 gwlarson 3.8 break;
1433     s_id = sAlloc_samp(s,&replaced);
1434     cnt++;
1435     if(cnt > S_MAX_SAMP(s))
1436     error(CONSISTENCY,"smAlloc_samp():unable to find free samp\n");
1437     }
1438     return(s_id);
1439 gwlarson 3.1 }
1440    
1441     int
1442 gwlarson 3.10 smReplace_samp(sm,c,dir,p,s_id,t_id,o_id,on,which)
1443     SM *sm;
1444     COLR c;
1445     FVECT dir,p;
1446     int s_id,t_id,o_id,on,which;
1447     {
1448     int tonemap,v_id;
1449     TRI *t,*tri;
1450    
1451     tri = SM_NTH_TRI(sm,t_id);
1452     v_id = T_NTH_V(tri,which);
1453    
1454     /* If it is a base id, need new sample */
1455     if(SM_BASE_ID(sm,v_id))
1456     {
1457     tonemap = (SM_TONE_MAP(sm) > s_id);
1458     sInit_samp(SM_SAMP(sm),s_id,c,dir,p,o_id,tonemap);
1459     SM_NTH_VERT(sm,s_id) = t_id;
1460     T_NTH_V(tri,which) = s_id;
1461     if(!(SM_BASE_ID(sm,T_NTH_V(tri,(which+1)%3)) ||
1462     SM_BASE_ID(sm,T_NTH_V(tri,(which+2)%3))))
1463     SM_CLR_NTH_T_BASE(sm,t_id);
1464     t_id = smTri_next_ccw_nbr(sm,tri,v_id);
1465     while(t_id != INVALID)
1466     {
1467     t = SM_NTH_TRI(sm,t_id);
1468     which = T_WHICH_V(t,v_id);
1469     T_NTH_V(t,which) = s_id;
1470     /* Check if still a base triangle */
1471     if(!(SM_BASE_ID(sm,T_NTH_V(t,(which+1)%3)) ||
1472     SM_BASE_ID(sm,T_NTH_V(t,(which+2)%3))))
1473     SM_CLR_NTH_T_BASE(sm,t_id);
1474     t_id = smTri_next_ccw_nbr(sm,t,v_id);
1475     }
1476     return(s_id);
1477     }
1478     else
1479     if(on == ON_V || !p)
1480     {
1481     tonemap = (SM_TONE_MAP(sm) > v_id);
1482     sInit_samp(SM_SAMP(sm),v_id,c,dir,p,o_id,tonemap);
1483     return(v_id);
1484     }
1485     else /* on == ON_P */
1486     {
1487     FVECT spt,npt;
1488     double d,d2;
1489    
1490     /* Need to check which sample is preferable */
1491     VSUB(spt,p,SM_VIEW_CENTER(sm));
1492     normalize(spt);
1493    
1494     VSUB(npt,SM_NTH_WV(sm,v_id),SM_VIEW_CENTER(sm));
1495     normalize(npt);
1496     d = fdir2diff(SM_NTH_W_DIR(sm,v_id), npt);
1497     d2 = 2. - 2.*DOT(dir,spt);
1498     /* The existing sample is a better sample:punt */
1499     if(d2 < d)
1500     {
1501     /* The new sample has better information- replace values */
1502     tonemap = (SM_TONE_MAP(sm) > v_id);
1503     sInit_samp(SM_SAMP(sm),v_id,c,dir,p,o_id,tonemap);
1504     }
1505     return(v_id);
1506     }
1507     }
1508    
1509     /* Add sample to the mesh:
1510    
1511     the sample can be a world space or directional point. If o_id !=INVALID,
1512     then we have an existing sample containing the appropriate information
1513     already converted into storage format.
1514     The spherical projection of the point/ray can:
1515     a) coincide with existing vertex
1516     i) conincides with existing sample
1517     ii) projects same as existing sample
1518     iii) hits a base vertex
1519     b) coincide with existing edge
1520     c) Fall in existing triangle
1521     */
1522     int
1523 gwlarson 3.9 smAdd_samp(sm,c,dir,p,o_id)
1524 gwlarson 3.1 SM *sm;
1525 gwlarson 3.9 COLR c;
1526     FVECT dir,p;
1527     int o_id;
1528 gwlarson 3.1 {
1529 gwlarson 3.10 int t_id,s_id,r_id,on,which,n_id;
1530 gwlarson 3.9 double d;
1531     FVECT wpt;
1532 gwlarson 3.8
1533 gwlarson 3.9 r_id = INVALID;
1534 gwlarson 3.10
1535     /* Must do this first-as may change mesh */
1536     s_id = smAlloc_samp(sm);
1537 gwlarson 3.9 /* If sample is a world space point */
1538     if(p)
1539     {
1540 gwlarson 3.10 t_id = smPointLocateTri(sm,p,&on,&which);
1541 gwlarson 3.8 if(t_id == INVALID)
1542 gwlarson 3.9 {
1543     #ifdef DEBUG
1544     eputs("smAddSamp(): unable to locate tri containing sample \n");
1545     #endif
1546 gwlarson 3.10 smUnalloc_samp(sm,s_id);
1547 gwlarson 3.9 return(INVALID);
1548     }
1549 gwlarson 3.10 /* If spherical projection coincides with existing sample: replace */
1550     if((on == ON_V || on == ON_P))
1551 gwlarson 3.1 {
1552 gwlarson 3.10 if((n_id = smReplace_samp(sm,c,dir,p,s_id,t_id,o_id,on,which))!= s_id)
1553     smUnalloc_samp(sm,s_id);
1554     return(n_id);
1555 gwlarson 3.9 }
1556 gwlarson 3.10 if((!(smTest_sample(sm,t_id,dir,p,&r_id))))
1557     {
1558     smUnalloc_samp(sm,s_id);
1559     return(INVALID);
1560     }
1561     /* If sample is being added in the middle of the sample array: tone
1562     map individually
1563     */
1564     /* Initialize sample */
1565     sInit_samp(SM_SAMP(sm),s_id,c,dir,p,o_id,(SM_TONE_MAP(sm)>s_id));
1566    
1567 gwlarson 3.9 }
1568     /* If sample is a direction vector */
1569     else
1570     {
1571 gwlarson 3.10 VADD(wpt,dir,SM_VIEW_CENTER(sm));
1572     t_id = smPointLocateTri(sm,wpt,&on,&which);
1573 gwlarson 3.9 if(t_id == INVALID)
1574     {
1575 gwlarson 3.1 #ifdef DEBUG
1576 gwlarson 3.9 eputs("smAddSamp(): unable to locate tri containing sample \n");
1577 gwlarson 3.1 #endif
1578 gwlarson 3.10 smUnalloc_samp(sm,s_id);
1579 gwlarson 3.9 return(INVALID);
1580     }
1581 gwlarson 3.10 if(on == ON_V || on == ON_P)
1582     {
1583     if((n_id = smReplace_samp(sm,c,wpt,NULL,s_id,t_id,o_id,on,which))!= s_id)
1584     smUnalloc_samp(sm,s_id);
1585     return(n_id);
1586     }
1587 gwlarson 3.9 /* Allocate space for a sample and initialize */
1588 gwlarson 3.10 sInit_samp(SM_SAMP(sm),s_id,c,wpt,NULL,o_id,(SM_TONE_MAP(sm)>s_id));
1589 gwlarson 3.1 }
1590 gwlarson 3.10 if(!SM_DIR_ID(sm,s_id))
1591 gwlarson 3.1 {
1592 gwlarson 3.9 /* If not a base or sky point, add distance from the
1593     viewcenter to average*/
1594 gwlarson 3.10 d = DIST(SM_NTH_WV(sm,s_id),SM_VIEW_CENTER(sm));
1595 gwlarson 3.9 smDist_sum += 1.0/d;
1596 gwlarson 3.1 }
1597 gwlarson 3.10 smInsert_samp(sm,s_id,t_id,on,which);
1598 gwlarson 3.5
1599 gwlarson 3.9 /* If new sample replaces existing one- remove that vertex now */
1600     if(r_id != INVALID)
1601     {
1602     smRemoveVertex(sm,r_id);
1603     sDelete_samp(SM_SAMP(sm),r_id);
1604     }
1605     return(s_id);
1606 gwlarson 3.1 }
1607    
1608     /*
1609     * int
1610     * smNewSamp(c, dir, p) : register new sample point and return index
1611     * COLR c; : pixel color (RGBE)
1612     * FVECT dir; : ray direction vector
1613     * FVECT p; : world intersection point
1614     *
1615     * Add new sample point to data structures, removing old values as necessary.
1616     * New sample representation will be output in next call to smUpdate().
1617     * If the point is a sky point: then v=NULL
1618     */
1619     int
1620     smNewSamp(c,dir,p)
1621     COLR c;
1622     FVECT dir;
1623     FVECT p;
1624     {
1625     int s_id;
1626 gwlarson 3.9
1627 gwlarson 3.10
1628 gwlarson 3.1 /* First check if this the first sample: if so initialize mesh */
1629 gwlarson 3.10
1630 gwlarson 3.1 if(SM_NUM_SAMP(smMesh) == 0)
1631 gwlarson 3.9 {
1632     smInit_sm(smMesh,odev.v.vp);
1633     sClear_all_flags(SM_SAMP(smMesh));
1634     }
1635     /* Add the sample to the mesh */
1636     s_id = smAdd_samp(smMesh,c,dir,p,INVALID);
1637 gwlarson 3.8
1638 gwlarson 3.1 return(s_id);
1639    
1640     }
1641     int
1642 gwlarson 3.8 smAdd_base_vertex(sm,v)
1643 gwlarson 3.1 SM *sm;
1644 gwlarson 3.8 FVECT v;
1645 gwlarson 3.1 {
1646     int id;
1647    
1648     /* First add coordinate to the sample array */
1649 gwlarson 3.8 id = sAdd_base_point(SM_SAMP(sm),v);
1650     if(id == INVALID)
1651     return(INVALID);
1652 gwlarson 3.1 /* Initialize triangle pointer to -1 */
1653     smClear_vert(sm,id);
1654     return(id);
1655     }
1656    
1657    
1658    
1659     /* Initialize a the point location DAG based on a 6 triangle tesselation
1660     of the unit sphere centered on the view center. The DAG structure
1661     contains 6 roots: one for each initial base triangle
1662     */
1663    
1664     int
1665     smCreate_base_mesh(sm,type)
1666     SM *sm;
1667     int type;
1668     {
1669 gwlarson 3.8 int i,s_id,tri_id,nbr_id;
1670 gwlarson 3.10 int p[SM_EXTRA_POINTS],ids[SM_BASE_TRIS];
1671 gwlarson 3.1 int v0_id,v1_id,v2_id;
1672 gwlarson 3.8 FVECT d,pt,cntr,v0,v1,v2;
1673 gwlarson 3.10 TRI *t;
1674    
1675 gwlarson 3.1 /* First insert the base vertices into the sample point array */
1676    
1677 gwlarson 3.10 for(i=0; i < SM_EXTRA_POINTS; i++)
1678 gwlarson 3.1 {
1679 gwlarson 3.10 VADD(cntr,icosa_verts[i],SM_VIEW_CENTER(sm));
1680 gwlarson 3.8 p[i] = smAdd_base_vertex(sm,cntr);
1681 gwlarson 3.1 }
1682     /* Create the base triangles */
1683 gwlarson 3.10 for(i=0;i < SM_BASE_TRIS; i++)
1684 gwlarson 3.1 {
1685 gwlarson 3.10 v0_id = p[icosa_tris[i][0]];
1686     v1_id = p[icosa_tris[i][1]];
1687     v2_id = p[icosa_tris[i][2]];
1688 gwlarson 3.8 ids[i] = smAdd_tri(sm, v0_id,v1_id,v2_id);
1689 gwlarson 3.10 /* Set neighbors */
1690     for(nbr_id=0; nbr_id < 3; nbr_id++)
1691     T_NTH_NBR(SM_NTH_TRI(sm,ids[i]),nbr_id) = icosa_nbrs[i][nbr_id];
1692    
1693 gwlarson 3.1 }
1694 gwlarson 3.10 for(i=0;i < SM_BASE_TRIS; i++)
1695 gwlarson 3.8 {
1696 gwlarson 3.10 t = SM_NTH_TRI(sm,ids[i]);
1697     smLocator_add_tri(sm,ids[i],T_NTH_V(t,0),T_NTH_V(t,1),T_NTH_V(t,2));
1698 gwlarson 3.8 }
1699     return(1);
1700 gwlarson 3.1
1701     }
1702    
1703    
1704     int
1705     smNext_tri_flag_set(sm,i,which,b)
1706     SM *sm;
1707     int i,which;
1708 gwlarson 3.5 int b;
1709 gwlarson 3.1 {
1710    
1711 gwlarson 3.8 for(; i < SM_NUM_TRI(sm);i++)
1712 gwlarson 3.1 {
1713    
1714 gwlarson 3.8 if(!T_IS_VALID(SM_NTH_TRI(sm,i)))
1715 gwlarson 3.5 continue;
1716 gwlarson 3.8 if(!SM_IS_NTH_T_FLAG(sm,i,which))
1717     continue;
1718 gwlarson 3.1 if(!b)
1719     break;
1720     if((b==1) && !SM_BG_TRI(sm,i))
1721     break;
1722     if((b==2) && SM_BG_TRI(sm,i))
1723     break;
1724     }
1725    
1726     return(i);
1727     }
1728    
1729    
1730     int
1731     smNext_valid_tri(sm,i)
1732     SM *sm;
1733     int i;
1734     {
1735    
1736 gwlarson 3.8 while( i < SM_NUM_TRI(sm) && !T_IS_VALID(SM_NTH_TRI(sm,i)))
1737 gwlarson 3.1 i++;
1738    
1739     return(i);
1740     }
1741    
1742    
1743    
1744 gwlarson 3.8 qtTri_from_id(t_id,v0,v1,v2)
1745 gwlarson 3.1 int t_id;
1746     FVECT v0,v1,v2;
1747     {
1748     TRI *t;
1749    
1750     t = SM_NTH_TRI(smMesh,t_id);
1751 gwlarson 3.8 if(!T_IS_VALID(t))
1752     return(0);
1753     VSUB(v0,SM_T_NTH_WV(smMesh,t,0),SM_VIEW_CENTER(smMesh));
1754     VSUB(v1,SM_T_NTH_WV(smMesh,t,1),SM_VIEW_CENTER(smMesh));
1755     VSUB(v2,SM_T_NTH_WV(smMesh,t,2),SM_VIEW_CENTER(smMesh));
1756     return(1);
1757 gwlarson 3.1 }
1758    
1759    
1760 gwlarson 3.8 smRebuild_mesh(sm,v)
1761 gwlarson 3.1 SM *sm;
1762 gwlarson 3.8 VIEW *v;
1763 gwlarson 3.1 {
1764 gwlarson 3.8 int i,j,cnt;
1765     FVECT p,ov,dir;
1766     double d;
1767 gwlarson 3.1
1768 gwlarson 3.8 #ifdef DEBUG
1769 gwlarson 3.10 fprintf(stderr,"smRebuild_mesh(): rebuilding....");
1770 gwlarson 3.8 #endif
1771 gwlarson 3.1 /* Clear the mesh- and rebuild using the current sample array */
1772 gwlarson 3.8 /* Remember the current number of samples */
1773     cnt = SM_NUM_SAMP(sm);
1774     /* Calculate the difference between the current and new viewpoint*/
1775     /* Will need to subtract this off of sky points */
1776 gwlarson 3.9 VCOPY(ov,SM_VIEW_CENTER(sm));
1777 gwlarson 3.8 /* Initialize the mesh to 0 samples and the base triangles */
1778 gwlarson 3.1
1779 gwlarson 3.8 /* Go through all samples and add in if in the new view frustum, and
1780     the dir is <= 30 degrees from new view
1781     */
1782 gwlarson 3.9 smInit_sm(sm,v->vp);
1783 gwlarson 3.8 for(i=0; i < cnt; i++)
1784 gwlarson 3.1 {
1785 gwlarson 3.8 /* First check if sample visible(conservative approx: if one of tris
1786     attached to sample is in frustum) */
1787     if(!S_IS_FLAG(i))
1788     continue;
1789    
1790     /* Next test if direction > 30 from current direction */
1791     if(SM_NTH_W_DIR(sm,i)!=-1)
1792     {
1793     VSUB(p,SM_NTH_WV(sm,i),v->vp);
1794     normalize(p);
1795 gwlarson 3.9 decodedir(dir,SM_NTH_W_DIR(sm,i));
1796     d = 2. - 2.*DOT(dir,p);
1797 gwlarson 3.8 if (d > MAXDIFF2)
1798     continue;
1799 gwlarson 3.9 VCOPY(p,SM_NTH_WV(sm,i));
1800     smAdd_samp(sm,NULL,dir,p,i);
1801 gwlarson 3.8 }
1802 gwlarson 3.9 else
1803     {
1804     VSUB(dir,SM_NTH_WV(sm,i),ov);
1805 gwlarson 3.10 VADD(SM_NTH_WV(sm,i),dir,SM_VIEW_CENTER(sm));
1806 gwlarson 3.9 smAdd_samp(sm,NULL,dir,NULL,i);
1807     }
1808    
1809 gwlarson 3.1 }
1810 gwlarson 3.10
1811 gwlarson 3.9 smNew_tri_cnt = SM_SAMPLE_TRIS(sm);
1812 gwlarson 3.8 #ifdef DEBUG
1813 gwlarson 3.10 fprintf(stderr,"smRebuild_mesh():done\n");
1814 gwlarson 3.8 #endif
1815 gwlarson 3.1 }
1816 gwlarson 3.5
1817     int
1818     intersect_tri_set(t_set,orig,dir,pt)
1819     OBJECT *t_set;
1820     FVECT orig,dir,pt;
1821     {
1822     OBJECT *optr;
1823 gwlarson 3.10 int i,t_id,id,base;
1824 gwlarson 3.5 int pid0,pid1,pid2;
1825     FVECT p0,p1,p2,p;
1826     TRI *t;
1827 gwlarson 3.10 double d,d1;
1828 gwlarson 3.5
1829     optr = QT_SET_PTR(t_set);
1830     for(i = QT_SET_CNT(t_set); i > 0; i--)
1831     {
1832     t_id = QT_SET_NEXT_ELEM(optr);
1833     t = SM_NTH_TRI(smMesh,t_id);
1834 gwlarson 3.8 if(!T_IS_VALID(t))
1835     continue;
1836 gwlarson 3.5 pid0 = T_NTH_V(t,0);
1837     pid1 = T_NTH_V(t,1);
1838     pid2 = T_NTH_V(t,2);
1839 gwlarson 3.10 if(base = SM_IS_NTH_T_BASE(smMesh,t_id))
1840     if(SM_BASE_ID(smMesh,pid0) && SM_BASE_ID(smMesh,pid1) &&
1841     SM_BASE_ID(smMesh,pid2))
1842     continue;
1843 gwlarson 3.5 VCOPY(p0,SM_NTH_WV(smMesh,pid0));
1844     VCOPY(p1,SM_NTH_WV(smMesh,pid1));
1845     VCOPY(p2,SM_NTH_WV(smMesh,pid2));
1846     if(ray_intersect_tri(orig,dir,p0,p1,p2,p))
1847     {
1848 gwlarson 3.10 if(!base)
1849     {
1850     d = DIST_SQ(p,p0);
1851     d1 = DIST_SQ(p,p1);
1852     if(d < d1)
1853     {
1854     d1 = DIST_SQ(p,p2);
1855     id = (d1 < d)?pid2:pid0;
1856     }
1857     else
1858     {
1859     d = DIST_SQ(p,p2);
1860     id = (d < d1)? pid2:pid1;
1861     }
1862     }
1863     else
1864     {
1865     if(SM_BASE_ID(smMesh,pid0))
1866     {
1867     if(SM_BASE_ID(smMesh,pid1))
1868     id = pid2;
1869     else
1870     if(SM_BASE_ID(smMesh,pid2))
1871     id = pid1;
1872     else
1873     {
1874     d = DIST_SQ(p,p1);
1875     d1 = DIST_SQ(p,p2);
1876     if(d < d1)
1877     id = pid1;
1878     else
1879     id = pid2;
1880     }
1881     }
1882     else
1883     {
1884     if(SM_BASE_ID(smMesh,pid1))
1885     {
1886     if(SM_BASE_ID(smMesh,pid2))
1887     id = pid0;
1888     else
1889     {
1890     d = DIST_SQ(p,p0);
1891     d1 = DIST_SQ(p,p2);
1892     if(d < d1)
1893     id = pid0;
1894     else
1895     id = pid2;
1896     }
1897     }
1898     else
1899     {
1900     d = DIST_SQ(p,p0);
1901     d1 = DIST_SQ(p,p1);
1902     if(d < d1)
1903     id = pid0;
1904     else
1905     id = pid1;
1906     }
1907     }
1908     }
1909 gwlarson 3.5 if(pt)
1910 gwlarson 3.10 VCOPY(pt,p);
1911     #ifdef TEST_DRIVER
1912 gwlarson 3.5 Pick_tri = t_id;
1913     Pick_samp = id;
1914     VCOPY(Pick_point[0],p);
1915     #endif
1916     return(id);
1917     }
1918     }
1919     return(-1);
1920     }
1921    
1922 gwlarson 3.8 /* OS is constrained to be <= QT_MAXCSET : if the set exceeds this, the
1923     results of check_set are conservative
1924     */
1925 gwlarson 3.10 int
1926     compare_ids(id1,id2)
1927     OBJECT *id1,*id2;
1928     {
1929     int d;
1930 gwlarson 3.8
1931 gwlarson 3.10 d = *id2-*id1;
1932    
1933     if(d > 0)
1934     return(-1);
1935     if(d < 0)
1936     return(1);
1937    
1938     return(0);
1939     }
1940    
1941     ray_trace_check_set(qt,argptr,fptr)
1942     QUADTREE qt;
1943     RT_ARGS *argptr;
1944 gwlarson 3.8 int *fptr;
1945 gwlarson 3.5 {
1946 gwlarson 3.8 OBJECT tset[QT_MAXSET+1],*tptr;
1947 gwlarson 3.5 double dt,t;
1948     int found;
1949 gwlarson 3.10 if(QT_IS_EMPTY(qt))
1950     return;
1951     if(QT_LEAF_IS_FLAG(qt))
1952     {
1953     QT_FLAG_SET_DONE(*fptr);
1954     #if DEBUG
1955     eputs("ray_trace_check_set():Already visited this node:aborting\n");
1956     #endif
1957     return;
1958     }
1959     else
1960     QT_LEAF_SET_FLAG(qt);
1961    
1962     tptr = qtqueryset(qt);
1963     if(QT_SET_CNT(tptr) > QT_MAXSET)
1964     tptr = (OBJECT *)malloc((QT_SET_CNT(tptr)+1)*sizeof(OBJECT));
1965     else
1966     tptr = tset;
1967     if(!tptr)
1968     goto memerr;
1969 gwlarson 3.5
1970 gwlarson 3.10 qtgetset(tptr,qt);
1971     /* Must sort */
1972     qsort((void *)(&(tptr[1])),tptr[0],sizeof(OBJECT),compare_ids);
1973     /* Check triangles in set against those seen so far(os):only
1974     check new triangles for intersection (t_set')
1975     */
1976     check_set_large(tptr,argptr->os);
1977    
1978     if(!QT_SET_CNT(tptr))
1979 gwlarson 3.8 return;
1980 gwlarson 3.10 found = intersect_tri_set(tptr,argptr->orig,argptr->dir,NULL);
1981     if(tptr != tset)
1982     free(tptr);
1983     if(found != INVALID)
1984     {
1985     argptr->t_id = found;
1986     QT_FLAG_SET_DONE(*fptr);
1987     }
1988     return;
1989 gwlarson 3.8 memerr:
1990     error(SYSTEM,"ray_trace_check_set():Unable to allocate memory");
1991 gwlarson 3.5 }
1992    
1993 gwlarson 3.8
1994     /*
1995     * int
1996     * smFindSamp(FVECT orig, FVECT dir)
1997     *
1998     * Find the closest sample to the given ray. Returns sample id, -1 on failure.
1999     * "dir" is assumed to be normalized
2000     */
2001    
2002 gwlarson 3.5 int
2003 gwlarson 3.6 smFindSamp(orig,dir)
2004 gwlarson 3.5 FVECT orig,dir;
2005     {
2006     FVECT b,p,o;
2007     OBJECT *ts;
2008     QUADTREE qt;
2009 gwlarson 3.10 int s_id,test;
2010 gwlarson 3.5 double d;
2011    
2012     /* r is the normalized vector from the view center to the current
2013     * ray point ( starting with "orig"). Find the cell that r falls in,
2014     * and test the ray against all triangles stored in the cell. If
2015     * the test fails, trace the projection of the ray across to the
2016     * next cell it intersects: iterate until either an intersection
2017     * is found, or the projection ray is // to the direction. The sample
2018     * corresponding to the triangle vertex closest to the intersection
2019     * point is returned.
2020     */
2021    
2022     /* First test if "orig" coincides with the View_center or if "dir" is
2023     parallel to r formed by projecting "orig" on the sphere. In
2024     either case, do a single test against the cell containing the
2025     intersection of "dir" and the sphere
2026     */
2027     /* orig will be updated-so preserve original value */
2028     if(!smMesh)
2029     return;
2030 gwlarson 3.9 #ifdef TEST_DRIVER
2031     Picking= TRUE;
2032     #endif
2033 gwlarson 3.10
2034     if(EQUAL_VEC3(orig,SM_VIEW_CENTER(smMesh)))
2035 gwlarson 3.5 {
2036 gwlarson 3.10 qt = smPointLocateCell(smMesh,dir);
2037     if(QT_IS_EMPTY(qt))
2038     goto Lerror;
2039     ts = qtqueryset(qt);
2040     s_id = intersect_tri_set(ts,orig,dir,p);
2041     #ifdef DEBUG_TEST_DRIVER
2042     VCOPY(Pick_point[0],p);
2043     #endif
2044 gwlarson 3.8 #ifdef DEBUG
2045 gwlarson 3.10 fprintf(stderr,"Simple pick returning %d\n",s_id);
2046 gwlarson 3.8 #endif
2047 gwlarson 3.10
2048     return(s_id);
2049     }
2050     d = point_on_sphere(b,orig,SM_VIEW_CENTER(smMesh));
2051     if(EQUAL_VEC3(b,dir))
2052     {
2053     qt = smPointLocateCell(smMesh,dir);
2054     if(QT_IS_EMPTY(qt))
2055     goto Lerror;
2056     ts = qtqueryset(qt);
2057     s_id = intersect_tri_set(ts,orig,dir,p);
2058 gwlarson 3.5 #ifdef DEBUG_TEST_DRIVER
2059 gwlarson 3.10 VCOPY(Pick_point[0],p);
2060 gwlarson 3.5 #endif
2061 gwlarson 3.10 #ifdef DEBUG
2062     fprintf(stderr,"Simple pick2 returning %d\n",s_id);
2063     #endif
2064     return(s_id);
2065 gwlarson 3.5 }
2066 gwlarson 3.10 if(OPP_EQUAL_VEC3(b,dir))
2067 gwlarson 3.5 {
2068 gwlarson 3.10 qt = smPointLocateCell(smMesh,orig);
2069     if(QT_IS_EMPTY(qt))
2070     goto Lerror;
2071     ts = qtqueryset(qt);
2072     s_id = intersect_tri_set(ts,orig,dir,p);
2073     #ifdef DEBUG_TEST_DRIVER
2074     VCOPY(Pick_point[0],p);
2075     #endif
2076     if(s_id != INVALID)
2077     {
2078     #ifdef DEBUG
2079     fprintf(stderr,"Front pick returning %d\n",s_id);
2080     #endif
2081     return(s_id);
2082     }
2083     qt = smPointLocateCell(smMesh,dir);
2084     if(QT_IS_EMPTY(qt))
2085     goto Lerror;
2086     ts = qtqueryset(qt);
2087     s_id = intersect_tri_set(ts,orig,dir,p);
2088     #ifdef DEBUG_TEST_DRIVER
2089     VCOPY(Pick_point[0],p);
2090     #endif
2091     #ifdef DEBUG
2092     fprintf(stderr,"Back pick returning %d\n",s_id);
2093     #endif
2094     return(s_id);
2095     }
2096     {
2097 gwlarson 3.8 OBJECT t_set[QT_MAXCSET + 1];
2098     RT_ARGS rt;
2099 gwlarson 3.10 FUNC func;
2100 gwlarson 3.8
2101 gwlarson 3.5 /* Test each of the root triangles against point id */
2102     QT_CLEAR_SET(t_set);
2103     VSUB(o,orig,SM_VIEW_CENTER(smMesh));
2104     ST_CLEAR_FLAGS(SM_LOCATOR(smMesh));
2105 gwlarson 3.8 rt.t_id = -1;
2106     rt.os = t_set;
2107     VCOPY(rt.orig,orig);
2108     VCOPY(rt.dir,dir);
2109 gwlarson 3.10 F_FUNC(func) = ray_trace_check_set;
2110     F_ARGS(func) = (int *)(&rt);
2111     stTrace_ray(SM_LOCATOR(smMesh),o,dir,func);
2112 gwlarson 3.8 s_id = rt.t_id;
2113 gwlarson 3.10 }
2114     #ifdef DEBUG
2115     fprintf(stderr,"Trace pick returning %d\n",s_id);
2116     #endif
2117 gwlarson 3.5 return(s_id);
2118 gwlarson 3.10
2119     Lerror:
2120     #ifdef DEBUG
2121     eputs("smFindSamp(): point not found");
2122     #endif
2123     return(INVALID);
2124    
2125 gwlarson 3.5 }
2126    
2127 gwlarson 3.10
2128     mark_active_tris(argptr,root,qt)
2129     int *argptr;
2130     int root;
2131     QUADTREE qt;
2132     {
2133     OBJECT *os,*optr;
2134     register int i,t_id;
2135     TRI *tri;
2136    
2137     if(QT_IS_EMPTY(qt) || QT_LEAF_IS_FLAG(qt))
2138     return;
2139    
2140     /* For each triangle in the set, set the which flag*/
2141     os = qtqueryset(qt);
2142    
2143     for (i = QT_SET_CNT(os), optr = QT_SET_PTR(os); i > 0; i--)
2144     {
2145     t_id = QT_SET_NEXT_ELEM(optr);
2146     /* Set the render flag */
2147     tri = SM_NTH_TRI(smMesh,t_id);
2148     if(!T_IS_VALID(tri) || SM_IS_NTH_T_BASE(smMesh,t_id))
2149     continue;
2150     SM_SET_NTH_T_ACTIVE(smMesh,t_id);
2151     /* Set the Active bits of the Vertices */
2152     S_SET_FLAG(T_NTH_V(tri,0));
2153     S_SET_FLAG(T_NTH_V(tri,1));
2154     S_SET_FLAG(T_NTH_V(tri,2));
2155     }
2156     }
2157    
2158    
2159     mark_tris_in_frustum(view)
2160     VIEW *view;
2161     {
2162     FVECT nr[4],far[4];
2163     FPEQ peq;
2164     int debug=0;
2165     FUNC f;
2166    
2167     F_FUNC(f) = mark_active_tris;
2168     F_ARGS(f) = NULL;
2169    
2170     /* Mark triangles in approx. view frustum as being active:set
2171     LRU counter: for use in discarding samples when out
2172     of space
2173     Radiance often has no far clipping plane: but driver will set
2174     dev_zmin,dev_zmax to satisfy OGL
2175     */
2176    
2177     /* First clear all the quadtree node and triangle active flags */
2178     qtClearAllFlags();
2179     smClear_flags(smMesh,T_ACTIVE_FLAG);
2180     /* Clear all of the active sample flags*/
2181     sClear_all_flags(SM_SAMP(smMesh));
2182    
2183    
2184     /* calculate the world space coordinates of the view frustum */
2185     calculate_view_frustum(view->vp,view->hvec,view->vvec,view->horiz,
2186     view->vert, dev_zmin,dev_zmax,nr,far);
2187    
2188     #ifdef TEST_DRIVER
2189     VCOPY(FrustumFar[0],far[0]);
2190     VCOPY(FrustumFar[1],far[1]);
2191     VCOPY(FrustumFar[2],far[2]);
2192     VCOPY(FrustumFar[3],far[3]);
2193     VCOPY(FrustumNear[0],nr[0]);
2194     VCOPY(FrustumNear[1],nr[1]);
2195     VCOPY(FrustumNear[2],nr[2]);
2196     VCOPY(FrustumNear[3],nr[3]);
2197     #endif
2198     /* Project the view frustum onto the spherical quadtree */
2199     /* For every cell intersected by the projection of the faces
2200    
2201     of the frustum: mark all triangles in the cell as ACTIVE-
2202     Also set the triangles LRU clock counter
2203     */
2204    
2205     if(EQUAL_VEC3(view->vp,SM_VIEW_CENTER(smMesh)))
2206     {/* Near face triangles */
2207    
2208     smLocator_apply(smMesh,nr[0],nr[2],nr[3],f);
2209     smLocator_apply(smMesh,nr[2],nr[0],nr[1],f);
2210     return;
2211     }
2212     /* Test the view against the planes: and swap orientation if inside:*/
2213     tri_plane_equation(nr[0],nr[2],nr[3], &peq,FALSE);
2214     if(PT_ON_PLANE(SM_VIEW_CENTER(smMesh),peq) < 0.0)
2215     {/* Near face triangles */
2216     smLocator_apply(smMesh,nr[3],nr[2],nr[0],f);
2217     smLocator_apply(smMesh,nr[1],nr[0],nr[2],f);
2218     }
2219     else
2220     {/* Near face triangles */
2221     smLocator_apply(smMesh,nr[0],nr[2],nr[3],f);
2222     smLocator_apply(smMesh,nr[2],nr[0],nr[1],f);
2223     }
2224     tri_plane_equation(nr[0],far[3],far[0], &peq,FALSE);
2225     if(PT_ON_PLANE(SM_VIEW_CENTER(smMesh),peq) < 0.0)
2226     { /* Right face triangles */
2227     smLocator_apply(smMesh,far[0],far[3],nr[0],f);
2228     smLocator_apply(smMesh,nr[3],nr[0],far[3],f);
2229     }
2230     else
2231     {/* Right face triangles */
2232     smLocator_apply(smMesh,nr[0],far[3],far[0],f);
2233     smLocator_apply(smMesh,far[3],nr[0],nr[3],f);
2234     }
2235    
2236     tri_plane_equation(nr[1],far[2],nr[2], &peq,FALSE);
2237     if(PT_ON_PLANE(SM_VIEW_CENTER(smMesh),peq) < 0.0)
2238     { /* Left face triangles */
2239     smLocator_apply(smMesh,nr[2],far[2],nr[1],f);
2240     smLocator_apply(smMesh,far[1],nr[1],far[2],f);
2241     }
2242     else
2243     { /* Left face triangles */
2244     smLocator_apply(smMesh,nr[1],far[2],nr[2],f);
2245     smLocator_apply(smMesh,far[2],nr[1],far[1],f);
2246     }
2247     tri_plane_equation(nr[0],far[0],nr[1], &peq,FALSE);
2248     if(PT_ON_PLANE(SM_VIEW_CENTER(smMesh),peq) < 0.0)
2249     {/* Top face triangles */
2250     smLocator_apply(smMesh,nr[1],far[0],nr[0],f);
2251     smLocator_apply(smMesh,far[1],far[0],nr[1],f);
2252     }
2253     else
2254     {/* Top face triangles */
2255     smLocator_apply(smMesh,nr[0],far[0],nr[1],f);
2256     smLocator_apply(smMesh,nr[1],far[0],far[1],f);
2257     }
2258     tri_plane_equation(nr[3],nr[2],far[3], &peq,FALSE);
2259     if(PT_ON_PLANE(SM_VIEW_CENTER(smMesh),peq) < 0.0)
2260     {/* Bottom face triangles */
2261     smLocator_apply(smMesh,far[3],nr[2],nr[3],f);
2262     smLocator_apply(smMesh,far[3],far[2],nr[2],f);
2263     }
2264     else
2265     { /* Bottom face triangles */
2266     smLocator_apply(smMesh,nr[3],nr[2],far[3],f);
2267     smLocator_apply(smMesh,nr[2],far[2],far[3],f);
2268     }
2269     tri_plane_equation(far[2],far[0],far[1], &peq,FALSE);
2270     if(PT_ON_PLANE(SM_VIEW_CENTER(smMesh),peq) < 0.0)
2271     {/* Far face triangles */
2272     smLocator_apply(smMesh,far[0],far[2],far[1],f);
2273     smLocator_apply(smMesh,far[2],far[0],far[3],f);
2274     }
2275     else
2276     {/* Far face triangles */
2277     smLocator_apply(smMesh,far[1],far[2],far[0],f);
2278     smLocator_apply(smMesh,far[3],far[0],far[2],f);
2279     }
2280    
2281     }
2282 gwlarson 3.5
2283    
2284    
2285    
2286    
2287    
2288    
2289    
2290 gwlarson 3.1
2291