Branch data Line data Source code
1 : : /*
2 : : *
3 : : *
4 : : * Copyright (C) 2004 Sandia Corporation. Under the terms of Contract DE-AC04-94AL85000
5 : : * with Sandia Corporation, the U.S. Government retains certain rights in this software.
6 : : *
7 : : * This file is part of facetbool--contact via [email protected]
8 : : *
9 : : * This library is free software; you can redistribute it and/or
10 : : * modify it under the terms of the GNU Lesser General Public
11 : : * License as published by the Free Software Foundation; either
12 : : * version 2.1 of the License, or (at your option) any later version.
13 : : *
14 : : * This library is distributed in the hope that it will be useful,
15 : : * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 : : * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 : : * Lesser General Public License for more details.
18 : : *
19 : : * You should have received a copy of the GNU Lesser General Public
20 : : * License along with this library; if not, write to the Free Software
21 : : * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
22 : : *
23 : : *
24 : : *
25 : : */
26 : :
27 : : #include <vector>
28 : : #include "FBIntersect.hpp"
29 : : #include "FBPolyhedron.hpp"
30 : : #include "FBRetriangulate.hpp"
31 : : #include "CubitMessage.hpp"
32 : : #include "GfxDebug.hpp"
33 : : #include "GeometryDefines.h"
34 : :
35 : 0 : FBIntersect::FBIntersect()
36 : : {
37 : 0 : poly1 = poly2 = 0;
38 : 0 : do_classify = false;
39 : 0 : do_edges_only = false;
40 : 0 : do_imprint = false;
41 : 0 : classify1 = classify2 = 0;
42 : 0 : body1_is_plane = body2_is_plane = false;
43 : 0 : f_c_indices1 = f_c_indices2 = 0;
44 : 0 : nothing_intersected = false;
45 : 0 : }
46 : :
47 : 0 : FBIntersect::~FBIntersect()
48 : : {
49 [ # # ][ # # ]: 0 : if ( poly1 ) delete poly1;
50 [ # # ][ # # ]: 0 : if ( poly2 ) delete poly2;
51 [ # # ][ # # ]: 0 : if ( classify1 ) delete classify1;
52 [ # # ][ # # ]: 0 : if ( classify2 ) delete classify2;
53 : 0 : }
54 : :
55 : 0 : CubitStatus FBIntersect::intersect(const std::vector<double>& Ticoords,
56 : : const std::vector<int>& Ticonnections,
57 : : const std::vector<double>& Tjcoords,
58 : : const std::vector<int>& Tjconnections,
59 : : std::vector<int>& duddedTiFacets,
60 : : std::vector<int>& duddedTjFacets,
61 : : std::vector<int>& newTiFacets,
62 : : std::vector<int>& newTjFacets,
63 : : std::vector<int>& newTiFacetsIndex,
64 : : std::vector<int>& newTjFacetsIndex,
65 : : std::vector<double>& newTiPoints,
66 : : std::vector<double>& newTjPoints,
67 : : std::vector<int>& edgesTi,
68 : : std::vector<int>& edgesTj
69 : : )
70 : : {
71 : : CubitStatus status;
72 : : bool boxes_intersect;
73 : :
74 : 0 : status = CUBIT_SUCCESS;
75 : :
76 [ # # ]: 0 : if ( Ticoords.size()%3 != 0 ) {
77 [ # # ][ # # ]: 0 : PRINT_ERROR("Bad coordinates for first part fed to FBIntersect.\n");
78 : 0 : return CUBIT_FAILURE;
79 : : }
80 [ # # ]: 0 : if ( Tjcoords.size()%3 != 0 ) {
81 [ # # ][ # # ]: 0 : PRINT_ERROR("Bad coordinates for second part fed to FBIntersect.\n");
82 : 0 : return CUBIT_FAILURE;
83 : : }
84 [ # # ]: 0 : if ( Ticonnections.size()%3 != 0 ) {
85 [ # # ][ # # ]: 0 : PRINT_ERROR("Bad connection list for first part fed to FBIntersect.\n");
86 : 0 : return CUBIT_FAILURE;
87 : : }
88 [ # # ]: 0 : if ( Tjconnections.size()%3 != 0 ) {
89 [ # # ][ # # ]: 0 : PRINT_ERROR("Bad connection list for second part fed to FBIntersect.\n");
90 : 0 : return CUBIT_FAILURE;
91 : : }
92 : :
93 [ # # ]: 0 : poly1 = new FBPolyhedron;
94 : 0 : status = poly1->makepoly(Ticoords,Ticonnections,f_c_indices1);
95 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
96 : : {
97 : 0 : return status;
98 : : }
99 : :
100 [ # # ]: 0 : poly2 = new FBPolyhedron;
101 : 0 : status = poly2->makepoly(Tjcoords,Tjconnections,f_c_indices2);
102 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
103 : : {
104 : 0 : return status;
105 : : }
106 : :
107 : 0 : boxes_intersect = poly1->boxintersection(poly2);
108 [ # # ]: 0 : if ( boxes_intersect == false ) {
109 : 0 : nothing_intersected = true;
110 : 0 : return status;
111 : : }
112 : :
113 : 0 : status = pair_intersect();
114 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
115 : : {
116 : 0 : return status;
117 : : }
118 : :
119 : 0 : poly1->putnewpoints(newTiPoints);
120 : 0 : poly2->putnewpoints(newTjPoints);
121 : :
122 : 0 : poly1->putedges(edgesTi);
123 : 0 : poly2->putedges(edgesTj);
124 : :
125 : 0 : status = poly1->retriangulate(newTiFacets, newTiFacetsIndex);
126 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
127 : : {
128 : 0 : return status;
129 : : }
130 : :
131 : 0 : status = poly2->retriangulate(newTjFacets, newTjFacetsIndex);
132 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
133 : : {
134 : 0 : return status;
135 : : }
136 : :
137 : : unsigned int i;
138 [ # # ]: 0 : for ( i = 0; i < poly1->tris.size(); i++ ) {
139 [ # # ][ # # ]: 0 : if ( poly1->tris[i]->dudded == true ) duddedTiFacets.push_back(i);
140 : : }
141 : :
142 [ # # ]: 0 : for ( i = 0; i < poly2->tris.size(); i++ ) {
143 [ # # ][ # # ]: 0 : if ( poly2->tris[i]->dudded == true ) duddedTjFacets.push_back(i);
144 : : }
145 : :
146 : : // newxxFacets is updated in FBRetriangulate::make_this_tri() after
147 : : // a new triangle has been made. Also newxxFacetsIndex is updated if
148 : : // any new facets were made. But this update points to the start of
149 : : // newxxFacets. Therefore we have to add on the last member to
150 : : // newxxFacetsIndex here.
151 [ # # ]: 0 : newTiFacetsIndex.push_back(newTiFacets.size());
152 [ # # ]: 0 : newTjFacetsIndex.push_back(newTjFacets.size());
153 : :
154 : : // If classification of the intersected surface parts is needed, we need to
155 : : // do some things.
156 : : bool other_is_planar;
157 [ # # ]: 0 : if ( do_classify == true ) {
158 : 0 : poly1->add_new_triangle_data();
159 : 0 : poly2->add_new_triangle_data();
160 : 0 : poly1->removeduddedtriangles();
161 : 0 : poly2->removeduddedtriangles();
162 [ # # ]: 0 : classify1 = new FBClassify;
163 : 0 : classify1->SetPoly(poly1, poly2);
164 : 0 : status = classify1->Group(1);
165 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
166 : : {
167 : 0 : return status;
168 : : }
169 : 0 : other_is_planar = body2_is_plane;
170 : 0 : status = classify1->CharacterizeGroups(1,other_is_planar);
171 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
172 : : {
173 : 0 : return status;
174 : : }
175 [ # # ]: 0 : classify2 = new FBClassify;
176 : 0 : classify2->SetPoly(poly1, poly2);
177 : 0 : status = classify2->Group(2);
178 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
179 : : {
180 : 0 : return status;
181 : : }
182 : 0 : other_is_planar = body1_is_plane;
183 : 0 : status = classify2->CharacterizeGroups(2,other_is_planar);
184 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
185 : : {
186 : 0 : return status;
187 : : }
188 : : }
189 : :
190 : 0 : return status;
191 : :
192 : : }
193 : :
194 : 0 : CubitStatus FBIntersect::intersect(const std::vector<double>& Ticoords,
195 : : const std::vector<int>& Ticonnections,
196 : : const std::vector<double>& Tjcoords,
197 : : const std::vector<int>& Tjconnections,
198 : : std::vector<int>& newTiFacets,
199 : : std::vector<int>& newTjFacets,
200 : : std::vector<int> *indices1,
201 : : std::vector<int> *indices2
202 : : )
203 : : {
204 : : CubitStatus status;
205 : : bool boxes_intersect;
206 : :
207 : 0 : status = CUBIT_SUCCESS;
208 : :
209 [ # # ][ # # ]: 0 : if ( Ticoords.size()%3 != 0 ) {
210 [ # # ][ # # ]: 0 : PRINT_ERROR("Bad coordinates for first part fed to FBIntersect.\n");
[ # # ][ # # ]
211 : 0 : return CUBIT_FAILURE;
212 : : }
213 [ # # ][ # # ]: 0 : if ( Tjcoords.size()%3 != 0 ) {
214 [ # # ][ # # ]: 0 : PRINT_ERROR("Bad coordinates for second part fed to FBIntersect.\n");
[ # # ][ # # ]
215 : 0 : return CUBIT_FAILURE;
216 : : }
217 [ # # ][ # # ]: 0 : if ( Ticonnections.size()%3 != 0 ) {
218 [ # # ][ # # ]: 0 : PRINT_ERROR("Bad connection list for first part fed to FBIntersect.\n");
[ # # ][ # # ]
219 : 0 : return CUBIT_FAILURE;
220 : : }
221 [ # # ][ # # ]: 0 : if ( Tjconnections.size()%3 != 0 ) {
222 [ # # ][ # # ]: 0 : PRINT_ERROR("Bad connection list for second part fed to FBIntersect.\n");
[ # # ][ # # ]
223 : 0 : return CUBIT_FAILURE;
224 : : }
225 : :
226 : 0 : f_c_indices1 = indices1;
227 : 0 : f_c_indices2 = indices2;
228 : :
229 [ # # ][ # # ]: 0 : poly1 = new FBPolyhedron;
230 [ # # ]: 0 : status = poly1->makepoly(Ticoords,Ticonnections,f_c_indices1);
231 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
232 : : {
233 : 0 : return status;
234 : : }
235 : :
236 [ # # ][ # # ]: 0 : poly2 = new FBPolyhedron;
237 [ # # ]: 0 : status = poly2->makepoly(Tjcoords,Tjconnections,f_c_indices2);
238 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
239 : : {
240 : 0 : return status;
241 : : }
242 : :
243 : 0 : double min_angle=0.0, max_angle=0.0;
244 : 0 : int mydebug = 0;
245 [ # # ]: 0 : if(mydebug){
246 [ # # ]: 0 : poly1->min_max_angles_in_polyhedron(min_angle, max_angle);
247 : :
248 [ # # ][ # # ]: 0 : PRINT_INFO("(0) Min angle in poly1 %f, max angle %f\n", min_angle,
[ # # ]
249 [ # # ]: 0 : max_angle);
250 [ # # ]: 0 : poly2->min_max_angles_in_polyhedron(min_angle, max_angle);
251 [ # # ][ # # ]: 0 : PRINT_INFO("(0) Min angle in poly2 %f, max angle %f\n", min_angle,
[ # # ]
252 [ # # ]: 0 : max_angle);
253 : : }
254 [ # # ]: 0 : boxes_intersect = poly1->boxintersection(poly2);
255 [ # # ]: 0 : if ( boxes_intersect == false ) {
256 : 0 : nothing_intersected = true;
257 : : //don't return because in the case of union we need to unite the
258 : : //bodies even if they do not intersect
259 : : //return status;
260 : : }
261 : : else{
262 [ # # ]: 0 : status = pair_intersect();
263 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
264 : : {
265 : 0 : return status;
266 : : }
267 : : }
268 : :
269 [ # # ]: 0 : if(mydebug){
270 [ # # ]: 0 : poly1->min_max_angles_in_polyhedron(min_angle, max_angle);
271 [ # # ][ # # ]: 0 : PRINT_INFO("(1) Min angle in poly1 %f, max angle %f\n", min_angle,
[ # # ]
272 [ # # ]: 0 : max_angle);
273 [ # # ]: 0 : poly2->min_max_angles_in_polyhedron(min_angle, max_angle);
274 [ # # ][ # # ]: 0 : PRINT_INFO("(1) Min angle in poly2 %f, max angle %f\n", min_angle,
[ # # ]
275 [ # # ]: 0 : max_angle);
276 : : }
277 : :
278 [ # # ]: 0 : status = poly1->retriangulate(newTiFacets);
279 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
280 : : {
281 : 0 : return status;
282 : : }
283 [ # # ]: 0 : status = poly2->retriangulate(newTjFacets);
284 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
285 : : {
286 : 0 : return status;
287 : : }
288 : :
289 [ # # ]: 0 : if(mydebug){
290 [ # # ]: 0 : GfxDebug::clear();
291 [ # # ]: 0 : poly1->debug_draw_boundary_edges(CUBIT_BLUE_INDEX);
292 [ # # ]: 0 : GfxDebug::mouse_xforms();
293 [ # # ]: 0 : GfxDebug::clear();
294 [ # # ]: 0 : poly2->debug_draw_boundary_edges(CUBIT_RED_INDEX);
295 [ # # ]: 0 : GfxDebug::mouse_xforms();
296 [ # # ]: 0 : poly1->min_max_angles_in_polyhedron(min_angle, max_angle);
297 : :
298 [ # # ][ # # ]: 0 : PRINT_INFO("(2) Min angle in poly1 %f, max angle %f\n", min_angle,
[ # # ]
299 [ # # ]: 0 : max_angle);
300 [ # # ]: 0 : poly2->min_max_angles_in_polyhedron(min_angle, max_angle);
301 [ # # ][ # # ]: 0 : PRINT_INFO("(2) Min angle in poly2 %f, max angle %f\n", min_angle,
[ # # ]
302 [ # # ]: 0 : max_angle);
303 : : }
304 : : // If classification of the intersected surface parts is needed, we need to
305 : : // do some things.
306 : : bool other_is_planar;
307 [ # # ]: 0 : if ( do_classify == true ) {
308 [ # # ]: 0 : poly1->add_new_triangle_data();
309 [ # # ]: 0 : poly2->add_new_triangle_data();
310 [ # # ]: 0 : poly1->removeduddedtriangles();
311 [ # # ]: 0 : poly2->removeduddedtriangles();
312 [ # # ][ # # ]: 0 : classify1 = new FBClassify;
313 [ # # ]: 0 : classify1->SetPoly(poly1, poly2);
314 [ # # ]: 0 : status = classify1->Group(1);
315 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
316 : : {
317 : 0 : return status;
318 : : }
319 : 0 : other_is_planar = body2_is_plane;
320 [ # # ]: 0 : status = classify1->CharacterizeGroups(1,other_is_planar);
321 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
322 : : {
323 : 0 : return status;
324 : : }
325 [ # # ][ # # ]: 0 : classify2 = new FBClassify;
326 [ # # ]: 0 : classify2->SetPoly(poly1, poly2);
327 [ # # ]: 0 : status = classify2->Group(2);
328 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
329 : : {
330 : 0 : return status;
331 : : }
332 : 0 : other_is_planar = body1_is_plane;
333 [ # # ]: 0 : status = classify2->CharacterizeGroups(2,other_is_planar);
334 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
335 : : {
336 : 0 : return status;
337 : : }
338 : : }
339 : :
340 : 0 : return status;
341 : :
342 : : }
343 : :
344 : 0 : CubitStatus FBIntersect::pair_intersect()
345 : : {
346 : : CubitStatus status;
347 : : unsigned int i, j, k;
348 : : int numboxesfound, *boxlist;
349 : : double dtemp;
350 : :
351 [ # # ][ # # ]: 0 : boxlist = new int[poly2->tris.size()];
[ # # ]
352 : :
353 : 0 : status = CUBIT_SUCCESS;
354 [ # # ][ # # ]: 0 : for ( i = 0; i < poly1->tris.size(); i++ ) {
355 [ # # ]: 0 : FB_Triangle *tri1 = poly1->tris[i];
356 [ # # ][ # # ]: 0 : if ( (tri1->boundingbox.xmax < poly2->polyxmin) ||
357 [ # # ]: 0 : (tri1->boundingbox.xmin > poly2->polyxmax) ||
358 [ # # ]: 0 : (tri1->boundingbox.ymax < poly2->polyymin) ||
359 [ # # ]: 0 : (tri1->boundingbox.ymin > poly2->polyymax) ||
360 [ # # ]: 0 : (tri1->boundingbox.zmax < poly2->polyzmin) ||
361 : 0 : (tri1->boundingbox.zmin > poly2->polyzmax) ) continue;
362 [ # # ]: 0 : poly2->kdtree->box_kdtree_intersect(tri1->boundingbox,&numboxesfound,boxlist);
363 : :
364 [ # # ]: 0 : for ( j = 0; j < (unsigned int)numboxesfound; j++ ) {
365 [ # # ]: 0 : FB_Triangle *tri2 = poly2->tris[boxlist[j]];
366 [ # # ][ # # ]: 0 : if ( (tri1->boundingbox.xmax < tri2->boundingbox.xmin) ||
367 [ # # ]: 0 : (tri1->boundingbox.xmin > tri2->boundingbox.xmax) ||
368 [ # # ]: 0 : (tri1->boundingbox.ymax < tri2->boundingbox.ymin) ||
369 [ # # ]: 0 : (tri1->boundingbox.ymin > tri2->boundingbox.ymax) ||
370 [ # # ]: 0 : (tri1->boundingbox.zmax < tri2->boundingbox.zmin ) ||
371 : 0 : (tri1->boundingbox.zmin > tri2->boundingbox.zmax) ) continue;
372 : :
373 : 0 : int mydebug = 0;
374 [ # # ]: 0 : if ( mydebug )
375 : : {
376 : : // Draw the 2 facets which are intersecting.
377 [ # # ]: 0 : GfxDebug::clear();
378 [ # # ]: 0 : poly1->debug_draw_boundary_edges(CUBIT_YELLOW_INDEX);
379 [ # # ]: 0 : poly2->debug_draw_boundary_edges(CUBIT_PINK_INDEX);
380 [ # # ]: 0 : poly1->debug_draw_fb_triangle(tri1);
381 [ # # ]: 0 : poly2->debug_draw_fb_triangle(tri2);
382 [ # # ]: 0 : GfxDebug::mouse_xforms();
383 : : }
384 : :
385 : : // Find the line of intersection of the two triangle planes.
386 : 0 : linecoeff[0] = tri1->c*tri2->b - tri1->b*tri2->c;
387 : 0 : linecoeff[1] = tri1->a*tri2->c - tri1->c*tri2->a;
388 : 0 : linecoeff[2] = tri1->b*tri2->a - tri1->a*tri2->b;
389 : :
390 [ # # ][ # # ]: 0 : if ( (fabs(linecoeff[0]) < EPSILON) &&
391 [ # # ]: 0 : (fabs(linecoeff[1]) < EPSILON) &&
392 : 0 : (fabs(linecoeff[2]) < EPSILON) ) {
393 [ # # ]: 0 : if ( do_imprint == true ) continue;
394 : : // Don't intersect nearly-coplanar triangles.
395 : :
396 : : // calculate the distance between the triangles. Just because they are coplanar
397 : : // does not mean we have to intersect. If the distance between them is larger than
398 : : // GEOMETRY_RESABS, then just continue on.
399 [ # # ]: 0 : double dist = poly1->verts[tri1->v2]->coord[0]*tri2->a +
400 [ # # ]: 0 : poly1->verts[tri1->v2]->coord[1]*tri2->b +
401 [ # # ]: 0 : poly1->verts[tri1->v2]->coord[2]*tri2->c + tri2->d;
402 [ # # ]: 0 : if ( dist > GEOMETRY_RESABS )
403 : : {
404 : 0 : continue;
405 : : }
406 : :
407 : : // coplanar triangles; tilt each vertex of the triangle up successively
408 : : // and then do the usual tri-tri intersection.
409 : : double ta, tb, tc, td;
410 : : double tx1[3],ty1[3],tz1[3];
411 : :
412 : : // triangle 1
413 : 0 : ta = tri1->a; tb = tri1->b; tc = tri1->c; td = tri1->d;
414 [ # # ]: 0 : for ( k = 0; k < 3; k++ ) {
415 [ # # ]: 0 : tx1[k] = poly1->verts[tri1->v0]->coord[k];
416 [ # # ]: 0 : ty1[k] = poly1->verts[tri1->v1]->coord[k];
417 [ # # ]: 0 : tz1[k] = poly1->verts[tri1->v2]->coord[k];
418 : : }
419 : :
420 [ # # ]: 0 : poly1->verts[tri1->v2]->coord[0] += ta;
421 [ # # ]: 0 : poly1->verts[tri1->v2]->coord[1] += tb;
422 [ # # ]: 0 : poly1->verts[tri1->v2]->coord[2] += tc;
423 : : // Compute new normal and linecoeff
424 [ # # ]: 0 : newplanecoefficients(poly1, tri1);
425 : 0 : linecoeff[0] = tri1->c*tri2->b - tri1->b*tri2->c;
426 : 0 : linecoeff[1] = tri1->a*tri2->c - tri1->c*tri2->a;
427 : 0 : linecoeff[2] = tri1->b*tri2->a - tri1->a*tri2->b;
428 : 0 : dtemp = sqrt(linecoeff[0]*linecoeff[0] +
429 : 0 : linecoeff[1]*linecoeff[1] +
430 : 0 : linecoeff[2]*linecoeff[2]);
431 : 0 : linecoeff[0] /= dtemp;
432 : 0 : linecoeff[1] /= dtemp;
433 : 0 : linecoeff[2] /= dtemp;
434 : :
435 [ # # ]: 0 : status = tri_tri_intersect(tri1,tri2);
436 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
437 : : {
438 : 0 : return status;
439 : : }
440 : :
441 : : // Restore values.
442 [ # # ]: 0 : poly1->verts[tri1->v2]->coord[0] = tz1[0];
443 [ # # ]: 0 : poly1->verts[tri1->v2]->coord[1] = tz1[1];
444 [ # # ]: 0 : poly1->verts[tri1->v2]->coord[2] = tz1[2];
445 : 0 : tri1->a = ta; tri1->b = tb; tri1->c = tc; tri1->d = td;
446 : :
447 [ # # ]: 0 : poly1->verts[tri1->v1]->coord[0] += ta;
448 [ # # ]: 0 : poly1->verts[tri1->v1]->coord[1] += tb;
449 [ # # ]: 0 : poly1->verts[tri1->v1]->coord[2] += tc;
450 : : // Compute new normal and linecoeff
451 [ # # ]: 0 : newplanecoefficients(poly1, tri1);
452 : 0 : linecoeff[0] = tri1->c*tri2->b - tri1->b*tri2->c;
453 : 0 : linecoeff[1] = tri1->a*tri2->c - tri1->c*tri2->a;
454 : 0 : linecoeff[2] = tri1->b*tri2->a - tri1->a*tri2->b;
455 : 0 : dtemp = sqrt(linecoeff[0]*linecoeff[0] +
456 : 0 : linecoeff[1]*linecoeff[1] +
457 : 0 : linecoeff[2]*linecoeff[2]);
458 : 0 : linecoeff[0] /= dtemp;
459 : 0 : linecoeff[1] /= dtemp;
460 : 0 : linecoeff[2] /= dtemp;
461 : :
462 [ # # ]: 0 : status = tri_tri_intersect(tri1,tri2);
463 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
464 : : {
465 : 0 : return status;
466 : : }
467 : :
468 : : // Restore values.
469 [ # # ]: 0 : poly1->verts[tri1->v1]->coord[0] = ty1[0];
470 [ # # ]: 0 : poly1->verts[tri1->v1]->coord[1] = ty1[1];
471 [ # # ]: 0 : poly1->verts[tri1->v1]->coord[2] = ty1[2];
472 : 0 : tri1->a = ta; tri1->b = tb; tri1->c = tc; tri1->d = td;
473 : :
474 [ # # ]: 0 : poly1->verts[tri1->v0]->coord[0] += ta;
475 [ # # ]: 0 : poly1->verts[tri1->v0]->coord[1] += tb;
476 [ # # ]: 0 : poly1->verts[tri1->v0]->coord[2] += tc;
477 : : // Compute new normal and linecoeff
478 [ # # ]: 0 : newplanecoefficients(poly1, tri1);
479 : 0 : linecoeff[0] = tri1->c*tri2->b - tri1->b*tri2->c;
480 : 0 : linecoeff[1] = tri1->a*tri2->c - tri1->c*tri2->a;
481 : 0 : linecoeff[2] = tri1->b*tri2->a - tri1->a*tri2->b;
482 : 0 : dtemp = sqrt(linecoeff[0]*linecoeff[0] +
483 : 0 : linecoeff[1]*linecoeff[1] +
484 : 0 : linecoeff[2]*linecoeff[2]);
485 : 0 : linecoeff[0] /= dtemp;
486 : 0 : linecoeff[1] /= dtemp;
487 : 0 : linecoeff[2] /= dtemp;
488 : :
489 [ # # ]: 0 : status = tri_tri_intersect(tri1,tri2);
490 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
491 : : {
492 : 0 : return status;
493 : : }
494 : :
495 : : // Restore values.
496 [ # # ]: 0 : poly1->verts[tri1->v0]->coord[0] = tx1[0];
497 [ # # ]: 0 : poly1->verts[tri1->v0]->coord[1] = tx1[1];
498 [ # # ]: 0 : poly1->verts[tri1->v0]->coord[2] = tx1[2];
499 : 0 : tri1->a = ta; tri1->b = tb; tri1->c = tc; tri1->d = td;
500 : :
501 : : // triangle 2
502 : 0 : ta = tri2->a; tb = tri2->b; tc = tri2->c; td = tri2->d;
503 [ # # ]: 0 : for ( k = 0; k < 3; k++ ) {
504 [ # # ]: 0 : tx1[k] = poly2->verts[tri2->v0]->coord[k];
505 [ # # ]: 0 : ty1[k] = poly2->verts[tri2->v1]->coord[k];
506 [ # # ]: 0 : tz1[k] = poly2->verts[tri2->v2]->coord[k];
507 : : }
508 : :
509 [ # # ]: 0 : poly2->verts[tri2->v2]->coord[0] += ta;
510 [ # # ]: 0 : poly2->verts[tri2->v2]->coord[1] += tb;
511 [ # # ]: 0 : poly2->verts[tri2->v2]->coord[2] += tc;
512 : : // Compute new normal and linecoeff
513 [ # # ]: 0 : newplanecoefficients(poly2, tri2);
514 : 0 : linecoeff[0] = tri1->c*tri2->b - tri1->b*tri2->c;
515 : 0 : linecoeff[1] = tri1->a*tri2->c - tri1->c*tri2->a;
516 : 0 : linecoeff[2] = tri1->b*tri2->a - tri1->a*tri2->b;
517 : 0 : dtemp = sqrt(linecoeff[0]*linecoeff[0] +
518 : 0 : linecoeff[1]*linecoeff[1] +
519 : 0 : linecoeff[2]*linecoeff[2]);
520 : 0 : linecoeff[0] /= dtemp;
521 : 0 : linecoeff[1] /= dtemp;
522 : 0 : linecoeff[2] /= dtemp;
523 : :
524 [ # # ]: 0 : status = tri_tri_intersect(tri1,tri2);
525 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
526 : : {
527 : 0 : return status;
528 : : }
529 : :
530 : : // Restore values.
531 [ # # ]: 0 : poly2->verts[tri2->v2]->coord[0] = tz1[0];
532 [ # # ]: 0 : poly2->verts[tri2->v2]->coord[1] = tz1[1];
533 [ # # ]: 0 : poly2->verts[tri2->v2]->coord[2] = tz1[2];
534 : 0 : tri2->a = ta; tri2->b = tb; tri2->c = tc; tri2->d = td;
535 : :
536 [ # # ]: 0 : poly2->verts[tri2->v1]->coord[0] += ta;
537 [ # # ]: 0 : poly2->verts[tri2->v1]->coord[1] += tb;
538 [ # # ]: 0 : poly2->verts[tri2->v1]->coord[2] += tc;
539 : : // Compute new normal and linecoeff
540 [ # # ]: 0 : newplanecoefficients(poly2, tri2);
541 : 0 : linecoeff[0] = tri1->c*tri2->b - tri1->b*tri2->c;
542 : 0 : linecoeff[1] = tri1->a*tri2->c - tri1->c*tri2->a;
543 : 0 : linecoeff[2] = tri1->b*tri2->a - tri1->a*tri2->b;
544 : 0 : dtemp = sqrt(linecoeff[0]*linecoeff[0] +
545 : 0 : linecoeff[1]*linecoeff[1] +
546 : 0 : linecoeff[2]*linecoeff[2]);
547 : 0 : linecoeff[0] /= dtemp;
548 : 0 : linecoeff[1] /= dtemp;
549 : 0 : linecoeff[2] /= dtemp;
550 : :
551 [ # # ]: 0 : status = tri_tri_intersect(tri1,tri2);
552 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
553 : : {
554 : 0 : return status;
555 : : }
556 : :
557 : : // Restore values.
558 [ # # ]: 0 : poly2->verts[tri2->v1]->coord[0] = ty1[0];
559 [ # # ]: 0 : poly2->verts[tri2->v1]->coord[1] = ty1[1];
560 [ # # ]: 0 : poly2->verts[tri2->v1]->coord[2] = ty1[2];
561 : 0 : tri2->a = ta; tri2->b = tb; tri2->c = tc; tri2->d = td;
562 : :
563 [ # # ]: 0 : poly2->verts[tri2->v0]->coord[0] += ta;
564 [ # # ]: 0 : poly2->verts[tri2->v0]->coord[1] += tb;
565 [ # # ]: 0 : poly2->verts[tri2->v0]->coord[2] += tc;
566 : : // Compute new normal and linecoeff
567 [ # # ]: 0 : newplanecoefficients(poly2, tri2);
568 : 0 : linecoeff[0] = tri1->c*tri2->b - tri1->b*tri2->c;
569 : 0 : linecoeff[1] = tri1->a*tri2->c - tri1->c*tri2->a;
570 : 0 : linecoeff[2] = tri1->b*tri2->a - tri1->a*tri2->b;
571 : 0 : dtemp = sqrt(linecoeff[0]*linecoeff[0] +
572 : 0 : linecoeff[1]*linecoeff[1] +
573 : 0 : linecoeff[2]*linecoeff[2]);
574 : 0 : linecoeff[0] /= dtemp;
575 : 0 : linecoeff[1] /= dtemp;
576 : 0 : linecoeff[2] /= dtemp;
577 : :
578 [ # # ]: 0 : status = tri_tri_intersect(tri1,tri2);
579 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
580 : : {
581 : 0 : return status;
582 : : }
583 : :
584 : : // Restore values.
585 [ # # ]: 0 : poly2->verts[tri2->v0]->coord[0] = tx1[0];
586 [ # # ]: 0 : poly2->verts[tri2->v0]->coord[1] = tx1[1];
587 [ # # ]: 0 : poly2->verts[tri2->v0]->coord[2] = tx1[2];
588 : 0 : tri2->a = ta; tri2->b = tb; tri2->c = tc; tri2->d = td;
589 : :
590 : 0 : continue;
591 : : }
592 : :
593 [ # # ]: 0 : if ( do_imprint == true ) {
594 : : // Don't intersect nearly-coplanar triangles.
595 [ # # ]: 0 : if ( fabs(tri1->a*tri2->a + tri1->b*tri2->b +tri1->c*tri2->c) > 0.8 ) continue;
596 : : }
597 : : double dtemp;
598 : 0 : dtemp = sqrt(linecoeff[0]*linecoeff[0] +
599 : 0 : linecoeff[1]*linecoeff[1] +
600 : 0 : linecoeff[2]*linecoeff[2]);
601 : 0 : linecoeff[0] /= dtemp;
602 : 0 : linecoeff[1] /= dtemp;
603 : 0 : linecoeff[2] /= dtemp;
604 : :
605 [ # # ]: 0 : status = tri_tri_intersect(tri1,tri2);
606 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
607 : : {
608 : 0 : return status;
609 : : }
610 : :
611 : : } // end of loop over poly2
612 : :
613 : : } // end of loop over poly1
614 : :
615 : : //delete entire array.
616 [ # # ]: 0 : delete [] boxlist;
617 : :
618 : 0 : return status;
619 : : }
620 : :
621 : 0 : CubitStatus FBIntersect::tri_tri_intersect(FB_Triangle *tri1,
622 : : FB_Triangle *tri2)
623 : : {
624 : : int ret1, ret2, i;
625 : : double xc10[3], xc11[3],xc12[3]; // coords for poly1 triangle
626 : : double xc20[3], xc21[3],xc22[3]; // coords for poly2 triangle
627 : : double d10, d11, d12, d20, d21, d22; // distance of vertex from plane of other triangle
628 : : double tt[4];
629 : : int edge_vert_type[4];
630 : : CubitStatus status;
631 : :
632 : 0 : int mydebug = 0;
633 : 0 : status = CUBIT_SUCCESS;
634 : :
635 : 0 : tt[0] = tt[1] = tt[2] = tt[3] = CUBIT_DBL_MAX;
636 : : // Is tri1 entirely on one side of tri2?
637 [ # # ]: 0 : for ( i = 0; i < 3; i++ ) {
638 [ # # ]: 0 : xc10[i] = poly1->verts[tri1->v0]->coord[i];
639 [ # # ]: 0 : xc11[i] = poly1->verts[tri1->v1]->coord[i];
640 [ # # ]: 0 : xc12[i] = poly1->verts[tri1->v2]->coord[i];
641 : : }
642 : :
643 : : // distance of each tri1 vert to plane of tri2
644 : 0 : d10 = xc10[0]*tri2->a + xc10[1]*tri2->b + xc10[2]*tri2->c + tri2->d;
645 : 0 : d11 = xc11[0]*tri2->a + xc11[1]*tri2->b + xc11[2]*tri2->c + tri2->d;
646 : 0 : d12 = xc12[0]*tri2->a + xc12[1]*tri2->b + xc12[2]*tri2->c + tri2->d;
647 [ # # ][ # # ]: 0 : if ( ( (d10 < -EPSILON2) && (d11 < -EPSILON2) && (d12 < -EPSILON2) ) ||
[ # # ][ # # ]
648 [ # # ][ # # ]: 0 : ( (d10 > EPSILON2) && (d11 > EPSILON2) && (d12 > EPSILON2) ) )
649 : 0 : return CUBIT_SUCCESS;
650 : : // Is tri2 entirely on one side of tri1?
651 [ # # ]: 0 : for ( i = 0; i < 3; i++ ) {
652 [ # # ]: 0 : xc20[i] = poly2->verts[tri2->v0]->coord[i];
653 [ # # ]: 0 : xc21[i] = poly2->verts[tri2->v1]->coord[i];
654 [ # # ]: 0 : xc22[i] = poly2->verts[tri2->v2]->coord[i];
655 : : }
656 : :
657 : : // distance of each tri2 vert to plane of tri1
658 : 0 : d20 = xc20[0]*tri1->a + xc20[1]*tri1->b + xc20[2]*tri1->c + tri1->d;
659 : 0 : d21 = xc21[0]*tri1->a + xc21[1]*tri1->b + xc21[2]*tri1->c + tri1->d;
660 : 0 : d22 = xc22[0]*tri1->a + xc22[1]*tri1->b + xc22[2]*tri1->c + tri1->d;
661 [ # # ][ # # ]: 0 : if ( ( (d20 < -EPSILON2) && (d21 < -EPSILON2) && (d22 < -EPSILON2) ) ||
[ # # ][ # # ]
662 [ # # ][ # # ]: 0 : ( (d20 > EPSILON2) && (d21 > EPSILON2) && (d22 > EPSILON2) ) )
663 : 0 : return CUBIT_SUCCESS;
664 : : // Get a point on the line of intersection to serve as a reference point.
665 : : double ta, tb, ts1, ts2, tdot11;
666 : 0 : ts1 = -tri1->d; ts2 = -tri2->d;
667 : :
668 : 0 : tdot11 = tri1->a*tri2->a + tri1->b*tri2->b + tri1->c*tri2->c;
669 : 0 : ta = (ts2*tdot11 - ts1)/(tdot11*tdot11 - 1);
670 : 0 : tb = (ts1*tdot11 - ts2)/(tdot11*tdot11 - 1);
671 : 0 : linept[0] = ta*tri1->a + tb*tri2->a;
672 : 0 : linept[1] = ta*tri1->b + tb*tri2->b;
673 : 0 : linept[2] = ta*tri1->c + tb*tri2->c;
674 : :
675 : : // There are several cases for the distances.
676 : : // ret1 holds the number of intersections of the triangle with the
677 : : // intersection line.
678 : : // Do tri1.
679 : : ret1 = get_intersectionline_parameter_values(d10,d11,d12,
680 : : xc10,xc11,xc12,
681 : : tt[0],tt[1],
682 [ # # ]: 0 : edge_vert_type[0],edge_vert_type[1]);
683 : : // Do tri2.
684 : : ret2 = get_intersectionline_parameter_values(d20,d21,d22,
685 : : xc20,xc21,xc22,
686 : : tt[2],tt[3],
687 [ # # ]: 0 : edge_vert_type[2],edge_vert_type[3]);
688 : : // If not two intersections for each triangle, no intersection edge exists.
689 [ # # ][ # # ]: 0 : if ( (ret1 == 2) && (ret2 == 2) )
690 : : {
691 [ # # ]: 0 : status = add_intersection_edges(tri1,tri2,tt,edge_vert_type);
692 [ # # ]: 0 : if ( status != CUBIT_SUCCESS )
693 : : {
694 : 0 : return status;
695 : : }
696 : : }
697 [ # # ]: 0 : if(mydebug){
698 : : double min_angle_1, max_angle_1;
699 : : double min_angle_2, max_angle_2;
700 [ # # ]: 0 : poly1->min_max_angles_in_fb_triangle(tri1, min_angle_1, max_angle_1);
701 [ # # ]: 0 : poly2->min_max_angles_in_fb_triangle(tri2, min_angle_2, max_angle_2);
702 : :
703 : :
704 [ # # ][ # # ]: 0 : if(min_angle_1 < 10 || min_angle_2 < 10){
705 [ # # ][ # # ]: 0 : PRINT_INFO(" Tri 1 min angle = %f\n",min_angle_1);
[ # # ][ # # ]
706 [ # # ][ # # ]: 0 : PRINT_INFO(" Tri 2 min angle = %f\n",min_angle_2);
[ # # ][ # # ]
707 : : }
708 : : }
709 : :
710 : 0 : return status;
711 : : }
712 : :
713 : 0 : CubitStatus FBIntersect::add_intersection_edges(FB_Triangle *tri1,
714 : : FB_Triangle *tri2,
715 : : double *tt,
716 : : int *edge_vert_type
717 : : )
718 : : {
719 : : int v10, v11, v20, v21;
720 : : bool ifoundit, exists;
721 : : FB_Edge *edge;
722 : : double x1pt, y1pt, z1pt, x2pt, y2pt, z2pt;
723 : : int edge1_vert0_type, edge1_vert1_type, edge2_vert0_type, edge2_vert1_type;
724 : :
725 : 0 : ifoundit = false;
726 : 0 : edge1_vert0_type = edge1_vert1_type = edge2_vert0_type = edge2_vert1_type = UNKNOWN;
727 : :
728 : : // Try to handle the epsilon cases by forcing tt[] values that are almost
729 : : // equal to be equal.
730 : :
731 [ # # ]: 0 : if ( fabs(tt[0]-tt[2]) < EPSILON ) tt[2] = tt[0];
732 [ # # ]: 0 : if ( fabs(tt[0]-tt[3]) < EPSILON ) tt[3] = tt[0];
733 [ # # ]: 0 : if ( fabs(tt[1]-tt[2]) < EPSILON ) tt[2] = tt[1];
734 [ # # ]: 0 : if ( fabs(tt[1]-tt[3]) < EPSILON ) tt[3] = tt[1];
735 : :
736 : : // cases 0 and 9, no overlap
737 : : // if ( (tt[1] < tt[2]) || (tt[3] < tt[0]) ) return CUBIT_SUCCESS;
738 [ # # ]: 0 : if ( tt[1] < tt[2] ) { return CUBIT_SUCCESS; }
739 : :
740 [ # # ]: 0 : if ( tt[3] < tt[0] ) { return CUBIT_SUCCESS; }
741 : :
742 : : // These next four cases are by far the most common forms of overlap,
743 : : // so check for them first.
744 : :
745 : : // The 1's (e.g.: 11111111) refer to the portion of the line of intersection
746 : : // that is in tri1; the 2's for tri2. The case numbers were arbirtarily assigned.
747 : :
748 : : // case 4
749 : : // 11111111
750 : : // 2222
751 [ # # ][ # # ]: 0 : if ( (tt[2] > tt[0]) && (tt[3] < tt[1]) ) {
752 : 0 : ifoundit = true;
753 [ # # ]: 0 : get_point_from_parameter(tt[2],&x1pt,&y1pt,&z1pt);
754 [ # # ]: 0 : get_point_from_parameter(tt[3],&x2pt,&y2pt,&z2pt);
755 : 0 : edge1_vert0_type = INTERIOR_VERT;
756 : 0 : edge1_vert1_type = INTERIOR_VERT;
757 : 0 : edge2_vert0_type = edge_vert_type[2];
758 : 0 : edge2_vert1_type = edge_vert_type[3];
759 : :
760 : : // case 2
761 : : // 11111111
762 : : // 22222222
763 [ # # ][ # # ]: 0 : } else if ( (tt[2] < tt[1]) && (tt[2] > tt[0]) && (tt[3] > tt[1]) ) {
[ # # ]
764 : 0 : ifoundit = true;
765 [ # # ]: 0 : get_point_from_parameter(tt[2],&x1pt,&y1pt,&z1pt);
766 [ # # ]: 0 : get_point_from_parameter(tt[1],&x2pt,&y2pt,&z2pt);
767 : 0 : edge1_vert0_type = INTERIOR_VERT;
768 : 0 : edge1_vert1_type = edge_vert_type[1];
769 : 0 : edge2_vert0_type = edge_vert_type[2];
770 : 0 : edge2_vert1_type = INTERIOR_VERT;
771 : :
772 : : // case 7
773 : : // 11111111
774 : : // 22222222
775 [ # # ][ # # ]: 0 : } else if ( (tt[0] < tt[3]) && (tt[0] > tt[2]) && (tt[1] > tt[3]) ) {
[ # # ]
776 : 0 : ifoundit = true;
777 [ # # ]: 0 : get_point_from_parameter(tt[0],&x1pt,&y1pt,&z1pt);
778 [ # # ]: 0 : get_point_from_parameter(tt[3],&x2pt,&y2pt,&z2pt);
779 : 0 : edge1_vert0_type = edge_vert_type[0];
780 : 0 : edge1_vert1_type = INTERIOR_VERT;
781 : 0 : edge2_vert0_type = INTERIOR_VERT;
782 : 0 : edge2_vert1_type = edge_vert_type[3];
783 : :
784 : : // case 12
785 : : // 1111
786 : : // 22222222
787 [ # # ][ # # ]: 0 : } else if ( (tt[0] > tt[2]) && (tt[1] < tt[3]) ) {
788 : 0 : ifoundit = true;
789 [ # # ]: 0 : get_point_from_parameter(tt[0],&x1pt,&y1pt,&z1pt);
790 [ # # ]: 0 : get_point_from_parameter(tt[1],&x2pt,&y2pt,&z2pt);
791 : 0 : edge1_vert0_type = edge_vert_type[0];
792 : 0 : edge1_vert1_type = edge_vert_type[1];
793 : 0 : edge2_vert0_type = INTERIOR_VERT;
794 : 0 : edge2_vert1_type = INTERIOR_VERT;
795 : :
796 : : // case 1
797 : : // 11111111
798 : : // 22222222
799 [ # # ]: 0 : } else if ( fabs(tt[1]-tt[2]) < EPSILON ) {
800 : 0 : ifoundit = true;
801 : : // return CUBIT_SUCCESS;
802 [ # # ]: 0 : get_point_from_parameter(tt[1],&x1pt,&y1pt,&z1pt);
803 [ # # ]: 0 : get_point_from_parameter(tt[1],&x2pt,&y2pt,&z2pt);
804 : 0 : edge1_vert0_type = edge_vert_type[1];
805 : 0 : edge1_vert1_type = edge_vert_type[1];
806 : 0 : edge2_vert0_type = edge_vert_type[2];
807 : 0 : edge2_vert1_type = edge_vert_type[2];
808 : :
809 : : // case 3
810 [ # # ]: 0 : } else if ( tt[0] < tt[2] ) {
811 : : // case 3
812 : : // 11111111
813 : : // 22222
814 [ # # ]: 0 : if ( tt[1] == tt[3] ) {
815 : 0 : ifoundit = true;
816 [ # # ]: 0 : get_point_from_parameter(tt[2],&x1pt,&y1pt,&z1pt);
817 [ # # ]: 0 : get_point_from_parameter(tt[3],&x2pt,&y2pt,&z2pt);
818 : 0 : edge1_vert0_type = INTERIOR_VERT;
819 : 0 : edge1_vert1_type = edge_vert_type[1];
820 : 0 : edge2_vert0_type = edge_vert_type[2];
821 : 0 : edge2_vert1_type = edge_vert_type[3];
822 : : }
823 : :
824 : : // case 5, 6, or 11
825 [ # # ]: 0 : } else if ( fabs(tt[0]-tt[2]) < EPSILON ) {
826 : :
827 : : // case 5
828 : : // 11111111
829 : : // 22222222
830 [ # # ]: 0 : if ( tt[1] == tt[3] ) {
831 : 0 : ifoundit = true;
832 [ # # ]: 0 : get_point_from_parameter(tt[0],&x1pt,&y1pt,&z1pt);
833 [ # # ]: 0 : get_point_from_parameter(tt[1],&x2pt,&y2pt,&z2pt);
834 : 0 : edge1_vert0_type = edge_vert_type[0];
835 : 0 : edge1_vert1_type = edge_vert_type[1];
836 : 0 : edge2_vert0_type = edge_vert_type[2];
837 : 0 : edge2_vert1_type = edge_vert_type[3];
838 : :
839 : : // case 6
840 : : // 11111111
841 : : // 2222
842 [ # # ]: 0 : } else if ( tt[1] > tt[3] ) {
843 : 0 : ifoundit = true;
844 [ # # ]: 0 : get_point_from_parameter(tt[2],&x1pt,&y1pt,&z1pt);
845 [ # # ]: 0 : get_point_from_parameter(tt[3],&x2pt,&y2pt,&z2pt);
846 : 0 : edge1_vert0_type = edge_vert_type[0];
847 : 0 : edge1_vert1_type = INTERIOR_VERT;
848 : 0 : edge2_vert0_type = edge_vert_type[2];
849 : 0 : edge2_vert1_type = edge_vert_type[3];
850 : :
851 : : // case 11
852 : : // 1111
853 : : // 22222222
854 : : } else {
855 : 0 : ifoundit = true;
856 [ # # ]: 0 : get_point_from_parameter(tt[0],&x1pt,&y1pt,&z1pt);
857 [ # # ]: 0 : get_point_from_parameter(tt[1],&x2pt,&y2pt,&z2pt);
858 : 0 : edge1_vert0_type = edge_vert_type[0];
859 : 0 : edge1_vert1_type = edge_vert_type[1];
860 : 0 : edge2_vert0_type = edge_vert_type[2];
861 : 0 : edge2_vert1_type = INTERIOR_VERT;
862 : : }
863 : :
864 : : // case 8 or 10
865 : : } else {
866 : :
867 : : // case 8
868 : : // 11111111
869 : : // 22222222
870 [ # # ]: 0 : if ( fabs(tt[0]-tt[3]) < EPSILON ) {
871 : 0 : ifoundit = true;
872 : : // return CUBIT_SUCCESS;
873 [ # # ]: 0 : get_point_from_parameter(tt[0],&x1pt,&y1pt,&z1pt);
874 [ # # ]: 0 : get_point_from_parameter(tt[0],&x2pt,&y2pt,&z2pt);
875 : 0 : edge1_vert0_type = edge_vert_type[0];
876 : 0 : edge1_vert1_type = edge_vert_type[0];
877 : 0 : edge2_vert0_type = edge_vert_type[3];
878 : 0 : edge2_vert1_type = edge_vert_type[3];
879 : :
880 : : // case 10
881 : : // 1111
882 : : // 22222222
883 [ # # ]: 0 : } else if ( fabs(tt[1]-tt[3]) < EPSILON ) {
884 : 0 : ifoundit = true;
885 [ # # ]: 0 : get_point_from_parameter(tt[0],&x1pt,&y1pt,&z1pt);
886 [ # # ]: 0 : get_point_from_parameter(tt[1],&x2pt,&y2pt,&z2pt);
887 : 0 : edge1_vert0_type = edge_vert_type[0];
888 : 0 : edge1_vert1_type = edge_vert_type[1];
889 : 0 : edge2_vert0_type = INTERIOR_VERT;
890 : 0 : edge2_vert1_type = edge_vert_type[3];
891 : : }
892 : :
893 : : }
894 : :
895 [ # # ]: 0 : if ( ifoundit == false ) {
896 [ # # ][ # # ]: 0 : PRINT_ERROR("unaccounted for case in add_intersection_edges: tt[] = %e %e %e %e\n",
[ # # ]
897 [ # # ]: 0 : tt[0],tt[1],tt[2],tt[3]);
898 : 0 : return CUBIT_FAILURE;
899 : : } else {
900 : 0 : tri1->dudded = true;
901 [ # # ]: 0 : v10 = poly1->addavertex(x1pt,y1pt,z1pt);
902 [ # # ]: 0 : v11 = poly1->addavertex(x2pt,y2pt,z2pt);
903 [ # # ]: 0 : if ( v10 != v11 ) {
904 [ # # ]: 0 : exists = poly1->edge_exists_in_tri(*tri1,v10,v11);
905 [ # # ]: 0 : if ( exists == false ) {
906 [ # # ]: 0 : if ( edge1_vert0_type == INTERIOR_VERT )
907 : : edge1_vert0_type = determine_edge_vert_type(edge_vert_type[0],
908 [ # # ]: 0 : edge_vert_type[1]);
909 [ # # ]: 0 : if ( edge1_vert1_type == INTERIOR_VERT )
910 : : edge1_vert1_type = determine_edge_vert_type(edge_vert_type[0],
911 [ # # ]: 0 : edge_vert_type[1]);
912 [ # # ][ # # ]: 0 : edge = new FB_Edge(v10,v11,edge1_vert0_type,edge1_vert1_type,true);
913 [ # # ]: 0 : tri1->edge_list.push_back(edge);
914 [ # # ][ # # ]: 0 : if ( poly1->edge_exists(v10,v11) == false )
915 [ # # ]: 0 : poly1->intersection_edges.push_back(edge);
916 : : }
917 : : } else {
918 : 0 : int edge_type = UNKNOWN;
919 : 0 : int vtype1 = UNKNOWN_VERT, vtype2 = UNKNOWN_VERT;
920 : 0 : int v_other1 = UNKNOWN_VERT, v_other2 = UNKNOWN_VERT;
921 : : // edge_type = UNKNOWN;
922 [ # # ]: 0 : if ( edge1_vert0_type == EDGE_2 ) {
923 : 0 : v_other1 = tri1->v0;
924 : 0 : v_other2 = tri1->v2;
925 : 0 : edge_type = EDGE_2;
926 : 0 : vtype1 = VERTEX_0;
927 : 0 : vtype2 = VERTEX_2;
928 [ # # ]: 0 : } else if ( edge1_vert0_type == EDGE_1 ) {
929 : 0 : v_other1 = tri1->v1;
930 : 0 : v_other2 = tri1->v2;
931 : 0 : edge_type = EDGE_1;
932 : 0 : vtype1 = VERTEX_1;
933 : 0 : vtype2 = VERTEX_2;
934 : :
935 [ # # ]: 0 : } else if ( edge1_vert0_type == EDGE_0 ) {
936 : 0 : v_other1 = tri1->v0;
937 : 0 : v_other2 = tri1->v1;
938 : 0 : edge_type = EDGE_0;
939 : 0 : vtype1 = VERTEX_0;
940 : 0 : vtype2 = VERTEX_1;
941 : :
942 : : }
943 [ # # ]: 0 : if ( edge_type != UNKNOWN ) {
944 [ # # ]: 0 : exists = poly1->edge_exists_in_tri(*tri1,v_other1,v10);
945 [ # # ]: 0 : if ( exists == false ) {
946 [ # # ][ # # ]: 0 : edge = new FB_Edge(v_other1,v10,vtype1,edge_type,false);
947 [ # # ]: 0 : tri1->edge_list.push_back(edge);
948 [ # # ][ # # ]: 0 : if ( poly1->edge_exists(v_other1,v10) == false )
949 [ # # ]: 0 : poly1->intersection_edges.push_back(edge);
950 : : }
951 [ # # ]: 0 : exists = poly1->edge_exists_in_tri(*tri1,v10,v_other2);
952 [ # # ]: 0 : if ( exists == false ) {
953 [ # # ][ # # ]: 0 : edge = new FB_Edge(v10,v_other2,edge_type,vtype2,false);
954 [ # # ]: 0 : tri1->edge_list.push_back(edge);
955 [ # # ][ # # ]: 0 : if ( poly1->edge_exists(v10,v_other2) == false )
956 [ # # ]: 0 : poly1->intersection_edges.push_back(edge);
957 : : }
958 : : }
959 : : }
960 : :
961 : 0 : tri2->dudded = true;
962 [ # # ]: 0 : v20 = poly2->addavertex(x1pt,y1pt,z1pt);
963 [ # # ]: 0 : v21 = poly2->addavertex(x2pt,y2pt,z2pt);
964 [ # # ]: 0 : if ( v20 != v21) {
965 [ # # ]: 0 : exists = poly2->edge_exists_in_tri(*tri2,v20,v21);
966 [ # # ]: 0 : if ( exists == false ) {
967 [ # # ]: 0 : if ( edge2_vert0_type == INTERIOR_VERT )
968 : 0 : edge2_vert0_type = determine_edge_vert_type(edge_vert_type[2],
969 [ # # ]: 0 : edge_vert_type[3]);
970 [ # # ]: 0 : if ( edge2_vert1_type == INTERIOR_VERT )
971 : 0 : edge2_vert1_type = determine_edge_vert_type(edge_vert_type[2],
972 [ # # ]: 0 : edge_vert_type[3]);
973 [ # # ][ # # ]: 0 : edge = new FB_Edge(v20,v21,edge2_vert0_type,edge2_vert1_type,true);
974 [ # # ]: 0 : tri2->edge_list.push_back(edge);
975 [ # # ][ # # ]: 0 : if ( poly2->edge_exists(v20,v21) == false )
976 [ # # ]: 0 : poly2->intersection_edges.push_back(edge);
977 : : }
978 : : } else {
979 : 0 : int edge_type = UNKNOWN;
980 : 0 : int vtype1 = UNKNOWN_VERT, vtype2 = UNKNOWN_VERT;
981 : 0 : int v_other1 = UNKNOWN_VERT, v_other2 = UNKNOWN_VERT;
982 : : // edge_type = UNKNOWN;
983 [ # # ]: 0 : if ( edge2_vert0_type == EDGE_2 ) {
984 : 0 : v_other1 = tri2->v0;
985 : 0 : v_other2 = tri2->v2;
986 : 0 : edge_type = EDGE_2;
987 : 0 : vtype1 = VERTEX_0;
988 : 0 : vtype2 = VERTEX_2;
989 [ # # ]: 0 : } else if ( edge2_vert0_type == EDGE_1 ) {
990 : 0 : v_other1 = tri2->v1;
991 : 0 : v_other2 = tri2->v2;
992 : 0 : edge_type = EDGE_1;
993 : 0 : vtype1 = VERTEX_1;
994 : 0 : vtype2 = VERTEX_2;
995 : :
996 [ # # ]: 0 : } else if ( edge2_vert0_type == EDGE_0 ) {
997 : 0 : v_other1 = tri2->v0;
998 : 0 : v_other2 = tri2->v1;
999 : 0 : edge_type = EDGE_0;
1000 : 0 : vtype1 = VERTEX_0;
1001 : 0 : vtype2 = VERTEX_1;
1002 : :
1003 : : }
1004 [ # # ]: 0 : if ( edge_type != UNKNOWN ) {
1005 [ # # ]: 0 : exists = poly2->edge_exists_in_tri(*tri2,v_other1,v20);
1006 [ # # ]: 0 : if ( exists == false ) {
1007 [ # # ][ # # ]: 0 : edge = new FB_Edge(v_other1,v20,vtype1,edge_type,false);
1008 [ # # ]: 0 : tri2->edge_list.push_back(edge);
1009 [ # # ][ # # ]: 0 : if ( poly2->edge_exists(v_other1,v20) == false )
1010 [ # # ]: 0 : poly2->intersection_edges.push_back(edge);
1011 : : }
1012 [ # # ]: 0 : exists = poly2->edge_exists_in_tri(*tri2,v20,v_other2);
1013 [ # # ]: 0 : if ( exists == false ) {
1014 [ # # ][ # # ]: 0 : edge = new FB_Edge(v20,v_other2,edge_type,vtype2,false);
1015 [ # # ]: 0 : tri2->edge_list.push_back(edge);
1016 [ # # ][ # # ]: 0 : if ( poly2->edge_exists(v20,v_other2) == false )
1017 [ # # ]: 0 : poly2->intersection_edges.push_back(edge);
1018 : : }
1019 : : }
1020 : : }
1021 : : }
1022 : :
1023 : 0 : return CUBIT_SUCCESS;
1024 : :
1025 : : }
1026 : :
1027 : 0 : void FBIntersect::get_point_from_parameter(double parameter,
1028 : : double *x, double *y, double *z)
1029 : : {
1030 : 0 : *x = linept[0] + linecoeff[0]*parameter;
1031 : 0 : *y = linept[1] + linecoeff[1]*parameter;
1032 : 0 : *z = linept[2] + linecoeff[2]*parameter;
1033 : :
1034 : 0 : }
1035 : :
1036 : 0 : double FBIntersect::get_distance_parameter(double *xc0,
1037 : : double *xc1,
1038 : : double d0, double d1)
1039 : : {
1040 : : double v1dot, v2dot;
1041 : :
1042 : : // dot the coordinates onto the line.
1043 : : // Assume that it has been checked already that d0 != d1.
1044 : :
1045 : 0 : v1dot = linecoeff[0]*(xc0[0] - linept[0]) +
1046 : 0 : linecoeff[1]*(xc0[1] - linept[1]) +
1047 : 0 : linecoeff[2]*(xc0[2] - linept[2]);
1048 : 0 : v2dot = linecoeff[0]*(xc1[0] - linept[0]) +
1049 : 0 : linecoeff[1]*(xc1[1] - linept[1]) +
1050 : 0 : linecoeff[2]*(xc1[2] - linept[2]);
1051 : :
1052 : 0 : return v1dot + (v2dot - v1dot)*d0/(d0-d1);
1053 : : }
1054 : :
1055 : :
1056 : 0 : double FBIntersect::get_distance_parameter_single(double *xc)
1057 : : {
1058 : : double coeff;
1059 : : int coord;
1060 : :
1061 [ # # ][ # # ]: 0 : if ( (fabs(linecoeff[0]) >= fabs(linecoeff[1])) &&
1062 : 0 : (fabs(linecoeff[0]) >= fabs(linecoeff[2])) ){
1063 : 0 : coord = 0;
1064 [ # # ]: 0 : } else if ( fabs(linecoeff[1]) >= fabs(linecoeff[2]) ) {
1065 : 0 : coord = 1;
1066 : : } else {
1067 : 0 : coord = 2;
1068 : : }
1069 : 0 : coeff = linecoeff[coord];
1070 : 0 : return (xc[coord] - linept[coord])/coeff;
1071 : : }
1072 : :
1073 : 0 : int FBIntersect::get_intersectionline_parameter_values(
1074 : : double d0,
1075 : : double d1,
1076 : : double d2,
1077 : : double *pt0,
1078 : : double *pt1,
1079 : : double *pt2,
1080 : : double& t0,
1081 : : double& t1,
1082 : : int& vert_type_0,
1083 : : int& vert_type_1)
1084 : : {
1085 : 0 : int ret = 0;
1086 : : // ret holds the nnumber of intersections of the line with the triangle.
1087 : :
1088 [ # # ]: 0 : if ( fabs(d0) < EPSILON ) {
1089 [ # # ]: 0 : if ( fabs(d1) < EPSILON ) {
1090 [ # # ]: 0 : if ( fabs(d2) < EPSILON ) {
1091 : : // coplanar
1092 : 0 : ret = 0;
1093 : : } else {
1094 : : // d0 = d1 = 0
1095 : 0 : t0 = get_distance_parameter_single(pt0);
1096 : 0 : t1 = get_distance_parameter_single(pt1);
1097 : 0 : vert_type_0 = VERTEX_0;
1098 : 0 : vert_type_1 = VERTEX_1;
1099 : 0 : ret = 2;
1100 : : }
1101 [ # # ]: 0 : } else if ( fabs(d2) < EPSILON ) {
1102 : : // d0 = d2 = 0
1103 : 0 : t0 = get_distance_parameter_single(pt0);
1104 : 0 : t1 = get_distance_parameter_single(pt2);
1105 : 0 : vert_type_0 = VERTEX_0;
1106 : 0 : vert_type_1 = VERTEX_2;
1107 : 0 : ret = 2;
1108 [ # # ]: 0 : } else if ( d1*d2 < 0.0 ) {
1109 : : // d0 = 0 and edge 12 crosses
1110 : 0 : t0 = get_distance_parameter_single(pt0);
1111 : 0 : t1 = get_distance_parameter(pt1,pt2,d1,d2);
1112 : 0 : vert_type_0 = VERTEX_0;
1113 : 0 : vert_type_1 = EDGE_1;
1114 : 0 : ret = 2;
1115 : : } else {
1116 : : // d0 = 0 and no edges cross
1117 : 0 : ret = 0;
1118 : : }
1119 [ # # ]: 0 : } else if ( fabs(d1) < EPSILON ) {
1120 [ # # ]: 0 : if ( fabs(d2) < EPSILON ) {
1121 : : // d1 = d2 = 0
1122 : 0 : t0 = get_distance_parameter_single(pt1);
1123 : 0 : t1 = get_distance_parameter_single(pt2);
1124 : 0 : vert_type_0 = VERTEX_1;
1125 : 0 : vert_type_1 = VERTEX_2;
1126 : 0 : ret = 2;
1127 [ # # ]: 0 : } else if ( d0*d2 < 0.0 ) {
1128 : : // d1 = 0 and edge 20 crosses
1129 : 0 : t0 = get_distance_parameter_single(pt1);
1130 : 0 : t1 = get_distance_parameter(pt0,pt2,d0,d2);
1131 : 0 : vert_type_0 = VERTEX_1;
1132 : 0 : vert_type_1 = EDGE_2;
1133 : 0 : ret = 2;
1134 : : } else {
1135 : : // d1 = 0 and no edges cross
1136 : 0 : ret = 0;
1137 : : }
1138 [ # # ]: 0 : } else if ( fabs(d2) < EPSILON ) {
1139 [ # # ]: 0 : if ( d0*d1 < 0.0 ) {
1140 : : // d2 = 0 and edge 01 crosses
1141 : 0 : t0 = get_distance_parameter_single(pt2);
1142 : 0 : t1 = get_distance_parameter(pt0,pt1,d0,d1);
1143 : 0 : vert_type_0 = VERTEX_2;
1144 : 0 : vert_type_1 = EDGE_0;
1145 : 0 : ret = 2;
1146 : : } else {
1147 : : // d2 = 0 and no edges cross
1148 : 0 : ret = 0;
1149 : : }
1150 [ # # ]: 0 : } else if ( d0*d1 < 0.0 ) {
1151 [ # # ]: 0 : if ( d0*d2 < 0.0 ) {
1152 : : // edges 01 and 02 cross
1153 : 0 : t0 = get_distance_parameter(pt0,pt1,d0,d1);
1154 : 0 : t1 = get_distance_parameter(pt0,pt2,d0,d2);
1155 : 0 : vert_type_0 = EDGE_0;
1156 : 0 : vert_type_1 = EDGE_2;
1157 : 0 : ret = 2;
1158 : : } else {
1159 : : // edges 01 and 12 cross
1160 : 0 : t0 = get_distance_parameter(pt0,pt1,d0,d1);
1161 : 0 : t1 = get_distance_parameter(pt1,pt2,d1,d2);
1162 : 0 : vert_type_0 = EDGE_0;
1163 : 0 : vert_type_1 = EDGE_1;
1164 : 0 : ret = 2;
1165 : : }
1166 : : } else {
1167 : : // edges 02 and 12 cross
1168 : 0 : t0 = get_distance_parameter(pt0,pt2,d0,d2);
1169 : 0 : t1 = get_distance_parameter(pt1,pt2,d1,d2);
1170 : 0 : vert_type_0 = EDGE_2;
1171 : 0 : vert_type_1 = EDGE_1;
1172 : 0 : ret = 2;
1173 : : }
1174 : :
1175 [ # # ]: 0 : if ( ret == 2 ) { // Sort them if necessary.
1176 [ # # ]: 0 : if ( t0 > t1 ) {
1177 : : int itemp;
1178 : : double dtemp;
1179 : 0 : dtemp = t0;
1180 : 0 : t0 = t1;
1181 : 0 : t1 = dtemp;
1182 : 0 : itemp = vert_type_0;
1183 : 0 : vert_type_0 = vert_type_1;
1184 : 0 : vert_type_1 = itemp;
1185 : : }
1186 : :
1187 : : }
1188 : :
1189 : 0 : return ret;
1190 : : }
1191 : :
1192 : 0 : void FBIntersect::set_classify_flag(bool value)
1193 : : {
1194 : 0 : do_classify = value;
1195 : 0 : }
1196 : :
1197 : 0 : CubitStatus FBIntersect::get_persistent_entity_info(bool *surfs_in,
1198 : : bool *curves_in,
1199 : : bool *surfs_out,
1200 : : bool *curves_out,
1201 : : const CubitFacetboolOp op,
1202 : : const int whichparent
1203 : : )
1204 : : {
1205 : : std::vector<int> *group, *groupcharacterization;
1206 [ # # ][ # # ]: 0 : std::vector<int>::iterator it, ig;
1207 : : int booltest1, booltest2, booltest1a, booltest2a, booltest_in, booltest_out;
1208 : : FBPolyhedron *poly;
1209 : :
1210 [ # # ]: 0 : if ( nothing_intersected == true ) return CUBIT_SUCCESS;
1211 [ # # ][ # # ]: 0 : if ( !poly1 || !poly2 ) {
1212 [ # # ][ # # ]: 0 : PRINT_ERROR("Error: Objects for Booleans must first be created.\n");
[ # # ][ # # ]
1213 : 0 : return CUBIT_FAILURE;
1214 : : }
1215 [ # # ][ # # ]: 0 : if ( !classify1 || !classify2 ) {
1216 [ # # ][ # # ]: 0 : PRINT_ERROR("Error: Objects for Booleans must first be classified.\n");
[ # # ][ # # ]
1217 : 0 : return CUBIT_FAILURE;
1218 : : }
1219 [ # # ][ # # ]: 0 : if ( (whichparent < 1) || (whichparent > 2) ) {
1220 [ # # ][ # # ]: 0 : PRINT_ERROR("Error: Requested nonexistent object.\n");
[ # # ][ # # ]
1221 : 0 : return CUBIT_FAILURE;
1222 : : }
1223 [ # # ]: 0 : if ( op == CUBIT_FB_UNION ) {
1224 : 0 : booltest1 = FB_ORIENTATION_INSIDE + FB_ORIENTATION_SAME;
1225 : 0 : booltest2 = FB_ORIENTATION_INSIDE;
1226 : 0 : booltest1a = FB_ORIENTATION_OUTSIDE + FB_ORIENTATION_SAME;
1227 : 0 : booltest2a = FB_ORIENTATION_OUTSIDE;
1228 [ # # ]: 0 : } else if ( op == CUBIT_FB_INTERSECTION ) {
1229 : 0 : booltest1 = FB_ORIENTATION_INSIDE + FB_ORIENTATION_SAME;
1230 : 0 : booltest2 = FB_ORIENTATION_INSIDE;
1231 : 0 : booltest1a = FB_ORIENTATION_OUTSIDE+ FB_ORIENTATION_OPPOSITE;
1232 : 0 : booltest2a = FB_ORIENTATION_INSIDE;
1233 [ # # ]: 0 : } else if ( op == CUBIT_FB_SUBTRACTION ) {
1234 : 0 : booltest1 = FB_ORIENTATION_OUTSIDE+ FB_ORIENTATION_OPPOSITE;
1235 : 0 : booltest2 = FB_ORIENTATION_INSIDE;
1236 : 0 : booltest1a = FB_ORIENTATION_INSIDE + FB_ORIENTATION_SAME;
1237 : 0 : booltest2a = FB_ORIENTATION_INSIDE;
1238 : : } else {
1239 [ # # ][ # # ]: 0 : PRINT_ERROR("Error: Unrecognized Boolean operation.\n");
[ # # ][ # # ]
1240 : 0 : return CUBIT_FAILURE;
1241 : : }
1242 : :
1243 [ # # ]: 0 : if ( whichparent == 1 ) {
1244 : 0 : booltest_in = booltest1;
1245 : 0 : booltest_out = booltest1a;
1246 [ # # ]: 0 : classify1->get_group(&group, &groupcharacterization);
1247 : 0 : poly = poly1;
1248 : : } else {
1249 : 0 : booltest_in = booltest2;
1250 : 0 : booltest_out = booltest2a;
1251 [ # # ]: 0 : classify2->get_group(&group, &groupcharacterization);
1252 : 0 : poly = poly2;
1253 : : }
1254 : : unsigned int i;
1255 : : int isurf, icurve0, icurve1, icurve2;
1256 : :
1257 [ # # ]: 0 : it = group->begin();
1258 [ # # ]: 0 : ig = groupcharacterization->begin();
1259 : :
1260 : 0 : i = 0;
1261 [ # # ][ # # ]: 0 : while ( it != group->end() ) {
[ # # ]
1262 [ # # ][ # # ]: 0 : if ( ( *(ig + *it) & booltest_in ) != 0 ) {
[ # # ][ # # ]
1263 : : // vtx = poly1->tris[i]->v0;
1264 [ # # ]: 0 : isurf = poly->tris[i]->cubitsurfaceindex;
1265 [ # # ]: 0 : if ( isurf > 0 )
1266 : 0 : surfs_in[isurf] = true;
1267 [ # # ]: 0 : icurve0 = poly->tris[i]->cubitedge0index;
1268 [ # # ]: 0 : if ( icurve0 > 0 )
1269 : 0 : curves_in[icurve0] = true;
1270 [ # # ]: 0 : icurve1 = poly->tris[i]->cubitedge1index;
1271 [ # # ]: 0 : if ( icurve1 > 0 )
1272 : 0 : curves_in[icurve1] = true;
1273 [ # # ]: 0 : icurve2 = poly->tris[i]->cubitedge2index;
1274 [ # # ]: 0 : if ( icurve2 > 0 )
1275 : 0 : curves_in[icurve2] = true;
1276 : : }
1277 [ # # ]: 0 : it++;
1278 : 0 : i++;
1279 : : }
1280 : :
1281 [ # # ]: 0 : it = group->begin();
1282 [ # # ]: 0 : ig = groupcharacterization->begin();
1283 : :
1284 : 0 : i = 0;
1285 [ # # ][ # # ]: 0 : while ( it != group->end() ) {
[ # # ]
1286 [ # # ][ # # ]: 0 : if ( ( *(ig + *it) & booltest_out ) != 0 ) {
[ # # ][ # # ]
1287 [ # # ]: 0 : isurf = poly->tris[i]->cubitsurfaceindex;
1288 [ # # ]: 0 : if ( isurf > 0 )
1289 : 0 : surfs_out[isurf] = true;
1290 [ # # ]: 0 : icurve0 = poly->tris[i]->cubitedge0index;
1291 [ # # ]: 0 : if ( icurve0 > 0 )
1292 : 0 : curves_out[icurve0] = true;
1293 [ # # ]: 0 : icurve1 = poly->tris[i]->cubitedge1index;
1294 [ # # ]: 0 : if ( icurve1 > 0 )
1295 : 0 : curves_out[icurve1] = true;
1296 [ # # ]: 0 : icurve2 = poly->tris[i]->cubitedge2index;
1297 [ # # ]: 0 : if ( icurve2 > 0 )
1298 : 0 : curves_out[icurve2] = true;
1299 : : }
1300 [ # # ]: 0 : it++;
1301 : 0 : i++;
1302 : : }
1303 : :
1304 : 0 : return CUBIT_SUCCESS;
1305 : :
1306 : : }
1307 : :
1308 : 0 : CubitStatus FBIntersect::update_surfs_and_curves(std::vector<double>& out_coords,
1309 : : std::vector<int>& out_connections,
1310 : : std::vector<int> *out_surf_index,
1311 : : std::vector<int> *out_curve_index,
1312 : : const int whichone)
1313 : : {
1314 : : unsigned int i;
1315 : : FBPolyhedron *poly;
1316 : :
1317 [ # # ]: 0 : if ( whichone == 1 ) poly = poly1;
1318 : 0 : else poly = poly2;
1319 [ # # ]: 0 : for ( i = 0; i < poly->verts.size(); i++ ) {
1320 : 0 : out_coords.push_back(poly->verts[i]->coord[0]);
1321 : 0 : out_coords.push_back(poly->verts[i]->coord[1]);
1322 : 0 : out_coords.push_back(poly->verts[i]->coord[2]);
1323 : : }
1324 [ # # ]: 0 : for ( i = 0; i < poly->tris.size(); i++ ) {
1325 [ # # ]: 0 : if ( poly->tris[i]->dudded == true ) continue;
1326 : 0 : out_connections.push_back(poly->tris[i]->v0);
1327 : 0 : out_connections.push_back(poly->tris[i]->v1);
1328 : 0 : out_connections.push_back(poly->tris[i]->v2);
1329 : 0 : out_surf_index->push_back(poly->tris[i]->cubitsurfaceindex);
1330 : 0 : out_curve_index->push_back(poly->tris[i]->cubitedge0index);
1331 : 0 : out_curve_index->push_back(poly->tris[i]->cubitedge1index);
1332 : 0 : out_curve_index->push_back(poly->tris[i]->cubitedge2index);
1333 : : }
1334 : 0 : return CUBIT_SUCCESS;
1335 : : }
1336 : :
1337 : 0 : CubitStatus FBIntersect::gather_by_boolean(std::vector<double>& out_coords,
1338 : : std::vector<int>& out_connections,
1339 : : std::vector<int> *out_surf_index,
1340 : : std::vector<int> *out_curve_index,
1341 : : std::vector<bool> *is_body_1,
1342 : : const CubitFacetboolOp op)
1343 : : {
1344 : : CubitStatus status;
1345 : : std::vector<int> *group1,
1346 : : *group2,
1347 : : *groupcharacterization1,
1348 : : *groupcharacterization2;
1349 [ # # ][ # # ]: 0 : std::vector<int>::iterator it, ig;
1350 : : int booltest1, booltest2;
1351 : :
1352 [ # # ][ # # ]: 0 : if ( nothing_intersected == true && op !=CUBIT_FB_UNION)
1353 : : {
1354 : 0 : return CUBIT_SUCCESS;
1355 : : }
1356 [ # # ][ # # ]: 0 : if ( !poly1 || !poly2 )
1357 : : {
1358 [ # # ][ # # ]: 0 : PRINT_ERROR("Error: Objects for Booleans must first be created.\n");
[ # # ][ # # ]
1359 : 0 : return CUBIT_FAILURE;
1360 : : }
1361 [ # # ][ # # ]: 0 : if ( !classify1 || !classify2 )
1362 : : {
1363 [ # # ][ # # ]: 0 : PRINT_ERROR("Error: Objects for Booleans must first be classified.\n");
[ # # ][ # # ]
1364 : 0 : return CUBIT_FAILURE;
1365 : : }
1366 : :
1367 [ # # ]: 0 : if ( op == CUBIT_FB_UNION )
1368 : : {
1369 : 0 : booltest1 = FB_ORIENTATION_OUTSIDE + FB_ORIENTATION_SAME;
1370 : 0 : booltest2 = FB_ORIENTATION_OUTSIDE;
1371 : : }
1372 [ # # ]: 0 : else if ( op == CUBIT_FB_INTERSECTION )
1373 : : {
1374 : 0 : booltest1 = FB_ORIENTATION_INSIDE + FB_ORIENTATION_SAME;
1375 : 0 : booltest2 = FB_ORIENTATION_INSIDE;
1376 : : }
1377 [ # # ]: 0 : else if ( op == CUBIT_FB_SUBTRACTION )
1378 : : {
1379 : 0 : booltest1 = FB_ORIENTATION_OUTSIDE+ FB_ORIENTATION_OPPOSITE;
1380 : 0 : booltest2 = FB_ORIENTATION_INSIDE;
1381 : : }
1382 : : else
1383 : : {
1384 [ # # ][ # # ]: 0 : PRINT_ERROR("Error: Unrecognized Boolean operation.\n");
[ # # ][ # # ]
1385 : 0 : return CUBIT_FAILURE;
1386 : : }
1387 : :
1388 [ # # ]: 0 : classify1->get_group(&group1, &groupcharacterization1);
1389 [ # # ]: 0 : classify2->get_group(&group2, &groupcharacterization2);
1390 : :
1391 [ # # ][ # # ]: 0 : IntegerHash *hashobj = new IntegerHash(NUMHASHBINS,20);
1392 : :
1393 [ # # ]: 0 : it = group1->begin();
1394 [ # # ]: 0 : ig = groupcharacterization1->begin();
1395 : : unsigned int i;
1396 : : int vtx, vertnum1, vertnum2, vertnum3, verts_sofar;
1397 : :
1398 : 0 : i = 0;
1399 : 0 : verts_sofar = 0;
1400 [ # # ][ # # ]: 0 : while ( it != group1->end() )
[ # # ]
1401 : : {
1402 [ # # ][ # # ]: 0 : if ( ( *(ig + *it) & booltest1 ) != 0 )
[ # # ][ # # ]
1403 : : {
1404 [ # # ]: 0 : vtx = poly1->tris[i]->v0;
1405 [ # # ]: 0 : vertnum1 = get_vertex(poly1, vtx, hashobj, out_coords, verts_sofar);
1406 : :
1407 [ # # ]: 0 : vtx = poly1->tris[i]->v1;
1408 [ # # ]: 0 : vertnum2 = get_vertex(poly1, vtx, hashobj, out_coords, verts_sofar);
1409 : :
1410 [ # # ]: 0 : vtx = poly1->tris[i]->v2;
1411 [ # # ]: 0 : vertnum3 = get_vertex(poly1, vtx, hashobj, out_coords, verts_sofar);
1412 : :
1413 [ # # ]: 0 : if ( store_connectivity( out_connections,
1414 : : vertnum1,
1415 : : vertnum2,
1416 [ # # ]: 0 : vertnum3 ) != CUBIT_SUCCESS )
1417 : : {
1418 : 0 : return CUBIT_FAILURE;
1419 : : }
1420 [ # # ]: 0 : if ( out_surf_index )
1421 : : {
1422 [ # # ][ # # ]: 0 : out_surf_index->push_back(poly1->tris[i]->cubitsurfaceindex);
1423 : : }
1424 [ # # ]: 0 : if ( out_curve_index )
1425 : : {
1426 [ # # ][ # # ]: 0 : out_curve_index->push_back(poly1->tris[i]->cubitedge0index);
1427 [ # # ][ # # ]: 0 : out_curve_index->push_back(poly1->tris[i]->cubitedge1index);
1428 [ # # ][ # # ]: 0 : out_curve_index->push_back(poly1->tris[i]->cubitedge2index);
1429 : : }
1430 [ # # ][ # # ]: 0 : if ( is_body_1 ) is_body_1->push_back(true);
1431 : : }
1432 [ # # ]: 0 : it++;
1433 : 0 : i++;
1434 : : }
1435 [ # # ]: 0 : it = group2->begin();
1436 [ # # ]: 0 : ig = groupcharacterization2->begin();
1437 : 0 : i = 0;
1438 [ # # ][ # # ]: 0 : while ( it != group2->end() )
[ # # ]
1439 : : {
1440 [ # # ][ # # ]: 0 : if ( ( *(ig + *it) & booltest2 ) != 0 )
[ # # ][ # # ]
1441 : : {
1442 [ # # ]: 0 : vtx = poly2->tris[i]->v0;
1443 [ # # ]: 0 : vertnum1 = get_vertex(poly2, vtx, hashobj, out_coords, verts_sofar);
1444 : :
1445 [ # # ]: 0 : vtx = poly2->tris[i]->v1;
1446 [ # # ]: 0 : vertnum2 = get_vertex(poly2, vtx, hashobj, out_coords, verts_sofar);
1447 : :
1448 [ # # ]: 0 : vtx = poly2->tris[i]->v2;
1449 [ # # ]: 0 : vertnum3 = get_vertex(poly2, vtx, hashobj, out_coords, verts_sofar);
1450 : :
1451 [ # # ]: 0 : if ( op == CUBIT_FB_SUBTRACTION ) // reverse the winding
1452 : : {
1453 [ # # ]: 0 : if ( store_connectivity( out_connections,
1454 : : vertnum1,
1455 : : vertnum3,
1456 [ # # ]: 0 : vertnum2 ) != CUBIT_SUCCESS )
1457 : : {
1458 : 0 : return CUBIT_FAILURE;
1459 : : }
1460 [ # # ]: 0 : if ( out_curve_index )
1461 : : {
1462 [ # # ][ # # ]: 0 : out_curve_index->push_back(poly2->tris[i]->cubitedge2index);
1463 [ # # ][ # # ]: 0 : out_curve_index->push_back(poly2->tris[i]->cubitedge1index);
1464 [ # # ][ # # ]: 0 : out_curve_index->push_back(poly2->tris[i]->cubitedge0index);
1465 : : }
1466 : : }
1467 : : else
1468 : : {
1469 [ # # ]: 0 : if ( store_connectivity( out_connections,
1470 : : vertnum1,
1471 : : vertnum2,
1472 [ # # ]: 0 : vertnum3 ) != CUBIT_SUCCESS )
1473 : : {
1474 : 0 : return CUBIT_FAILURE;
1475 : : }
1476 [ # # ]: 0 : if ( out_curve_index )
1477 : : {
1478 [ # # ][ # # ]: 0 : out_curve_index->push_back(poly2->tris[i]->cubitedge0index);
1479 [ # # ][ # # ]: 0 : out_curve_index->push_back(poly2->tris[i]->cubitedge1index);
1480 [ # # ][ # # ]: 0 : out_curve_index->push_back(poly2->tris[i]->cubitedge2index);
1481 : : }
1482 : : }
1483 [ # # ][ # # ]: 0 : if ( out_surf_index ) out_surf_index->push_back(poly2->tris[i]->cubitsurfaceindex);
[ # # ]
1484 [ # # ][ # # ]: 0 : if ( is_body_1 ) is_body_1->push_back(false);
1485 : : }
1486 [ # # ]: 0 : it++;
1487 : 0 : i++;
1488 : : }
1489 : : // GfxDebug::clear();
1490 : : // poly1->debug_draw_boundary_edges(CUBIT_MAGENTA_INDEX);
1491 : : // GfxDebug::mouse_xforms();
1492 : : // GfxDebug::clear();
1493 : : // poly2->debug_draw_boundary_edges(CUBIT_GREEN_INDEX);
1494 [ # # ][ # # ]: 0 : delete hashobj;
1495 : :
1496 : 0 : status = CUBIT_SUCCESS;
1497 : 0 : return status;
1498 : : }
1499 : :
1500 : 0 : CubitStatus FBIntersect::store_connectivity
1501 : : (
1502 : : std::vector<int>& out_connections,
1503 : : int vertnum1,
1504 : : int vertnum2,
1505 : : int vertnum3
1506 : : )
1507 : : {
1508 [ # # ][ # # ]: 0 : if ( vertnum1 == vertnum2 ||
1509 [ # # ]: 0 : vertnum2 == vertnum3 ||
1510 : 0 : vertnum1 == vertnum3 )
1511 : : {
1512 [ # # ][ # # ]: 0 : PRINT_ERROR( "Cannot continue without generating a degenerate facet.\n" );
1513 : 0 : return CUBIT_FAILURE;
1514 : : }
1515 : 0 : out_connections.push_back(vertnum1);
1516 : 0 : out_connections.push_back(vertnum2);
1517 : 0 : out_connections.push_back(vertnum3);
1518 : 0 : return CUBIT_SUCCESS;
1519 : : }
1520 : :
1521 : 0 : int FBIntersect::get_vertex(FBPolyhedron *poly, int vtx,
1522 : : IntegerHash *hashobj,
1523 : : std::vector<double>& out_coords,
1524 : : int &num_sofar)
1525 : : {
1526 : : double xx, yy, zz, xval, yval, zval;
1527 : : int i, hashvalue, *hasharrayptr, hasharraysize, hptr, ifoundit;
1528 : :
1529 [ # # ]: 0 : xx = poly->verts[vtx]->coord[0];
1530 [ # # ]: 0 : yy = poly->verts[vtx]->coord[1];
1531 [ # # ]: 0 : zz = poly->verts[vtx]->coord[2];
1532 [ # # ]: 0 : hashvalue = makeahashvaluefrom_coord(xx,yy,zz);
1533 [ # # ]: 0 : hasharrayptr = hashobj->getHashBin(hashvalue,&hasharraysize);
1534 : 0 : ifoundit = -1;
1535 [ # # ]: 0 : for ( i = 0; i < hasharraysize; i++ ) {
1536 : 0 : hptr = hasharrayptr[i];
1537 [ # # ]: 0 : xval = out_coords[3*hptr];
1538 [ # # ]: 0 : yval = out_coords[3*hptr+1];
1539 [ # # ]: 0 : zval = out_coords[3*hptr+2];
1540 [ # # ][ # # ]: 0 : if ( ( fabs(xval-xx) < EPSILON ) &&
1541 [ # # ]: 0 : ( fabs(yval-yy) < EPSILON ) &&
1542 : 0 : ( fabs(zval-zz) < EPSILON ) ) {
1543 : 0 : ifoundit = hasharrayptr[i];
1544 : 0 : break;
1545 : : }
1546 : : }
1547 [ # # ]: 0 : if ( ifoundit == -1 ) {
1548 : 0 : ifoundit = num_sofar;
1549 [ # # ]: 0 : hashobj->addtoHashList(hashvalue,num_sofar);
1550 [ # # ]: 0 : out_coords.push_back(xx);
1551 [ # # ]: 0 : out_coords.push_back(yy);
1552 [ # # ]: 0 : out_coords.push_back(zz);
1553 : 0 : num_sofar++;
1554 : : }
1555 : :
1556 : 0 : return ifoundit;
1557 : : }
1558 : :
1559 : 0 : int FBIntersect::makeahashvaluefrom_coord(double x, double y, double z)
1560 : : {
1561 : : double mantissasum;
1562 : :
1563 [ # # ]: 0 : if ( fabs(x) < 1.e-3 ) x = 0.0;
1564 [ # # ]: 0 : if ( fabs(y) < 1.e-3 ) y = 0.0;
1565 [ # # ]: 0 : if ( fabs(z) < 1.e-3 ) z = 0.0;
1566 : 0 : mantissasum = (int)(10000.0*fabs(x) + 0.5) +
1567 : 0 : (int)(10000.0*fabs(y) + 0.5) +
1568 : 0 : (int)(10000.0*fabs(z) + 0.5);
1569 : :
1570 : 0 : return (int)(mantissasum) % NUMHASHBINS;
1571 : : }
1572 : :
1573 : :
1574 : 0 : void FBIntersect::set_body1_planar()
1575 : : {
1576 : 0 : body1_is_plane = true;
1577 : 0 : }
1578 : :
1579 : 0 : void FBIntersect::set_body2_planar()
1580 : : {
1581 : 0 : body2_is_plane = true;
1582 : 0 : }
1583 : :
1584 : 0 : void FBIntersect::set_imprint()
1585 : : {
1586 : 0 : do_imprint = true;
1587 : 0 : }
1588 : :
1589 : 0 : void FBIntersect::newplanecoefficients(FBPolyhedron *poly, FB_Triangle *tri)
1590 : : {
1591 : : FB_Coord *mycoord;
1592 : : double x1, x2, x3, y1, y2, y3, z1, z2, z3, e1x, e1y, e1z, e2x, e2y, e2z;
1593 : : double a, b, c, d, dtemp;
1594 : :
1595 : 0 : mycoord = poly->verts[tri->v0];
1596 : 0 : x1 = mycoord->coord[0];
1597 : 0 : y1 = mycoord->coord[1];
1598 : 0 : z1 = mycoord->coord[2];
1599 : 0 : mycoord = poly->verts[tri->v1];
1600 : 0 : x2 = mycoord->coord[0];
1601 : 0 : y2 = mycoord->coord[1];
1602 : 0 : z2 = mycoord->coord[2];
1603 : 0 : mycoord = poly->verts[tri->v2];
1604 : 0 : x3 = mycoord->coord[0];
1605 : 0 : y3 = mycoord->coord[1];
1606 : 0 : z3 = mycoord->coord[2];
1607 : 0 : e1x = x1 - x2; e1y = y1 - y2; e1z = z1 - z2;
1608 : 0 : e2x = x3 - x2; e2y = y3 - y2; e2z = z3 - z2;
1609 : 0 : a = e1z*e2y - e2z*e1y;
1610 : 0 : b = e1x*e2z - e2x*e1z;
1611 : 0 : c = e1y*e2x - e2y*e1x;
1612 : 0 : dtemp = sqrt(a*a + b*b + c*c);
1613 [ # # ]: 0 : if ( dtemp > EPSILON2 ) {
1614 : 0 : a /= dtemp;
1615 : 0 : b /= dtemp;
1616 : 0 : c /= dtemp;
1617 : : } else {
1618 [ # # ][ # # ]: 0 : PRINT_WARNING("small-area triangle\n");
1619 : : }
1620 : 0 : d = -(a*x1 + b*y1 + c*z1);
1621 : 0 : tri->a = a; tri->b = b; tri->c = c; tri->d = d;
1622 : :
1623 [ + - ][ + - ]: 6540 : }
|
