gts/split.c

   1 /* GTS - Library for the manipulation of triangulated surfaces
   2  * Copyright (C) 1999 Stéphane Popinet
   3  *
   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.
   8  *
   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.
  13  *
  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.
  18  */
  19
  20 #include <stdlib.h>
  21 #include <string.h>
  22 #include "gts.h"
  23
  24 #define DYNAMIC_SPLIT
  25 #define NEW
  26
  27 /* #define DEBUG
  28    #define DEBUG_HEXPAND
  29    #define DEBUG_EXPAND */
  30
  31 struct _GtsSplitCFace {
  32   GtsFace * f;
  33   GtsTriangle ** a1, ** a2;
  34 };
  35
  36 typedef struct _CFace      CFace;
  37 typedef struct _CFaceClass CFaceClass;
  38
  39 struct _CFace {
  40   GtsObject object;
  41
  42   GtsSplit * parent_split;
  43   GtsTriangle * t;
  44   guint flags;
  45 };
  46 /* the size of the CFace structure must be smaller or equal to the size
  47    of the GtsFace structure as both structures use the same memory location */
  48
  49 struct _CFaceClass {
  50   GtsObjectClass parent_class;
  51 };
  52
  53 #define IS_CFACE(obj) (gts_object_is_from_class (obj, cface_class ()))
  54 #define CFACE(obj)    ((CFace *) obj)
  55 #define CFACE_ORIENTATION(cf) ((cf)->flags & 0x1)
  56 #define CFACE_ORIENTATION_DIRECT(cf) ((cf)->flags |= 0x1)
  57 #define CFACE_VVS(cf)                ((cf)->flags & 0x2)
  58 #define CFACE_VVS_DIRECT(cf)         ((cf)->flags |= 0x2)
  59 #define CFACE_E1                     0x4
  60 #define CFACE_E2                     0x8
  61 #define CFACE_KEEP_VVS               0x10
  62
  63 #define ROTATE_ORIENT(e, e1, e2, e3)  { if (e1 == e) { e1 = e2; e2 = e3; }\
  64                                  else if (e2 == e) { e2 = e1; e1 = e3; }\
  65                                  else g_assert (e3 == e); }
  66 #define SEGMENT_USE_VERTEX(s, v) ((s)->v1 == v || (s)->v2 == v)
  67 #define TRIANGLE_REPLACE_EDGE(t, e, with) { if ((t)->e1 == e)\
  68                                               (t)->e1 = with;\
  69                                             else if ((t)->e2 == e)\
  70                                               (t)->e2 = with;\
  71                                             else {\
  72                                               g_assert ((t)->e3 == e);\
  73                                               (t)->e3 = with;\
  74                                             }\
  75                                           }
  76
  77 #define HEAP_INSERT_OBJECT(h, e) (GTS_OBJECT (e)->reserved =\
  78                                   gts_eheap_insert (h, e))
  79 #define HEAP_REMOVE_OBJECT(h, e) (gts_eheap_remove (h, GTS_OBJECT (e)->reserved),\
  80                                    GTS_OBJECT (e)->reserved = NULL)
  81
  82 static GtsObjectClass * cface_class (void)
  83 {
  84   static GtsObjectClass * klass = NULL;
  85
  86   if (klass == NULL) {
  87     GtsObjectClassInfo cface_info = {
  88       "GtsCFace",
  89       sizeof (CFace),
  90       sizeof (CFaceClass),
  91       (GtsObjectClassInitFunc) NULL,
  92       (GtsObjectInitFunc) NULL,
  93       (GtsArgSetFunc) NULL,
  94       (GtsArgGetFunc) NULL
  95     };
  96     klass = gts_object_class_new (gts_object_class (), &cface_info);
  97     g_assert (sizeof (CFace) <= sizeof (GtsFace));
  98   }
  99
 100   return klass;
 101 }
 102
 103 /* Replace @e with @with for all the triangles using @e but @f.
 104    Destroys @e and removes it from @heap (if not %NULL).
 105    Returns a triangle using e different from f or %NULL. */
 106 static GtsTriangle * replace_edge_collapse (GtsEdge * e,
 107                                             GtsEdge * with,
 108                                             CFace * cf,
 109                                             GtsEHeap * heap
 110 #ifdef DYNAMIC_SPLIT
 111                                             , GtsTriangle *** a1
 112 #endif
 113 #ifdef NEW
 114                                             , guint edge_flag
 115 #endif
 116                                             )
 117 {
 118   GSList * i;
 119   GtsTriangle * rt = NULL;
 120 #ifdef DYNAMIC_SPLIT
 121   guint size;
 122   GtsTriangle ** a;
 123 #endif
 124
 125 #ifdef NEW
 126   i = e->triangles;
 127   e->triangles = NULL;
 128   size = g_slist_length (i)*sizeof (GtsTriangle *);
 129   *a1 = a = g_malloc (size > 0 ? size : sizeof (GtsTriangle *));
 130   while (i) {
 131     GtsTriangle * t = i->data;
 132     GSList * next = i->next;
 133     if (t != ((GtsTriangle *) cf)) {
 134       if (IS_CFACE (t)) {
 135         i->next = e->triangles;
 136         e->triangles = i;
 137         /* set the edge given by edge_flag (CFACE_E1 or CFACE_E2) */
 138         GTS_OBJECT (t)->reserved = GUINT_TO_POINTER (edge_flag);
 139         cf->flags |= CFACE_KEEP_VVS;
 140       }
 141       else {
 142         TRIANGLE_REPLACE_EDGE (t, e, with);
 143         i->next = with->triangles;
 144         with->triangles = i;
 145         rt = t;
 146         *(a++) = t;
 147       }
 148     }
 149     else
 150       g_slist_free_1 (i);
 151     i = next;
 152   }
 153   *a = NULL;
 154   if (!e->triangles) {
 155     if (heap)
 156       HEAP_REMOVE_OBJECT (heap, e);
 157     gts_object_destroy (GTS_OBJECT (e));
 158   }
 159 #else /* not NEW */
 160   i = e->triangles;
 161 #ifdef DYNAMIC_SPLIT
 162   size = g_slist_length (i)*sizeof (GtsTriangle *);
 163   *a1 = a = g_malloc (size > 0 ? size : sizeof (GtsTriangle *));
 164 #endif
 165   while (i) {
 166     GtsTriangle * t = i->data;
 167     GSList * next = i->next;
 168     if (t != ((GtsTriangle *) cf)) {
 169       TRIANGLE_REPLACE_EDGE (t, e, with);
 170       i->next = with->triangles;
 171       with->triangles = i;
 172       rt = t;
 173 #ifdef DYNAMIC_SPLIT
 174       *(a++) = t;
 175 #endif
 176     }
 177     else
 178       g_slist_free_1 (i);
 179     i = next;
 180   }
 181 #ifdef DYNAMIC_SPLIT
 182   *a = NULL;
 183 #endif
 184   if (heap)
 185     HEAP_REMOVE_OBJECT (heap, e);
 186   e->triangles = NULL;
 187   gts_object_destroy (GTS_OBJECT (e));
 188 #endif /* NEW */
 189
 190   return rt;
 191 }
 192
 193 static CFace * cface_new (GtsFace * f,
 194                           GtsEdge * e,
 195                           GtsVertex * v1,
 196                           GtsVertex * v2,
 197                           GtsSplit * vs,
 198                           GtsEHeap * heap,
 199                           GtsEdgeClass * klass
 200 #ifdef DYNAMIC_SPLIT
 201                           , GtsSplitCFace * scf
 202 #endif
 203                           )
 204 {
 205   CFace * cf;
 206   GtsVertex * v;
 207   GtsEdge * e1, * e2, * e3, * vvs;
 208   GSList * i;
 209   GtsTriangle * t, * t1 = NULL, * t2 = NULL;
 210   guint flags;
 211
 212   g_return_val_if_fail (f != NULL, NULL);
 213 #ifndef NEW
 214   g_return_val_if_fail (GTS_IS_FACE (f), NULL);
 215 #endif
 216   g_return_val_if_fail (e != NULL, NULL);
 217   g_return_val_if_fail (vs != NULL, NULL);
 218
 219   t = ((GtsTriangle *) f);
 220   if (heap)
 221     g_return_val_if_fail (!gts_triangle_is_duplicate (t), NULL);
 222
 223 #ifdef NEW
 224   /* get CFACE_E1 and CFACE_E2 info */
 225   flags = GPOINTER_TO_UINT (GTS_OBJECT (f)->reserved);
 226 #endif
 227   GTS_OBJECT_SET_FLAGS (f, GTS_DESTROYED);
 228
 229   i = f->surfaces;
 230   while (i) {
 231     GSList * next = i->next;
 232     gts_surface_remove_face (i->data, f);
 233     i = next;
 234   }
 235   g_slist_free (f->surfaces);
 236
 237   e1 = t->e1; e2 = t->e2; e3 = t->e3;
 238   ROTATE_ORIENT (e, e1, e2, e3);
 239
 240   cf = (CFace *) f;
 241 #ifndef NEW
 242   GTS_OBJECT (cf)->klass = cface_class ();
 243 #else
 244   cf->flags = flags;
 245 #endif
 246   gts_object_init (GTS_OBJECT (cf), cface_class ());
 247   cf->parent_split = vs;
 248
 249   if (SEGMENT_USE_VERTEX (GTS_SEGMENT (e1), v2)) {
 250     CFACE_ORIENTATION_DIRECT (cf); /* v1->v2->v */
 251     e3 = e1; e1 = e2; e2 = e3;
 252   }
 253   v = GTS_SEGMENT (e1)->v1 == v1 ?
 254     GTS_SEGMENT (e1)->v2 : GTS_SEGMENT (e1)->v1;
 255 #ifdef NEW
 256   if ((cf->flags & CFACE_E1) || (cf->flags & CFACE_E2)) {
 257     vvs = GTS_EDGE (gts_vertices_are_connected (vs->v, v));
 258     g_assert (vvs != NULL);
 259   } else
 260 #endif
 261   vvs = gts_edge_new (klass, v, vs->v);
 262
 263   t1 = replace_edge_collapse (e1, vvs, cf, heap
 264 #ifdef DYNAMIC_SPLIT
 265                               , &scf->a1
 266 #endif
 267 #ifdef NEW
 268                               , CFACE_E1
 269 #endif
 270                               );
 271   t2 = replace_edge_collapse (e2, vvs, cf, heap
 272 #ifdef DYNAMIC_SPLIT
 273                               , &scf->a2
 274 #endif
 275 #ifdef NEW
 276                               , CFACE_E2
 277 #endif
 278                               );
 279   t = cf->t = t1 ? t1 : t2;
 280   g_assert (t);
 281
 282   /* set up flags necessary to find vvs */
 283   if (t->e1 == vvs) e2 = t->e2;
 284   else if (t->e2 == vvs) e2 = t->e3;
 285   else {
 286     g_assert (t->e3 == vvs);
 287     e2 = t->e1;
 288   }
 289   if (SEGMENT_USE_VERTEX (GTS_SEGMENT (e2), v))
 290     CFACE_VVS_DIRECT (cf);
 291
 292   return cf;
 293 }
 294
 295 static void find_vvs (GtsVertex * vs,
 296                       GtsTriangle * t,
 297                       GtsVertex ** v, GtsEdge ** vvs,
 298                       gboolean orientation)
 299 {
 300   GtsEdge * e1 = t->e1, * e2 = t->e2, * e3 = t->e3, * tmp;
 301
 302   if (SEGMENT_USE_VERTEX (GTS_SEGMENT (e2), vs)) {
 303     tmp = e1; e1 = e2; e2 = e3; e3 = tmp;
 304   }
 305   else if (SEGMENT_USE_VERTEX (GTS_SEGMENT (e3), vs)) {
 306     tmp = e1; e1 = e3; e3 = e2; e2 = tmp;
 307   }
 308   else
 309     g_assert (SEGMENT_USE_VERTEX (GTS_SEGMENT (e1), vs));
 310   if (SEGMENT_USE_VERTEX (GTS_SEGMENT (e2), vs) ||
 311       !gts_segments_touch (GTS_SEGMENT (e1), GTS_SEGMENT (e2))) {
 312     tmp = e1; e1 = e2; e2 = e3; e3 = tmp;
 313     g_assert (gts_segments_touch (GTS_SEGMENT (e1), GTS_SEGMENT (e2)));
 314   }
 315
 316   *vvs = orientation ? e1 : e3;
 317
 318   if (GTS_SEGMENT (*vvs)->v1 != vs) {
 319     g_assert (GTS_SEGMENT (*vvs)->v2 == vs);
 320     *v = GTS_SEGMENT (*vvs)->v1;
 321   }
 322   else
 323     *v = GTS_SEGMENT (*vvs)->v2;
 324 }
 325
 326 static void replace_edge_expand (GtsEdge * e,
 327                                  GtsEdge * with,
 328                                  GtsTriangle ** a,
 329                                  GtsVertex * v)
 330 {
 331   GtsTriangle ** i = a, * t;
 332
 333   while ((t = *(i++))) {
 334 #ifdef DEBUG_EXPAND
 335     g_assert (!IS_CFACE (t));
 336     fprintf (stderr, "replacing %p->%d: e: %p->%d with: %p->%d\n",
 337              t, id (t), e, id (e), with, id (with));
 338 #endif
 339     TRIANGLE_REPLACE_EDGE (t, e, with);
 340     with->triangles = g_slist_prepend (with->triangles, t);
 341     if (GTS_OBJECT (t)->reserved) {
 342       /* apart from the triangles having e as an edge, t is the only
 343          triangle using v */
 344       g_assert (GTS_OBJECT (t)->reserved == v);
 345       GTS_OBJECT (t)->reserved = NULL;
 346     }
 347     else
 348       GTS_OBJECT (t)->reserved = v;
 349   }
 350 }
 351
 352 static void cface_expand (CFace * cf,
 353                           GtsTriangle ** a1,
 354                           GtsTriangle ** a2,
 355                           GtsEdge * e,
 356                           GtsVertex * v1,
 357                           GtsVertex * v2,
 358                           GtsVertex * vs,
 359                           GtsEdgeClass * klass)
 360 {
 361   GtsVertex * v;
 362   GtsEdge * e1, * e2, * vvs;
 363   gboolean orientation;
 364   guint flags;
 365
 366   g_return_if_fail (cf != NULL);
 367   g_return_if_fail (IS_CFACE (cf));
 368   g_return_if_fail (e != NULL);
 369   g_return_if_fail (vs != NULL);
 370
 371   flags = cf->flags;
 372   orientation = CFACE_ORIENTATION (cf);
 373
 374   find_vvs (vs, cf->t, &v, &vvs, CFACE_VVS (cf));
 375
 376 #ifdef NEW
 377   if (flags & CFACE_E1)
 378     e1 = GTS_EDGE (gts_vertices_are_connected (v1, v));
 379   else
 380     e1 = gts_edge_new (klass, v, v1);
 381   if (flags & CFACE_E2)
 382     e2 = GTS_EDGE (gts_vertices_are_connected (v2, v));
 383   else
 384     e2 = gts_edge_new (klass, v, v2);
 385 #else
 386   e1 = gts_edge_new (v, v1);
 387   e2 = gts_edge_new (v, v2);
 388 #endif
 389
 390   replace_edge_expand (vvs, e1, a1, v1);
 391   replace_edge_expand (vvs, e2, a2, v2);
 392
 393 #ifdef NEW
 394   if (!(flags & CFACE_KEEP_VVS)) {
 395     g_slist_free (vvs->triangles);
 396     vvs->triangles = NULL;
 397     gts_object_destroy (GTS_OBJECT (vvs));
 398   }
 399 #else
 400   g_slist_free (vvs->triangles);
 401   vvs->triangles = NULL;
 402   gts_object_destroy (GTS_OBJECT (vvs));
 403 #endif
 404
 405   /* gts_face_new : because I am "creating" a face */
 406   GTS_OBJECT (cf)->klass = GTS_OBJECT_CLASS (gts_face_class ());
 407   gts_object_init (GTS_OBJECT (cf), GTS_OBJECT (cf)->klass);
 408
 409   if (orientation)
 410     gts_triangle_set (GTS_TRIANGLE (cf), e, e2, e1);
 411   else
 412     gts_triangle_set (GTS_TRIANGLE (cf), e, e1, e2);
 413 }
 414
 415 static void split_destroy (GtsObject * object)
 416 {
 417   GtsSplit * vs = GTS_SPLIT (object);
 418   guint i = vs->ncf;
 419   GtsSplitCFace * cf = vs->cfaces;
 420
 421   while (i--) {
 422     if (IS_CFACE (cf->f))
 423       gts_object_destroy (GTS_OBJECT (cf->f));
 424     g_free (cf->a1);
 425     g_free (cf->a2);
 426     cf++;
 427   }
 428   g_free (vs->cfaces);
 429
 430   if (!gts_allow_floating_vertices && vs->v && vs->v->segments == NULL)
 431     gts_object_destroy (GTS_OBJECT (vs->v));
 432
 433   (* GTS_OBJECT_CLASS (gts_split_class ())->parent_class->destroy) (object);
 434 }
 435
 436 static void split_class_init (GtsObjectClass * klass)
 437 {
 438   klass->destroy = split_destroy;
 439 }
 440
 441 static void split_init (GtsSplit * split)
 442 {
 443   split->v1 = split->v2 = NULL;
 444   split->v = NULL;
 445   split->cfaces = NULL;
 446   split->ncf = 0;
 447 }
 448
 449 /**
 450  * gts_split_class:
 451  *
 452  * Returns: the #GtsSplitClass.
 453  */
 454 GtsSplitClass * gts_split_class (void)
 455 {
 456   static GtsSplitClass * klass = NULL;
 457
 458   if (klass == NULL) {
 459     GtsObjectClassInfo split_info = {
 460       "GtsSplit",
 461       sizeof (GtsSplit),
 462       sizeof (GtsSplitClass),
 463       (GtsObjectClassInitFunc) split_class_init,
 464       (GtsObjectInitFunc) split_init,
 465       (GtsArgSetFunc) NULL,
 466       (GtsArgGetFunc) NULL
 467     };
 468     klass = gts_object_class_new (gts_object_class (),
 469                                   &split_info);
 470   }
 471
 472   return klass;
 473 }
 474
 475 #ifdef DEBUG
 476 static gboolean edge_collapse_is_valid (GtsEdge * e)
 477 {
 478   GSList * i;
 479
 480   g_return_val_if_fail (e != NULL, FALSE);
 481
 482   if (gts_segment_is_duplicate (GTS_SEGMENT (e))) {
 483     g_warning ("collapsing duplicate edge");
 484     return FALSE;
 485   }
 486
 487   i = GTS_SEGMENT (e)->v1->segments;
 488   while (i) {
 489     GtsEdge * e1 = i->data;
 490     if (e1 != e && GTS_IS_EDGE (e1)) {
 491       GtsEdge * e2 = NULL;
 492       GSList * j = GTS_SEGMENT (e1)->v1 == GTS_SEGMENT (e)->v1 ?
 493         GTS_SEGMENT (e1)->v2->segments : GTS_SEGMENT (e1)->v1->segments;
 494       while (j && !e2) {
 495         GtsEdge * e1 = j->data;
 496         if (GTS_IS_EDGE (e1) &&
 497             (GTS_SEGMENT (e1)->v1 == GTS_SEGMENT (e)->v2 ||
 498              GTS_SEGMENT (e1)->v2 == GTS_SEGMENT (e)->v2))
 499           e2 = e1;
 500         j = j->next;
 501       }
 502       if (e2 && !gts_triangle_use_edges (e, e1, e2)) {
 503         g_warning ("collapsing empty triangle");
 504         return FALSE;
 505       }
 506     }
 507     i = i->next;
 508   }
 509
 510   if (gts_edge_is_boundary (e, NULL)) {
 511     GtsTriangle * t = e->triangles->data;
 512     if (gts_edge_is_boundary (t->e1, NULL) &&
 513         gts_edge_is_boundary (t->e2, NULL) &&
 514         gts_edge_is_boundary (t->e3, NULL)) {
 515       g_warning ("collapsing single triangle");
 516       return FALSE;
 517     }
 518   }
 519   else {
 520     if (gts_vertex_is_boundary (GTS_SEGMENT (e)->v1, NULL) &&
 521         gts_vertex_is_boundary (GTS_SEGMENT (e)->v2, NULL)) {
 522       g_warning ("collapsing two sides of a strip");
 523       return FALSE;
 524     }
 525     if (gts_edge_belongs_to_tetrahedron (e)) {
 526       g_warning ("collapsing tetrahedron");
 527       return FALSE;
 528     }
 529   }
 530
 531   return TRUE;
 532 }
 533 #endif /* DEBUG */
 534
 535 /* Not currently used.  May be useful for some debug code */
 536 #ifdef DEBUG
 537 static void print_split (GtsSplit * vs, FILE * fptr)
 538 {
 539   guint j;
 540   GtsSplitCFace * cf;
 541
 542   g_return_if_fail (vs != NULL);
 543   g_return_if_fail (fptr != NULL);
 544
 545   fprintf (fptr, "%p: v: %p v1: %p v2: %p ncf: %u cfaces: %p\n",
 546            vs, vs->v, vs->v1, vs->v2, vs->ncf, vs->cfaces);
 547   cf = vs->cfaces;
 548   j = vs->ncf;
 549   while (j--) {
 550     fprintf (stderr, "  f: %p a1: %p a2: %p\n",
 551              cf->f, cf->a1, cf->a2);
 552     cf++;
 553   }
 554 }
 555 #endif
 556
 557 /**
 558  * gts_split_collapse:
 559  * @vs: a #GtsSplit.
 560  * @klass: a #GtsEdgeClass.
 561  * @heap: a #GtsEHeap or %NULL.
 562  *
 563  * Collapses the vertex split @vs. Any new edge created during the process will
 564  * be of class @klass. If heap is not %NULL, the new edges will be inserted
 565  * into it and the destroyed edges will be removed from it.
 566  */
 567 void gts_split_collapse (GtsSplit * vs,
 568                          GtsEdgeClass * klass,
 569                          GtsEHeap * heap)
 570 {
 571   GtsEdge * e;
 572   GtsVertex * v, * v1, * v2;
 573   GSList * i, * end;
 574 #ifdef DYNAMIC_SPLIT
 575   GtsSplitCFace * cf;
 576   guint j;
 577 #endif
 578 #ifdef DEBUG
 579   gboolean invalid = FALSE;
 580   static guint ninvalid = 0;
 581 #endif
 582
 583   g_return_if_fail (vs != NULL);
 584   g_return_if_fail (klass != NULL);
 585
 586   v = vs->v;
 587
 588   g_return_if_fail (v->segments == NULL);
 589
 590   /* we don't want to destroy vertices */
 591   gts_allow_floating_vertices = TRUE;
 592
 593   v1 = GTS_SPLIT_V1 (vs);
 594   v2 = GTS_SPLIT_V2 (vs);
 595   e = GTS_EDGE (gts_vertices_are_connected (v1, v2));
 596   g_assert (e != NULL);
 597
 598 #ifdef DEBUG
 599   fprintf (stderr, "collapsing %p: v1: %p v2: %p v: %p\n", vs, v1, v2, v);
 600   if (!edge_collapse_is_valid (e)) {
 601     char fname[80];
 602     FILE * fptr;
 603     GSList * triangles, * i;
 604
 605     g_warning ("invalid edge collapse");
 606     invalid = TRUE;
 607     sprintf (fname, "invalid.%d", ninvalid);
 608     fptr = fopen (fname, "wt");
 609     gts_write_segment (GTS_SEGMENT (e), GTS_POINT (v), fptr);
 610     triangles = gts_vertex_triangles (v1, NULL);
 611     triangles = gts_vertex_triangles (v2, triangles);
 612     i = triangles;
 613     while (i) {
 614       gts_write_triangle (i->data, GTS_POINT (v), fptr);
 615       i = i->next;
 616     }
 617     g_slist_free (triangles);
 618     fclose (fptr);
 619   }
 620 #endif
 621
 622   i = e->triangles;
 623 #ifdef DYNAMIC_SPLIT
 624   cf = vs->cfaces;
 625   j = vs->ncf;
 626   while (j--) {
 627     g_free (cf->a1);
 628     g_free (cf->a2);
 629     cf++;
 630   }
 631   g_free (vs->cfaces);
 632
 633   vs->ncf = g_slist_length (i);
 634   g_assert (vs->ncf > 0);
 635   cf = vs->cfaces = g_malloc (vs->ncf*sizeof (GtsSplitCFace));
 636 #endif /* DYNAMIC_SPLIT */
 637 #ifdef NEW
 638   while (i) {
 639     cf->f = i->data;
 640     g_assert (GTS_IS_FACE (cf->f));
 641     GTS_OBJECT (cf->f)->klass = GTS_OBJECT_CLASS (cface_class ());
 642     cf++;
 643     i = i->next;
 644   }
 645   i = e->triangles;
 646   cf = vs->cfaces;
 647   while (i) {
 648     cface_new (i->data, e, v1, v2, vs, heap, klass, cf);
 649 #ifdef DEBUG
 650     fprintf (stderr, "cface: %p->%d t: %p->%d a1: ",
 651              cf->f, id (cf->f), CFACE (cf->f)->t, id (CFACE (cf->f)->t));
 652     {
 653       GtsTriangle * t, ** a;
 654       a = cf->a1;
 655       while ((t = *(a++)))
 656         fprintf (stderr, "%p->%d ", t, id (t));
 657       fprintf (stderr, "a2: ");
 658       a = cf->a2;
 659       while ((t = *(a++)))
 660         fprintf (stderr, "%p->%d ", t, id (t));
 661       fprintf (stderr, "\n");
 662     }
 663 #endif
 664     cf++;
 665     i = i->next;
 666   }
 667 #else /* not NEW */
 668   while (i) {
 669     cface_new (i->data, e, v1, v2, vs, heap
 670 #ifdef DYNAMIC_SPLIT
 671                , cf
 672 #endif /* DYNAMIC_SPLIT */
 673                );
 674 #ifdef DYNAMIC_SPLIT
 675     cf->f = i->data;
 676     cf++;
 677 #endif /* DYNAMIC_SPLIT */
 678     i = i->next;
 679   }
 680 #endif /* NEW */
 681   g_slist_free (e->triangles);
 682   e->triangles = NULL;
 683   gts_object_destroy (GTS_OBJECT (e));
 684
 685   gts_allow_floating_vertices = FALSE;
 686
 687   end = NULL;
 688   i = v1->segments;
 689   while (i) {
 690     GtsSegment * s = i->data;
 691     if (s->v1 == v1)
 692       s->v1 = v;
 693     else
 694       s->v2 = v;
 695     end = i;
 696     i = i->next;
 697   }
 698   if (end) {
 699     end->next = v->segments;
 700     v->segments = v1->segments;
 701     v1->segments = NULL;
 702   }
 703
 704   end = NULL;
 705   i = v2->segments;
 706   while (i) {
 707     GtsSegment * s = i->data;
 708     if (s->v1 == v2)
 709       s->v1 = v;
 710     else
 711       s->v2 = v;
 712     end = i;
 713     i = i->next;
 714   }
 715   if (end) {
 716     end->next = v->segments;
 717     v->segments = v2->segments;
 718     v2->segments = NULL;
 719   }
 720
 721 #ifdef DEBUG
 722   if (invalid) {
 723     char fname[80];
 724     FILE * fptr;
 725     GSList * triangles, * i;
 726     GtsSurface * surface = NULL;
 727
 728     sprintf (fname, "invalid_after.%d", ninvalid);
 729     fptr = fopen (fname, "wt");
 730     triangles = gts_vertex_triangles (v, NULL);
 731     i = triangles;
 732     while (i) {
 733       GtsTriangle * t = i->data;
 734       fprintf (stderr, "checking %p->%d\n", t, id (t));
 735       g_assert (GTS_IS_FACE (t));
 736       gts_write_triangle (t, GTS_POINT (v), fptr);
 737       surface = GTS_FACE (t)->surfaces->data;
 738       if (gts_triangle_is_duplicate (t))
 739         fprintf (stderr, "%p->%d is duplicate\n", t, id (t));
 740       if (gts_segment_is_duplicate (GTS_SEGMENT (t->e1)))
 741         fprintf (stderr, "e1 of %p->%d is duplicate\n", t, id (t));
 742       if (gts_segment_is_duplicate (GTS_SEGMENT (t->e2)))
 743         fprintf (stderr, "e2 of %p->%d is duplicate\n", t, id (t));
 744       if (gts_segment_is_duplicate (GTS_SEGMENT (t->e3)))
 745         fprintf (stderr, "e3 of %p->%d is duplicate\n", t, id (t));
 746       i = i->next;
 747     }
 748     fclose (fptr);
 749     g_slist_free (triangles);
 750 #if 0
 751     gts_split_expand (vs, surface);
 752
 753     sprintf (fname, "invalid_after_after.%d", ninvalid);
 754     fptr = fopen (fname, "wt");
 755     triangles = gts_vertex_triangles (v1, NULL);
 756     triangles = gts_vertex_triangles (v2, triangles);
 757     i = triangles;
 758     while (i) {
 759       GtsTriangle * t = i->data;
 760       gts_write_triangle (t, GTS_POINT (v), fptr);
 761       surface = GTS_FACE (t)->surfaces->data;
 762       if (gts_triangle_is_duplicate (t))
 763         fprintf (stderr, "%p->%d is duplicate\n", t, id (t));
 764       if (gts_segment_is_duplicate (GTS_SEGMENT (t->e1)))
 765         fprintf (stderr, "e1 of %p->%d is duplicate\n", t, id (t));
 766       if (gts_segment_is_duplicate (GTS_SEGMENT (t->e2)))
 767         fprintf (stderr, "e2 of %p->%d is duplicate\n", t, id (t));
 768       if (gts_segment_is_duplicate (GTS_SEGMENT (t->e3)))
 769         fprintf (stderr, "e3 of %p->%d is duplicate\n", t, id (t));
 770       i = i->next;
 771     }
 772     fclose (fptr);
 773     g_slist_free (triangles);
 774
 775     exit (1);
 776 #endif
 777     ninvalid++;
 778   }
 779 #endif
 780 }
 781
 782 /**
 783  * gts_split_expand:
 784  * @vs: a #GtsSplit.
 785  * @s: a #GtsSurface.
 786  * @klass: a #GtsEdgeClass.
 787  *
 788  * Expands the vertex split @vs adding the newly created faces to @s. Any
 789  * new edge will be of class @klass.
 790  */
 791 void gts_split_expand (GtsSplit * vs,
 792                        GtsSurface * s,
 793                        GtsEdgeClass * klass)
 794 {
 795   GSList * i;
 796   GtsEdge * e;
 797   GtsVertex * v, * v1, * v2;
 798   gboolean changed = FALSE;
 799   GtsSplitCFace * cf;
 800   guint j;
 801
 802   g_return_if_fail (vs != NULL);
 803   g_return_if_fail (s != NULL);
 804   g_return_if_fail (klass != NULL);
 805
 806   /* we don't want to destroy vertices */
 807   gts_allow_floating_vertices = TRUE;
 808
 809   v1 = GTS_SPLIT_V1 (vs);
 810   v2 = GTS_SPLIT_V2 (vs);
 811   v = vs->v;
 812 #ifdef DEBUG_EXPAND
 813   fprintf (stderr, "expanding %p->%d: v1: %p->%d v2: %p->%d v: %p->%d\n",
 814            vs, id (vs), v1, id (v1), v2, id (v2), v, id (v));
 815 #endif
 816   e = gts_edge_new (klass, v1, v2);
 817   cf = vs->cfaces;
 818   j = vs->ncf;
 819   while (j--) {
 820     cface_expand (CFACE (cf->f), cf->a1, cf->a2, e, v1, v2, v, klass);
 821     gts_surface_add_face (s, cf->f);
 822     cf++;
 823   }
 824
 825   gts_allow_floating_vertices = FALSE;
 826
 827   /* this part is described by figure "expand.fig" */
 828   i = v->segments;
 829   while (i) {
 830     GtsEdge * e1 = i->data;
 831     GtsVertex * with = NULL;
 832     GSList * j = e1->triangles, * next = i->next;
 833     // fprintf (stderr, "e1: %p->%d\n", e1, id (e1));
 834     while (j && !with) {
 835       with = GTS_OBJECT (j->data)->reserved;
 836       j = j->next;
 837     }
 838     if (with) {
 839       j = e1->triangles;
 840       while (j) {
 841         GtsTriangle * t = j->data;
 842         if (GTS_OBJECT (t)->reserved) {
 843           g_assert (GTS_OBJECT (t)->reserved == with);
 844           GTS_OBJECT (t)->reserved = NULL;
 845         }
 846         else
 847           GTS_OBJECT (t)->reserved = with;
 848         j = j->next;
 849       }
 850       if (GTS_SEGMENT (e1)->v1 == v)
 851         GTS_SEGMENT (e1)->v1 = with;
 852       else
 853         GTS_SEGMENT (e1)->v2 = with;
 854
 855       v->segments = g_slist_remove_link (v->segments, i);
 856       i->next = with->segments;
 857       with->segments = i;
 858       changed = TRUE;
 859     }
 860     if (next)
 861       i = next;
 862     else {
 863       /* check for infinite loop (the crossed out case in
 864          figure "expand.fig") */
 865       g_assert (changed);
 866       changed = FALSE;
 867       i = v->segments;
 868     }
 869   }
 870 }
 871
 872 #ifndef DYNAMIC_SPLIT
 873 static void cface_neighbors (GtsSplitCFace * cf,
 874                              GtsEdge * e,
 875                              GtsVertex * v1,
 876                              GtsVertex * v2)
 877 {
 878   GtsTriangle * t = GTS_TRIANGLE (cf->f), ** a;
 879   GtsEdge * e1 = t->e1, * e2 = t->e2, * e3 = t->e3;
 880   GSList * i;
 881   guint size;
 882
 883   ROTATE_ORIENT (e, e1, e2, e3);
 884   if (SEGMENT_USE_VERTEX (GTS_SEGMENT (e1), v2)) {
 885     e3 = e1; e1 = e2; e2 = e3;
 886   }
 887
 888   i = e1->triangles;
 889   size = g_slist_length (i)*sizeof (GtsTriangle *);
 890   a = cf->a1 = g_malloc (size > 0 ? size : sizeof (GtsTriangle *));
 891   while (i) {
 892     if (i->data != t)
 893       *(a++) = i->data;
 894     i = i->next;
 895   }
 896   *a = NULL;
 897
 898   i = e2->triangles;
 899   size = g_slist_length (i)*sizeof (GtsTriangle *);
 900   a = cf->a2 = g_malloc (size > 0 ? size : sizeof (GtsTriangle *));
 901   while (i) {
 902     if (i->data != t)
 903       *(a++) = i->data;
 904     i = i->next;
 905   }
 906   *a = NULL;
 907 }
 908 #endif /*ifndef DYNAMIC_SPLIT */
 909
 910 /**
 911  * gts_split_new:
 912  * @klass: a #GtsSplitClass.
 913  * @v: a #GtsVertex.
 914  * @o1: either a #GtsVertex or a #GtsSplit.
 915  * @o2: either a #GtsVertex or a #GtsSplit.
 916  *
 917  * Creates a new #GtsSplit which would collapse @o1 and @o2 into @v. The
 918  * collapse itself is not performed.
 919  *
 920  * Returns: the new #GtsSplit.
 921  */
 922 GtsSplit * gts_split_new (GtsSplitClass * klass,
 923                           GtsVertex * v,
 924                           GtsObject * o1,
 925                           GtsObject * o2)
 926 {
 927   GtsSplit * vs;
 928 #ifndef DYNAMIC_SPLIT
 929   GtsVertex * v1, * v2;
 930   GtsEdge * e;
 931   GSList * i;
 932   GtsSplitCFace * cf;
 933 #endif
 934
 935   g_return_val_if_fail (klass != NULL, NULL);
 936   g_return_val_if_fail (v != NULL, NULL);
 937   g_return_val_if_fail (GTS_IS_SPLIT (o1) || GTS_IS_VERTEX (o1), NULL);
 938   g_return_val_if_fail (GTS_IS_SPLIT (o2) || GTS_IS_VERTEX (o2), NULL);
 939
 940   vs = GTS_SPLIT (gts_object_new (GTS_OBJECT_CLASS (klass)));
 941   vs->v = v;
 942   vs->v1 = o1;
 943   vs->v2 = o2;
 944 #ifdef DYNAMIC_SPLIT
 945   vs->ncf = 0;
 946   vs->cfaces = NULL;
 947 #else
 948   v1 = GTS_SPLIT_V1 (vs);
 949   v2 = GTS_SPLIT_V2 (vs);
 950   e = GTS_EDGE (gts_vertices_are_connected (v1, v2));
 951   g_assert (e != NULL);
 952   i = e->triangles;
 953   vs->ncf = g_slist_length (i);
 954   g_assert (vs->ncf > 0);
 955   cf = vs->cfaces = g_malloc (vs->ncf*sizeof (GtsSplitCFace));
 956   while (i) {
 957     cf->f = i->data;
 958     cface_neighbors (cf, e, v1, v2);
 959     i = i->next;
 960     cf++;
 961   }
 962 #endif
 963
 964   return vs;
 965 }
 966
 967 static gboolean
 968 split_traverse_pre_order (GtsSplit *           vs,
 969                           GtsSplitTraverseFunc func,
 970                           gpointer             data)
 971 {
 972   if (func (vs, data))
 973     return TRUE;
 974   if (GTS_IS_SPLIT (vs->v1) &&
 975       split_traverse_pre_order (GTS_SPLIT (vs->v1), func, data))
 976     return TRUE;
 977   if (GTS_IS_SPLIT (vs->v2) &&
 978       split_traverse_pre_order (GTS_SPLIT (vs->v2), func, data))
 979     return TRUE;
 980   return FALSE;
 981 }
 982
 983 static gboolean
 984 split_depth_traverse_pre_order (GtsSplit *             vs,
 985                                 guint                  depth,
 986                                 GtsSplitTraverseFunc   func,
 987                                 gpointer               data)
 988 {
 989   if (func (vs, data))
 990       return TRUE;
 991
 992   depth--;
 993   if (!depth)
 994     return FALSE;
 995
 996   if (GTS_IS_SPLIT (vs->v1) &&
 997       split_depth_traverse_pre_order (GTS_SPLIT (vs->v1), depth, func, data))
 998     return TRUE;
 999   if (GTS_IS_SPLIT (vs->v2) &&
1000       split_depth_traverse_pre_order (GTS_SPLIT (vs->v2), depth, func, data))
1001     return TRUE;
1002   return FALSE;
1003 }
1004
1005 static gboolean
1006 split_traverse_post_order (GtsSplit *           vs,
1007                            GtsSplitTraverseFunc func,
1008                            gpointer             data)
1009 {
1010   if (GTS_IS_SPLIT (vs->v1) &&
1011       split_traverse_post_order (GTS_SPLIT (vs->v1), func, data))
1012     return TRUE;
1013   if (GTS_IS_SPLIT (vs->v2) &&
1014       split_traverse_post_order (GTS_SPLIT (vs->v2), func, data))
1015     return TRUE;
1016   if (func (vs, data))
1017     return TRUE;
1018   return FALSE;
1019 }
1020
1021 static gboolean
1022 split_depth_traverse_post_order (GtsSplit *           vs,
1023                                  guint                depth,
1024                                  GtsSplitTraverseFunc func,
1025                                  gpointer             data)
1026 {
1027   depth--;
1028   if (depth) {
1029     if (GTS_IS_SPLIT (vs->v1) &&
1030         split_depth_traverse_post_order (GTS_SPLIT (vs->v1),
1031                                          depth, func, data))
1032       return TRUE;
1033     if (GTS_IS_SPLIT (vs->v2) &&
1034         split_depth_traverse_post_order (GTS_SPLIT (vs->v2),
1035                                          depth, func, data))
1036       return TRUE;
1037   }
1038   if (func (vs, data))
1039     return TRUE;
1040   return FALSE;
1041 }
1042
1043 /**
1044  * gts_split_traverse:
1045  * @root: the #GtsSplit to start the traversal from.
1046  * @order: the order in which nodes are visited - G_PRE_ORDER or G_POST_ORDER.
1047  * @depth: the maximum depth of the traversal. Nodes below this depth
1048  * will not be visited. If depth is -1 all nodes in the tree are
1049  * visited. If depth is 1, only the root is visited. If depth is 2,
1050  * the root and its children are visited. And so on.
1051  * @func: the function to call for each visited #GtsHSplit.
1052  * @data: user data to pass to the function.
1053  *
1054  * Traverses the #GtsSplit tree having @root as root. Calls @func for each
1055  * #GtsSplit of the tree in the order specified by @order. If order is set
1056  * to G_PRE_ORDER @func is called for the #GtsSplit then its children, if order
1057  * is set to G_POST_ORDER @func is called for the children and then for the
1058  * #GtsSplit.
1059  */
1060 void gts_split_traverse (GtsSplit *           root,
1061                          GTraverseType        order,
1062                          gint                 depth,
1063                          GtsSplitTraverseFunc func,
1064                          gpointer             data)
1065 {
1066   g_return_if_fail (root != NULL);
1067   g_return_if_fail (func != NULL);
1068   g_return_if_fail (order < G_LEVEL_ORDER);
1069   g_return_if_fail (depth == -1 || depth > 0);
1070
1071   switch (order) {
1072   case G_PRE_ORDER:
1073     if (depth < 0)
1074       split_traverse_pre_order (root, func, data);
1075     else
1076       split_depth_traverse_pre_order (root, depth, func, data);
1077     break;
1078   case G_POST_ORDER:
1079     if (depth < 0)
1080       split_traverse_post_order (root, func, data);
1081     else
1082       split_depth_traverse_post_order (root, depth, func, data);
1083     break;
1084   default:
1085     g_assert_not_reached ();
1086   }
1087 }
1088
1089 /**
1090  * gts_split_height:
1091  * @root: a #GtsSplit.
1092  *
1093  * Returns: the maximum height of the vertex split tree having @root as root.
1094  */
1095 guint gts_split_height (GtsSplit * root)
1096 {
1097   guint height = 0, tmp_height;
1098
1099   g_return_val_if_fail (root != NULL, 0);
1100
1101   if (GTS_IS_SPLIT (root->v1)) {
1102     tmp_height = gts_split_height (GTS_SPLIT (root->v1));
1103     if (tmp_height > height)
1104       height = tmp_height;
1105   }
1106   if (GTS_IS_SPLIT (root->v2)) {
1107     tmp_height = gts_split_height (GTS_SPLIT (root->v2));
1108     if (tmp_height > height)
1109       height = tmp_height;
1110   }
1111
1112   return height + 1;
1113 }
1114
1115 #ifndef DYNAMIC_SPLIT
1116 static gboolean list_array_are_identical (GSList * list,
1117                                           gpointer * array,
1118                                           gpointer excluded)
1119 {
1120   while (list) {
1121     gpointer data = list->data;
1122     if (data != excluded) {
1123       gboolean found = FALSE;
1124       gpointer * a = array;
1125
1126       while (!found && *a)
1127         if (*(a++) == data)
1128           found = TRUE;
1129       if (!found)
1130         return FALSE;
1131     }
1132     list = list->next;
1133   }
1134   return TRUE;
1135 }
1136 #endif /* ifndef DYNAMIC_SPLIT */
1137
1138 #ifndef NEW
1139 gboolean gts_split_is_collapsable (GtsSplit * vs)
1140 {
1141   guint i;
1142   GtsSplitCFace * cf;
1143   GtsVertex * v1, * v2;
1144   GtsEdge * e;
1145
1146   g_return_val_if_fail (vs != NULL, FALSE);
1147
1148   v1 = GTS_SPLIT_V1 (vs);
1149   v2 = GTS_SPLIT_V2 (vs);
1150   g_return_val_if_fail ((e = GTS_EDGE (gts_vertices_are_connected (v1, v2))),
1151                         FALSE);
1152
1153 #ifdef DYNAMIC_SPLIT
1154   if (!gts_edge_collapse_is_valid (e))
1155     return FALSE;
1156 #else
1157   i = vs->ncf;
1158   cf = vs->cfaces;
1159   while (i--) {
1160     GtsTriangle * t = GTS_TRIANGLE (cf->f);
1161     GtsEdge * e1 = t->e1, * e2 = t->e2, * e3 = t->e3;
1162
1163     ROTATE_ORIENT (e, e1, e2, e3);
1164     if (SEGMENT_USE_VERTEX (GTS_SEGMENT (e1), v2)) {
1165       e3 = e1; e1 = e2; e2 = e3;
1166     }
1167
1168     if (!list_array_are_identical (e1->triangles, (gpointer *) cf->a1, t))
1169       return FALSE;
1170     if (!list_array_are_identical (e2->triangles, (gpointer *) cf->a2, t))
1171       return FALSE;
1172
1173     cf++;
1174   }
1175 #endif
1176   return TRUE;
1177 }
1178 #endif /* not NEW */
1179
1180 #ifdef DEBUG_HEXPAND
1181 static guint expand_level = 0;
1182
1183 static void expand_indent (FILE * fptr)
1184 {
1185   guint i = expand_level;
1186   while (i--)
1187     fputc (' ', fptr);
1188 }
1189 #endif
1190
1191 /**
1192  * gts_hsplit_force_expand:
1193  * @hs: a #GtsHSplit.
1194  * @hsurface: a #GtsHSurface.
1195  *
1196  * Forces the expansion of @hs by first expanding all its dependencies not
1197  * already expanded.
1198  */
1199 void gts_hsplit_force_expand (GtsHSplit * hs,
1200                               GtsHSurface * hsurface)
1201 {
1202   guint i;
1203   GtsSplitCFace * cf;
1204
1205   g_return_if_fail (hs != NULL);
1206   g_return_if_fail (hsurface != NULL);
1207   g_return_if_fail (hs->nchild == 0);
1208
1209 #ifdef DEBUG_HEXPAND
1210   expand_level += 2;
1211 #endif
1212
1213   if (hs->parent && hs->parent->nchild == 0) {
1214 #ifdef DEBUG_HEXPAND
1215     expand_indent (stderr);
1216     fprintf (stderr, "expand parent %p\n", hs->parent);
1217 #endif
1218     gts_hsplit_force_expand (hs->parent, hsurface);
1219   }
1220
1221   i = GTS_SPLIT (hs)->ncf;
1222   cf = GTS_SPLIT (hs)->cfaces;
1223   while (i--) {
1224     GtsTriangle ** j, * t;
1225
1226     j = cf->a1;
1227     while ((t = *(j++)))
1228       if (IS_CFACE (t)) {
1229 #ifdef DEBUG_HEXPAND
1230         expand_indent (stderr);
1231         fprintf (stderr, "expand a1: cf->f: %p t: %p parent_split: %p\n",
1232                  cf->f,
1233                  t,
1234                  GTS_HSPLIT (CFACE (t)->parent_split));
1235 #endif
1236         gts_hsplit_force_expand (GTS_HSPLIT (CFACE (t)->parent_split),
1237                                  hsurface);
1238 #ifdef DEBUG_HEXPAND
1239         g_assert (!IS_CFACE (t));
1240 #endif
1241       }
1242     j = cf->a2;
1243     while ((t = *(j++)))
1244       if (IS_CFACE (t)) {
1245 #ifdef DEBUG_HEXPAND
1246         expand_indent (stderr);
1247         fprintf (stderr, "expand a2: cf->f: %p t: %p parent_split: %p\n",
1248                  cf->f,
1249                  t,
1250                  GTS_HSPLIT (CFACE (t)->parent_split));
1251 #endif
1252         gts_hsplit_force_expand (GTS_HSPLIT (CFACE (t)->parent_split),
1253                                  hsurface);
1254       }
1255     cf++;
1256   }
1257
1258   gts_hsplit_expand (hs, hsurface);
1259
1260 #ifdef DEBUG_HEXPAND
1261   expand_level -= 2;
1262   expand_indent (stderr);
1263   fprintf (stderr, "%p expanded\n", hs);
1264 #endif
1265 }
1266
1267 static void index_object (GtsObject * o, guint * n)
1268 {
1269   o->reserved = GUINT_TO_POINTER ((*n)++);
1270 }
1271
1272 static void index_face (GtsFace * f, gpointer * data)
1273 {
1274   guint * nf = data[1];
1275
1276   g_hash_table_insert (data[0], f, GUINT_TO_POINTER ((*nf)++));
1277 }
1278
1279 /**
1280  * gts_psurface_write:
1281  * @ps: a #GtsPSurface.
1282  * @fptr: a file pointer.
1283  *
1284  * Writes to @fptr a GTS progressive surface description.
1285  */
1286 void gts_psurface_write (GtsPSurface * ps, FILE * fptr)
1287 {
1288   guint nv = 1;
1289   guint nf = 1;
1290   GHashTable * hash;
1291   gpointer data[2];
1292
1293   g_return_if_fail (ps != NULL);
1294   g_return_if_fail (fptr != NULL);
1295   g_return_if_fail (GTS_PSURFACE_IS_CLOSED (ps));
1296
1297   while (gts_psurface_remove_vertex (ps))
1298     ;
1299
1300   GTS_POINT_CLASS (ps->s->vertex_class)->binary = FALSE;
1301   gts_surface_write (ps->s, fptr);
1302
1303   gts_surface_foreach_vertex (ps->s, (GtsFunc) index_object, &nv);
1304   hash = g_hash_table_new (NULL, NULL);
1305   data[0] = hash;
1306   data[1] = &nf;
1307   gts_surface_foreach_face (ps->s, (GtsFunc) index_face, data);
1308
1309   fprintf (fptr, "%u\n", ps->split->len);
1310   while (ps->pos) {
1311     GtsSplit * vs = g_ptr_array_index (ps->split, --ps->pos);
1312     GtsSplitCFace * scf = vs->cfaces;
1313     GtsVertex * v1, * v2;
1314     guint i = vs->ncf;
1315
1316     fprintf (fptr, "%u %u",
1317              GPOINTER_TO_UINT (GTS_OBJECT (vs->v)->reserved),
1318              vs->ncf);
1319     if (GTS_OBJECT (vs)->klass->write)
1320       (*GTS_OBJECT (vs)->klass->write) (GTS_OBJECT (vs), fptr);
1321     fputc ('\n', fptr);
1322
1323     v1 = GTS_IS_SPLIT (vs->v1) ? GTS_SPLIT (vs->v1)->v : GTS_VERTEX (vs->v1);
1324     GTS_OBJECT (v1)->reserved = GUINT_TO_POINTER (nv++);
1325     v2 = GTS_IS_SPLIT (vs->v2) ? GTS_SPLIT (vs->v2)->v : GTS_VERTEX (vs->v2);
1326     GTS_OBJECT (v2)->reserved = GUINT_TO_POINTER (nv++);
1327
1328     (*GTS_OBJECT (v1)->klass->write) (GTS_OBJECT (v1), fptr);
1329     fputc ('\n', fptr);
1330
1331     (*GTS_OBJECT (v2)->klass->write) (GTS_OBJECT (v2), fptr);
1332     fputc ('\n', fptr);
1333
1334     while (i--) {
1335       CFace * cf = CFACE (scf->f);
1336       GtsTriangle ** a, * t;
1337
1338       fprintf (fptr, "%u %u",
1339                GPOINTER_TO_UINT (g_hash_table_lookup (hash, cf->t)),
1340                cf->flags);
1341       if (GTS_OBJECT_CLASS (ps->s->face_class)->write)
1342         (*GTS_OBJECT_CLASS (ps->s->face_class)->write) (GTS_OBJECT (cf), fptr);
1343       fputc ('\n', fptr);
1344
1345       a = scf->a1;
1346       while ((t = *(a++)))
1347         fprintf (fptr, "%u ",
1348                  GPOINTER_TO_UINT (g_hash_table_lookup (hash, t)));
1349       fprintf (fptr, "\n");
1350
1351       a = scf->a2;
1352       while ((t = *(a++)))
1353         fprintf (fptr, "%u ",
1354                  GPOINTER_TO_UINT (g_hash_table_lookup (hash, t)));
1355       fprintf (fptr, "\n");
1356
1357       g_hash_table_insert (hash, cf, GUINT_TO_POINTER (nf++));
1358
1359       scf++;
1360     }
1361
1362     gts_split_expand (vs, ps->s, ps->s->edge_class);
1363   }
1364
1365   gts_surface_foreach_vertex (ps->s,
1366                               (GtsFunc) gts_object_reset_reserved, NULL);
1367   g_hash_table_destroy (hash);
1368 }
1369
1370 static guint surface_read (GtsSurface * surface,
1371                            GtsFile * f,
1372                            GPtrArray * vertices,
1373                            GPtrArray * faces)
1374 {
1375   GtsEdge ** edges;
1376   guint n, nv, ne, nf;
1377
1378   g_return_val_if_fail (surface != NULL, 1);
1379   g_return_val_if_fail (f != NULL, 1);
1380
1381   if (f->type != GTS_INT) {
1382     gts_file_error (f, "expecting an integer (number of vertices)");
1383     return f->line;
1384   }
1385   nv = atoi (f->token->str);
1386
1387   gts_file_next_token (f);
1388   if (f->type != GTS_INT) {
1389     gts_file_error (f, "expecting an integer (number of edges)");
1390     return f->line;
1391   }
1392   ne = atoi (f->token->str);
1393
1394   gts_file_next_token (f);
1395   if (f->type != GTS_INT) {
1396     gts_file_error (f, "expecting an integer (number of faces)");
1397     return f->line;
1398   }
1399   nf = atoi (f->token->str);
1400
1401   gts_file_next_token (f);
1402   if (f->type == GTS_STRING) {
1403     if (f->type != GTS_STRING) {
1404       gts_file_error (f, "expecting a string (GtsSurfaceClass)");
1405       return f->line;
1406     }
1407     gts_file_next_token (f);
1408     if (f->type != GTS_STRING) {
1409       gts_file_error (f, "expecting a string (GtsFaceClass)");
1410       return f->line;
1411     }
1412     gts_file_next_token (f);
1413     if (f->type != GTS_STRING) {
1414       gts_file_error (f, "expecting a string (GtsEdgeClass)");
1415       return f->line;
1416     }
1417     gts_file_next_token (f);
1418     if (f->type != GTS_STRING) {
1419       gts_file_error (f, "expecting a string (GtsVertexClass)");
1420       return f->line;
1421     }
1422     if (!strcmp (f->token->str, "GtsVertexBinary"))
1423       GTS_POINT_CLASS (surface->vertex_class)->binary = TRUE;
1424     else
1425       gts_file_first_token_after (f, '\n');
1426   }
1427   else
1428     gts_file_first_token_after (f, '\n');
1429
1430   g_ptr_array_set_size (vertices, nv);
1431   g_ptr_array_set_size (faces, nf);
1432   /* allocate nv + 1 just in case nv == 0 */
1433   edges = g_malloc ((ne + 1)*sizeof (GtsEdge *));
1434
1435   n = 0;
1436   while (n < nv && f->type != GTS_ERROR) {
1437     GtsObject * new_vertex =
1438       gts_object_new (GTS_OBJECT_CLASS (surface->vertex_class));
1439
1440     (* GTS_OBJECT_CLASS (surface->vertex_class)->read) (&new_vertex, f);
1441     if (f->type != GTS_ERROR) {
1442       if (!GTS_POINT_CLASS (surface->vertex_class)->binary)
1443         gts_file_first_token_after (f, '\n');
1444       g_ptr_array_index (vertices, n++) = new_vertex;
1445     }
1446     else
1447       gts_object_destroy (new_vertex);
1448   }
1449   if (f->type == GTS_ERROR)
1450     nv = n;
1451   if (GTS_POINT_CLASS (surface->vertex_class)->binary)
1452     gts_file_first_token_after (f, '\n');
1453
1454   n = 0;
1455   while (n < ne && f->type != GTS_ERROR) {
1456     guint p1, p2;
1457
1458     if (f->type != GTS_INT)
1459       gts_file_error (f, "expecting an integer (first vertex index)");
1460     else {
1461       p1 = atoi (f->token->str);
1462       if (p1 == 0 || p1 > nv)
1463         gts_file_error (f, "vertex index `%d' is out of range `[1,%d]'",
1464                         p1, nv);
1465       else {
1466         gts_file_next_token (f);
1467         if (f->type != GTS_INT)
1468           gts_file_error (f, "expecting an integer (second vertex index)");
1469         else {
1470           p2 = atoi (f->token->str);
1471           if (p2 == 0 || p2 > nv)
1472             gts_file_error (f, "vertex index `%d' is out of range `[1,%d]'",
1473                             p2, nv);
1474           else {
1475             GtsEdge * new_edge =
1476               gts_edge_new (surface->edge_class,
1477                             g_ptr_array_index (vertices, p1 - 1),
1478                             g_ptr_array_index (vertices, p2 - 1));
1479
1480             gts_file_next_token (f);
1481             if (f->type != '\n')
1482               if (GTS_OBJECT_CLASS (surface->edge_class)->read)
1483                 (*GTS_OBJECT_CLASS (surface->edge_class)->read)
1484                   ((GtsObject **) &new_edge, f);
1485             gts_file_first_token_after (f, '\n');
1486             edges[n++] = new_edge;
1487           }
1488         }
1489       }
1490     }
1491   }
1492   if (f->type == GTS_ERROR)
1493     ne = n;
1494
1495   n = 0;
1496   while (n < nf && f->type != GTS_ERROR) {
1497     guint s1, s2, s3;
1498
1499     if (f->type != GTS_INT)
1500       gts_file_error (f, "expecting an integer (first edge index)");
1501     else {
1502       s1 = atoi (f->token->str);
1503       if (s1 == 0 || s1 > ne)
1504         gts_file_error (f, "edge index `%d' is out of range `[1,%d]'",
1505                         s1, ne);
1506       else {
1507         gts_file_next_token (f);
1508         if (f->type != GTS_INT)
1509           gts_file_error (f, "expecting an integer (second edge index)");
1510         else {
1511           s2 = atoi (f->token->str);
1512           if (s2 == 0 || s2 > ne)
1513             gts_file_error (f, "edge index `%d' is out of range `[1,%d]'",
1514                             s2, ne);
1515           else {
1516             gts_file_next_token (f);
1517             if (f->type != GTS_INT)
1518               gts_file_error (f, "expecting an integer (third edge index)");
1519             else {
1520               s3 = atoi (f->token->str);
1521               if (s3 == 0 || s3 > ne)
1522                 gts_file_error (f, "edge index `%d' is out of range `[1,%d]'",
1523                                 s3, ne);
1524               else {
1525                 GtsFace * new_face = gts_face_new (surface->face_class,
1526                                                    edges[s1 - 1],
1527                                                    edges[s2 - 1],
1528                                                    edges[s3 - 1]);
1529
1530                 gts_file_next_token (f);
1531                 if (f->type != '\n')
1532                   if (GTS_OBJECT_CLASS (surface->face_class)->read)
1533                     (*GTS_OBJECT_CLASS (surface->face_class)->read)
1534                       ((GtsObject **) &new_face, f);
1535                 gts_file_first_token_after (f, '\n');
1536                 gts_surface_add_face (surface, new_face);
1537                 g_ptr_array_index (faces, n++) = new_face;
1538               }
1539             }
1540           }
1541         }
1542       }
1543     }
1544   }
1545
1546   g_free (edges);
1547
1548   if (f->type == GTS_ERROR) {
1549     gts_allow_floating_vertices = TRUE;
1550     while (nv)
1551       gts_object_destroy (GTS_OBJECT (g_ptr_array_index (vertices, nv-- - 1)));
1552     gts_allow_floating_vertices = FALSE;
1553     return f->line;
1554   }
1555
1556   return 0;
1557 }
1558
1559 /**
1560  * gts_psurface_open:
1561  * @klass: a #GtsPSurfaceClass.
1562  * @s: a #GtsSurface.
1563  * @split_class: a #GtsSplitClass to use for the #GtsSplit.
1564  * @f: a #GtsFile.
1565  *
1566  * Creates a new #GtsPSurface prepared for input from the file @f
1567  * containing a valid GTS representation of a progressive surface. The initial
1568  * shape of the progressive surface is loaded into @s.
1569  *
1570  * Before being usable as such this progressive surface must be closed using
1571  * gts_psurface_close(). While open however, the functions
1572  * gts_psurface_get_vertex_number(), gts_psurface_min_vertex_number() and
1573  * gts_psurface_max_vertex_number() can still be used.
1574  *
1575  * Returns: a new #GtsPSurface or %NULL if there was a format error while
1576  * reading the file, in which case @f contains information about the error.
1577  */
1578 GtsPSurface * gts_psurface_open (GtsPSurfaceClass * klass,
1579                                  GtsSurface * s,
1580                                  GtsSplitClass * split_class,
1581                                  GtsFile * f)
1582 {
1583   GtsPSurface * ps;
1584
1585   g_return_val_if_fail (klass != NULL, NULL);
1586   g_return_val_if_fail (s != NULL, NULL);
1587   g_return_val_if_fail (split_class != NULL, NULL);
1588   g_return_val_if_fail (f != NULL, NULL);
1589
1590   ps = GTS_PSURFACE (gts_object_new (GTS_OBJECT_CLASS (klass)));
1591   ps->s = s;
1592   ps->split_class = split_class;
1593
1594   ps->vertices = g_ptr_array_new ();
1595   ps->faces = g_ptr_array_new ();
1596
1597   if (surface_read (s, f, ps->vertices, ps->faces)) {
1598     ps->s = NULL;
1599     gts_object_destroy (GTS_OBJECT (ps));
1600     return NULL;
1601   }
1602
1603   ps->min = gts_surface_vertex_number (ps->s);
1604   ps->pos = 0;
1605
1606   if (f->type == GTS_INT) {
1607     gint ns = atoi (f->token->str);
1608
1609     if (ns > 0) {
1610       g_ptr_array_set_size (ps->split, ns);
1611       gts_file_first_token_after (f, '\n');
1612     }
1613   }
1614
1615   return ps;
1616 }
1617
1618 /**
1619  * gts_psurface_read_vertex:
1620  * @ps: a #GtsPSurface prealably created with gts_psurface_open().
1621  * @fp: a #GtsFile.
1622  *
1623  * Reads in one vertex split operation from @fp and performs the expansion.
1624  *
1625  * If an error occurs while reading the file, the @error field of @fp is set.
1626  *
1627  * Returns: the newly created #GtsSplit or %NULL if no vertex split could be
1628  * read from @fp.
1629  */
1630 GtsSplit * gts_psurface_read_vertex (GtsPSurface * ps, GtsFile * fp)
1631 {
1632   guint nv, ncf;
1633   GtsSplit * vs, * parent;
1634   GtsSplitCFace * scf;
1635
1636   g_return_val_if_fail (ps != NULL, NULL);
1637   g_return_val_if_fail (fp != NULL, NULL);
1638   g_return_val_if_fail (!GTS_PSURFACE_IS_CLOSED (ps), NULL);
1639
1640   if (ps->pos >= ps->split->len)
1641     return NULL;
1642
1643   if (fp->type == GTS_NONE)
1644     return NULL;
1645   if (fp->type != GTS_INT) {
1646     gts_file_error (fp, "expecting an integer (vertex index)");
1647     return NULL;
1648   }
1649   nv = atoi (fp->token->str);
1650   if (nv == 0 || nv > ps->vertices->len) {
1651     gts_file_error (fp, "vertex index `%d' is out of range `[1,%d]'",
1652                     nv, ps->vertices->len);
1653     return NULL;
1654   }
1655
1656   gts_file_next_token (fp);
1657   if (fp->type != GTS_INT) {
1658     gts_file_error (fp, "expecting an integer (ncf)");
1659     return NULL;
1660   }
1661   ncf = atoi (fp->token->str);
1662
1663   vs = GTS_SPLIT (gts_object_new (GTS_OBJECT_CLASS (ps->split_class)));
1664
1665   vs->v = g_ptr_array_index (ps->vertices, nv - 1);
1666   vs->v1 = vs->v2 = NULL;
1667   vs->cfaces = NULL;
1668   vs->ncf = 0;
1669
1670   gts_file_next_token (fp);
1671   if (fp->type != '\n')
1672     if (GTS_OBJECT (vs)->klass->read)
1673       (* GTS_OBJECT (vs)->klass->read) ((GtsObject **) &vs, fp);
1674   gts_file_first_token_after (fp, '\n');
1675
1676   if (fp->type != GTS_ERROR) {
1677     vs->v1 = gts_object_new (GTS_OBJECT_CLASS (ps->s->vertex_class));
1678     (* GTS_OBJECT_CLASS (ps->s->vertex_class)->read) (&(vs->v1), fp);
1679     if (fp->type != GTS_ERROR) {
1680       vs->v1->reserved = vs;
1681       g_ptr_array_add (ps->vertices, vs->v1);
1682
1683       gts_file_first_token_after (fp, '\n');
1684
1685       vs->v2 = gts_object_new (GTS_OBJECT_CLASS (ps->s->vertex_class));
1686       (*GTS_OBJECT_CLASS (ps->s->vertex_class)->read) (&(vs->v2), fp);
1687       if (fp->type != GTS_ERROR) {
1688         vs->v2->reserved = vs;
1689         g_ptr_array_add (ps->vertices, vs->v2);
1690         gts_file_first_token_after (fp, '\n');
1691       }
1692     }
1693   }
1694
1695   if (fp->type != GTS_ERROR) {
1696     scf = vs->cfaces = g_malloc (sizeof (GtsSplitCFace)*ncf);
1697     while (fp->type != GTS_ERROR && ncf--) {
1698       guint it, flags;
1699       GtsFace * f;
1700       CFace * cf;
1701       GPtrArray * a;
1702
1703       if (fp->type != GTS_INT)
1704         gts_file_error (fp, "expecting an integer (face index)");
1705       else {
1706         it = atoi (fp->token->str);
1707         if (it == 0 || it > ps->faces->len)
1708           gts_file_error (fp, "face index `%d' is out of range `[1,%d]'",
1709                           it, ps->faces->len);
1710         else {
1711           gts_file_next_token (fp);
1712           if (fp->type != GTS_INT)
1713             gts_file_error (fp, "expecting an integer (flags)");
1714           else {
1715             flags = atoi (fp->token->str);
1716             f =
1717               GTS_FACE (gts_object_new (GTS_OBJECT_CLASS (ps->s->face_class)));
1718
1719             gts_file_next_token (fp);
1720             if (fp->type != '\n')
1721               if (GTS_OBJECT (f)->klass->read)
1722                 (*GTS_OBJECT (f)->klass->read) ((GtsObject **) &f, fp);
1723             gts_file_first_token_after (fp, '\n');
1724             if (fp->type != GTS_ERROR) {
1725               scf->f = f;
1726
1727               cf = (CFace *) f;
1728               GTS_OBJECT (cf)->klass = GTS_OBJECT_CLASS (cface_class ());
1729               cf->parent_split = vs;
1730               cf->t = g_ptr_array_index (ps->faces, it - 1);
1731               cf->flags = flags;
1732
1733               a = g_ptr_array_new ();
1734               do {
1735                 if (fp->type != GTS_INT)
1736                   gts_file_error (fp, "expecting an integer (face index)");
1737                 else {
1738                   it = atoi (fp->token->str);
1739                   if (it > ps->faces->len)
1740                     gts_file_error (fp,
1741                                     "face index `%d' is out of range `[1,%d]'",
1742                                     it, ps->faces->len);
1743                   else {
1744                     g_ptr_array_add (a, g_ptr_array_index (ps->faces,
1745                                                            it - 1));
1746                     gts_file_next_token (fp);
1747                   }
1748                 }
1749               } while (fp->type != GTS_ERROR && fp->type != '\n');
1750               gts_file_first_token_after (fp, '\n');
1751               g_ptr_array_add (a, NULL);
1752               scf->a1 = (GtsTriangle **) a->pdata;
1753               g_ptr_array_free (a, FALSE);
1754
1755               if (fp->type != GTS_ERROR) {
1756                 a = g_ptr_array_new ();
1757                 do {
1758                   if (fp->type != GTS_INT)
1759                     gts_file_error (fp, "expecting an integer (face index)");
1760                   else {
1761                     it = atoi (fp->token->str);
1762                     if (it > ps->faces->len)
1763                       gts_file_error (fp,
1764                                    "face index `%d' is out of range `[1,%d]'",
1765                                       it, ps->faces->len);
1766                     else {
1767                       g_ptr_array_add (a, g_ptr_array_index (ps->faces,
1768                                                              it - 1));
1769                       gts_file_next_token (fp);
1770                     }
1771                   }
1772                 } while (fp->type != GTS_ERROR && fp->type != '\n');
1773                 gts_file_first_token_after (fp, '\n');
1774                 g_ptr_array_add (a, NULL);
1775                 scf->a2 = (GtsTriangle **) a->pdata;
1776                 g_ptr_array_free (a, FALSE);
1777
1778                 g_ptr_array_add (ps->faces, f);
1779
1780                 vs->ncf++;
1781                 scf++;
1782               }
1783             }
1784           }
1785         }
1786       }
1787     }
1788   }
1789
1790   if (fp->type != GTS_ERROR) {
1791     if ((parent = GTS_OBJECT (vs->v)->reserved)) {
1792       GTS_OBJECT (vs->v)->reserved = NULL;
1793       if (parent->v1 == GTS_OBJECT (vs->v))
1794         parent->v1 = GTS_OBJECT (vs);
1795       else {
1796         g_assert (parent->v2 == GTS_OBJECT (vs->v));
1797         parent->v2 = GTS_OBJECT (vs);
1798       }
1799     }
1800     g_ptr_array_index (ps->split, ps->pos++) = vs;
1801     gts_split_expand (vs, ps->s, ps->s->edge_class);
1802
1803     return vs;
1804   }
1805
1806   if (vs->v1) gts_object_destroy (vs->v1);
1807   if (vs->v2) gts_object_destroy (vs->v2);
1808   gts_object_destroy (GTS_OBJECT (vs));
1809
1810   return NULL;
1811 }
1812
1813 /**
1814  * gts_psurface_close:
1815  * @ps: a #GtsPSurface prealably created with gts_psurface_open().
1816  *
1817  * Closes a progressive surface.
1818  */
1819 void gts_psurface_close (GtsPSurface * ps)
1820 {
1821   g_return_if_fail (ps != NULL);
1822   g_return_if_fail (!GTS_PSURFACE_IS_CLOSED (ps));
1823
1824   g_ptr_array_free (ps->vertices, TRUE);
1825   g_ptr_array_free (ps->faces, TRUE);
1826   ps->faces = ps->vertices = NULL;
1827
1828   gts_surface_foreach_vertex (ps->s,
1829                               (GtsFunc) gts_object_reset_reserved, NULL);
1830   if (ps->pos > 0)
1831     g_ptr_array_set_size (ps->split, ps->pos);
1832   if (ps->split->len > 1) {
1833     guint i, half = ps->split->len/2, n = ps->split->len - 1;
1834
1835     for (i = 0; i < half; i++) {
1836       gpointer p1 = g_ptr_array_index (ps->split, i);
1837       gpointer p2 = g_ptr_array_index (ps->split, n - i);
1838       g_ptr_array_index (ps->split, n - i) = p1;
1839       g_ptr_array_index (ps->split, i) = p2;
1840     }
1841   }
1842   ps->pos = 0;
1843 }