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 <math.h>
28 : : #include <algorithm>
29 : : #include <vector>
30 : : #include "CubitMessage.hpp"
31 : : #include "FBTiler.hpp"
32 : :
33 : 0 : FBTiler::FBTiler(std::vector<FB_Coord *>& my_verts, int pd, int sd, int sequence,
34 : : double a, double b, double c,
35 [ # # ]: 0 : std::vector<int> *tri_list)
36 : : {
37 : :
38 [ # # ]: 0 : verts = my_verts;
39 : 0 : p_dir = pd;
40 : 0 : s_dir = sd;
41 : 0 : xnorm = a;
42 : 0 : ynorm = b;
43 : 0 : znorm = c;
44 : 0 : parent = sequence;
45 : 0 : my_tri_list = tri_list;
46 : 0 : }
47 : :
48 [ # # ]: 0 : FBTiler::~FBTiler()
49 : : {
50 : :
51 : 0 : }
52 : :
53 : 0 : CubitStatus FBTiler::Tile_Poly(std::vector<int> *coordlist)
54 : : {
55 : : int v0, v1, v2, v3;
56 [ # # ][ # # ]: 0 : std::vector<int>::iterator it, itmin, itmax;
[ # # ]
57 : : double min_p, min_s, max_p, max_s, test_p, test_s;
58 : : CubitStatus status;
59 : :
60 : 0 : status = CUBIT_SUCCESS;
61 : :
62 [ # # ]: 0 : it = coordlist->begin();
63 : :
64 [ # # ][ # # ]: 0 : if ( coordlist->size() == 3 ) {
65 [ # # ][ # # ]: 0 : v1 = *it++; v2 = *it++; v3 = *it++;
[ # # ][ # # ]
[ # # ][ # # ]
66 [ # # ]: 0 : add_triangle(v1, v2, v3);
67 : 0 : return status;
68 : : }
69 : :
70 : : // Get location in list of min and max coord in p_dir.
71 : 0 : min_p = min_s = CUBIT_DBL_MAX;
72 : 0 : max_p = max_s = -CUBIT_DBL_MAX;
73 [ # # ][ # # ]: 0 : while ( it != coordlist->end() ) {
[ # # ]
74 [ # # ]: 0 : v0 = *it;
75 [ # # ]: 0 : test_p = verts[v0]->coord[p_dir];
76 [ # # ]: 0 : test_s = verts[v0]->coord[s_dir];
77 : : // if ( test_p < min_p ) {
78 [ # # ]: 0 : if ( (test_p - min_p) <= -EPSILON ) {
79 : 0 : itmin = it;
80 : 0 : min_p = test_p;
81 [ # # ]: 0 : if ( test_s < min_s )
82 : 0 : min_s = test_s;
83 : : }
84 : : // if ( test_p > max_p ) {
85 [ # # ]: 0 : if ( (test_p - max_p) >= EPSILON ) {
86 : 0 : itmax = it;
87 : 0 : max_p = test_p;
88 [ # # ]: 0 : if ( test_s > max_s )
89 : 0 : max_s = test_s;
90 : : }
91 [ # # ]: 0 : if ( fabs(test_p - min_p) < EPSILON ) {
92 [ # # ]: 0 : if ( test_s < min_s ) {
93 : 0 : itmin = it;
94 : 0 : min_s = test_s;
95 : : }
96 : : }
97 [ # # ]: 0 : if ( fabs(test_p - max_p) < EPSILON ) {
98 [ # # ]: 0 : if ( test_s > max_s ) {
99 : 0 : itmax = it;
100 : 0 : max_s = test_s;
101 : : }
102 : : }
103 : :
104 [ # # ]: 0 : it++;
105 : : }
106 : :
107 : : // Now we've got to start at the min coord and make a single sorted list
108 : : // of the coords, from min to max, noting whether the coord is on the
109 : : // LEFT or RIGHT chain, or -- for the top and bottom coords -- on BOTH chains.
110 : :
111 [ # # ]: 0 : sortedchainlist.clear();
112 : :
113 [ # # ][ # # ]: 0 : std::vector<int>::iterator itleft, itright, itlistbegin, itlistend;
[ # # ][ # # ]
114 : 0 : it = itleft = itright = itmin;
115 [ # # ]: 0 : itlistbegin = coordlist->begin();
116 [ # # ]: 0 : itlistend = coordlist->end();
117 : :
118 [ # # ]: 0 : decrement_list_ptr(itleft,itlistbegin,itlistend);
119 [ # # ]: 0 : increment_list_ptr(itright,itlistbegin,itlistend);
120 : :
121 [ # # ]: 0 : v0 = *itmin;
122 : : int v0left, v0right;
123 : : double leftpdircoord, rightpdircoord;
124 : : FBTilerChainVert *cvert;
125 [ # # ][ # # ]: 0 : cvert = new FBTilerChainVert(v0,BOTH);
126 [ # # ]: 0 : sortedchainlist.push_back(cvert);
127 [ # # ][ # # ]: 0 : v0left = *itleft; v0right = *itright;
128 [ # # ]: 0 : leftpdircoord = verts[v0left]->coord[p_dir];
129 [ # # ]: 0 : rightpdircoord = verts[v0right]->coord[p_dir];
130 : :
131 : : while(1) {
132 [ # # ]: 0 : while ( (leftpdircoord - rightpdircoord) > -EPSILON2 ) {
133 [ # # ][ # # ]: 0 : if ( itright == itmax ) break;
134 : 0 : v0 = v0right;
135 [ # # ][ # # ]: 0 : cvert = new FBTilerChainVert(v0,LEFT);
136 [ # # ]: 0 : sortedchainlist.push_back(cvert);
137 [ # # ]: 0 : increment_list_ptr(itright,itlistbegin,itlistend);
138 [ # # ]: 0 : v0right = *itright;
139 [ # # ]: 0 : rightpdircoord = verts[v0right]->coord[p_dir];
140 : :
141 : : }
142 [ # # ]: 0 : while ( (rightpdircoord - leftpdircoord) > -EPSILON2 ) {
143 [ # # ][ # # ]: 0 : if ( itleft == itmax ) break;
144 : 0 : v0 = v0left;
145 [ # # ][ # # ]: 0 : cvert = new FBTilerChainVert(v0,RIGHT);
146 [ # # ]: 0 : sortedchainlist.push_back(cvert);
147 [ # # ]: 0 : decrement_list_ptr(itleft,itlistbegin,itlistend);
148 [ # # ]: 0 : v0left = *itleft;
149 [ # # ]: 0 : leftpdircoord = verts[v0left]->coord[p_dir];
150 : :
151 : : }
152 [ # # ][ # # ]: 0 : if ( (itright == itmax) || (itleft == itmax) ) break;
[ # # ][ # # ]
[ # # ]
153 : : }
154 : :
155 : : // In case we have some horizontal edges (wrt p_dir) left on the list,
156 : : // the next two while statements will get them.
157 [ # # ][ # # ]: 0 : while ( itright != itmax ) {
158 : 0 : v0 = v0right;
159 [ # # ][ # # ]: 0 : cvert = new FBTilerChainVert(v0,LEFT);
160 [ # # ]: 0 : sortedchainlist.push_back(cvert);
161 [ # # ]: 0 : increment_list_ptr(itright,itlistbegin,itlistend);
162 [ # # ]: 0 : v0right = *itright;
163 : : }
164 [ # # ][ # # ]: 0 : while ( itleft != itmax ) {
165 : 0 : v0 = v0left;
166 [ # # ][ # # ]: 0 : cvert = new FBTilerChainVert(v0,RIGHT);
167 [ # # ]: 0 : sortedchainlist.push_back(cvert);
168 [ # # ]: 0 : decrement_list_ptr(itleft,itlistbegin,itlistend);
169 [ # # ]: 0 : v0left = *itleft;
170 : : }
171 : :
172 : : // Put the top coord on the list.
173 : :
174 [ # # ]: 0 : v0 = *itmax;
175 [ # # ][ # # ]: 0 : cvert = new FBTilerChainVert(v0,BOTH);
176 [ # # ]: 0 : sortedchainlist.push_back(cvert);
177 [ # # ]: 0 : status = retriangulate(coordlist);
178 [ # # ][ # # ]: 0 : std::vector<FBTilerChainVert*>::iterator itv, itvlast;
179 : :
180 [ # # ]: 0 : itvlast = sortedchainlist.end();
181 [ # # ]: 0 : itv = sortedchainlist.begin();
182 [ # # ][ # # ]: 0 : while ( itv != itvlast ) {
183 [ # # ][ # # ]: 0 : delete *itv;
[ # # ]
184 [ # # ]: 0 : itv++;
185 : : }
186 [ # # ]: 0 : sortedchainlist.clear();
187 : 0 : return status;
188 : : }
189 : :
190 : 0 : int FBTiler::add_triangle(int v1, int v2, int v3)
191 : : {
192 : : int status;
193 : :
194 : 0 : status = 1;
195 : :
196 : : // See if the winding is OK by comparing the normal to the parent's normal.
197 : : double xn, yn, zn;
198 : : double x1, y1, z1, x2, y2, z2, x3, y3, z3;
199 : : double ux, uy, uz, vx, vy, vz, product;
200 : :
201 : 0 : x1 = verts[v1]->coord[0]; y1 = verts[v1]->coord[1]; z1 = verts[v1]->coord[2];
202 : 0 : x2 = verts[v2]->coord[0]; y2 = verts[v2]->coord[1]; z2 = verts[v2]->coord[2];
203 : 0 : x3 = verts[v3]->coord[0]; y3 = verts[v3]->coord[1]; z3 = verts[v3]->coord[2];
204 : 0 : ux = x2 - x1; uy = y2 - y1; uz = z2 - z1;
205 : 0 : vx = x3 - x1; vy = y3 - y1; vz = z3 - z1;
206 : 0 : xn = uy*vz - uz*vy; yn = uz*vx - ux*vz; zn = ux*vy - uy*vx;
207 : 0 : product = xnorm*xn + ynorm*yn + znorm*zn;
208 [ # # ]: 0 : if ( product < 0.0 ) { // Must reverse the connections.
209 : : int itemp;
210 : 0 : itemp = v1; v1 = v2; v2 = itemp;
211 : : }
212 : 0 : my_tri_list->push_back(v1);
213 : 0 : my_tri_list->push_back(v2);
214 : 0 : my_tri_list->push_back(v3);
215 : :
216 : 0 : return status;
217 : :
218 : : }
219 : :
220 : 0 : void FBTiler::dud_from_coord_list(int val, std::vector<int> *ivec)
221 : : {
222 [ # # ]: 0 : std::vector<int>::iterator itv;
223 : : bool ifoundit;
224 : : // Also decrements ivec->size() if val was found.
225 : :
226 : 0 : ifoundit = false;
227 [ # # ]: 0 : itv = ivec->begin();
228 [ # # ][ # # ]: 0 : do {
229 [ # # ]: 0 : if ( ifoundit == false ) {
230 [ # # ][ # # ]: 0 : if ( *itv == val ) {
231 : 0 : ifoundit = true;
232 : 0 : break;
233 : : }
234 : : }
235 [ # # ]: 0 : itv++;
236 [ # # ]: 0 : } while ( itv != ivec->end() );
237 [ # # ]: 0 : if ( ifoundit == true ) {
238 [ # # ]: 0 : ivec->erase(itv);
239 : : }
240 : :
241 : 0 : }
242 : :
243 : 0 : bool FBTiler::reflex_angle(int v0, int v1, int v2, std::vector<int> *coordlist)
244 : : {
245 : : double v0x, v0y, v1x, v1y, v2x, v2y, xbary, ybary;
246 : :
247 [ # # ]: 0 : v0x = verts[v0]->coord[s_dir];
248 [ # # ]: 0 : v0y = verts[v0]->coord[p_dir];
249 [ # # ]: 0 : v1x = verts[v1]->coord[s_dir];
250 [ # # ]: 0 : v1y = verts[v1]->coord[p_dir];
251 [ # # ]: 0 : v2x = verts[v2]->coord[s_dir];
252 [ # # ]: 0 : v2y = verts[v2]->coord[p_dir];
253 : : // Are the three points colinear?
254 : :
255 [ # # ]: 0 : if ( fabs((v1x-v0x)*(v2y-v0y) - (v2x-v0x)*(v1y-v0y)) < EPSILON2 )
256 : 0 : return true;
257 : :
258 : : // Get coordinates of barycenter of the triangle formed by
259 : : // v0, v1, and v2
260 : :
261 : 0 : xbary = (v0x + v1x + v2x)/3.;
262 : 0 : ybary = (v0y + v1y + v2y)/3.;
263 : :
264 : : // Do a point-in-polygon test on this point and the polygon formed
265 : : // by the chain. If point is in polygon, angle is not reflex;
266 : : // otherwise, angle is reflex. Point in polygon algorithm from
267 : : // Schneider and Eberly, Geometric Tools for COmputer Graphics,
268 : : // Sec. 13.3.
269 : :
270 : 0 : bool inside = false;
271 [ # # ]: 0 : std::vector<int>::iterator itv;
272 : : int u0, u1;
273 : :
274 [ # # ]: 0 : itv = coordlist->begin();
275 [ # # ][ # # ]: 0 : while ( itv != coordlist->end() ) {
[ # # ]
276 [ # # ][ # # ]: 0 : u0 = *itv++;
277 [ # # ][ # # ]: 0 : if ( itv == coordlist->end() )
[ # # ]
278 [ # # ][ # # ]: 0 : u1 = *(coordlist->begin());
279 : : else
280 [ # # ]: 0 : u1 = *itv;
281 [ # # ][ # # ]: 0 : if ( ybary < verts[u1]->coord[p_dir] ) {
282 : : // u1 is above the hroizontal ray in s_dir from the barycenter.
283 [ # # ][ # # ]: 0 : if ( verts[u0]->coord[p_dir] <= ybary ) {
284 : : // u0 is on or below the ray.
285 [ # # ][ # # ]: 0 : if ( ( (ybary - verts[u0]->coord[p_dir])*
286 [ # # ][ # # ]: 0 : (verts[u1]->coord[s_dir] -verts[ u0]->coord[s_dir]) ) >
287 [ # # ]: 0 : ( (xbary - verts[u0]->coord[s_dir])*
288 [ # # ][ # # ]: 0 : (verts[u1]->coord[p_dir] - verts[u0]->coord[p_dir]) ) )
289 : 0 : inside = !inside;
290 : : }
291 [ # # ][ # # ]: 0 : } else if ( ybary < verts[u0]->coord[p_dir] ) {
292 : : // U1 is on or below the ray; u0 is above the ray.
293 [ # # ][ # # ]: 0 : if ( ( (ybary - verts[u0]->coord[p_dir])*
294 [ # # ][ # # ]: 0 : (verts[u1]->coord[s_dir] - verts[u0]->coord[s_dir]) ) <
295 [ # # ]: 0 : ( (xbary - verts[u0]->coord[s_dir])*
296 [ # # ][ # # ]: 0 : (verts[u1]->coord[p_dir] - verts[u0]->coord[p_dir]) ) )
297 : 0 : inside = !inside;
298 : : }
299 : : }
300 : :
301 : 0 : return !inside;
302 : :
303 : : }
304 : :
305 : 0 : CubitStatus FBTiler::retriangulate(std::vector<int> *coordlist)
306 : : {
307 : : int stacktopadjacency, this_adjacency, adjacency_case;
308 [ # # ]: 0 : std::vector<FBTilerChainVert*> vstack;
309 [ # # ][ # # ]: 0 : std::vector<FBTilerChainVert*>::iterator itv, itvlast;
310 : : int v0, v1, v2;
311 : : unsigned int i, numtris, totaltris;
312 : :
313 : : // Put the first two verts on vstack.
314 [ # # ]: 0 : itv = sortedchainlist.begin();
315 [ # # ][ # # ]: 0 : vstack.push_back(*itv++);
[ # # ]
316 [ # # ][ # # ]: 0 : vstack.push_back(*itv);
317 [ # # ]: 0 : stacktopadjacency = (*itv)->whichchain;
318 : 0 : itvlast = itv;
319 : 0 : numtris = 0;
320 [ # # ]: 0 : totaltris = sortedchainlist.size() - 2;
321 : :
322 : : while(1) {
323 : :
324 [ # # ]: 0 : if ( numtris == totaltris ) break;
325 [ # # ][ # # ]: 0 : *itv++;
326 [ # # ][ # # ]: 0 : if ( itv == sortedchainlist.end() ) {
[ # # ]
327 [ # # ][ # # ]: 0 : PRINT_ERROR("Tiler error in FacetBoolean: premature end of list.\n");
[ # # ][ # # ]
328 : 0 : return CUBIT_FAILURE;
329 : : }
330 : :
331 [ # # ]: 0 : this_adjacency = (*itv)->whichchain;
332 [ # # ]: 0 : v0 = (*itv)->v0;
333 : :
334 [ # # ]: 0 : adjacency_case = get_adjacency(stacktopadjacency,this_adjacency);
335 : :
336 [ # # # ]: 0 : switch ( adjacency_case ) {
337 : :
338 : : case 1:
339 : : {
340 [ # # ]: 0 : v1 = vstack[0]->v0;
341 [ # # ][ # # ]: 0 : for ( i = 1; i < vstack.size(); i++ ) {
342 [ # # ]: 0 : v2 = vstack[i]->v0;
343 [ # # ]: 0 : add_triangle(v0,v1,v2);
344 [ # # ]: 0 : dud_from_coord_list(v1,coordlist);
345 : 0 : numtris += 1;
346 : 0 : v1 = v2;
347 : : }
348 [ # # ]: 0 : vstack.clear();
349 [ # # ][ # # ]: 0 : vstack.push_back(*itvlast);
350 [ # # ][ # # ]: 0 : vstack.push_back(*itv);
351 [ # # ]: 0 : stacktopadjacency = (*itv)->whichchain;
352 : 0 : itvlast = itv;
353 : :
354 : 0 : break;
355 : : }
356 : :
357 : : case 2:
358 : : {
359 [ # # ][ # # ]: 0 : while ( vstack.size() > 1 ) {
360 [ # # ]: 0 : int size = vstack.size() - 1;
361 [ # # ]: 0 : v1 = vstack[size]->v0;
362 : : // int v1chain = vstack[size]->whichchain;
363 [ # # ]: 0 : v2 = vstack[size-1]->v0;
364 : : // bool isreflex = reflex_angle(v0,v1,v2,v1chain);
365 [ # # ]: 0 : bool isreflex = reflex_angle(v0,v1,v2,coordlist);
366 [ # # ]: 0 : if ( isreflex == true ) break;
367 [ # # ]: 0 : add_triangle(v0,v1,v2);
368 [ # # ]: 0 : dud_from_coord_list(v1,coordlist);
369 : 0 : numtris += 1;
370 [ # # ]: 0 : vstack.pop_back();
371 : :
372 : : }
373 : :
374 [ # # ][ # # ]: 0 : vstack.push_back(*itv);
375 [ # # ]: 0 : stacktopadjacency = (*itv)->whichchain;
376 : 0 : itvlast = itv;
377 : :
378 : :
379 : 0 : break;
380 : : }
381 : :
382 : : }
383 : :
384 : : } // endwhile
385 : :
386 [ # # ]: 0 : return CUBIT_SUCCESS;
387 : :
388 : : }
389 : :
390 : :
|
