1 |
|
---|
2 |
|
---|
3 |
|
---|
4 |
|
---|
5 |
|
---|
6 |
|
---|
7 |
|
---|
8 |
|
---|
9 |
|
---|
10 |
|
---|
11 |
|
---|
12 |
|
---|
13 |
|
---|
14 |
|
---|
15 |
|
---|
16 |
|
---|
17 |
|
---|
18 |
|
---|
19 |
|
---|
20 |
|
---|
21 |
|
---|
22 |
|
---|
23 |
|
---|
24 |
|
---|
25 |
|
---|
26 |
|
---|
27 |
#include "CubitMessage.hpp" |
---|
28 |
#include "FBClassify.hpp" |
---|
29 |
#include "IntegerHash.hpp" |
---|
30 |
#include "GfxDebug.hpp" |
---|
31 |
#include <stack> |
---|
32 |
|
---|
33 |
//const double EPSILON_CLASSIFY = 1.e-12; |
---|
34 |
const double EPSILON_CLASSIFY = 1.e-10; |
---|
35 |
FBClassify::FBClassify() |
---|
36 |
{ |
---|
37 |
number_of_groups = 0; |
---|
38 |
polya = polyb = 0; |
---|
39 |
} |
---|
40 |
|
---|
41 |
FBClassify::~FBClassify() |
---|
42 |
{ |
---|
43 |
|
---|
44 |
} |
---|
45 |
|
---|
46 |
void FBClassify::SetPoly(FBPolyhedron *poly1, FBPolyhedron *poly2) |
---|
47 |
{ |
---|
48 |
polya = poly1; |
---|
49 |
polyb = poly2; |
---|
50 |
} |
---|
51 |
|
---|
52 |
CubitStatus FBClassify::Group(int which) |
---|
53 |
{ |
---|
54 |
unsigned int i; |
---|
55 |
int hashvalue, v0, v1, v2, tv0, tv1, tv2; |
---|
56 |
|
---|
57 |
int *hasharrayptr, hasharraysize, j, itri; |
---|
58 |
FBPolyhedron *poly; |
---|
59 |
const int primes[] = { 307, 601, 1009, 3001, 6007, 10007, 30011, 60013, 100003 }; |
---|
60 |
int numhashbins; |
---|
61 |
|
---|
62 |
if ( which == 1 ) poly = polya; |
---|
63 |
else poly = polyb; |
---|
64 |
|
---|
65 |
if ( poly == 0 ) { |
---|
66 |
PRINT_ERROR("ERROR: Polyhedral object undefined in FBClassify::Group\n"); |
---|
67 |
return CUBIT_FAILURE; |
---|
68 |
} |
---|
69 |
|
---|
70 |
i = poly->verts.size(); |
---|
71 |
if ( i < 3000 ) numhashbins = primes[0]; |
---|
72 |
else if ( i < 6000 ) numhashbins = primes[1]; |
---|
73 |
else if ( i < 10000 ) numhashbins = primes[2]; |
---|
74 |
else if ( i < 30000 ) numhashbins = primes[3]; |
---|
75 |
else if ( i < 60000 ) numhashbins = primes[4]; |
---|
76 |
else if ( i < 100000 ) numhashbins = primes[5]; |
---|
77 |
else if ( i < 300000 ) numhashbins = primes[6]; |
---|
78 |
|
---|
79 |
else numhashbins = primes[7]; |
---|
80 |
|
---|
81 |
IntegerHash *hash = new IntegerHash(numhashbins,50); |
---|
82 |
e0 = new int[poly->tris.size()]; |
---|
83 |
e1 = new int[poly->tris.size()]; |
---|
84 |
e2 = new int[poly->tris.size()]; |
---|
85 |
|
---|
86 |
for ( i = 0; i < poly->tris.size(); i++ ) { |
---|
87 |
e0[i] = e1[i] = e2[i]= NO_EDGE_NBR; |
---|
88 |
group.push_back(UNSET); |
---|
89 |
v0 = poly->tris[i]->v0; |
---|
90 |
v1 = poly->tris[i]->v1; |
---|
91 |
v2 = poly->tris[i]->v2; |
---|
92 |
|
---|
93 |
hash->addtoHashList((v0+v1)%numhashbins,i); |
---|
94 |
hash->addtoHashList((v1+v2)%numhashbins,i); |
---|
95 |
hash->addtoHashList((v2+v0)%numhashbins,i); |
---|
96 |
} |
---|
97 |
|
---|
98 |
|
---|
99 |
|
---|
100 |
|
---|
101 |
|
---|
102 |
|
---|
103 |
|
---|
104 |
|
---|
105 |
|
---|
106 |
|
---|
107 |
|
---|
108 |
|
---|
109 |
|
---|
110 |
|
---|
111 |
|
---|
112 |
for ( i = 0; i < poly->tris.size(); i++ ) { |
---|
113 |
v0 = poly->tris[i]->v0; |
---|
114 |
v1 = poly->tris[i]->v1; |
---|
115 |
v2 = poly->tris[i]->v2; |
---|
116 |
hashvalue = (v0+v1)%numhashbins; |
---|
117 |
hasharrayptr = hash->getHashBin(hashvalue,&hasharraysize); |
---|
118 |
for ( j = 0; j < hasharraysize; j++ ) { |
---|
119 |
|
---|
120 |
|
---|
121 |
itri = hasharrayptr[j]; |
---|
122 |
if ( (unsigned int)itri == i ) continue; |
---|
123 |
tv0 = poly->tris[itri]->v0; |
---|
124 |
tv1 = poly->tris[itri]->v1; |
---|
125 |
tv2 = poly->tris[itri]->v2; |
---|
126 |
if ( ((v1 == tv0) && (v0 == tv1)) || ((v0 == tv0) && (v1 == tv1)) ) { |
---|
127 |
e0[i] = itri; |
---|
128 |
} else if ( ((v1 == tv1) && (v0 == tv2)) || ((v0 == tv1) && (v1 == tv2)) ) { |
---|
129 |
e0[i] = itri; |
---|
130 |
} else if ( ((v1 == tv2) && (v0 == tv0)) || ((v0 == tv2) && (v1 == tv0)) ) { |
---|
131 |
e0[i] = itri; |
---|
132 |
} |
---|
133 |
} |
---|
134 |
hashvalue = (v1+v2)%numhashbins; |
---|
135 |
hasharrayptr = hash->getHashBin(hashvalue,&hasharraysize); |
---|
136 |
for ( j = 0; j < hasharraysize; j++ ) { |
---|
137 |
itri = hasharrayptr[j]; |
---|
138 |
if ( (unsigned int)itri == i ) continue; |
---|
139 |
tv0 = poly->tris[itri]->v0; |
---|
140 |
tv1 = poly->tris[itri]->v1; |
---|
141 |
tv2 = poly->tris[itri]->v2; |
---|
142 |
if ( ((v1 == tv0) && (v2 == tv1)) || ((v2 == tv0) && (v1 == tv1)) ) { |
---|
143 |
e1[i] = itri; |
---|
144 |
} else if ( ((v1 == tv1) && (v2 == tv2)) || ((v2 == tv1) && (v1 == tv2)) ) { |
---|
145 |
e1[i] = itri; |
---|
146 |
} else if ( ((v1 == tv2) && (v2 == tv0)) || ((v2 == tv2) && (v1 == tv0)) ) { |
---|
147 |
e1[i] = itri; |
---|
148 |
} |
---|
149 |
} |
---|
150 |
hashvalue = (v2+v0)%numhashbins; |
---|
151 |
hasharrayptr = hash->getHashBin(hashvalue,&hasharraysize); |
---|
152 |
for ( j = 0; j < hasharraysize; j++ ) { |
---|
153 |
itri = hasharrayptr[j]; |
---|
154 |
if ( (unsigned int)itri == i ) continue; |
---|
155 |
tv0 = poly->tris[itri]->v0; |
---|
156 |
tv1 = poly->tris[itri]->v1; |
---|
157 |
tv2 = poly->tris[itri]->v2; |
---|
158 |
if ( ((v0 == tv0) && (v2 == tv1)) || ((v2 == tv0) && (v0 == tv1)) ) { |
---|
159 |
e2[i] = itri; |
---|
160 |
} else if ( ((v0 == tv1) && (v2 == tv2)) || ((v2 == tv1) && (v0 == tv2)) ) { |
---|
161 |
e2[i] = itri; |
---|
162 |
} else if ( ((v0 == tv2) && (v2 == tv0)) || ((v2 == tv2) && (v0 == tv0)) ) { |
---|
163 |
e2[i] = itri; |
---|
164 |
} |
---|
165 |
} |
---|
166 |
} |
---|
167 |
|
---|
168 |
|
---|
169 |
for ( i = 0; i < poly->intersection_edges.size(); i++ ) { |
---|
170 |
v0 = poly->intersection_edges[i]->v0; |
---|
171 |
v1 = poly->intersection_edges[i]->v1; |
---|
172 |
hashvalue = (v0+v1)%numhashbins; |
---|
173 |
hasharrayptr = hash->getHashBin(hashvalue,&hasharraysize); |
---|
174 |
for ( j = 0; j < hasharraysize; j++ ) { |
---|
175 |
itri = hasharrayptr[j]; |
---|
176 |
tv0 = poly->tris[itri]->v0; |
---|
177 |
tv1 = poly->tris[itri]->v1; |
---|
178 |
tv2 = poly->tris[itri]->v2; |
---|
179 |
if ( ((v0 == tv0) && (v1 == tv1)) || ((v0 == tv1) && (v1 == tv0)) ) { |
---|
180 |
e0[itri] = NO_EDGE_NBR; |
---|
181 |
}; |
---|
182 |
if ( ((v0 == tv0) && (v1 == tv2)) || ((v0 == tv2) && (v1 == tv0)) ) { |
---|
183 |
e2[itri] = NO_EDGE_NBR; |
---|
184 |
}; |
---|
185 |
if ( ((v0 == tv1) && (v1 == tv2)) || ((v0 == tv2) && (v1 == tv1)) ) { |
---|
186 |
e1[itri] = NO_EDGE_NBR; |
---|
187 |
} |
---|
188 |
} |
---|
189 |
} |
---|
190 |
|
---|
191 |
|
---|
192 |
for ( i = 0; i < poly->tris.size(); i++ ) { |
---|
193 |
if ( group[i] == UNSET ) { |
---|
194 |
fill_group(i,number_of_groups++); |
---|
195 |
} |
---|
196 |
} |
---|
197 |
|
---|
198 |
delete[] e0; delete[] e1; delete[] e2; |
---|
199 |
delete hash; |
---|
200 |
|
---|
201 |
return CUBIT_SUCCESS; |
---|
202 |
} |
---|
203 |
|
---|
204 |
void FBClassify::fill_group(int itri, int ngroup) |
---|
205 |
{ |
---|
206 |
std::stack<int> vstack; |
---|
207 |
int ktri; |
---|
208 |
|
---|
209 |
group[itri] = ngroup; |
---|
210 |
|
---|
211 |
if ( (e0[itri] != NO_EDGE_NBR) && (group[e0[itri]] == UNSET) ) |
---|
212 |
vstack.push(e0[itri]); |
---|
213 |
if ( (e1[itri] != NO_EDGE_NBR) && (group[e1[itri]] == UNSET) ) |
---|
214 |
vstack.push(e1[itri]); |
---|
215 |
if ( (e2[itri] != NO_EDGE_NBR) && (group[e2[itri]] == UNSET) ) |
---|
216 |
vstack.push(e2[itri]); |
---|
217 |
while (vstack.size() > 0 ) { |
---|
218 |
ktri = vstack.top(); |
---|
219 |
vstack.pop(); |
---|
220 |
group[ktri] = ngroup; |
---|
221 |
if ( (e0[ktri] != NO_EDGE_NBR) && (group[e0[ktri]] == UNSET) ) |
---|
222 |
vstack.push(e0[ktri]); |
---|
223 |
if ( (e1[ktri] != NO_EDGE_NBR) && (group[e1[ktri]] == UNSET) ) |
---|
224 |
vstack.push(e1[ktri]); |
---|
225 |
if ( (e2[ktri] != NO_EDGE_NBR) && (group[e2[ktri]] == UNSET) ) |
---|
226 |
vstack.push(e2[ktri]); |
---|
227 |
} |
---|
228 |
} |
---|
229 |
|
---|
230 |
CubitStatus FBClassify::CharacterizeGroups(int which, bool other_is_planar) |
---|
231 |
{ |
---|
232 |
int numberdone; |
---|
233 |
unsigned int i; |
---|
234 |
bool ifoundit; |
---|
235 |
int itri = -1, type; |
---|
236 |
FBPolyhedron *poly; |
---|
237 |
|
---|
238 |
if ( which == 1 ) poly = polya; |
---|
239 |
else poly = polyb; |
---|
240 |
|
---|
241 |
if ( poly == 0 ) { |
---|
242 |
PRINT_ERROR("ERROR: Polyhedral object undefined in FBClassify::Group\n"); |
---|
243 |
return CUBIT_FAILURE; |
---|
244 |
} |
---|
245 |
|
---|
246 |
numberdone = 0; |
---|
247 |
i = 0; |
---|
248 |
while ( numberdone < number_of_groups ) { |
---|
249 |
ifoundit = false; |
---|
250 |
while ( i < poly->tris.size() ) { |
---|
251 |
if ( group[i] == numberdone ) { |
---|
252 |
ifoundit = true; |
---|
253 |
itri = (int)i; |
---|
254 |
break; |
---|
255 |
} |
---|
256 |
i++; |
---|
257 |
} |
---|
258 |
if ( ifoundit == true ) { |
---|
259 |
if ( other_is_planar == false ) |
---|
260 |
type = classify(itri, which); |
---|
261 |
else |
---|
262 |
type = classify_against_plane(itri, which); |
---|
263 |
group_characterization.push_back(type); |
---|
264 |
numberdone++; |
---|
265 |
} else { |
---|
266 |
PRINT_ERROR("Error in FBClassify::CharacterizeGroups\n"); |
---|
267 |
return CUBIT_FAILURE; |
---|
268 |
} |
---|
269 |
} |
---|
270 |
return CUBIT_SUCCESS; |
---|
271 |
} |
---|
272 |
|
---|
273 |
int FBClassify::classify_against_plane(int itri, int which) |
---|
274 |
{ |
---|
275 |
FBPolyhedron *polyref, *polyobj; |
---|
276 |
int type, v0, v1, v2; |
---|
277 |
double xbary, ybary, zbary, a, b, c; |
---|
278 |
; |
---|
279 |
|
---|
280 |
type = FB_ORIENTATION_UNDEFINED; |
---|
281 |
if ( which == 1 ) { |
---|
282 |
polyobj = polyb; |
---|
283 |
polyref = polya; |
---|
284 |
} else if ( which == 2 ) { |
---|
285 |
polyobj = polya; |
---|
286 |
polyref = polya; |
---|
287 |
} else { |
---|
288 |
PRINT_ERROR("ERROR in FBClassify::classify\n"); |
---|
289 |
return type; |
---|
290 |
} |
---|
291 |
v0 = polyref->tris[itri]->v0; |
---|
292 |
v1 = polyref->tris[itri]->v1; |
---|
293 |
v2 = polyref->tris[itri]->v2; |
---|
294 |
|
---|
295 |
xbary = ( polyref->verts[v0]->coord[0] + |
---|
296 |
polyref->verts[v1]->coord[0] + |
---|
297 |
polyref->verts[v2]->coord[0] )/3.; |
---|
298 |
ybary = ( polyref->verts[v0]->coord[1] + |
---|
299 |
polyref->verts[v1]->coord[1] + |
---|
300 |
polyref->verts[v2]->coord[1] )/3.; |
---|
301 |
zbary = ( polyref->verts[v0]->coord[2] + |
---|
302 |
polyref->verts[v1]->coord[2] + |
---|
303 |
polyref->verts[v2]->coord[2] )/3.; |
---|
304 |
a = polyref->tris[itri]->a; |
---|
305 |
b = polyref->tris[itri]->b; |
---|
306 |
c = polyref->tris[itri]->c; |
---|
307 |
|
---|
308 |
|
---|
309 |
|
---|
310 |
|
---|
311 |
|
---|
312 |
double obj_tri_a, obj_tri_b, obj_tri_c, obj_tri_d, dotprod, disttoplane; |
---|
313 |
|
---|
314 |
obj_tri_a = polyobj->tris[0]->a; |
---|
315 |
obj_tri_b = polyobj->tris[0]->b; |
---|
316 |
obj_tri_c = polyobj->tris[0]->c; |
---|
317 |
obj_tri_d = polyobj->tris[0]->d; |
---|
318 |
|
---|
319 |
disttoplane = obj_tri_a*xbary + obj_tri_b*ybary + obj_tri_c*zbary + obj_tri_d; |
---|
320 |
if ( disttoplane > EPSILON ) return FB_ORIENTATION_OUTSIDE; |
---|
321 |
else if ( disttoplane < -EPSILON ) return FB_ORIENTATION_INSIDE; |
---|
322 |
|
---|
323 |
dotprod = obj_tri_a*a + obj_tri_b*b + obj_tri_c*c; |
---|
324 |
if ( dotprod > 0. ) return FB_ORIENTATION_SAME; |
---|
325 |
else return FB_ORIENTATION_OPPOSITE; |
---|
326 |
|
---|
327 |
return type; |
---|
328 |
} |
---|
329 |
|
---|
330 |
|
---|
331 |
|
---|
332 |
|
---|
333 |
|
---|
334 |
|
---|
335 |
|
---|
336 |
int FBClassify::classify(int itri, int which) |
---|
337 |
{ |
---|
338 |
|
---|
339 |
|
---|
340 |
|
---|
341 |
FBPolyhedron *polyref, *polyobj; |
---|
342 |
|
---|
343 |
int type, v0, v1, v2; |
---|
344 |
double xbary, ybary, zbary, a, b, c; |
---|
345 |
|
---|
346 |
type = FB_ORIENTATION_UNDEFINED; |
---|
347 |
if ( which == 1 ) { |
---|
348 |
polyobj = polyb; |
---|
349 |
polyref = polya; |
---|
350 |
} else if ( which == 2 ) { |
---|
351 |
polyobj = polya; |
---|
352 |
polyref = polyb; |
---|
353 |
} else { |
---|
354 |
PRINT_ERROR("ERROR in FBClassify::classify\n"); |
---|
355 |
return type; |
---|
356 |
} |
---|
357 |
int mydebug = 0; |
---|
358 |
if(mydebug) |
---|
359 |
polyref->debug_draw_fb_triangle(polyref->tris[itri]); |
---|
360 |
|
---|
361 |
v0 = polyref->tris[itri]->v0; |
---|
362 |
v1 = polyref->tris[itri]->v1; |
---|
363 |
v2 = polyref->tris[itri]->v2; |
---|
364 |
|
---|
365 |
xbary = ( polyref->verts[v0]->coord[0] + |
---|
366 |
polyref->verts[v1]->coord[0] + |
---|
367 |
polyref->verts[v2]->coord[0] )/3.; |
---|
368 |
ybary = ( polyref->verts[v0]->coord[1] + |
---|
369 |
polyref->verts[v1]->coord[1] + |
---|
370 |
polyref->verts[v2]->coord[1] )/3.; |
---|
371 |
zbary = ( polyref->verts[v0]->coord[2] + |
---|
372 |
polyref->verts[v1]->coord[2] + |
---|
373 |
polyref->verts[v2]->coord[2] )/3.; |
---|
374 |
a = polyref->tris[itri]->a; |
---|
375 |
b = polyref->tris[itri]->b; |
---|
376 |
c = polyref->tris[itri]->c; |
---|
377 |
|
---|
378 |
unsigned int i, num_perturb; |
---|
379 |
double obj_tri_a, obj_tri_b, obj_tri_c, obj_tri_d, dotprod; |
---|
380 |
double distance_to_other_sqr, closest_distance_to_plane, t; |
---|
381 |
double closest_distance_to_other_sqr; |
---|
382 |
double xint, yint, zint, distance_to_plane, closest_dotproduct; |
---|
383 |
bool perturb, done, foundone; |
---|
384 |
double other_xbar, other_ybar, other_zbar; |
---|
385 |
int other_tri_0, other_tri_1, other_tri_2; |
---|
386 |
|
---|
387 |
perturb = false; |
---|
388 |
num_perturb = 0; |
---|
389 |
done = false; |
---|
390 |
|
---|
391 |
while ( (done == false) && (num_perturb < 20) ) { |
---|
392 |
closest_dotproduct = -CUBIT_DBL_MAX + 1.; |
---|
393 |
closest_distance_to_plane = CUBIT_DBL_MAX; |
---|
394 |
closest_distance_to_other_sqr = CUBIT_DBL_MAX; |
---|
395 |
foundone = false; |
---|
396 |
for ( i = 0; i < polyobj->tris.size(); i++ ) { |
---|
397 |
obj_tri_a = polyobj->tris[i]->a; |
---|
398 |
obj_tri_b = polyobj->tris[i]->b; |
---|
399 |
obj_tri_c = polyobj->tris[i]->c; |
---|
400 |
obj_tri_d = polyobj->tris[i]->d; |
---|
401 |
|
---|
402 |
dotprod = obj_tri_a*a + obj_tri_b*b + obj_tri_c*c; |
---|
403 |
|
---|
404 |
distance_to_plane = (obj_tri_a*xbary + obj_tri_b*ybary + |
---|
405 |
obj_tri_c*zbary + obj_tri_d); |
---|
406 |
|
---|
407 |
|
---|
408 |
if ( fabs(dotprod) < EPSILON_CLASSIFY ) { |
---|
409 |
|
---|
410 |
if ( fabs(distance_to_plane) < EPSILON_CLASSIFY ) { |
---|
411 |
|
---|
412 |
perturb = true; |
---|
413 |
num_perturb += 1; |
---|
414 |
break; |
---|
415 |
} |
---|
416 |
continue; |
---|
417 |
} |
---|
418 |
|
---|
419 |
t =-(distance_to_plane)/dotprod; |
---|
420 |
if ( t < -EPSILON_CLASSIFY ) continue; |
---|
421 |
xint = xbary + a*t; |
---|
422 |
yint = ybary + b*t; |
---|
423 |
zint = zbary + c*t; |
---|
424 |
|
---|
425 |
|
---|
426 |
|
---|
427 |
|
---|
428 |
if ( (polyobj->tris[i]->boundingbox.xmin - EPSILON > xint) || |
---|
429 |
(polyobj->tris[i]->boundingbox.xmax + EPSILON < xint) || |
---|
430 |
(polyobj->tris[i]->boundingbox.ymin - EPSILON > yint) || |
---|
431 |
(polyobj->tris[i]->boundingbox.ymax + EPSILON < yint) || |
---|
432 |
(polyobj->tris[i]->boundingbox.zmin - EPSILON > zint) || |
---|
433 |
(polyobj->tris[i]->boundingbox.zmax + EPSILON < zint) ) |
---|
434 |
continue; |
---|
435 |
|
---|
436 |
|
---|
437 |
|
---|
438 |
double xp1, yp1, xp2, yp2, xp3, yp3, ptx, pty; |
---|
439 |
int retval; |
---|
440 |
|
---|
441 |
if ( (fabs(obj_tri_b) >= fabs(obj_tri_a)) && |
---|
442 |
(fabs(obj_tri_b) >= fabs(obj_tri_c)) ) { |
---|
443 |
xp1 = polyobj->verts[polyobj->tris[i]->v0]->coord[0]; |
---|
444 |
yp1 = polyobj->verts[polyobj->tris[i]->v0]->coord[2]; |
---|
445 |
xp2 = polyobj->verts[polyobj->tris[i]->v1]->coord[0]; |
---|
446 |
yp2 = polyobj->verts[polyobj->tris[i]->v1]->coord[2]; |
---|
447 |
xp3 = polyobj->verts[polyobj->tris[i]->v2]->coord[0]; |
---|
448 |
yp3 = polyobj->verts[polyobj->tris[i]->v2]->coord[2]; |
---|
449 |
ptx = xint; |
---|
450 |
pty = zint; |
---|
451 |
} else if ( fabs(obj_tri_a) >= fabs(obj_tri_c) ) { |
---|
452 |
xp1 = polyobj->verts[polyobj->tris[i]->v0]->coord[1]; |
---|
453 |
yp1 = polyobj->verts[polyobj->tris[i]->v0]->coord[2]; |
---|
454 |
xp2 = polyobj->verts[polyobj->tris[i]->v1]->coord[1]; |
---|
455 |
yp2 = polyobj->verts[polyobj->tris[i]->v1]->coord[2]; |
---|
456 |
xp3 = polyobj->verts[polyobj->tris[i]->v2]->coord[1]; |
---|
457 |
yp3 = polyobj->verts[polyobj->tris[i]->v2]->coord[2]; |
---|
458 |
ptx = yint; |
---|
459 |
pty = zint; |
---|
460 |
} else { |
---|
461 |
xp1 = polyobj->verts[polyobj->tris[i]->v0]->coord[0]; |
---|
462 |
yp1 = polyobj->verts[polyobj->tris[i]->v0]->coord[1]; |
---|
463 |
xp2 = polyobj->verts[polyobj->tris[i]->v1]->coord[0]; |
---|
464 |
yp2 = polyobj->verts[polyobj->tris[i]->v1]->coord[1]; |
---|
465 |
xp3 = polyobj->verts[polyobj->tris[i]->v2]->coord[0]; |
---|
466 |
yp3 = polyobj->verts[polyobj->tris[i]->v2]->coord[1]; |
---|
467 |
ptx = xint; |
---|
468 |
pty = yint; |
---|
469 |
} |
---|
470 |
retval = pt_in_tri_2d(ptx,pty,xp1,yp1,xp2,yp2,xp3,yp3); |
---|
471 |
if ( (retval == FB_ORIENTATION_INSIDE) || |
---|
472 |
(retval == FB_ORIENTATION_ON) ) { |
---|
473 |
|
---|
474 |
other_tri_0 = polyobj->tris[i]->v0; |
---|
475 |
other_tri_1 = polyobj->tris[i]->v1; |
---|
476 |
other_tri_2 = polyobj->tris[i]->v2; |
---|
477 |
other_xbar = ( polyobj->verts[other_tri_0]->coord[0] + |
---|
478 |
polyobj->verts[other_tri_1]->coord[0] + |
---|
479 |
polyobj->verts[other_tri_2]->coord[0] )/3.; |
---|
480 |
other_ybar = ( polyobj->verts[other_tri_0]->coord[1] + |
---|
481 |
polyobj->verts[other_tri_1]->coord[1] + |
---|
482 |
polyobj->verts[other_tri_2]->coord[1] )/3.; |
---|
483 |
other_zbar = ( polyobj->verts[other_tri_0]->coord[2] + |
---|
484 |
polyobj->verts[other_tri_1]->coord[2] + |
---|
485 |
polyobj->verts[other_tri_2]->coord[2] )/3.; |
---|
486 |
|
---|
487 |
|
---|
488 |
distance_to_other_sqr = ( (xbary-other_xbar)*(xbary-other_xbar) + |
---|
489 |
(ybary-other_ybar)*(ybary-other_ybar) + |
---|
490 |
(zbary-other_zbar)*(zbary-other_zbar) ); |
---|
491 |
|
---|
492 |
if(closest_distance_to_other_sqr > distance_to_other_sqr){ |
---|
493 |
|
---|
494 |
|
---|
495 |
foundone = true; |
---|
496 |
closest_distance_to_other_sqr = distance_to_other_sqr; |
---|
497 |
closest_dotproduct = dotprod; |
---|
498 |
if(mydebug){ |
---|
499 |
polyobj->debug_draw_fb_triangle(polyobj->tris[i]); |
---|
500 |
GfxDebug::mouse_xforms(); |
---|
501 |
} |
---|
502 |
closest_distance_to_plane = distance_to_plane; |
---|
503 |
|
---|
504 |
if ( fabs(closest_distance_to_plane) < EPSILON_CLASSIFY ) |
---|
505 |
break; |
---|
506 |
} |
---|
507 |
} |
---|
508 |
} |
---|
509 |
if ( perturb == false ) done = true; |
---|
510 |
else { |
---|
511 |
|
---|
512 |
perturb_the_ray(xbary, ybary, zbary); |
---|
513 |
perturb = false; |
---|
514 |
} |
---|
515 |
} |
---|
516 |
|
---|
517 |
|
---|
518 |
|
---|
519 |
if ( (fabs(closest_distance_to_plane) < EPSILON_CLASSIFY ) ){ |
---|
520 |
if ( closest_dotproduct > 0.0 ) |
---|
521 |
type = FB_ORIENTATION_SAME; |
---|
522 |
else |
---|
523 |
type = FB_ORIENTATION_OPPOSITE; |
---|
524 |
} |
---|
525 |
|
---|
526 |
|
---|
527 |
|
---|
528 |
|
---|
529 |
else { |
---|
530 |
if ( foundone == false ) |
---|
531 |
type = FB_ORIENTATION_OUTSIDE; |
---|
532 |
else if ( closest_dotproduct > 0.0 ) |
---|
533 |
type = FB_ORIENTATION_INSIDE; |
---|
534 |
else |
---|
535 |
type = FB_ORIENTATION_OUTSIDE; |
---|
536 |
} |
---|
537 |
if(mydebug) |
---|
538 |
GfxDebug::display_all(); |
---|
539 |
return type; |
---|
540 |
} |
---|
541 |
|
---|
542 |
void FBClassify::perturb_the_ray(double &xbary, double &ybary, double &zbary) |
---|
543 |
{ |
---|
544 |
xbary += 1.e-4*(double(rand())/(RAND_MAX+1.0)-0.5); |
---|
545 |
ybary += 1.e-4*(double(rand())/(RAND_MAX+1.0)-0.5); |
---|
546 |
zbary += 1.e-4*(double(rand())/(RAND_MAX+1.0)-0.5); |
---|
547 |
} |
---|
548 |
|
---|
549 |
int FBClassify::pt_in_tri_2d(double xpt, double ypt, |
---|
550 |
double x0, double y0, |
---|
551 |
double x1, double y1, |
---|
552 |
double x2, double y2) |
---|
553 |
{ |
---|
554 |
|
---|
555 |
|
---|
556 |
|
---|
557 |
|
---|
558 |
double c0, c1, c2; |
---|
559 |
double e0x, e1x, e2x, e0y, e1y, e2y; |
---|
560 |
double n0x, n1x, n2x, n0y, n1y, n2y; |
---|
561 |
double denom0, denom1, denom2; |
---|
562 |
int result; |
---|
563 |
|
---|
564 |
e0x = x1 - x0; e0y = y1 - y0; |
---|
565 |
e1x = x2 - x1; e1y = y2 - y1; |
---|
566 |
e2x = x0 - x2; e2y = y0 - y2; |
---|
567 |
n0x = e0y; n0y = -e0x; |
---|
568 |
n1x = e1y; n1y = -e1x; |
---|
569 |
n2x = e2y; n2y = -e2x; |
---|
570 |
denom0 = n1x*e0x + n1y*e0y; |
---|
571 |
if ( fabs(denom0) < EPSILON_CLASSIFY ) { |
---|
572 |
PRINT_ERROR("Failure in pt_in_tri_2d; needle-thin triangle encountered.\n"); |
---|
573 |
return FB_ORIENTATION_UNDEFINED; |
---|
574 |
} |
---|
575 |
denom1 = n2x*e1x + n2y*e1y; |
---|
576 |
if ( fabs(denom1) < EPSILON_CLASSIFY ) { |
---|
577 |
PRINT_ERROR("Failure in pt_in_tri_2d; needle-thin triangle encountered.\n"); |
---|
578 |
return FB_ORIENTATION_UNDEFINED; |
---|
579 |
} |
---|
580 |
denom2 = n0x*e2x + n0y*e2y; |
---|
581 |
if ( fabs(denom2) < EPSILON_CLASSIFY ) { |
---|
582 |
PRINT_ERROR("Failure in pt_in_tri_2d; needle-thin triangle encountered.\n"); |
---|
583 |
return FB_ORIENTATION_UNDEFINED; |
---|
584 |
} |
---|
585 |
|
---|
586 |
c0 = -( n1x*(xpt-x1) + n1y*(ypt-y1) )/denom0; |
---|
587 |
c1 = -( n2x*(xpt-x2) + n2y*(ypt-y2) )/denom1; |
---|
588 |
c2 = -( n0x*(xpt-x0) + n0y*(ypt-y0) )/denom2; |
---|
589 |
|
---|
590 |
if ( (c0 > 0.0) && (c1 > 0.0) && (c2 > 0.0) ) result = FB_ORIENTATION_INSIDE; |
---|
591 |
else if ( (c0 < 0.0) || (c1 < 0.0) || (c2 < 0.0) ) result = FB_ORIENTATION_OUTSIDE; |
---|
592 |
else result = FB_ORIENTATION_ON; |
---|
593 |
|
---|
594 |
return result; |
---|
595 |
|
---|
596 |
} |
---|
597 |
|
---|
598 |
void FBClassify::get_group(std::vector<int> **this_group, |
---|
599 |
std::vector<int> **this_group_characterization) |
---|
600 |
{ |
---|
601 |
*this_group = &group; |
---|
602 |
*this_group_characterization = &group_characterization; |
---|
603 |
} |
---|