Branch data Line data Source code
1 : : #include "moab/CartVect.hpp"
2 : : #include "moab/BSPTreePoly.hpp"
3 : : #include <assert.h>
4 : : #include <stdlib.h>
5 : : #include <set>
6 : :
7 : : #undef DEBUG_IDS
8 : :
9 : : namespace moab
10 : : {
11 : :
12 : : struct BSPTreePoly::Vertex : public CartVect
13 : : {
14 : 670 : Vertex( const CartVect& v )
15 : 670 : : CartVect( v ), usePtr( 0 ), markVal( 0 )
16 : : #ifdef DEBUG_IDS
17 : : ,
18 : : id( nextID++ )
19 : : #endif
20 : : {
21 : 670 : }
22 : 670 : ~Vertex()
23 : : {
24 [ - + ]: 670 : assert( !usePtr );
25 : 670 : }
26 : : BSPTreePoly::VertexUse* usePtr;
27 : : int markVal;
28 : : #ifdef DEBUG_IDS
29 : : int id;
30 : : static int nextID;
31 : : #endif
32 : : };
33 : :
34 : : struct BSPTreePoly::VertexUse
35 : : {
36 : : VertexUse( Edge* edge, Vertex* vtx );
37 : : ~VertexUse();
38 : :
39 : : void set_vertex( BSPTreePoly::Vertex*& vtx_ptr );
40 : :
41 : : BSPTreePoly::VertexUse *nextPtr, *prevPtr;
42 : : BSPTreePoly::Vertex* vtxPtr;
43 : : BSPTreePoly::Edge* edgePtr;
44 : : };
45 : :
46 : : struct BSPTreePoly::EdgeUse
47 : : {
48 : : EdgeUse( Edge* edge );
49 : : EdgeUse( Edge* edge, Face* face );
50 : : ~EdgeUse();
51 : :
52 : : BSPTreePoly::EdgeUse *prevPtr, *nextPtr;
53 : : BSPTreePoly::Edge* edgePtr;
54 : : BSPTreePoly::Face* facePtr;
55 : :
56 : : inline BSPTreePoly::Vertex* start() const;
57 : : inline BSPTreePoly::Vertex* end() const;
58 : : int sense() const;
59 : :
60 : : void insert_after( BSPTreePoly::EdgeUse* prev );
61 : : void insert_before( BSPTreePoly::EdgeUse* next );
62 : : };
63 : :
64 : : struct BSPTreePoly::Edge
65 : : {
66 : : BSPTreePoly::VertexUse *startPtr, *endPtr;
67 : : BSPTreePoly::EdgeUse *forwardPtr, *reversePtr;
68 : : #ifdef DEBUG_IDS
69 : : int id;
70 : : static int nextID;
71 : : #endif
72 : :
73 : 1232 : Edge( Vertex* vstart, Vertex* vend )
74 : 1232 : : forwardPtr( 0 ), reversePtr( 0 )
75 : : #ifdef DEBUG_IDS
76 : : ,
77 : : id( nextID++ )
78 : : #endif
79 : : {
80 [ + - ]: 1232 : startPtr = new VertexUse( this, vstart );
81 [ + - ]: 1232 : endPtr = new VertexUse( this, vend );
82 : 1232 : }
83 : :
84 : : ~Edge();
85 : :
86 : 14576 : BSPTreePoly::Vertex* start() const
87 : : {
88 : 14576 : return startPtr->vtxPtr;
89 : : }
90 : 14595 : BSPTreePoly::Vertex* end() const
91 : : {
92 : 14595 : return endPtr->vtxPtr;
93 : : }
94 : :
95 : : BSPTreePoly::Face* forward() const
96 : : {
97 : : return forwardPtr ? forwardPtr->facePtr : 0;
98 : : }
99 : : BSPTreePoly::Face* reverse() const
100 : : {
101 : : return reversePtr ? reversePtr->facePtr : 0;
102 : : }
103 : :
104 : 351 : BSPTreePoly::VertexUse* use( BSPTreePoly::Vertex* vtx ) const
105 : : {
106 [ + + ][ + - ]: 351 : return ( vtx == startPtr->vtxPtr ) ? startPtr : ( vtx == endPtr->vtxPtr ) ? endPtr : 0;
107 : : }
108 : : BSPTreePoly::Edge* next( BSPTreePoly::Vertex* about ) const
109 : : {
110 : : return use( about )->nextPtr->edgePtr;
111 : : }
112 : : BSPTreePoly::Edge* prev( BSPTreePoly::Vertex* about ) const
113 : : {
114 : : return use( about )->prevPtr->edgePtr;
115 : : }
116 : :
117 : : BSPTreePoly::EdgeUse* use( BSPTreePoly::Face* face ) const
118 : : {
119 : : return ( face == forwardPtr->facePtr ) ? forwardPtr : ( face == reversePtr->facePtr ) ? reversePtr : 0;
120 : : }
121 : : BSPTreePoly::Edge* next( BSPTreePoly::Face* about ) const
122 : : {
123 : : return use( about )->nextPtr->edgePtr;
124 : : }
125 : : BSPTreePoly::Edge* prev( BSPTreePoly::Face* about ) const
126 : : {
127 : : return use( about )->prevPtr->edgePtr;
128 : : }
129 : :
130 : : BSPTreePoly::VertexUse* other( BSPTreePoly::VertexUse* vuse ) const
131 : : {
132 : : return vuse == startPtr ? endPtr : vuse == endPtr ? startPtr : 0;
133 : : }
134 : 389 : BSPTreePoly::EdgeUse* other( BSPTreePoly::EdgeUse* vuse ) const
135 : : {
136 [ + + ][ + - ]: 389 : return vuse == forwardPtr ? reversePtr : vuse == reversePtr ? forwardPtr : 0;
137 : : }
138 : : BSPTreePoly::Vertex* other( BSPTreePoly::Vertex* vtx ) const
139 : : {
140 : : return vtx == startPtr->vtxPtr ? endPtr->vtxPtr : vtx == endPtr->vtxPtr ? startPtr->vtxPtr : 0;
141 : : }
142 : : BSPTreePoly::Vertex* common( BSPTreePoly::Edge* eother ) const
143 : : {
144 : : return start() == eother->start() || start() == eother->end()
145 : : ? start()
146 : : : end() == eother->start() || end() == eother->end() ? end() : 0;
147 : : }
148 : :
149 : : int sense( BSPTreePoly::Face* face ) const;
150 : :
151 : : void remove_from_vertex( BSPTreePoly::Vertex*& vtx_ptr );
152 : : void remove_from_face( BSPTreePoly::Face*& face_ptr );
153 : : void add_to_vertex( BSPTreePoly::Vertex* vtx_ptr );
154 : : };
155 : :
156 : : struct BSPTreePoly::Face
157 : : {
158 : 269 : Face( Face* next )
159 : 269 : : usePtr( 0 ), nextPtr( next )
160 : : #ifdef DEBUG_IDS
161 : : ,
162 : : id( nextID++ )
163 : : #endif
164 : : {
165 : 269 : }
166 : 450 : Face()
167 : 450 : : usePtr( 0 ), nextPtr( 0 )
168 : : #ifdef DEBUG_IDS
169 : : ,
170 : : id( nextID++ )
171 : : #endif
172 : : {
173 : 450 : }
174 : : ~Face();
175 : : BSPTreePoly::EdgeUse* usePtr;
176 : : BSPTreePoly::Face* nextPtr;
177 : : #ifdef DEBUG_IDS
178 : : int id;
179 : : static int nextID;
180 : : #endif
181 : : double signed_volume() const;
182 : : };
183 : :
184 : : #ifdef DEBUG_IDS
185 : : int BSPTreePoly::Vertex::nextID = 1;
186 : : int BSPTreePoly::Edge::nextID = 1;
187 : : int BSPTreePoly::Face::nextID = 1;
188 : : #endif
189 : 3 : void BSPTreePoly::reset_debug_ids()
190 : : {
191 : : #ifdef DEBUG_IDS
192 : : BSPTreePoly::Vertex::nextID = 1;
193 : : BSPTreePoly::Edge::nextID = 1;
194 : : BSPTreePoly::Face::nextID = 1;
195 : : #endif
196 : 3 : }
197 : :
198 : : // static void merge_edges( BSPTreePoly::Edge* keep_edge,
199 : : // BSPTreePoly::Edge* dead_edge );
200 : :
201 : : static BSPTreePoly::Edge* split_edge( BSPTreePoly::Vertex*& new_vtx, BSPTreePoly::Edge* into_edge );
202 : :
203 : 2464 : BSPTreePoly::VertexUse::VertexUse( BSPTreePoly::Edge* edge, BSPTreePoly::Vertex* vtx ) : vtxPtr( vtx ), edgePtr( edge )
204 : : {
205 [ + + ]: 2464 : if( !vtx->usePtr )
206 : : {
207 : 670 : vtx->usePtr = prevPtr = nextPtr = this;
208 : 670 : return;
209 : : }
210 : :
211 : 1794 : nextPtr = vtx->usePtr;
212 : 1794 : prevPtr = nextPtr->prevPtr;
213 [ - + ]: 1794 : assert( prevPtr->nextPtr == nextPtr );
214 : 1794 : nextPtr->prevPtr = this;
215 : 1794 : prevPtr->nextPtr = this;
216 : : }
217 : :
218 : 2464 : BSPTreePoly::VertexUse::~VertexUse()
219 : : {
220 [ + + ]: 2464 : if( nextPtr == this )
221 : : {
222 [ - + ]: 670 : assert( prevPtr == this );
223 [ - + ]: 670 : assert( vtxPtr->usePtr == this );
224 : 670 : vtxPtr->usePtr = 0;
225 [ + - ]: 670 : delete vtxPtr;
226 : : }
227 [ + + ]: 1794 : else if( vtxPtr->usePtr == this )
228 : 730 : vtxPtr->usePtr = nextPtr;
229 : :
230 : 2464 : nextPtr->prevPtr = prevPtr;
231 : 2464 : prevPtr->nextPtr = nextPtr;
232 : 2464 : nextPtr = prevPtr = 0;
233 : 2464 : }
234 : :
235 : 446 : void BSPTreePoly::VertexUse::set_vertex( BSPTreePoly::Vertex*& vtx )
236 : : {
237 [ + - ]: 446 : if( vtxPtr )
238 : : {
239 [ - + ]: 446 : if( nextPtr == prevPtr )
240 : : {
241 [ # # ]: 0 : assert( nextPtr == this );
242 : 0 : vtxPtr->usePtr = 0;
243 [ # # ]: 0 : delete vtx;
244 : 0 : vtx = 0;
245 : : }
246 : : else
247 : : {
248 : 446 : nextPtr->prevPtr = prevPtr;
249 : 446 : prevPtr->nextPtr = nextPtr;
250 [ + + ]: 446 : if( vtxPtr->usePtr == this ) vtxPtr->usePtr = nextPtr;
251 : : }
252 : : }
253 : :
254 [ + - ]: 446 : if( vtx )
255 : : {
256 : 446 : vtxPtr = vtx;
257 : 446 : nextPtr = vtxPtr->usePtr->nextPtr;
258 : 446 : prevPtr = vtxPtr->usePtr;
259 : 446 : nextPtr->prevPtr = this;
260 : 446 : vtxPtr->usePtr->nextPtr = this;
261 : : }
262 : 446 : }
263 : :
264 : 2647 : BSPTreePoly::EdgeUse::EdgeUse( BSPTreePoly::Edge* edge ) : prevPtr( 0 ), nextPtr( 0 ), edgePtr( edge ), facePtr( 0 ) {}
265 : :
266 : 269 : BSPTreePoly::EdgeUse::EdgeUse( BSPTreePoly::Edge* edge, BSPTreePoly::Face* face ) : edgePtr( edge ), facePtr( face )
267 : : {
268 [ - + ]: 269 : assert( !face->usePtr );
269 : 269 : face->usePtr = prevPtr = nextPtr = this;
270 : :
271 [ - + ]: 269 : if( !face->usePtr )
272 : : {
273 : 0 : face->usePtr = prevPtr = nextPtr = this;
274 : 0 : return;
275 : : }
276 : :
277 : 269 : nextPtr = face->usePtr;
278 : 269 : prevPtr = nextPtr->prevPtr;
279 [ - + ]: 269 : assert( prevPtr->nextPtr == nextPtr );
280 : 269 : nextPtr->prevPtr = this;
281 : 269 : prevPtr->nextPtr = this;
282 : : }
283 : :
284 : 1301 : void BSPTreePoly::EdgeUse::insert_after( BSPTreePoly::EdgeUse* prev )
285 : : {
286 : : // shouldn't already be in a face
287 [ - + ]: 1301 : assert( !facePtr );
288 : : // adjacent edges should share vertices
289 [ - + ]: 1301 : assert( start() == prev->end() );
290 : :
291 : 1301 : facePtr = prev->facePtr;
292 : 1301 : nextPtr = prev->nextPtr;
293 : 1301 : prevPtr = prev;
294 : 1301 : nextPtr->prevPtr = this;
295 : 1301 : prevPtr->nextPtr = this;
296 : 1301 : }
297 : :
298 : 446 : void BSPTreePoly::EdgeUse::insert_before( BSPTreePoly::EdgeUse* next )
299 : : {
300 : : // shouldn't already be in a face
301 [ - + ]: 446 : assert( !facePtr );
302 : : // adjacent edges should share vertices
303 [ - + ]: 446 : assert( end() == next->start() );
304 : :
305 : 446 : facePtr = next->facePtr;
306 : 446 : prevPtr = next->prevPtr;
307 : 446 : nextPtr = next;
308 : 446 : nextPtr->prevPtr = this;
309 : 446 : prevPtr->nextPtr = this;
310 : 446 : }
311 : :
312 : 2916 : BSPTreePoly::EdgeUse::~EdgeUse()
313 : : {
314 [ + - ][ + + ]: 2916 : if( facePtr->usePtr == this ) facePtr->usePtr = ( nextPtr == this ) ? 0 : nextPtr;
315 : :
316 [ + + ]: 2916 : if( edgePtr->forwardPtr == this ) edgePtr->forwardPtr = 0;
317 [ + + ]: 2916 : if( edgePtr->reversePtr == this ) edgePtr->reversePtr = 0;
318 : :
319 [ + + ][ + + ]: 2916 : if( !edgePtr->forwardPtr && !edgePtr->reversePtr ) delete edgePtr;
[ + - ]
320 : :
321 : 2916 : nextPtr->prevPtr = prevPtr;
322 : 2916 : prevPtr->nextPtr = nextPtr;
323 : 2916 : nextPtr = prevPtr = 0;
324 : 2916 : }
325 : :
326 : 0 : int BSPTreePoly::EdgeUse::sense() const
327 : : {
328 [ # # ]: 0 : if( edgePtr->forwardPtr == this )
329 : 0 : return 1;
330 [ # # ]: 0 : else if( edgePtr->reversePtr == this )
331 : 0 : return -1;
332 : : else
333 : 0 : return 0;
334 : : }
335 : :
336 : 7020 : BSPTreePoly::Vertex* BSPTreePoly::EdgeUse::start() const
337 : : {
338 [ + + ]: 7020 : if( edgePtr->forwardPtr == this )
339 : 3708 : return edgePtr->start();
340 [ + - ]: 3312 : else if( edgePtr->reversePtr == this )
341 : 3312 : return edgePtr->end();
342 : : else
343 : 0 : return 0;
344 : : }
345 : :
346 : 7505 : BSPTreePoly::Vertex* BSPTreePoly::EdgeUse::end() const
347 : : {
348 [ + + ]: 7505 : if( edgePtr->forwardPtr == this )
349 : 3737 : return edgePtr->end();
350 [ + - ]: 3768 : else if( edgePtr->reversePtr == this )
351 : 3768 : return edgePtr->start();
352 : : else
353 : 0 : return 0;
354 : : }
355 : :
356 : 1232 : BSPTreePoly::Edge::~Edge()
357 : : {
358 [ + - ]: 1232 : delete startPtr;
359 [ + - ]: 1232 : delete endPtr;
360 [ - + ]: 1232 : delete forwardPtr;
361 [ - + ]: 1232 : delete reversePtr;
362 : 1232 : }
363 : :
364 : 0 : int BSPTreePoly::Edge::sense( BSPTreePoly::Face* face ) const
365 : : {
366 [ # # ][ # # ]: 0 : if( forwardPtr && forwardPtr->facePtr == face )
367 : 0 : return 1;
368 [ # # ][ # # ]: 0 : else if( reversePtr && reversePtr->facePtr == face )
369 : 0 : return -1;
370 : : else
371 : 0 : return 0;
372 : : }
373 : :
374 : 719 : BSPTreePoly::Face::~Face()
375 : : {
376 : 719 : BSPTreePoly::EdgeUse* nextEdgeUsePtr = usePtr;
377 [ + + ]: 3635 : while( nextEdgeUsePtr )
378 : : {
379 [ + - ]: 2916 : delete nextEdgeUsePtr; // This is tricky: ~EdgeUse() might change the value of usePtr
380 [ + + ][ + - ]: 2916 : if( usePtr && usePtr != nextEdgeUsePtr )
381 : 2197 : nextEdgeUsePtr = usePtr;
382 : : else
383 : 719 : nextEdgeUsePtr = 0;
384 : : }
385 : 719 : usePtr = 0;
386 : 719 : }
387 : :
388 : 55 : void BSPTreePoly::clear()
389 : : {
390 [ + + ]: 224 : while( faceList )
391 : : {
392 : 169 : Face* face = faceList;
393 : 169 : faceList = faceList->nextPtr;
394 [ + - ]: 169 : delete face;
395 : : }
396 : 55 : }
397 : :
398 : 28 : ErrorCode BSPTreePoly::set( const CartVect hex_corners[8] )
399 : : {
400 [ + - ]: 28 : clear();
401 : :
402 : : Vertex* vertices[8];
403 [ + + ]: 252 : for( int i = 0; i < 8; ++i )
404 [ + - ][ + - ]: 224 : vertices[i] = new Vertex( hex_corners[i] );
405 : :
406 : : Edge* edges[12];
407 : : #ifdef DEBUG_IDS
408 : : int start_id = Edge::nextID;
409 : : #endif
410 [ + + ]: 140 : for( int i = 0; i < 4; ++i )
411 : : {
412 : 112 : int j = ( i + 1 ) % 4;
413 [ + - ][ + - ]: 112 : edges[i] = new Edge( vertices[i], vertices[j] );
414 [ + - ][ + - ]: 112 : edges[i + 4] = new Edge( vertices[i], vertices[i + 4] );
415 [ + - ][ + - ]: 112 : edges[i + 8] = new Edge( vertices[i + 4], vertices[j + 4] );
416 : : }
417 : : #ifdef DEBUG_IDS
418 : : for( int i = 0; i < 12; ++i )
419 : : edges[i]->id = start_id++;
420 : : #endif
421 : :
422 : : static const int face_conn[6][4] = { { 0, 5, -8, -4 }, { 1, 6, -9, -5 }, { 2, 7, -10, -6 },
423 : : { 3, 4, -11, -7 }, { -3, -2, -1, -12 }, { 8, 9, 10, 11 } };
424 [ + + ]: 196 : for( int i = 0; i < 6; ++i )
425 : : {
426 [ + - ][ + - ]: 168 : faceList = new Face( faceList );
427 : 168 : EdgeUse* prev = 0;
428 [ + + ]: 840 : for( int j = 0; j < 4; ++j )
429 : : {
430 : 672 : int e = face_conn[i][j];
431 [ + + ]: 672 : if( e < 0 )
432 : : {
433 : 336 : e = ( -e ) % 12;
434 [ - + ]: 336 : assert( !edges[e]->reversePtr );
435 [ + + ][ + - ]: 336 : if( !prev ) { edges[e]->reversePtr = new EdgeUse( edges[e], faceList ); }
[ + - ]
436 : : else
437 : : {
438 [ + - ][ + - ]: 308 : edges[e]->reversePtr = new EdgeUse( edges[e] );
439 [ + - ]: 308 : edges[e]->reversePtr->insert_after( prev );
440 : : }
441 : 336 : prev = edges[e]->reversePtr;
442 : : }
443 : : else
444 : : {
445 [ - + ]: 336 : assert( !edges[e]->forwardPtr );
446 [ + + ][ + - ]: 336 : if( !prev ) { edges[e]->forwardPtr = new EdgeUse( edges[e], faceList ); }
[ + - ]
447 : : else
448 : : {
449 [ + - ][ + - ]: 196 : edges[e]->forwardPtr = new EdgeUse( edges[e] );
450 [ + - ]: 196 : edges[e]->forwardPtr->insert_after( prev );
451 : : }
452 : 336 : prev = edges[e]->forwardPtr;
453 : : }
454 : : }
455 : : }
456 : :
457 : 28 : return MB_SUCCESS;
458 : : }
459 : :
460 : 22 : void BSPTreePoly::get_faces( std::vector< const Face* >& face_list ) const
461 : : {
462 : 22 : face_list.clear();
463 [ + + ]: 162 : for( Face* face = faceList; face; face = face->nextPtr )
464 [ + - ]: 140 : face_list.push_back( face );
465 : 22 : }
466 : :
467 : 70 : void BSPTreePoly::get_vertices( const Face* face, std::vector< CartVect >& vertices ) const
468 : : {
469 : 70 : vertices.clear();
470 [ + - ][ - + ]: 70 : if( !face || !face->usePtr ) return;
471 : :
472 : 70 : EdgeUse* coedge = face->usePtr;
473 [ + + ]: 286 : do
474 : : {
475 : 286 : vertices.push_back( *coedge->end() );
476 : 286 : coedge = coedge->nextPtr;
477 : 286 : } while( coedge != face->usePtr );
478 : : }
479 : :
480 : 150 : double BSPTreePoly::Face::signed_volume() const
481 : : {
482 [ + - ]: 150 : CartVect sum( 0.0 );
483 [ + - ]: 150 : const CartVect* base = usePtr->start();
484 [ + - ][ + - ]: 150 : CartVect d1 = ( *usePtr->end() - *base );
485 [ + + ]: 598 : for( EdgeUse* coedge = usePtr->nextPtr; coedge != usePtr; coedge = coedge->nextPtr )
486 : : {
487 [ + - ][ + - ]: 448 : CartVect d2 = ( *coedge->end() - *base );
488 [ + - ][ + - ]: 448 : sum += d1 * d2;
489 : 448 : d1 = d2;
490 : : }
491 [ + - ]: 150 : return ( 1.0 / 6.0 ) * ( sum % *base );
492 : : }
493 : :
494 : 25 : double BSPTreePoly::volume() const
495 : : {
496 : 25 : double result = 0;
497 [ + + ]: 175 : for( Face* ptr = faceList; ptr; ptr = ptr->nextPtr )
498 : 150 : result += ptr->signed_volume();
499 : 25 : return result;
500 : : }
501 : :
502 : 102 : void BSPTreePoly::set_vertex_marks( int value )
503 : : {
504 [ + + ]: 750 : for( Face* face = faceList; face; face = face->nextPtr )
505 : : {
506 : 648 : EdgeUse* edge = face->usePtr;
507 [ + - ]: 2664 : do
508 : : {
509 : 2664 : edge->edgePtr->start()->markVal = value;
510 : 2664 : edge->edgePtr->end()->markVal = value;
511 : 2664 : edge = edge->nextPtr;
512 [ + + ]: 2664 : } while( edge && edge != face->usePtr );
513 : : }
514 : 102 : }
515 : : /*
516 : : static void merge_edges( BSPTreePoly::Edge* keep_edge,
517 : : BSPTreePoly::Edge* dead_edge )
518 : : {
519 : : // edges must share a vertex
520 : : BSPTreePoly::Vertex* dead_vtx = keep_edge->common(dead_edge);
521 : : assert(dead_vtx);
522 : : // vertex may have only two adjacent edges
523 : : BSPTreePoly::VertexUse* dead_vtxuse = dead_edge->use(dead_vtx);
524 : : assert(dead_vtxuse);
525 : : BSPTreePoly::VertexUse* keep_vtxuse = dead_vtxuse->nextPtr;
526 : : assert(keep_vtxuse);
527 : : assert(keep_vtxuse->edgePtr == keep_edge);
528 : : assert(keep_vtxuse->nextPtr == dead_vtxuse);
529 : : assert(keep_vtxuse->prevPtr == dead_vtxuse);
530 : : assert(dead_vtxuse->prevPtr == keep_vtxuse);
531 : :
532 : : // kept edge now ends with the kept vertex on the dead edge
533 : : keep_vtxuse->set_vertex( dead_edge->other(dead_vtx) );
534 : :
535 : : // destructors should take care of everything else
536 : : // (including removing dead edge from face loops)
537 : : delete dead_edge;
538 : : }
539 : : */
540 : 446 : static BSPTreePoly::Edge* split_edge( BSPTreePoly::Vertex*& new_vtx, BSPTreePoly::Edge* into_edge )
541 : : {
542 : : // split edge, creating new edge
543 [ + - ]: 446 : BSPTreePoly::Edge* new_edge = new BSPTreePoly::Edge( new_vtx, into_edge->end() );
544 : 446 : into_edge->endPtr->set_vertex( new_vtx ); // This call might delete new_vtx
545 : :
546 : : // update coedge loops in faces
547 [ + - ]: 446 : if( into_edge->forwardPtr )
548 : : {
549 [ + - ]: 446 : new_edge->forwardPtr = new BSPTreePoly::EdgeUse( new_edge );
550 : 446 : new_edge->forwardPtr->insert_after( into_edge->forwardPtr );
551 : : }
552 [ + - ]: 446 : if( into_edge->reversePtr )
553 : : {
554 [ + - ]: 446 : new_edge->reversePtr = new BSPTreePoly::EdgeUse( new_edge );
555 : 446 : new_edge->reversePtr->insert_before( into_edge->reversePtr );
556 : : }
557 : :
558 : 446 : return new_edge;
559 : : }
560 : :
561 : 450 : static BSPTreePoly::Face* split_face( BSPTreePoly::EdgeUse* start, BSPTreePoly::EdgeUse* end )
562 : : {
563 : 450 : BSPTreePoly::Face* face = start->facePtr;
564 [ - + ]: 450 : assert( face == end->facePtr );
565 [ + - ]: 450 : BSPTreePoly::Face* new_face = new BSPTreePoly::Face;
566 : 450 : BSPTreePoly::EdgeUse* keep_start = start->prevPtr;
567 : 450 : BSPTreePoly::EdgeUse* keep_end = end->nextPtr;
568 [ + + ]: 1764 : for( BSPTreePoly::EdgeUse* ptr = start; ptr != keep_end; ptr = ptr->nextPtr )
569 : : {
570 [ - + ]: 1314 : if( face->usePtr == ptr ) face->usePtr = keep_start;
571 : 1314 : ptr->facePtr = new_face;
572 : : }
573 : 450 : new_face->usePtr = start;
574 [ + - ]: 450 : BSPTreePoly::Edge* edge = new BSPTreePoly::Edge( start->start(), end->end() );
575 [ + - ]: 450 : edge->forwardPtr = new BSPTreePoly::EdgeUse( edge );
576 [ + - ]: 450 : edge->reversePtr = new BSPTreePoly::EdgeUse( edge );
577 : :
578 : 450 : edge->forwardPtr->facePtr = face;
579 : 450 : edge->forwardPtr->prevPtr = keep_start;
580 : 450 : keep_start->nextPtr = edge->forwardPtr;
581 : 450 : edge->forwardPtr->nextPtr = keep_end;
582 : 450 : keep_end->prevPtr = edge->forwardPtr;
583 : :
584 : 450 : edge->reversePtr->facePtr = new_face;
585 : 450 : edge->reversePtr->nextPtr = start;
586 : 450 : start->prevPtr = edge->reversePtr;
587 : 450 : edge->reversePtr->prevPtr = end;
588 : 450 : end->nextPtr = edge->reversePtr;
589 : :
590 : 450 : return new_face;
591 : : }
592 : :
593 : 102 : bool BSPTreePoly::cut_polyhedron( const CartVect& plane_normal, double plane_coeff )
594 : : {
595 : 102 : const double EPSILON = 1e-6; // points this close are considered coincident
596 : :
597 : : // scale epsilon rather than normalizing normal vector
598 : 102 : const double epsilon = EPSILON * ( plane_normal % plane_normal );
599 : :
600 : : // Classify all points above/below plane and destroy any faces
601 : : // that have no vertices below the plane.
602 : 102 : const int UNKNOWN = 0;
603 : 102 : const int ABOVE = 1;
604 : 102 : const int ON = 2;
605 : 102 : const int BELOW = 3;
606 : 102 : int num_above = 0;
607 : 102 : set_vertex_marks( UNKNOWN );
608 : :
609 : : // Classify all points above/below plane and
610 : : // split any edge that intersect the plane.
611 [ + + ]: 750 : for( Face* face = faceList; face; face = face->nextPtr )
612 : : {
613 : 648 : EdgeUse* edge = face->usePtr;
614 : :
615 [ + - ]: 3332 : do
616 : : {
617 : 3332 : Vertex* start = edge->edgePtr->start();
618 : 3332 : Vertex* end = edge->edgePtr->end();
619 : :
620 [ + + ]: 3332 : if( !start->markVal )
621 : : {
622 : 454 : double d = plane_normal % *start + plane_coeff;
623 [ + + ]: 454 : if( d * d <= epsilon )
624 : 1 : start->markVal = ON;
625 [ + + ]: 453 : else if( d < 0.0 )
626 : 230 : start->markVal = BELOW;
627 : : else
628 : : {
629 : 223 : start->markVal = ABOVE;
630 : 454 : ++num_above;
631 : : }
632 : : }
633 : :
634 [ + + ]: 3332 : if( !end->markVal )
635 : : {
636 : 434 : double d = plane_normal % *end + plane_coeff;
637 [ + + ]: 434 : if( d * d <= epsilon )
638 : 9 : end->markVal = ON;
639 [ + + ]: 425 : else if( d < 0.0 )
640 : 203 : end->markVal = BELOW;
641 : : else
642 : : {
643 : 222 : end->markVal = ABOVE;
644 : 434 : ++num_above;
645 : : }
646 : : }
647 : :
648 [ + + ][ + + ]: 3332 : if( ( end->markVal == ABOVE && start->markVal == BELOW ) ||
[ + + ]
649 [ + + ]: 1383 : ( end->markVal == BELOW && start->markVal == ABOVE ) )
650 : : {
651 [ + - ]: 446 : CartVect dir = *end - *start;
652 [ + - ][ + - ]: 446 : double t = -( plane_normal % *start + plane_coeff ) / ( dir % plane_normal );
653 [ + - ][ + - ]: 446 : Vertex* new_vtx = new Vertex( *start + t * dir );
[ + - ][ + - ]
654 : 446 : new_vtx->markVal = ON;
655 [ + - ]: 446 : split_edge( new_vtx, edge->edgePtr ); // This call might delete new_vtx
656 : 446 : end = new_vtx;
657 : : }
658 : :
659 : 3332 : edge = edge->nextPtr;
660 [ + + ]: 3332 : } while( edge && edge != face->usePtr );
661 : : }
662 : :
663 [ + + ]: 102 : if( !num_above ) return false;
664 : :
665 : : // Split faces
666 [ + + ]: 743 : for( Face* face = faceList; face; face = face->nextPtr )
667 : : {
668 : 642 : EdgeUse* edge = face->usePtr;
669 : :
670 : 642 : EdgeUse *split_start = 0, *split_end = 0, *other_split = 0;
671 [ + - ]: 3532 : do
672 : : {
673 [ + + ][ + + ]: 3532 : if( edge->end()->markVal == ON && edge->start()->markVal != ON )
[ + + ]
674 : : {
675 [ + + ]: 906 : if( !split_start )
676 : 456 : split_start = edge->nextPtr;
677 [ + - ]: 450 : else if( !split_end )
678 : 450 : split_end = edge;
679 : : else
680 : 906 : other_split = edge;
681 : : }
682 : :
683 : 3532 : edge = edge->nextPtr;
684 [ + + ]: 3532 : } while( edge && edge != face->usePtr );
685 : :
686 : : // If two vertices are on plane (but not every vertex)
687 : : // then split the face
688 [ + + ][ + - ]: 642 : if( split_end && !other_split )
689 : : {
690 [ - + ]: 450 : assert( split_start );
691 : 450 : Face* new_face = split_face( split_start, split_end );
692 : 450 : new_face->nextPtr = faceList;
693 : 450 : faceList = new_face;
694 : : }
695 : : }
696 : :
697 : : // Destroy all faces that are above the plane
698 : 101 : Face** lptr = &faceList;
699 [ + + ]: 1193 : while( *lptr )
700 : : {
701 : 1092 : EdgeUse* edge = ( *lptr )->usePtr;
702 : 1092 : bool some_above = false;
703 [ + - ]: 2441 : do
704 : : {
705 [ + + ]: 2991 : if( edge->start()->markVal == ABOVE )
706 : : {
707 : 550 : some_above = true;
708 : 550 : break;
709 : : }
710 : 2441 : edge = edge->nextPtr;
711 [ + + ]: 2441 : } while( edge && edge != ( *lptr )->usePtr );
712 : :
713 [ + + ]: 1092 : if( some_above )
714 : : {
715 : 550 : Face* dead = *lptr;
716 : 550 : *lptr = ( *lptr )->nextPtr;
717 [ + - ]: 550 : delete dead;
718 : : }
719 : : else
720 : : {
721 : 542 : lptr = &( ( *lptr )->nextPtr );
722 : : }
723 : : }
724 : :
725 : : // Construct a new face in the cut plane
726 : :
727 : : // First find an edge to start at
728 : 101 : Edge* edge_ptr = 0;
729 [ + - ][ + + ]: 203 : for( Face* face = faceList; face && !edge_ptr; face = face->nextPtr )
730 : : {
731 : 102 : EdgeUse* co_edge = face->usePtr;
732 [ + - ]: 288 : do
733 : : {
734 [ + + ]: 389 : if( 0 == co_edge->edgePtr->other( co_edge ) )
735 : : {
736 : 101 : edge_ptr = co_edge->edgePtr;
737 : 101 : break;
738 : : }
739 : 288 : co_edge = co_edge->nextPtr;
740 [ + + ]: 288 : } while( co_edge && co_edge != face->usePtr );
741 : : }
742 [ - + ]: 101 : if( !edge_ptr ) return false;
743 : :
744 : : // Constuct new face and first CoEdge
745 [ + - ]: 101 : faceList = new Face( faceList );
746 : : Vertex *next_vtx, *start_vtx;
747 : : EdgeUse* prev_coedge;
748 [ + + ]: 101 : if( edge_ptr->forwardPtr )
749 : : {
750 : 12 : next_vtx = edge_ptr->start();
751 : 12 : start_vtx = edge_ptr->end();
752 [ - + ]: 12 : assert( !edge_ptr->reversePtr );
753 [ + - ]: 12 : prev_coedge = edge_ptr->reversePtr = new EdgeUse( edge_ptr, faceList );
754 : : }
755 : : else
756 : : {
757 : 89 : next_vtx = edge_ptr->end();
758 : 89 : start_vtx = edge_ptr->start();
759 [ + - ]: 89 : prev_coedge = edge_ptr->forwardPtr = new EdgeUse( edge_ptr, faceList );
760 : : }
761 : :
762 : : // Construct coedges until loop is closed
763 [ + + ]: 452 : while( next_vtx != start_vtx )
764 : : {
765 : : // find next edge adjacent to vertex with only one adjacent face
766 : 351 : VertexUse* this_use = edge_ptr->use( next_vtx );
767 : 351 : VertexUse* use = this_use->nextPtr;
768 [ + - ]: 549 : while( use != this_use )
769 : : {
770 [ + + ]: 549 : if( use->edgePtr->forwardPtr == 0 )
771 : : {
772 : 124 : edge_ptr = use->edgePtr;
773 [ - + ]: 124 : assert( edge_ptr->start() == next_vtx );
774 : 124 : next_vtx = edge_ptr->end();
775 [ + - ]: 124 : edge_ptr->forwardPtr = new EdgeUse( edge_ptr );
776 : 124 : edge_ptr->forwardPtr->insert_after( prev_coedge );
777 : 124 : prev_coedge = edge_ptr->forwardPtr;
778 : 124 : break;
779 : : }
780 [ + + ]: 425 : else if( use->edgePtr->reversePtr == 0 )
781 : : {
782 : 227 : edge_ptr = use->edgePtr;
783 [ - + ]: 227 : assert( edge_ptr->end() == next_vtx );
784 : 227 : next_vtx = edge_ptr->start();
785 [ + - ]: 227 : edge_ptr->reversePtr = new EdgeUse( edge_ptr );
786 : 227 : edge_ptr->reversePtr->insert_after( prev_coedge );
787 : 227 : prev_coedge = edge_ptr->reversePtr;
788 : 227 : break;
789 : : }
790 : :
791 : 198 : use = use->nextPtr;
792 [ - + ]: 198 : assert( use != this_use ); // failed to close loop!
793 : : }
794 : : }
795 : :
796 : 102 : return true;
797 : : }
798 : :
799 : 27 : bool BSPTreePoly::is_valid() const
800 : : {
801 [ + - ]: 27 : std::set< Face* > list_faces;
802 : :
803 : 27 : int i = 0;
804 [ + + ]: 190 : for( Face* ptr = faceList; ptr; ptr = ptr->nextPtr )
805 : : {
806 [ - + ]: 163 : if( ++i > 10000 ) return false;
807 [ + - ][ - + ]: 163 : if( !list_faces.insert( ptr ).second ) return false;
808 : : }
809 : :
810 [ + - ]: 54 : std::set< Vertex* > vertices;
811 [ + + ]: 190 : for( Face* face = faceList; face; face = face->nextPtr )
812 : : {
813 : 163 : i = 0;
814 : 163 : EdgeUse* coedge = face->usePtr;
815 [ + + ]: 652 : do
816 : : {
817 [ - + ]: 652 : if( ++i > 10000 ) return false;
818 : :
819 [ - + ]: 652 : if( coedge->facePtr != face ) return false;
820 : :
821 : 652 : Edge* edge = coedge->edgePtr;
822 [ + - ][ - + ]: 652 : if( !edge->startPtr || !edge->endPtr ) return false;
823 : :
824 [ + - ][ + - ]: 652 : vertices.insert( edge->start() );
825 [ + - ][ + - ]: 652 : vertices.insert( edge->end() );
826 : :
827 : : EdgeUse* other;
828 [ + + ]: 652 : if( edge->forwardPtr == coedge )
829 : 326 : other = edge->reversePtr;
830 [ - + ]: 326 : else if( edge->reversePtr != coedge )
831 : 0 : return false;
832 : : else
833 : 326 : other = edge->forwardPtr;
834 [ - + ]: 652 : if( !other ) return false;
835 [ + - ][ + - ]: 652 : if( list_faces.find( other->facePtr ) == list_faces.end() ) return false;
[ - + ]
836 : :
837 : 652 : EdgeUse* next = coedge->nextPtr;
838 [ - + ]: 652 : if( next->prevPtr != coedge ) return false;
839 [ + - ][ + - ]: 652 : if( coedge->end() != next->start() ) return false;
[ - + ]
840 : :
841 : 652 : coedge = next;
842 : 652 : } while( coedge != face->usePtr );
843 : : }
844 : :
845 [ + - ][ + - ]: 244 : for( std::set< Vertex* >::iterator j = vertices.begin(); j != vertices.end(); ++j )
[ + + ]
846 : : {
847 [ + - ]: 217 : Vertex* vtx = *j;
848 : :
849 : 217 : i = 0;
850 : 217 : VertexUse* use = vtx->usePtr;
851 [ + + ]: 652 : do
852 : : {
853 [ - + ]: 652 : if( ++i > 10000 ) return false;
854 : :
855 [ - + ]: 652 : if( use->vtxPtr != vtx ) return false;
856 : :
857 : 652 : Edge* edge = use->edgePtr;
858 [ - + ]: 652 : if( !edge ) return false;
859 [ + + ][ - + ]: 652 : if( edge->startPtr != use && edge->endPtr != use ) return false;
860 : :
861 : 652 : VertexUse* next = use->nextPtr;
862 [ - + ]: 652 : if( next->prevPtr != use ) return false;
863 : :
864 : 652 : use = next;
865 : 652 : } while( use != vtx->usePtr );
866 : : }
867 : :
868 : 54 : return true;
869 : : }
870 : :
871 : 20 : bool BSPTreePoly::is_point_contained( const CartVect& point ) const
872 : : {
873 [ - + ]: 20 : if( !faceList ) // empty (zero-dimension) polyhedron
874 : 0 : return false;
875 : :
876 : 20 : const double EPSILON = 1e-6;
877 : : // Test that point is below the plane of each face
878 : : // NOTE: This will NOT work for polyhedra w/ concavities
879 [ + + ]: 140 : for( Face* face = faceList; face; face = face->nextPtr )
880 : : {
881 : : Vertex *pt1, *pt2, *pt3;
882 [ + - ]: 120 : pt1 = face->usePtr->start();
883 [ + - ]: 120 : pt2 = face->usePtr->end();
884 [ + - ]: 120 : pt3 = face->usePtr->nextPtr->end();
885 : :
886 [ - + ]: 120 : if( pt3 == pt1 ) // degenerate
887 : 0 : continue;
888 : :
889 [ + - ][ + - ]: 120 : CartVect norm = ( *pt3 - *pt2 ) * ( *pt1 - *pt2 );
[ + - ]
890 [ + - ]: 120 : double coeff = -( norm % *pt2 );
891 [ + - ][ - + ]: 120 : if( ( norm % point + coeff ) > EPSILON ) // if above plane, with some -epsilon
892 : 120 : return false;
893 : : }
894 : :
895 : 20 : return true;
896 : : }
897 : :
898 : : } // namespace moab
|