1 /* GTS - Library for the manipulation of triangulated surfaces
2 * Copyright (C) 1999--2002 Stéphane Popinet
4 * This library is free software; you can redistribute it and/or
5 * modify it under the terms of the GNU Library General Public
6 * License as published by the Free Software Foundation; either
7 * version 2 of the License, or (at your option) any later version.
9 * This library is distributed in the hope that it will be useful,
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 * Library General Public License for more details.
14 * You should have received a copy of the GNU Library General Public
15 * License along with this library; if not, write to the
16 * Free Software Foundation, Inc., 59 Temple Place - Suite 330,
17 * Boston, MA 02111-1307, USA.
24 /*#define DEBUG_BOOLEAN*/
25 /*#define CHECK_ORIENTED*/
28 # include "gts-private.h"
31 static void surface_inter_destroy (GtsObject
* object
)
33 GtsSurfaceInter
* si
= GTS_SURFACE_INTER (object
);
35 gts_object_destroy (GTS_OBJECT (si
->s1
));
36 gts_object_destroy (GTS_OBJECT (si
->s2
));
37 g_slist_free (si
->edges
);
39 (* GTS_OBJECT_CLASS (gts_surface_inter_class ())->parent_class
->destroy
)
43 static void surface_inter_class_init (GtsObjectClass
* klass
)
45 klass
->destroy
= surface_inter_destroy
;
48 static void surface_inter_init (GtsSurfaceInter
* si
)
50 si
->s1
= si
->s2
= NULL
;
55 * gts_surface_inter_class:
57 * Returns: the #GtsSurfaceInterClass.
59 GtsSurfaceInterClass
* gts_surface_inter_class (void)
61 static GtsSurfaceInterClass
* klass
= NULL
;
64 GtsObjectClassInfo surface_inter_info
= {
66 sizeof (GtsSurfaceInter
),
67 sizeof (GtsSurfaceInterClass
),
68 (GtsObjectClassInitFunc
) surface_inter_class_init
,
69 (GtsObjectInitFunc
) surface_inter_init
,
73 klass
= gts_object_class_new (gts_object_class (), &surface_inter_info
);
79 /* EdgeInter: Header */
81 typedef struct _EdgeInter EdgeInter
;
86 GtsTriangle
* t1
, * t2
;
89 #define EDGE_INTER(obj) GTS_OBJECT_CAST (obj,\
92 #define IS_EDGE_INTER(obj) (gts_object_is_from_class (obj,\
95 static GtsEdgeClass
* edge_inter_class (void);
96 static EdgeInter
* edge_inter_new (GtsVertex
* v1
, GtsVertex
* v2
,
97 GtsTriangle
* t1
, GtsTriangle
* t2
);
99 /* EdgeInter: Object */
101 static GtsEdgeClass
* edge_inter_class (void)
103 static GtsEdgeClass
* klass
= NULL
;
106 GtsObjectClassInfo edge_inter_info
= {
109 sizeof (GtsEdgeClass
),
110 (GtsObjectClassInitFunc
) NULL
,
111 (GtsObjectInitFunc
) NULL
,
112 (GtsArgSetFunc
) NULL
,
115 klass
= gts_object_class_new (GTS_OBJECT_CLASS (gts_constraint_class ()),
122 static EdgeInter
* edge_inter_new (GtsVertex
* v1
, GtsVertex
* v2
,
123 GtsTriangle
* t1
, GtsTriangle
* t2
)
127 object
= EDGE_INTER (gts_edge_new (GTS_EDGE_CLASS (edge_inter_class ()),
136 static void write_surface_graph (GtsSurface
* s
, FILE * fp
)
140 static void add_to_list (gpointer data
, GSList
** l
) {
141 *l
= g_slist_prepend (*l
, data
);
144 gts_surface_foreach_vertex (s
, (GtsFunc
) gts_object_reset_reserved
, NULL
);
145 gts_surface_foreach_edge (s
, (GtsFunc
) gts_object_reset_reserved
, NULL
);
146 gts_surface_foreach_edge (s
, (GtsFunc
) add_to_list
, &l
);
147 g
= gts_segments_graph_new (gts_graph_class (), l
);
148 gts_graph_write_dot (g
, fp
);
149 gts_object_destroy (GTS_OBJECT (g
));
154 static GtsPoint
* segment_triangle_intersection (GtsSegment
* s
,
156 GtsPointClass
* klass
)
158 GtsPoint
* A
, * B
, * C
, * D
, * E
;
159 gint ABCE
, ABCD
, ADCE
, ABDE
, BCDE
;
160 GtsEdge
* AB
, * BC
, * CA
;
163 g_return_val_if_fail (s
!= NULL
, NULL
);
164 g_return_val_if_fail (t
!= NULL
, NULL
);
165 g_return_val_if_fail (klass
!= NULL
, NULL
);
167 gts_triangle_vertices_edges (t
, NULL
,
172 D
= GTS_POINT (s
->v1
);
173 E
= GTS_POINT (s
->v2
);
175 ABCE
= gts_point_orientation_3d_sos (A
, B
, C
, E
);
176 ABCD
= gts_point_orientation_3d_sos (A
, B
, C
, D
);
177 if (ABCE
< 0 || ABCD
> 0) {
181 tmpp
= E
; E
= D
; D
= tmpp
;
182 tmp
= ABCE
; ABCE
= ABCD
; ABCD
= tmp
;
184 if (ABCE
< 0 || ABCD
> 0)
186 ADCE
= gts_point_orientation_3d_sos (A
, D
, C
, E
);
189 ABDE
= gts_point_orientation_3d_sos (A
, B
, D
, E
);
192 BCDE
= gts_point_orientation_3d_sos (B
, C
, D
, E
);
195 a
= gts_point_orientation_3d (A
, B
, C
, E
);
196 b
= gts_point_orientation_3d (A
, B
, C
, D
);
199 return gts_point_new (klass
,
200 E
->x
+ c
*(D
->x
- E
->x
),
201 E
->y
+ c
*(D
->y
- E
->y
),
202 E
->z
+ c
*(D
->z
- E
->z
));
204 /* D and E are contained within ABC */
207 "segment: %p:%s triangle: %p:%s intersection\n"
208 "D and E contained in ABC\n",
209 s
, GTS_NEDGE (s
)->name
, t
, GTS_NFACE (t
)->name
);
212 return gts_point_new (klass
,
218 static gint
triangle_triangle_orientation (GtsPoint
* p1
,
219 GtsPoint
* p2
, GtsPoint
* p3
,
220 GtsPoint
* p4
, GtsPoint
* p5
,
223 gint o4
= 0, o5
= 0, o6
= 0;
225 if (p4
!= p1
&& p4
!= p2
&& p4
!= p3
)
226 o4
= gts_point_orientation_3d_sos (p1
, p2
, p3
, p4
);
227 if (p5
!= p1
&& p5
!= p2
&& p5
!= p3
)
228 o5
= gts_point_orientation_3d_sos (p1
, p2
, p3
, p5
);
231 if (p6
!= p1
&& p6
!= p2
&& p6
!= p3
)
232 o6
= gts_point_orientation_3d_sos (p1
, p2
, p3
, p6
);
233 if (o4
*o6
< 0 || o5
*o6
< 0)
241 static gint
triangle_point_orientation (GtsTriangle
* t1
,
246 GtsPoint
* p1
= GTS_POINT (GTS_SEGMENT (t1
->e1
)->v1
);
247 GtsPoint
* p2
= GTS_POINT (GTS_SEGMENT (t1
->e1
)->v2
);
248 GtsPoint
* p3
= GTS_POINT (gts_triangle_vertex (t1
));
249 GtsPoint
* p4
= GTS_POINT (GTS_SEGMENT (t2
->e1
)->v1
);
250 GtsPoint
* p5
= GTS_POINT (GTS_SEGMENT (t2
->e1
)->v2
);
251 GtsPoint
* p6
= GTS_POINT (gts_triangle_vertex (t2
));
252 gint o
= triangle_triangle_orientation (p1
, p2
, p3
, p4
, p5
, p6
);
256 o
= triangle_triangle_orientation (p4
, p5
, p6
, p1
, p2
, p3
);
258 gint o2
= gts_point_orientation_3d_sos (p4
, p5
, p6
, p
);
265 static void add_edge_inter (GtsEdge
* e
,
269 GtsVertex
* ev1
= GTS_SEGMENT (e
)->v1
, * ev2
= GTS_SEGMENT (e
)->v2
;
270 GList
* i
= GTS_OBJECT (e
)->reserved
;
272 GTS_OBJECT (v
)->reserved
= t
;
274 GTS_OBJECT (e
)->reserved
= g_list_prepend (NULL
, v
);
276 fprintf (stderr
, "add_edge_inter: inserting %p (%p,%p)\n", v
, e
, t
);
280 GtsPoint
* p1
= GTS_POINT (GTS_SEGMENT (t
->e1
)->v1
);
281 GtsPoint
* p2
= GTS_POINT (GTS_SEGMENT (t
->e1
)->v2
);
282 GtsPoint
* p3
= GTS_POINT (gts_triangle_vertex (t
));
283 gint o1
, oref
= gts_point_orientation_3d_sos (p1
, p2
, p3
, GTS_POINT (ev1
));
287 gint o2
= triangle_point_orientation (t
, GTS_OBJECT (i
->data
)->reserved
,
288 oref
, GTS_POINT (ev1
));
292 g_warning ("add_edge_inter: safe sign evaluation failed\n");
294 o2
= gts_point_orientation_3d_sos (p1
, p2
, p3
, i
->data
);
304 GList
* n
= g_list_prepend (NULL
, v
);
311 GTS_OBJECT (e
)->reserved
= n
;
316 g_assert (o1
*gts_point_orientation_3d_sos (p1
, p2
, p3
, GTS_POINT (ev2
))
318 GTS_OBJECT (e
)->reserved
= g_list_append (GTS_OBJECT (e
)->reserved
, v
);
322 "add_edge_inter: inserting %p (%p,%p) between %p and %p\n",
324 i
= GTS_OBJECT (e
)->reserved
;
326 fprintf (stderr
, " %p", i
->data
);
329 fprintf (stderr
, "\n");
334 static GtsVertex
* intersects (GtsEdge
* e
,
338 GList
* i
= GTS_OBJECT (e
)->reserved
;
342 if (GTS_OBJECT (i
->data
)->reserved
== t
)
347 v
= GTS_VERTEX (segment_triangle_intersection (GTS_SEGMENT (e
), t
,
348 GTS_POINT_CLASS (s
->vertex_class
)));
351 if (GTS_IS_NVERTEX (v
) && GTS_IS_NEDGE (e
) && GTS_IS_NFACE (t
) &&
352 GTS_NVERTEX (v
)->name
[0] == '\0')
353 g_snprintf (GTS_NVERTEX (v
)->name
, GTS_NAME_LENGTH
, "%s|%s",
354 GTS_NEDGE (e
)->name
, GTS_NFACE (t
)->name
);
356 if (s
->vertex_class
->intersection_attributes
)
357 (*s
->vertex_class
->intersection_attributes
)
358 (v
, GTS_OBJECT (e
), GTS_OBJECT (t
));
359 add_edge_inter (e
, t
, v
);
364 /* see figure misc/orientation.fig */
365 static gint
intersection_orientation (GtsTriangle
* t1
,
369 GtsVertex
* v1
, * v2
, * v3
;
371 GtsVertex
* v4
, * v5
, * v6
;
373 gts_triangle_vertices_edges (t1
, e
, &v1
, &v2
, &v3
, &e
, &e2
, &e3
);
374 gts_triangle_vertices (t2
, &v4
, &v5
, &v6
);
376 return gts_point_orientation_3d_sos (GTS_POINT (v4
),
382 #define UPDATE_ORIENTATION if (o > 0) { vi2 = v; /* e2 = e; */ } else { vi2 = vi1;\
387 static void intersect_edges (GtsBBox
* bb1
, GtsBBox
* bb2
,
388 GtsSurfaceInter
* si
)
390 GtsSurface
* s1
= GTS_OBJECT (si
->s1
)->reserved
;
391 GtsTriangle
* t1
= GTS_TRIANGLE (bb1
->bounded
);
392 GtsTriangle
* t2
= GTS_TRIANGLE (bb2
->bounded
);
393 GtsVertex
* v
, * vi1
= NULL
, * vi2
= NULL
;
394 //GtsEdge * e1 = NULL, * e2 = NULL, * e;
396 vi1
= intersects (t2
->e1
, t1
, s1
);
398 v
= intersects (t2
->e2
, t1
, s1
);
405 gint o
= intersection_orientation (t2
, t2
->e2
, t1
);
409 v
= intersects (t2
->e3
, t1
, s1
);
416 gint o
= intersection_orientation (t2
, t2
->e3
, t1
);
421 v
= intersects (t1
->e1
, t2
, s1
);
428 gint o
= - intersection_orientation (t1
, t1
->e1
, t2
);
433 v
= intersects (t1
->e2
, t2
, s1
);
440 gint o
= - intersection_orientation (t1
, t1
->e2
, t2
);
445 v
= intersects (t1
->e3
, t2
, s1
);
452 gint o
= - intersection_orientation (t1
, t1
->e3
, t2
);
457 g_assert ((!vi1
&& !vi2
) || (vi1
&& vi2
));
459 GtsEdge
* e
= GTS_EDGE (edge_inter_new (vi1
, vi2
, t1
, t2
));
462 fprintf (stderr
, "creating constraint %p: %p->%p: %p/%p\n",
463 e
, vi1
, vi2
, t1
, t2
);
465 gts_surface_add_face (si
->s1
, GTS_FACE (t1
));
466 gts_surface_add_face (si
->s2
, GTS_FACE (t2
));
467 si
->edges
= g_slist_prepend (si
->edges
, e
);
468 GTS_OBJECT (t1
)->reserved
= g_slist_prepend (GTS_OBJECT (t1
)->reserved
, e
);
469 GTS_OBJECT (t2
)->reserved
= g_slist_prepend (GTS_OBJECT (t2
)->reserved
, e
);
473 static GtsSurfaceInter
* surface_inter_new (GtsSurfaceInterClass
* klass
,
479 GtsSurfaceInter
* si
;
481 si
= GTS_SURFACE_INTER (gts_object_new (GTS_OBJECT_CLASS (klass
)));
482 si
->s1
= gts_surface_new (gts_surface_class (),
486 GTS_OBJECT (si
->s1
)->reserved
= s1
;
487 si
->s2
= gts_surface_new (gts_surface_class (),
491 GTS_OBJECT (si
->s2
)->reserved
= s2
;
492 gts_bb_tree_traverse_overlapping (faces_tree1
, faces_tree2
,
493 (GtsBBTreeTraverseFunc
) intersect_edges
,
499 static void free_slist (GtsObject
* o
)
501 g_slist_free (o
->reserved
);
505 static void free_glist (GtsObject
* o
)
507 g_list_foreach (o
->reserved
, (GFunc
) gts_object_reset_reserved
, NULL
);
508 g_list_free (o
->reserved
);
513 * gts_surface_intersection:
514 * @s1: a #GtsSurface.
515 * @s2: a #GtsSurface.
516 * @faces_tree1: a bounding box tree (see gts_bb_tree_new()) for
518 * @faces_tree2: a bounding box tree for the faces of @s2.
520 * Returns: a list of #GtsEdge defining the curve intersection of the
523 GSList
* gts_surface_intersection (GtsSurface
* s1
,
528 GtsSurfaceInter
* si
;
531 g_return_val_if_fail (s1
!= NULL
, NULL
);
532 g_return_val_if_fail (s2
!= NULL
, NULL
);
533 g_return_val_if_fail (faces_tree1
!= NULL
, NULL
);
534 g_return_val_if_fail (faces_tree2
!= NULL
, NULL
);
536 si
= surface_inter_new (gts_surface_inter_class (),
537 s1
, s2
, faces_tree1
, faces_tree2
);
539 gts_surface_foreach_face (si
->s1
, (GtsFunc
) free_slist
, NULL
);
540 gts_surface_foreach_face (si
->s2
, (GtsFunc
) free_slist
, NULL
);
541 gts_surface_foreach_edge (si
->s1
, (GtsFunc
) free_glist
, NULL
);
542 gts_surface_foreach_edge (si
->s2
, (GtsFunc
) free_glist
, NULL
);
545 gts_object_destroy (GTS_OBJECT (si
));
551 INTERIOR
= 1 << (GTS_USER_FLAG
),
552 RELEVANT
= 1 << (GTS_USER_FLAG
+ 1)
555 #define IS_SET(s, f) ((GTS_OBJECT_FLAGS (s) & (f)) != 0)
556 #define SET(s, f) (GTS_OBJECT_FLAGS (s) |= (f))
557 #define UNSET(s, f) (GTS_OBJECT_FLAGS (s) &= ~(f))
558 #define NEXT(s) (GTS_OBJECT (s)->reserved)
561 static void print_segment (GtsSegment
* s
)
563 fprintf (stderr
, "%p: %s->%s ", s
,
564 GTS_NVERTEX (s
->v1
)->name
,
565 GTS_NVERTEX (s
->v2
)->name
);
567 GtsSegment
* next
= NEXT (s
);
569 fprintf (stderr
, "next %p: %s->%s\n", next
,
570 GTS_NVERTEX (next
->v1
)->name
,
571 GTS_NVERTEX (next
->v2
)->name
);
574 fprintf (stderr
, "next NULL\n");
577 static void write_nodes (GSList
* i
, GHashTable
* hash
, guint
* nn
,
581 GtsSegment
* s
= i
->data
;
583 if (!g_hash_table_lookup (hash
, s
->v1
)) {
584 fprintf (fp
, " %u [ label = \"%p\" ];\n", *nn
, s
->v1
);
585 g_hash_table_insert (hash
, s
->v1
, GUINT_TO_POINTER ((*nn
)++));
587 if (!g_hash_table_lookup (hash
, s
->v2
)) {
588 fprintf (fp
, " %u [ label = \"%p\" ];\n", *nn
, s
->v2
);
589 g_hash_table_insert (hash
, s
->v2
, GUINT_TO_POINTER ((*nn
)++));
595 static void write_edges (GSList
* i
, GHashTable
* hash
,
596 GtsSurface
* surface
,
600 GtsSegment
* s
= i
->data
;
602 fprintf (fp
, " %u -> %u [ label = \"%p:%d\" ];\n",
603 GPOINTER_TO_UINT (g_hash_table_lookup (hash
, s
->v1
)),
604 GPOINTER_TO_UINT (g_hash_table_lookup (hash
, s
->v2
)),
606 gts_edge_face_number (GTS_EDGE (s
), surface
));
611 static void write_graph (GSList
* boundary
, GSList
* interior
,
612 GtsSurface
* surface
,
615 GHashTable
* hash
= g_hash_table_new (NULL
, NULL
);
618 fprintf (fp
, "digraph oriented_curve {\n");
619 write_nodes (boundary
, hash
, &nn
, fp
);
620 write_nodes (interior
, hash
, &nn
, fp
);
621 write_edges (boundary
, hash
, surface
, fp
);
622 fprintf (fp
, " edge [ color = red ];\n");
623 write_edges (interior
, hash
, surface
, fp
);
625 g_hash_table_destroy (hash
);
628 static void write_graph1 (GtsSegment
* start
, GSList
* i
,
629 GtsSurface
* surface
,
632 GSList
* boundary
= NULL
, * interior
= NULL
;
633 GtsSegment
* s
= start
;
636 boundary
= g_slist_prepend (boundary
, s
);
638 } while (s
!= start
);
640 if (IS_SET (i
->data
, INTERIOR
))
641 interior
= g_slist_prepend (interior
, i
->data
);
644 write_graph (boundary
, interior
, surface
, fp
);
645 g_slist_free (boundary
);
646 g_slist_free (interior
);
649 static void print_loop (GtsSegment
* start
, FILE * fp
)
651 GtsSegment
* s
= start
;
654 fprintf (fp
, " %p: %p:%s -> %p:%s\n",
656 s
->v1
, GTS_NVERTEX (s
->v1
)->name
,
657 s
->v2
, GTS_NVERTEX (s
->v2
)->name
);
659 } while (s
!= start
&& s
!= NULL
);
662 static void draw_vector (GtsPoint
* p1
, GtsPoint
* p2
, FILE * fp
)
664 gdouble x
= p2
->x
- p1
->x
;
665 gdouble y
= p2
->y
- p1
->y
;
666 gdouble z
= p2
->z
- p1
->z
;
668 fprintf (fp
, "VECT 1 3 0 3 0 %g %g %g %g %g %g %g %g %g\n",
669 p1
->x
+ x
- (x
- y
/2.)/5.,
670 p1
->y
+ y
- (x
/2. + y
)/5.,
671 p1
->z
+ z
- (x
/2. + z
)/5.,
675 p1
->x
+ x
- (x
+ y
/2.)/5.,
676 p1
->y
+ y
+ (x
/2. - y
)/5.,
677 p1
->z
+ z
+ (x
/2. - z
)/5.);
678 fprintf (fp
, "VECT 1 2 0 2 0 %g %g %g %g %g %g\n",
685 static void draw_vector1 (GtsPoint
* p1
, GtsPoint
* p2
, GtsPoint
* o
,
688 gdouble x1
= o
->x
+ 0.9*(p1
->x
- o
->x
);
689 gdouble y1
= o
->y
+ 0.9*(p1
->y
- o
->y
);
690 gdouble z1
= o
->z
+ 0.9*(p1
->z
- o
->z
);
691 gdouble x2
= o
->x
+ 0.9*(p2
->x
- o
->x
);
692 gdouble y2
= o
->y
+ 0.9*(p2
->y
- o
->y
);
693 gdouble z2
= o
->z
+ 0.9*(p2
->z
- o
->z
);
698 fprintf (fp
, "VECT 1 3 0 3 0 %g %g %g %g %g %g %g %g %g\n",
699 x1
+ x
- (x
- y
/2.)/5.,
700 y1
+ y
- (x
/2. + y
)/5.,
701 z1
+ z
- (x
/2. + z
)/5.,
705 x1
+ x
- (x
+ y
/2.)/5.,
706 y1
+ y
+ (x
/2. - y
)/5.,
707 z1
+ z
+ (x
/2. - z
)/5.);
708 fprintf (fp
, "VECT 1 2 0 2 0 %g %g %g %g %g %g\n",
715 static void write_segments (GSList
* boundary
, GSList
* interior
,
718 GSList
* i
= boundary
;
720 fprintf (fp
, "LIST {\n");
722 GSList
* inext
= i
->next
;
723 GtsSegment
* s
= i
->data
;
724 GtsSegment
* next
= inext
? inext
->data
: boundary
->data
;
725 GtsVertex
* v1
, * v2
;
727 if (s
->v1
!= next
->v1
&& s
->v1
!= next
->v2
) {
735 draw_vector (GTS_POINT (v1
), GTS_POINT (v2
), fp
);
740 GtsSegment
* s
= i
->data
;
742 draw_vector (GTS_POINT (s
->v1
), GTS_POINT (s
->v2
), fp
);
748 static void write_loops (GSList
* i
, FILE * fp
)
753 GtsSegment
* start
= i
->data
, * s
;
757 fprintf (fp
, "(geometry \"loop%d\" = LIST {\n", nl
++);
759 os
.x
= os
.y
= os
.z
= 0.;
762 GtsSegment
* next
= NEXT (s
);
765 if (s
->v1
!= next
->v1
&& s
->v1
!= next
->v2
)
766 p
= GTS_POINT (s
->v1
);
768 p
= GTS_POINT (s
->v2
);
769 os
.x
+= p
->x
; os
.y
+= p
->y
; os
.z
+= p
->z
; n
++;
771 } while (s
!= start
);
772 os
.x
/= n
; os
.y
/= n
; os
.z
/= n
;
776 GtsSegment
* next
= NEXT (s
);
778 if (s
->v1
!= next
->v1
&& s
->v1
!= next
->v2
)
779 draw_vector1 (GTS_POINT (s
->v1
), GTS_POINT (s
->v2
), &os
, fp
);
781 draw_vector1 (GTS_POINT (s
->v2
), GTS_POINT (s
->v1
), &os
, fp
);
783 } while (s
!= start
);
785 fprintf (fp
, "})\n");
791 #define NAME(v) (GTS_IS_NVERTEX (v) ? GTS_NVERTEX (v)->name : "")
794 static GtsSegment
* prev_flag (GtsSegment
* s
, CurveFlag flag
)
796 GSList
* i
= s
->v1
->segments
;
799 if (i
->data
!= s
&& IS_SET (i
->data
, flag
))
806 static GtsSegment
* next_flag (GtsSegment
* s
, CurveFlag flag
)
808 GSList
* i
= s
->v2
->segments
;
811 if (i
->data
!= s
&& IS_SET (i
->data
, flag
))
818 static GtsSegment
* next_interior (GtsVertex
* v
)
820 GSList
* i
= v
->segments
;
823 GtsSegment
* s
= i
->data
;
825 if (s
->v1
== v
&& IS_SET (s
, INTERIOR
))
832 static GtsSegment
* prev_interior (GtsVertex
* v
)
834 GSList
* i
= v
->segments
;
837 GtsSegment
* s
= i
->data
;
839 if (s
->v2
== v
&& IS_SET (s
, INTERIOR
))
846 static GtsSegment
* reverse (GtsSegment
* start
,
850 GtsSegment
* s
= start
, * prev
= NULL
, * rprev
= NULL
;
851 GtsSegment
* rstart
= NULL
, * rstart1
= NULL
;
856 g_assert (IS_EDGE_INTER (s
));
857 rs
= GTS_SEGMENT (edge_inter_new (s
->v2
, s
->v1
,
858 EDGE_INTER (s
)->t1
, EDGE_INTER (s
)->t2
));
862 else if (rstart1
== NULL
)
870 } while (s
!= NULL
&& s
!= start
);
872 NEXT (rstart
) = rprev
;
876 NEXT (rstart
) = start
;
883 static GSList
* interior_loops (GSList
* interior
)
885 GSList
* i
= interior
;
886 GSList
* loops
= NULL
;
890 GtsSegment
* s
= i
->data
;
892 if (IS_SET (s
, RELEVANT
)) {
893 GtsSegment
* start
= s
, * end
;
896 GtsSegment
* next
= next_flag (s
, INTERIOR
);
901 } while (s
!= NULL
&& s
!= start
);
904 loops
= g_slist_prepend (loops
, start
);
906 GtsSegment
* next
, * prev
;
909 s
= prev_flag (start
, INTERIOR
);
914 s
= prev_flag (s
, INTERIOR
);
916 next
= next_flag (end
, RELEVANT
);
917 prev
= prev_flag (start
, RELEVANT
);
919 SET (start
->v1
, INTERIOR
);
921 SET (end
->v2
, INTERIOR
);
922 if (next
== NULL
&& prev
== NULL
)
923 loops
= g_slist_prepend (loops
, start
);
925 reverse (start
, TRUE
, &isloop
);
933 #define ORIENTATION(p1,p2,p3,o) (gts_point_orientation_3d (p1, p2, o, p3))
934 #define ORIENTATION_SOS(p1,p2,p3,o) (gts_point_orientation_3d_sos (p1, p2, o, p3))
936 #define ORIENTED_VERTICES(s,next,w1,w2) {\
937 if ((s)->v1 == (next)->v1 || (s)->v1 == (next)->v2) {\
948 static GtsSegment
* segment_intersects (GtsPoint
* p1
, GtsPoint
* p2
,
953 GtsSegment
* s
= i
->data
;
954 GtsPoint
* p3
= GTS_POINT (s
->v1
);
955 GtsPoint
* p4
= GTS_POINT (s
->v2
);
957 if (p3
!= p1
&& p3
!= p2
&& p4
!= p1
&& p4
!= p2
) {
958 gdouble o1
= ORIENTATION (p3
, p4
, p1
, o
);
959 gdouble o2
= ORIENTATION (p3
, p4
, p2
, o
);
961 if ((o1
< 0. && o2
> 0.) || (o1
> 0. && o2
< 0.)) {
962 o1
= ORIENTATION (p1
, p2
, p3
, o
);
963 o2
= ORIENTATION (p1
, p2
, p4
, o
);
965 if ((o1
<= 0. && o2
>= 0.) || (o1
>= 0. && o2
<= 0.))
974 static GtsSegment
* segment_intersects (GtsPoint
* p1
, GtsPoint
* p2
,
979 GtsSegment
* s
= i
->data
;
980 GtsPoint
* p3
= GTS_POINT (s
->v1
);
981 GtsPoint
* p4
= GTS_POINT (s
->v2
);
983 if (p3
!= p1
&& p3
!= p2
&& p4
!= p1
&& p4
!= p2
) {
984 gint o1
= ORIENTATION_SOS (p3
, p4
, p1
, o
);
985 gint o2
= ORIENTATION_SOS (p3
, p4
, p2
, o
);
988 o1
= ORIENTATION_SOS (p1
, p2
, p3
, o
);
989 o2
= ORIENTATION_SOS (p1
, p2
, p4
, o
);
1001 static gboolean
is_inside_wedge (GtsSegment
* s1
, GtsSegment
* s2
,
1002 GtsPoint
* p
, GtsPoint
* o
)
1004 GtsVertex
* v1
, * v2
, * v3
;
1006 ORIENTED_VERTICES (s1
, s2
, v1
, v2
);
1007 v3
= s2
->v1
!= v2
? s2
->v1
: s2
->v2
;
1009 if (ORIENTATION (GTS_POINT (v1
), GTS_POINT (v2
),
1010 GTS_POINT (v3
), o
) >= 0.) {
1011 if (ORIENTATION (GTS_POINT (v1
), GTS_POINT (v2
), p
, o
) <= 0. ||
1012 ORIENTATION (GTS_POINT (v2
), GTS_POINT (v3
), p
, o
) <= 0.)
1015 else if (ORIENTATION (GTS_POINT (v1
), GTS_POINT (v2
), p
, o
) <= 0. &&
1016 ORIENTATION (GTS_POINT (v2
), GTS_POINT (v3
), p
, o
) <= 0.)
1021 static GtsSegment
* connection (GtsPoint
* p
,
1027 GtsSegment
* start
= bloops
->data
, * s
= start
;
1030 GtsSegment
* next
= NEXT (s
);
1031 GtsVertex
* v2
= s
->v1
== next
->v1
|| s
->v1
== next
->v2
? s
->v1
: s
->v2
;
1033 if (is_inside_wedge (s
, next
, p
, o
) &&
1034 !segment_intersects (p
, GTS_POINT (v2
), interior
, o
))
1037 } while (s
!= start
);
1038 bloops
= bloops
->next
;
1043 static gdouble
loop_orientation (GtsSegment
* start
,
1044 GtsPoint
* p
, GtsPoint
* o
)
1046 GtsSegment
* s
= start
;
1050 GtsSegment
* next
= NEXT (s
);
1051 GtsVertex
* v1
, * v2
;
1053 ORIENTED_VERTICES (s
, next
, v1
, v2
);
1054 or += ORIENTATION (p
, GTS_POINT (v1
), GTS_POINT (v2
), o
);
1056 } while (s
!= start
);
1059 fprintf (stderr
, "loop orientation: %g\n", or);
1065 static void connect_interior_loop (GtsSegment
* start
,
1068 GtsSurface
* surface
,
1071 GtsSegment
* s
= start
, * c
= NULL
, * next
, * s1
, * rs1
, * rs
;
1072 GtsVertex
* v
, * cv
;
1076 if (!(c
= connection (GTS_POINT (s
->v2
), *interior
, *bloops
, o
)))
1078 } while (s
!= start
&& !c
);
1081 v
= c
->v1
== next
->v1
|| c
->v1
== next
->v2
? c
->v1
: c
->v2
;
1084 fprintf (stderr
, "connecting %p:%s with %p:%s\n",
1085 cv
, NAME (cv
), v
, NAME (v
));
1086 fprintf (stderr
, " c: %p: %p:%s %p:%s\n", c
,
1087 c
->v1
, NAME (c
->v1
),
1088 c
->v2
, NAME (c
->v2
));
1089 fprintf (stderr
, " next: %p: %p:%s %p:%s\n", next
,
1090 next
->v1
, NAME (next
->v1
),
1091 next
->v2
, NAME (next
->v2
));
1093 rs
= reverse (s
, FALSE
, &isloop
);
1095 if (loop_orientation (rs
, GTS_POINT (v
), o
) < 0.) {
1096 GtsSegment
* tmp
= s
;
1100 *bloops
= g_slist_prepend (*bloops
, rs
);
1102 s1
= GTS_SEGMENT (gts_edge_new (surface
->edge_class
, v
, cv
));
1103 rs1
= GTS_SEGMENT (gts_edge_new (surface
->edge_class
, cv
, v
));
1106 *interior
= g_slist_prepend (*interior
, s1
);
1107 NEXT (s1
) = NEXT (s
);
1111 static GSList
* boundary_loops (GSList
* boundary
)
1113 GSList
* i
= boundary
;
1114 GtsSegment
* start
= i
->data
;
1115 GSList
* loops
= NULL
;
1118 GtsSegment
* s
= i
->data
;
1119 GSList
* inext
= i
->next
;
1120 GtsSegment
* next
= inext
? inext
->data
: start
;
1121 GtsVertex
* v
= s
->v1
== next
->v1
|| s
->v1
== next
->v2
? s
->v1
: s
->v2
;
1123 if (IS_SET (v
, INTERIOR
)) {
1124 GtsSegment
* intprev
= prev_interior (v
);
1126 NEXT (intprev
) = next
;
1127 NEXT (s
) = next_interior (v
);
1128 UNSET (v
, INTERIOR
);
1139 if (IS_SET (start
, RELEVANT
)) {
1140 GtsSegment
* s
= start
;
1143 UNSET (s
, RELEVANT
);
1144 UNSET (s
, INTERIOR
);
1146 } while (s
!= start
);
1147 loops
= g_slist_prepend (loops
, start
);
1155 typedef struct _Ear Ear
;
1158 GtsVertex
* v1
, * v2
, * v3
;
1159 GtsSegment
* s1
, * s2
, * s3
;
1162 static gboolean
point_in_wedge (GtsPoint
* p1
, GtsPoint
* p2
, GtsPoint
* p3
,
1163 GtsPoint
* p
, gboolean closed
, GtsPoint
* o
)
1167 if (p
== p2
|| p
== p3
)
1169 o1
= ORIENTATION (p1
, p2
, p
, o
);
1170 if ((closed
&& o1
< 0.) || (!closed
&& o1
<= 0.)) return FALSE
;
1171 o1
= ORIENTATION (p3
, p1
, p
, o
);
1172 if ((closed
&& o1
< 0.) || (!closed
&& o1
<= 0.)) return FALSE
;
1177 static gboolean
segment_intersects1 (GtsPoint
* p1
, GtsPoint
* p2
,
1178 GtsPoint
* p3
, GtsPoint
* p4
,
1179 gboolean closed
, GtsPoint
* o
)
1181 gdouble o1
= ORIENTATION (p3
, p4
, p1
, o
);
1182 gdouble o2
= ORIENTATION (p3
, p4
, p2
, o
);
1185 if ((closed
&& ((o1
> 0. && o2
> 0.) || (o1
< 0. && o2
< 0.))) ||
1186 (!closed
&& ((o1
>= 0. && o2
>= 0.) || (o1
<= 0. && o2
<= 0.))))
1188 o3
= ORIENTATION (p1
, p2
, p3
, o
);
1189 o4
= ORIENTATION (p1
, p2
, p4
, o
);
1190 if ((o3
> 0. && o4
> 0.) || (o3
< 0. && o4
< 0.))
1192 if (closed
) return TRUE
;
1193 if ((o3
== 0. && o4
> 0.) || (o4
== 0. && o3
> 0.))
1198 static gboolean
segment_intersects1 (GtsPoint
* p1
, GtsPoint
* p2
,
1199 GtsPoint
* p3
, GtsPoint
* p4
,
1200 gboolean closed
, GtsPoint
* o
)
1204 o1
= ORIENTATION_SOS (p3
, p4
, p1
, o
);
1205 o2
= ORIENTATION_SOS (p3
, p4
, p2
, o
);
1208 o1
= ORIENTATION_SOS (p1
, p2
, p3
, o
);
1209 o2
= ORIENTATION_SOS (p1
, p2
, p4
, o
);
1216 static GtsSegment
* triangle_intersects_segments (GtsPoint
* p1
,
1223 GtsSegment
* s
= start
;
1226 GtsPoint
* p4
= GTS_POINT (s
->v1
);
1227 GtsPoint
* p5
= GTS_POINT (s
->v2
);
1230 if (point_in_wedge (p1
, p2
, p3
, p5
, closed
, o
))
1233 else if (p4
== p2
) {
1234 if (point_in_wedge (p2
, p3
, p1
, p5
, closed
, o
))
1237 else if (p4
== p3
) {
1238 if (point_in_wedge (p3
, p1
, p2
, p5
, closed
, o
))
1241 else if (p5
== p1
) {
1242 if (point_in_wedge (p1
, p2
, p3
, p4
, closed
, o
))
1245 else if (p5
== p2
) {
1246 if (point_in_wedge (p2
, p3
, p1
, p4
, closed
, o
))
1249 else if (p5
== p3
) {
1250 if (point_in_wedge (p3
, p1
, p2
, p4
, closed
, o
))
1253 else if (segment_intersects1 (p1
, p2
, p4
, p5
, closed
, o
) ||
1254 segment_intersects1 (p2
, p3
, p4
, p5
, closed
, o
) ||
1255 segment_intersects1 (p3
, p1
, p4
, p5
, closed
, o
))
1258 } while (s
!= start
);
1262 static gboolean
new_ear (GtsSegment
* s
,
1273 g_return_val_if_fail (e
->s2
, FALSE
);
1274 g_return_val_if_fail (e
->s2
!= e
->s1
, FALSE
);
1276 ORIENTED_VERTICES (e
->s1
, e
->s2
, e
->v1
, e
->v2
);
1277 e
->v3
= e
->s2
->v1
!= e
->v2
? e
->s2
->v1
: e
->s2
->v2
;
1280 e
->s3
= NEXT (e
->s2
);
1281 if (gts_segment_connect (e
->s3
, e
->v1
, e
->v3
)) {
1282 if (NEXT (e
->s3
) != e
->s1
)
1285 else if (gts_vertices_are_connected (e
->v1
, e
->v3
))
1289 or = ORIENTATION (GTS_POINT (e
->v1
), GTS_POINT (e
->v2
), GTS_POINT (e
->v3
),o
);
1293 triangle_intersects_segments (GTS_POINT (e
->v1
), GTS_POINT (e
->v2
),
1294 GTS_POINT (e
->v3
), TRUE
, start
, o
))
1300 triangle_intersects_segments (GTS_POINT (e
->v1
), GTS_POINT (e
->v2
),
1301 GTS_POINT (e
->v3
), FALSE
, start
, o
)))
1306 triangle_intersects_segments (GTS_POINT (e
->v1
), GTS_POINT (e
->v2
),
1307 GTS_POINT (e
->v3
), FALSE
, start
, o
)) ||
1309 triangle_intersects_segments (GTS_POINT (e
->v2
), GTS_POINT (e
->v1
),
1310 GTS_POINT (e
->v3
), FALSE
, start
, o
)))
1320 fprintf (stderr
, "or: %g\n", or);
1322 g_assert (or > -1e-6);
1326 static void triangulate_loop (GtsSegment
* start
,
1327 GtsSurface
* surface
,
1330 GtsSegment
* prev
= start
, * s
;
1337 while (NEXT (s
) != s
) {
1338 GtsSegment
* next
= NEXT (s
);
1342 fprintf (stderr
, "prev: %p s: %p next: %p\n", prev
, s
, next
);
1345 if (!new_ear (s
, &e
, start
, sloppy
, o
)) {
1349 fprintf (stderr
, "sloppy: %u\n", sloppy
);
1358 if (!GTS_IS_EDGE (e
.s3
))
1359 e
.s3
= GTS_SEGMENT (gts_edge_new (surface
->edge_class
, e
.v1
, e
.v3
));
1360 f
= gts_face_new (surface
->face_class
,
1361 GTS_EDGE (e
.s1
), GTS_EDGE (e
.s2
), GTS_EDGE (e
.s3
));
1362 gts_surface_add_face (surface
, f
);
1363 UNSET (e
.s1
, RELEVANT
);
1364 UNSET (e
.s1
, INTERIOR
);
1365 UNSET (e
.s2
, RELEVANT
);
1366 UNSET (e
.s2
, INTERIOR
);
1368 NEXT (e
.s3
) = NEXT (e
.s2
);
1369 NEXT (e
.s1
) = NEXT (e
.s2
) = NULL
;
1378 fprintf (stderr
, " t.%u: (%p:%s,%p:%s,%p:%s)\n",
1383 sprintf (name
, "/tmp/t.%u", nt
++);
1384 fp
= fopen (name
, "wt");
1385 // gts_surface_write (surface, fp);
1386 gts_write_triangle (GTS_TRIANGLE (f
), NULL
, fp
);
1387 // write_graph1 (start, interior, surface, fp);
1389 print_loop (start
, stderr
);
1394 UNSET (s
, RELEVANT
);
1395 UNSET (s
, INTERIOR
);
1399 #ifdef CHECK_ORIENTED
1400 static void check_object (GtsObject
* o
)
1402 g_assert (o
->reserved
== NULL
);
1403 g_assert (o
->flags
== 0);
1406 static void check_boundary (GtsEdge
* e
, GtsSurface
* s
)
1408 check_object (GTS_OBJECT (e
));
1409 check_object (GTS_OBJECT (GTS_SEGMENT (e
)->v1
));
1410 check_object (GTS_OBJECT (GTS_SEGMENT (e
)->v2
));
1411 g_assert (gts_edge_face_number (e
, s
) == 1);
1414 static void check_interior (GtsEdge
* e
, GtsSurface
* s
)
1417 check_object (GTS_OBJECT (e
));
1418 check_object (GTS_OBJECT (GTS_SEGMENT (e
)->v1
));
1419 check_object (GTS_OBJECT (GTS_SEGMENT (e
)->v2
));
1421 n
= gts_edge_face_number (e
, s
);
1424 gts_surface_print_stats (s
, stderr
);
1429 static void check_boundary_interior_triangulation (GSList
* boundary
,
1431 GtsSurface
* surface
)
1433 g_slist_foreach (boundary
, (GFunc
) check_boundary
, surface
);
1434 g_slist_foreach (interior
, (GFunc
) check_interior
, surface
);
1436 #endif /*ifdef CHECK_ORIENTED */
1438 static void merge_duplicate (GtsEdge
* e
)
1440 GtsEdge
* dup
= gts_edge_is_duplicate (e
);
1443 gts_edge_replace (dup
, e
);
1444 gts_object_destroy (GTS_OBJECT (dup
));
1447 static void triangulate_boundary_interior (GSList
* boundary
,
1452 GSList
* iloops
, * bloops
, * i
;
1456 SET (i
->data
, RELEVANT
);
1461 SET (i
->data
, RELEVANT
);
1462 SET (i
->data
, INTERIOR
);
1466 iloops
= interior_loops (interior
);
1467 bloops
= boundary_loops (boundary
);
1472 fprintf (stderr
, "--- interior loop ---\n");
1473 print_loop (i
->data
, stderr
);
1475 connect_interior_loop (i
->data
, &interior
, &bloops
, s
, o
);
1481 FILE * fp
= fopen ("/tmp/bloops", "w");
1482 write_loops (bloops
, fp
);
1490 fprintf (stderr
, "--- boundary loop ---\n");
1491 print_loop (i
->data
, stderr
);
1493 triangulate_loop (i
->data
, s
, o
);
1497 g_slist_foreach (interior
, (GFunc
) merge_duplicate
, NULL
);
1498 g_slist_free (iloops
);
1499 g_slist_free (bloops
);
1501 #ifdef CHECK_ORIENTED
1502 check_boundary_interior_triangulation (boundary
, interior
, s
);
1503 #endif /* CHECK_ORIENTED */
1506 static void create_edges (GtsSegment
* s
, GtsSurface
* surface
)
1508 if (GTS_OBJECT (s
)->reserved
) {
1509 GList
* i
= GTS_OBJECT (s
)->reserved
;
1510 GtsVertex
* v1
= i
->data
;
1512 GTS_OBJECT (s
)->reserved
= g_list_prepend (i
,
1513 gts_edge_new (surface
->edge_class
, s
->v1
, v1
));
1515 GList
* next
= i
->next
;
1516 GtsVertex
* v2
= next
? next
->data
: s
->v2
;
1518 GTS_OBJECT (i
->data
)->reserved
= NULL
;
1519 i
->data
= gts_edge_new (surface
->edge_class
, v1
, v2
);
1526 static void add_boundary (GtsSegment
* s
, GtsSegment
* next
,
1529 if (GTS_OBJECT (s
)->reserved
== NULL
)
1530 *boundary
= g_slist_prepend (*boundary
, s
);
1532 if (s
->v2
== next
->v2
|| s
->v2
== next
->v1
) {
1533 GList
* i
= g_list_last (GTS_OBJECT (s
)->reserved
);
1536 *boundary
= g_slist_prepend (*boundary
, i
->data
);
1541 GList
* i
= GTS_OBJECT (s
)->reserved
;
1544 *boundary
= g_slist_prepend (*boundary
, i
->data
);
1551 static void triangulate_face (GtsTriangle
* t
, GtsSurface
* surface
)
1553 GSList
* interior
= GTS_OBJECT (t
)->reserved
;
1554 GSList
* boundary
= NULL
;
1555 GtsSurface
* s
= gts_surface_new (gts_surface_class (),
1556 surface
->face_class
,
1557 surface
->edge_class
,
1558 surface
->vertex_class
);
1560 GtsPoint
* p
= GTS_POINT (GTS_SEGMENT (t
->e1
)->v1
);
1563 GTS_OBJECT (t
)->reserved
= NULL
;
1564 gts_triangle_normal (t
, &x
, &y
, &z
);
1565 g_assert (x
!= 0. || y
!= 0. || z
!= 0.);
1566 o
= gts_point_new (gts_point_class (), p
->x
+ x
, p
->y
+ y
, p
->z
+ z
);
1567 add_boundary (GTS_SEGMENT (t
->e3
), GTS_SEGMENT (t
->e1
), &boundary
);
1568 add_boundary (GTS_SEGMENT (t
->e2
), GTS_SEGMENT (t
->e3
), &boundary
);
1569 add_boundary (GTS_SEGMENT (t
->e1
), GTS_SEGMENT (t
->e2
), &boundary
);
1572 static guint nt
= 0;
1576 fprintf (stderr
, "%u: triangulating %p\n", nt
, t
);
1578 fprintf (stderr
, "tintin!!!!\n");
1579 sprintf (name
, "/tmp/oc.%u", nt
++);
1580 fp
= fopen (name
, "wt");
1581 // write_graph (boundary, interior, s, fp);
1582 write_segments (boundary
, interior
, fp
);
1586 triangulate_boundary_interior (boundary
, interior
, s
, o
);
1587 g_slist_free (interior
);
1588 g_slist_free (boundary
);
1589 if (GTS_OBJECT (t
)->klass
->attributes
)
1590 gts_surface_foreach_face (s
, (GtsFunc
) gts_object_attributes
, t
);
1591 gts_surface_merge (surface
, s
);
1592 gts_object_destroy (GTS_OBJECT (s
));
1593 gts_object_destroy (GTS_OBJECT (o
));
1596 static void free_edge_list (GtsObject
* o
)
1598 g_list_free (o
->reserved
);
1603 * gts_surface_inter_new:
1604 * @klass: a #GtsSurfaceInterClass.
1605 * @s1: a #GtsSurface.
1606 * @s2: a #GtsSurface.
1607 * @faces_tree1: a bounding box tree (see gts_bb_tree_new()) for
1609 * @faces_tree2: a bounding box tree for the faces of @s2.
1610 * @is_open1: whether @s1 is an "open" surface.
1611 * @is_open2: whether @s2 is an "open" surface.
1613 * When triangulating the cut faces, the new faces inherit the
1614 * attributes of these original faces through their attributes()
1617 * Returns: a new #GtsSurfaceInter describing the intersection of @s1
1620 GtsSurfaceInter
* gts_surface_inter_new (GtsSurfaceInterClass
* klass
,
1623 GNode
* faces_tree1
,
1624 GNode
* faces_tree2
,
1628 GtsSurfaceInter
* si
;
1631 g_return_val_if_fail (klass
!= NULL
, NULL
);
1632 g_return_val_if_fail (s1
!= NULL
, NULL
);
1633 g_return_val_if_fail (s2
!= NULL
, NULL
);
1634 g_return_val_if_fail (faces_tree1
!= NULL
, NULL
);
1635 g_return_val_if_fail (faces_tree2
!= NULL
, NULL
);
1637 si
= surface_inter_new (klass
, s1
, s2
, faces_tree1
, faces_tree2
);
1639 gts_surface_foreach_edge (si
->s1
, (GtsFunc
) create_edges
, si
->s1
);
1640 gts_surface_foreach_edge (si
->s2
, (GtsFunc
) create_edges
, si
->s2
);
1643 fprintf (stderr
, "====== triangulating s1 ======\n");
1645 s
= gts_surface_new (gts_surface_class (),
1649 gts_surface_foreach_face (si
->s1
, (GtsFunc
) triangulate_face
, s
);
1650 gts_surface_foreach_edge (si
->s1
, (GtsFunc
) free_edge_list
, NULL
);
1651 gts_object_destroy (GTS_OBJECT (si
->s1
));
1653 GTS_OBJECT (si
->s1
)->reserved
= s1
;
1656 fprintf (stderr
, "====== triangulating s2 ======\n");
1658 s
= gts_surface_new (gts_surface_class (),
1662 gts_surface_foreach_face (si
->s2
, (GtsFunc
) triangulate_face
, s
);
1663 gts_surface_foreach_edge (si
->s2
, (GtsFunc
) free_edge_list
, NULL
);
1664 gts_object_destroy (GTS_OBJECT (si
->s2
));
1666 GTS_OBJECT (si
->s2
)->reserved
= s2
;
1671 static void check_surface_edge (GtsEdge
* e
, gpointer
* data
)
1673 gboolean
* ok
= data
[0];
1674 GtsSurface
* s
= data
[1];
1675 GtsSurface
* bs
= GTS_OBJECT (s
)->reserved
;
1676 guint nf
= gts_edge_face_number (e
, s
);
1678 if (nf
< 1 || nf
> 2) {
1680 g_return_if_fail (nf
>= 1 && nf
<= 2);
1682 if (nf
== 1 && gts_edge_face_number (e
, bs
) == 0) {
1684 g_return_if_fail (gts_edge_face_number (e
, bs
) > 0);
1688 static void mark_edge (GtsObject
* o
, gpointer data
)
1693 static gint
triangle_orientation (GtsTriangle
* t
, GtsEdge
* e
)
1695 GtsSegment
* s
= GTS_SEGMENT (t
->e1
== e
? t
->e2
1700 GtsVertex
* v2
= GTS_SEGMENT (e
)->v2
;
1702 if (s
->v1
== v2
|| s
->v2
== v2
)
1707 static gboolean
check_orientation (GtsEdge
* e
, GtsSurface
* s
)
1709 GtsTriangle
* t1
= NULL
, * t2
= NULL
;
1710 GSList
* i
= e
->triangles
;
1711 gint o1
= 0, o2
= 0;
1714 if (GTS_IS_FACE (i
->data
) &&
1715 gts_face_has_parent_surface (i
->data
, s
)) {
1718 o1
= triangle_orientation (t1
, e
);
1720 else if (t2
== NULL
) {
1722 o2
= triangle_orientation (t2
, e
);
1723 g_return_val_if_fail (o1
*o2
< 0, FALSE
);
1726 g_assert_not_reached ();
1730 g_return_val_if_fail (t1
&& t2
, FALSE
);
1734 static void check_edge (GtsSegment
* s
, gpointer
* data
)
1736 gboolean
* ok
= data
[0];
1737 GtsSurfaceInter
* si
= data
[1];
1738 gboolean
* closed
= data
[2];
1742 j
= s
->v1
->segments
;
1744 GtsSegment
* s1
= j
->data
;
1746 if (s1
!= s
&& GTS_OBJECT (s1
)->reserved
== si
) {
1747 if (s1
->v2
!= s
->v1
)
1753 j
= s
->v2
->segments
;
1755 GtsSegment
* s1
= j
->data
;
1757 if (s1
!= s
&& GTS_OBJECT (s1
)->reserved
== si
) {
1758 if (s1
->v1
!= s
->v2
)
1767 if (!check_orientation (GTS_EDGE (s
), si
->s1
))
1769 if (!check_orientation (GTS_EDGE (s
), si
->s2
))
1774 * gts_surface_inter_check:
1775 * @si: a #GtsSurfaceInter.
1776 * @closed: is set to %TRUE if @si->edges is a closed curve, %FALSE
1779 * Returns: %TRUE if the curve described by @si is an orientable
1780 * manifold, %FALSE otherwise.
1782 gboolean
gts_surface_inter_check (GtsSurfaceInter
* si
,
1788 g_return_val_if_fail (si
!= NULL
, FALSE
);
1789 g_return_val_if_fail (closed
!= NULL
, FALSE
);
1791 *closed
= si
->edges
? TRUE
: FALSE
;
1793 /* mark edges as used by si */
1794 g_slist_foreach (si
->edges
, (GFunc
) mark_edge
, si
);
1799 g_slist_foreach (si
->edges
, (GFunc
) check_edge
, data
);
1800 g_slist_foreach (si
->edges
, (GFunc
) gts_object_reset_reserved
, NULL
);
1802 /* check connectivity of the faces of @si */
1808 gts_surface_foreach_edge (si
->s1
, (GtsFunc
) check_surface_edge
, data
);
1810 gts_surface_foreach_edge (si
->s2
, (GtsFunc
) check_surface_edge
, data
);
1816 /* Given @e and @f returns a #GtsFace compatible with @f and belonging to
1818 static GtsFace
* next_compatible_face (GtsEdge
* e
,
1823 GSList
* i
= e
->triangles
;
1824 GtsFace
* f2
= NULL
, * f3
= NULL
;
1827 GtsFace
* f1
= i
->data
;
1829 if (f1
!= f
&& GTS_IS_FACE (f1
)) {
1830 if (gts_face_has_parent_surface (f1
, s1
))
1832 if (gts_face_has_parent_surface (f1
, s2
)) {
1833 if (f2
== NULL
) f2
= f1
;
1834 else if (f3
== NULL
) f3
= f1
;
1835 else g_assert_not_reached (); /* s2 is a non-manifold surface */
1841 if (gts_edge_is_boundary (e
, s2
))
1845 g_assert (gts_face_has_parent_surface (f
, s1
));
1846 if (gts_triangles_are_compatible (GTS_TRIANGLE (f
), GTS_TRIANGLE (f2
), e
))
1851 static void walk_faces (GtsEdge
* e
, GtsFace
* f
,
1856 GtsFifo
* faces
= gts_fifo_new ();
1857 GtsFifo
* edges
= gts_fifo_new ();
1859 gts_fifo_push (faces
, f
);
1860 gts_fifo_push (edges
, e
);
1861 while ((f
= gts_fifo_pop (faces
)) && (e
= gts_fifo_pop (edges
))) {
1862 if (!GTS_OBJECT (f
)->reserved
) {
1863 GtsTriangle
* t
= GTS_TRIANGLE (f
);
1866 gts_surface_add_face (s
, f
);
1867 GTS_OBJECT (f
)->reserved
= s
;
1868 if (t
->e1
!= e
&& !GTS_OBJECT (t
->e1
)->reserved
&&
1869 (f1
= next_compatible_face (t
->e1
, f
, s1
, s2
))) {
1870 gts_fifo_push (faces
, f1
);
1871 gts_fifo_push (edges
, t
->e1
);
1873 if (t
->e2
!= e
&& !GTS_OBJECT (t
->e2
)->reserved
&&
1874 (f1
= next_compatible_face (t
->e2
, f
, s1
, s2
))) {
1875 gts_fifo_push (faces
, f1
);
1876 gts_fifo_push (edges
, t
->e2
);
1878 if (t
->e3
!= e
&& !GTS_OBJECT (t
->e3
)->reserved
&&
1879 (f1
= next_compatible_face (t
->e3
, f
, s1
, s2
))) {
1880 gts_fifo_push (faces
, f1
);
1881 gts_fifo_push (edges
, t
->e3
);
1885 gts_fifo_destroy (faces
);
1886 gts_fifo_destroy (edges
);
1890 * gts_surface_inter_boolean:
1891 * @si: a #GtsSurfaceInter.
1892 * @surface: a #GtsSurface.
1893 * @op: a #GtsBooleanOperation.
1895 * Adds to @surface the part of the surface described by @si and @op.
1897 void gts_surface_inter_boolean (GtsSurfaceInter
* si
,
1898 GtsSurface
* surface
,
1899 GtsBooleanOperation op
)
1901 GtsSurface
* s
= NULL
;
1905 g_return_if_fail (si
!= NULL
);
1906 g_return_if_fail (surface
!= NULL
);
1909 case GTS_1_OUT_2
: s
= si
->s1
; orient
= 1; break;
1910 case GTS_1_IN_2
: s
= si
->s1
; orient
= -1; break;
1911 case GTS_2_OUT_1
: s
= si
->s2
; orient
= -1; break;
1912 case GTS_2_IN_1
: s
= si
->s2
; orient
= 1; break;
1913 default: g_assert_not_reached ();
1916 /* mark edges as belonging to intersection */
1917 g_slist_foreach (si
->edges
, (GFunc
) mark_edge
, si
);
1921 GtsEdge
* e
= i
->data
;
1922 GSList
* j
= e
->triangles
;
1925 if (gts_face_has_parent_surface (j
->data
, s
) &&
1926 orient
*triangle_orientation (j
->data
, e
) > 0) {
1927 #ifdef DEBUG_BOOLEAN
1928 GtsFace
* boundary
= gts_edge_is_boundary (e
, surface
);
1930 g_assert (!boundary
|| boundary
== j
->data
);
1931 #endif /* DEBUG_BOOLEAN */
1932 walk_faces (e
, j
->data
, s
, GTS_OBJECT (s
)->reserved
, surface
);
1939 g_slist_foreach (si
->edges
, (GFunc
) gts_object_reset_reserved
, NULL
);
1940 gts_surface_foreach_face (surface
,
1941 (GtsFunc
) gts_object_reset_reserved
, NULL
);
1944 static void self_intersecting (GtsBBox
* bb1
, GtsBBox
* bb2
,
1947 GtsTriangle
* t1
= bb1
->bounded
;
1948 GtsTriangle
* t2
= bb2
->bounded
;
1951 GtsSegment
* s1
= GTS_SEGMENT (t1
->e1
);
1952 GtsSegment
* s2
= GTS_SEGMENT (t1
->e2
);
1953 GtsSegment
* s3
= GTS_SEGMENT (t1
->e3
);
1954 GtsSegment
* s4
= GTS_SEGMENT (t2
->e1
);
1955 GtsSegment
* s5
= GTS_SEGMENT (t2
->e2
);
1956 GtsSegment
* s6
= GTS_SEGMENT (t2
->e3
);
1959 if ((!gts_segments_touch (s4
, s1
) &&
1960 !gts_segments_touch (s4
, s2
) &&
1961 !gts_segments_touch (s4
, s3
) &&
1962 (pi
= segment_triangle_intersection (s4
, t1
, gts_point_class ()))
1964 (!gts_segments_touch (s5
, s1
) &&
1965 !gts_segments_touch (s5
, s2
) &&
1966 !gts_segments_touch (s5
, s3
) &&
1967 (pi
= segment_triangle_intersection (s5
, t1
, gts_point_class ()))
1969 (!gts_segments_touch (s6
, s1
) &&
1970 !gts_segments_touch (s6
, s2
) &&
1971 !gts_segments_touch (s6
, s3
) &&
1972 (pi
= segment_triangle_intersection (s6
, t1
, gts_point_class ()))
1974 GtsBBTreeTraverseFunc func
= d
[0];
1975 gpointer data
= d
[1];
1976 gboolean
* self_inter
= d
[2];
1978 gts_object_destroy (GTS_OBJECT (pi
));
1980 (* func
) (bb1
, bb2
, data
);
1986 * gts_surface_foreach_intersecting_face:
1987 * @s: a #GtsSurface.
1988 * @func: a #GtsBBTreeTraverseFunc.
1989 * @data: user data to pass to @func.
1991 * Calls @func for each intersecting pair of faces of @s.
1993 * Returns: %TRUE if @func was called at least once, %FALSE otherwise.
1995 gboolean
gts_surface_foreach_intersecting_face (GtsSurface
* s
,
1996 GtsBBTreeTraverseFunc func
,
2001 gboolean self_inter
= FALSE
;
2003 g_return_val_if_fail (s
!= NULL
, FALSE
);
2004 g_return_val_if_fail (func
!= NULL
, FALSE
);
2006 tree
= gts_bb_tree_surface (s
);
2010 gts_bb_tree_traverse_overlapping (tree
, tree
,
2011 (GtsBBTreeTraverseFunc
) self_intersecting
,
2013 gts_bb_tree_destroy (tree
, TRUE
);
2018 static void add_intersecting (GtsBBox
* bb1
, GtsBBox
* bb2
,
2019 GtsSurface
* intersected
)
2021 gts_surface_add_face (intersected
, bb1
->bounded
);
2022 gts_surface_add_face (intersected
, bb2
->bounded
);
2026 * gts_surface_is_self_intersecting:
2027 * @s: a #GtsSurface.
2029 * Returns: a new #GtsSurface containing the faces of @s which are
2030 * self-intersecting or %NULL if no faces of @s are self-intersecting.
2032 GtsSurface
* gts_surface_is_self_intersecting (GtsSurface
* s
)
2034 GtsSurface
* intersected
;
2036 g_return_val_if_fail (s
!= NULL
, NULL
);
2038 intersected
= gts_surface_new (GTS_SURFACE_CLASS (GTS_OBJECT (s
)->klass
),
2042 if (!gts_surface_foreach_intersecting_face (s
,
2043 (GtsBBTreeTraverseFunc
) add_intersecting
, intersected
)) {
2044 gts_object_destroy (GTS_OBJECT (intersected
));