ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/sm.c
Revision: 3.16
Committed: Sat Feb 22 02:07:25 2003 UTC (21 years, 7 months ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad3R5
Changes since 3.15: +186 -431 lines
Log Message: 
Changes and check-in for 3.5 release
Includes new source files and modifications not recorded for many years
See ray/doc/notes/ReleaseNotes for notes between 3.1 and 3.5 release

File Contents

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