PostGIS  2.2.7dev-r@@SVN_REVISION@@
lwgeom_topo.c
Go to the documentation of this file.
1 /**********************************************************************
2  *
3  * PostGIS - Spatial Types for PostgreSQL
4  * http://postgis.net
5  *
6  * Copyright (C) 2015 Sandro Santilli <strk@keybit.net>
7  *
8  * This is free software; you can redistribute and/or modify it under
9  * the terms of the GNU General Public Licence. See the COPYING file.
10  *
11  **********************************************************************
12  *
13  * Topology extension for liblwgeom.
14  * Initially funded by Tuscany Region (Italy) - SITA (CIG: 60351023B8)
15  *
16  **********************************************************************/
17 
18 
19 #include "../postgis_config.h"
20 
21 /*#define POSTGIS_DEBUG_LEVEL 1*/
22 #include "lwgeom_log.h"
23 
24 #include "liblwgeom_internal.h"
26 #include "lwgeom_geos.h"
27 
28 #include <stdio.h>
29 #include <inttypes.h> /* for PRId64 */
30 #include <errno.h>
31 #include <math.h>
32 
33 #ifdef WIN32
34 # define LWTFMT_ELEMID "lld"
35 #else
36 # define LWTFMT_ELEMID PRId64
37 #endif
38 
39 /* TODO: move this to lwgeom_log.h */
40 #define LWDEBUGG(level, geom, msg) \
41  if (POSTGIS_DEBUG_LEVEL >= level) \
42  do { \
43  size_t sz; \
44  char *wkt1 = lwgeom_to_wkt(geom, WKT_EXTENDED, 15, &sz); \
45  /* char *wkt1 = lwgeom_to_hexwkb(geom, WKT_EXTENDED, &sz); */ \
46  LWDEBUGF(level, msg ": %s", wkt1); \
47  lwfree(wkt1); \
48  } while (0);
49 
50 
51 /*********************************************************************
52  *
53  * Backend iface
54  *
55  ********************************************************************/
56 
58 {
59  LWT_BE_IFACE *iface = lwalloc(sizeof(LWT_BE_IFACE));
60  iface->data = data;
61  iface->cb = NULL;
62  return iface;
63 }
64 
66  const LWT_BE_CALLBACKS* cb)
67 {
68  iface->cb = cb;
69 }
70 
72 {
73  lwfree(iface);
74 }
75 
76 /*********************************************************************
77  *
78  * Backend wrappers
79  *
80  ********************************************************************/
81 
82 #define CHECKCB(be, method) do { \
83  if ( ! (be)->cb || ! (be)->cb->method ) \
84  lwerror("Callback " # method " not registered by backend"); \
85 } while (0)
86 
87 #define CB0(be, method) \
88  CHECKCB(be, method);\
89  return (be)->cb->method((be)->data)
90 
91 #define CB1(be, method, a1) \
92  CHECKCB(be, method);\
93  return (be)->cb->method((be)->data, a1)
94 
95 #define CBT0(to, method) \
96  CHECKCB((to)->be_iface, method);\
97  return (to)->be_iface->cb->method((to)->be_topo)
98 
99 #define CBT1(to, method, a1) \
100  CHECKCB((to)->be_iface, method);\
101  return (to)->be_iface->cb->method((to)->be_topo, a1)
102 
103 #define CBT2(to, method, a1, a2) \
104  CHECKCB((to)->be_iface, method);\
105  return (to)->be_iface->cb->method((to)->be_topo, a1, a2)
106 
107 #define CBT3(to, method, a1, a2, a3) \
108  CHECKCB((to)->be_iface, method);\
109  return (to)->be_iface->cb->method((to)->be_topo, a1, a2, a3)
110 
111 #define CBT4(to, method, a1, a2, a3, a4) \
112  CHECKCB((to)->be_iface, method);\
113  return (to)->be_iface->cb->method((to)->be_topo, a1, a2, a3, a4)
114 
115 #define CBT5(to, method, a1, a2, a3, a4, a5) \
116  CHECKCB((to)->be_iface, method);\
117  return (to)->be_iface->cb->method((to)->be_topo, a1, a2, a3, a4, a5)
118 
119 #define CBT6(to, method, a1, a2, a3, a4, a5, a6) \
120  CHECKCB((to)->be_iface, method);\
121  return (to)->be_iface->cb->method((to)->be_topo, a1, a2, a3, a4, a5, a6)
122 
123 const char *
125 {
126  CB0(be, lastErrorMessage);
127 }
128 
131 {
132  CB1(be, loadTopologyByName, name);
133 }
134 
135 static int
137 {
138  CBT0(topo, topoGetSRID);
139 }
140 
141 static double
143 {
144  CBT0(topo, topoGetPrecision);
145 }
146 
147 static int
149 {
150  CBT0(topo, topoHasZ);
151 }
152 
153 int
155 {
156  CBT0(topo, freeTopology);
157 }
158 
161  int* numelems, int fields)
162 {
163  CBT3(topo, getNodeById, ids, numelems, fields);
164 }
165 
168  double dist, int* numelems, int fields,
169  int limit)
170 {
171  CBT5(topo, getNodeWithinDistance2D, pt, dist, numelems, fields, limit);
172 }
173 
174 static LWT_ISO_NODE*
176  const GBOX* box, int* numelems, int fields,
177  int limit )
178 {
179  CBT4(topo, getNodeWithinBox2D, box, numelems, fields, limit);
180 }
181 
182 static LWT_ISO_EDGE*
184  const GBOX* box, int* numelems, int fields,
185  int limit )
186 {
187  CBT4(topo, getEdgeWithinBox2D, box, numelems, fields, limit);
188 }
189 
190 static LWT_ISO_FACE*
192  const GBOX* box, int* numelems, int fields,
193  int limit )
194 {
195  CBT4(topo, getFaceWithinBox2D, box, numelems, fields, limit);
196 }
197 
198 int
199 lwt_be_insertNodes(LWT_TOPOLOGY* topo, LWT_ISO_NODE* node, int numelems)
200 {
201  CBT2(topo, insertNodes, node, numelems);
202 }
203 
204 static int
205 lwt_be_insertFaces(LWT_TOPOLOGY* topo, LWT_ISO_FACE* face, int numelems)
206 {
207  CBT2(topo, insertFaces, face, numelems);
208 }
209 
210 static int
211 lwt_be_deleteFacesById(const LWT_TOPOLOGY* topo, const LWT_ELEMID* ids, int numelems)
212 {
213  CBT2(topo, deleteFacesById, ids, numelems);
214 }
215 
216 static int
217 lwt_be_deleteNodesById(const LWT_TOPOLOGY* topo, const LWT_ELEMID* ids, int numelems)
218 {
219  CBT2(topo, deleteNodesById, ids, numelems);
220 }
221 
224 {
225  CBT0(topo, getNextEdgeId);
226 }
227 
230  int* numelems, int fields)
231 {
232  CBT3(topo, getEdgeById, ids, numelems, fields);
233 }
234 
235 static LWT_ISO_FACE*
237  int* numelems, int fields)
238 {
239  CBT3(topo, getFaceById, ids, numelems, fields);
240 }
241 
242 static LWT_ISO_EDGE*
244  int* numelems, int fields)
245 {
246  CBT3(topo, getEdgeByNode, ids, numelems, fields);
247 }
248 
249 static LWT_ISO_EDGE*
251  int* numelems, int fields, const GBOX *box)
252 {
253  CBT4(topo, getEdgeByFace, ids, numelems, fields, box);
254 }
255 
256 static LWT_ISO_NODE*
258  int* numelems, int fields, const GBOX *box)
259 {
260  CBT4(topo, getNodeByFace, ids, numelems, fields, box);
261 }
262 
265  double dist, int* numelems, int fields,
266  int limit)
267 {
268  CBT5(topo, getEdgeWithinDistance2D, pt, dist, numelems, fields, limit);
269 }
270 
271 int
272 lwt_be_insertEdges(LWT_TOPOLOGY* topo, LWT_ISO_EDGE* edge, int numelems)
273 {
274  CBT2(topo, insertEdges, edge, numelems);
275 }
276 
277 int
279  const LWT_ISO_EDGE* sel_edge, int sel_fields,
280  const LWT_ISO_EDGE* upd_edge, int upd_fields,
281  const LWT_ISO_EDGE* exc_edge, int exc_fields
282 )
283 {
284  CBT6(topo, updateEdges, sel_edge, sel_fields,
285  upd_edge, upd_fields,
286  exc_edge, exc_fields);
287 }
288 
289 static int
291  const LWT_ISO_NODE* sel_node, int sel_fields,
292  const LWT_ISO_NODE* upd_node, int upd_fields,
293  const LWT_ISO_NODE* exc_node, int exc_fields
294 )
295 {
296  CBT6(topo, updateNodes, sel_node, sel_fields,
297  upd_node, upd_fields,
298  exc_node, exc_fields);
299 }
300 
301 static int
303  const LWT_ISO_FACE* faces, int numfaces
304 )
305 {
306  CBT2(topo, updateFacesById, faces, numfaces);
307 }
308 
309 static int
311  const LWT_ISO_EDGE* edges, int numedges, int upd_fields
312 )
313 {
314  CBT3(topo, updateEdgesById, edges, numedges, upd_fields);
315 }
316 
317 static int
319  const LWT_ISO_NODE* nodes, int numnodes, int upd_fields
320 )
321 {
322  CBT3(topo, updateNodesById, nodes, numnodes, upd_fields);
323 }
324 
325 int
327  const LWT_ISO_EDGE* sel_edge, int sel_fields
328 )
329 {
330  CBT2(topo, deleteEdges, sel_edge, sel_fields);
331 }
332 
335 {
336  CBT1(topo, getFaceContainingPoint, pt);
337 }
338 
339 
340 int
342 {
343  CBT3(topo, updateTopoGeomEdgeSplit, split_edge, new_edge1, new_edge2);
344 }
345 
346 static int
348  LWT_ELEMID new_face1, LWT_ELEMID new_face2)
349 {
350  CBT3(topo, updateTopoGeomFaceSplit, split_face, new_face1, new_face2);
351 }
352 
353 static int
355  LWT_ELEMID face_left, LWT_ELEMID face_right)
356 {
357  CBT3(topo, checkTopoGeomRemEdge, edge_id, face_left, face_right);
358 }
359 
360 static int
362  LWT_ELEMID eid1, LWT_ELEMID eid2)
363 {
364  CBT3(topo, checkTopoGeomRemNode, node_id, eid1, eid2);
365 }
366 
367 static int
369  LWT_ELEMID face1, LWT_ELEMID face2,
370  LWT_ELEMID newface)
371 {
372  CBT3(topo, updateTopoGeomFaceHeal, face1, face2, newface);
373 }
374 
375 static int
377  LWT_ELEMID edge1, LWT_ELEMID edge2,
378  LWT_ELEMID newedge)
379 {
380  CBT3(topo, updateTopoGeomEdgeHeal, edge1, edge2, newedge);
381 }
382 
383 static LWT_ELEMID*
385  LWT_ELEMID edge, int *numedges, int limit )
386 {
387  CBT3(topo, getRingEdges, edge, numedges, limit);
388 }
389 
390 
391 /* wrappers of backend wrappers... */
392 
393 int
395 {
396  int exists = 0;
397  lwt_be_getNodeWithinDistance2D(topo, pt, 0, &exists, 0, -1);
398  if ( exists == -1 ) {
399  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
400  return 0;
401  }
402  return exists;
403 }
404 
405 int
407 {
408  int exists = 0;
409  lwt_be_getEdgeWithinDistance2D(topo, pt, 0, &exists, 0, -1);
410  if ( exists == -1 ) {
411  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
412  return 0;
413  }
414  return exists;
415 }
416 
417 /************************************************************************
418  *
419  * Utility functions
420  *
421  ************************************************************************/
422 
423 static LWGEOM *
424 _lwt_toposnap(LWGEOM *src, LWGEOM *tgt, double tol)
425 {
426  LWGEOM *tmp = src;
427  LWGEOM *tmp2;
428  int changed;
429  int iterations = 0;
430 
431  int maxiterations = lwgeom_count_vertices(tgt);
432 
433  /* GEOS snapping can be unstable */
434  /* See https://trac.osgeo.org/geos/ticket/760 */
435  do {
436  LWGEOM *tmp3;
437  tmp2 = lwgeom_snap(tmp, tgt, tol);
438  ++iterations;
439  changed = ( lwgeom_count_vertices(tmp2) != lwgeom_count_vertices(tmp) );
440 #if GEOS_NUMERIC_VERSION < 30309
441  /* Up to GEOS-3.3.8, snapping could duplicate points */
442  if ( changed ) {
443  tmp3 = lwgeom_remove_repeated_points( tmp2, 0 );
444  lwgeom_free(tmp2);
445  tmp2 = tmp3;
446  changed = ( lwgeom_count_vertices(tmp2) != lwgeom_count_vertices(tmp) );
447  }
448 #endif /* GEOS_NUMERIC_VERSION < 30309 */
449  LWDEBUGF(2, "After iteration %d, geometry changed ? %d (%d vs %d vertices)", iterations, changed, lwgeom_count_vertices(tmp2), lwgeom_count_vertices(tmp));
450  if ( tmp != src ) lwgeom_free(tmp);
451  tmp = tmp2;
452  } while ( changed && iterations <= maxiterations );
453 
454  LWDEBUGF(1, "It took %d/%d iterations to properly snap",
455  iterations, maxiterations);
456 
457  return tmp;
458 }
459 
460 static void
461 _lwt_release_faces(LWT_ISO_FACE *faces, int num_faces)
462 {
463  int i;
464  for ( i=0; i<num_faces; ++i ) {
465  if ( faces[i].mbr ) lwfree(faces[i].mbr);
466  }
467  lwfree(faces);
468 }
469 
470 static void
471 _lwt_release_edges(LWT_ISO_EDGE *edges, int num_edges)
472 {
473  int i;
474  for ( i=0; i<num_edges; ++i ) {
475  if ( edges[i].geom ) lwline_free(edges[i].geom);
476  }
477  lwfree(edges);
478 }
479 
480 static void
481 _lwt_release_nodes(LWT_ISO_NODE *nodes, int num_nodes)
482 {
483  int i;
484  for ( i=0; i<num_nodes; ++i ) {
485  if ( nodes[i].geom ) lwpoint_free(nodes[i].geom);
486  }
487  lwfree(nodes);
488 }
489 
490 /************************************************************************
491  *
492  * API implementation
493  *
494  ************************************************************************/
495 
496 LWT_TOPOLOGY *
497 lwt_LoadTopology( LWT_BE_IFACE *iface, const char *name )
498 {
499  LWT_BE_TOPOLOGY* be_topo;
500  LWT_TOPOLOGY* topo;
501 
502  be_topo = lwt_be_loadTopologyByName(iface, name);
503  if ( ! be_topo ) {
504  //lwerror("Could not load topology from backend: %s",
505  lwerror("%s", lwt_be_lastErrorMessage(iface));
506  return NULL;
507  }
508  topo = lwalloc(sizeof(LWT_TOPOLOGY));
509  topo->be_iface = iface;
510  topo->be_topo = be_topo;
511  topo->srid = lwt_be_topoGetSRID(topo);
512  topo->hasZ = lwt_be_topoHasZ(topo);
513  topo->precision = lwt_be_topoGetPrecision(topo);
514 
515  return topo;
516 }
517 
518 void
520 {
521  if ( ! lwt_be_freeTopology(topo) ) {
522  lwnotice("Could not release backend topology memory: %s",
524  }
525  lwfree(topo);
526 }
527 
530  LWPOINT* pt, int skipISOChecks )
531 {
532  LWT_ELEMID foundInFace = -1;
533 
534  if ( ! skipISOChecks )
535  {
536  if ( lwt_be_ExistsCoincidentNode(topo, pt) ) /*x*/
537  {
538  lwerror("SQL/MM Spatial exception - coincident node");
539  return -1;
540  }
541  if ( lwt_be_ExistsEdgeIntersectingPoint(topo, pt) ) /*x*/
542  {
543  lwerror("SQL/MM Spatial exception - edge crosses node.");
544  return -1;
545  }
546  }
547 
548  if ( face == -1 || ! skipISOChecks )
549  {
550  foundInFace = lwt_be_getFaceContainingPoint(topo, pt); /*x*/
551  if ( foundInFace == -2 ) {
552  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
553  return -1;
554  }
555  if ( foundInFace == -1 ) foundInFace = 0;
556  }
557 
558  if ( face == -1 ) {
559  face = foundInFace;
560  }
561  else if ( ! skipISOChecks && foundInFace != face ) {
562 #if 0
563  lwerror("SQL/MM Spatial exception - within face %d (not %d)",
564  foundInFace, face);
565 #else
566  lwerror("SQL/MM Spatial exception - not within face");
567 #endif
568  return -1;
569  }
570 
571  LWT_ISO_NODE node;
572  node.node_id = -1;
573  node.containing_face = face;
574  node.geom = pt;
575  if ( ! lwt_be_insertNodes(topo, &node, 1) )
576  {
577  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
578  return -1;
579  }
580 
581  return node.node_id;
582 }
583 
584 /* Check that an edge does not cross an existing node or edge
585  *
586  * @param myself the id of an edge to skip, if any
587  * (for ChangeEdgeGeom). Can use 0 for none.
588  *
589  * Return -1 on cross or error, 0 if everything is fine.
590  * Note that before returning -1, lwerror is invoked...
591  */
592 static int
594  LWT_ELEMID start_node, LWT_ELEMID end_node,
595  const LWLINE *geom, LWT_ELEMID myself )
596 {
597  int i, num_nodes, num_edges;
598  LWT_ISO_EDGE *edges;
599  LWT_ISO_NODE *nodes;
600  const GBOX *edgebox;
601  GEOSGeometry *edgegg;
602  const GEOSPreparedGeometry* prepared_edge;
603 
604  initGEOS(lwnotice, lwgeom_geos_error);
605 
606  edgegg = LWGEOM2GEOS( lwline_as_lwgeom(geom), 0);
607  if ( ! edgegg ) {
608  lwerror("Could not convert edge geometry to GEOS: %s", lwgeom_geos_errmsg);
609  return -1;
610  }
611  prepared_edge = GEOSPrepare( edgegg );
612  if ( ! prepared_edge ) {
613  lwerror("Could not prepare edge geometry: %s", lwgeom_geos_errmsg);
614  return -1;
615  }
616  edgebox = lwgeom_get_bbox( lwline_as_lwgeom(geom) );
617 
618  /* loop over each node within the edge's gbox */
619  nodes = lwt_be_getNodeWithinBox2D( topo, edgebox, &num_nodes,
620  LWT_COL_NODE_ALL, 0 );
621  LWDEBUGF(1, "lwt_be_getNodeWithinBox2D returned %d nodes", num_nodes);
622  if ( num_nodes == -1 ) {
623  GEOSPreparedGeom_destroy(prepared_edge);
624  GEOSGeom_destroy(edgegg);
625  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
626  return -1;
627  }
628  for ( i=0; i<num_nodes; ++i )
629  {
630  LWT_ISO_NODE* node = &(nodes[i]);
631  GEOSGeometry *nodegg;
632  int contains;
633  if ( node->node_id == start_node ) continue;
634  if ( node->node_id == end_node ) continue;
635  /* check if the edge contains this node (not on boundary) */
636  nodegg = LWGEOM2GEOS( lwpoint_as_lwgeom(node->geom) , 0);
637  /* ST_RelateMatch(rec.relate, 'T********') */
638  contains = GEOSPreparedContains( prepared_edge, nodegg );
639  GEOSGeom_destroy(nodegg);
640  if (contains == 2)
641  {
642  GEOSPreparedGeom_destroy(prepared_edge);
643  GEOSGeom_destroy(edgegg);
644  _lwt_release_nodes(nodes, num_nodes);
645  lwerror("GEOS exception on PreparedContains: %s", lwgeom_geos_errmsg);
646  return -1;
647  }
648  if ( contains )
649  {
650  GEOSPreparedGeom_destroy(prepared_edge);
651  GEOSGeom_destroy(edgegg);
652  _lwt_release_nodes(nodes, num_nodes);
653  lwerror("SQL/MM Spatial exception - geometry crosses a node");
654  return -1;
655  }
656  }
657  if ( nodes ) _lwt_release_nodes(nodes, num_nodes);
658  /* may be NULL if num_nodes == 0 */
659 
660  /* loop over each edge within the edge's gbox */
661  edges = lwt_be_getEdgeWithinBox2D( topo, edgebox, &num_edges, LWT_COL_EDGE_ALL, 0 );
662  LWDEBUGF(1, "lwt_be_getEdgeWithinBox2D returned %d edges", num_edges);
663  if ( num_edges == -1 ) {
664  GEOSPreparedGeom_destroy(prepared_edge);
665  GEOSGeom_destroy(edgegg);
666  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
667  return -1;
668  }
669  for ( i=0; i<num_edges; ++i )
670  {
671  LWT_ISO_EDGE* edge = &(edges[i]);
672  LWT_ELEMID edge_id = edge->edge_id;
673  GEOSGeometry *eegg;
674  char *relate;
675  int match;
676 
677  if ( edge_id == myself ) continue;
678 
679  if ( ! edge->geom ) {
680  _lwt_release_edges(edges, num_edges);
681  lwerror("Edge %d has NULL geometry!", edge_id);
682  return -1;
683  }
684 
685  eegg = LWGEOM2GEOS( lwline_as_lwgeom(edge->geom), 0 );
686  if ( ! eegg ) {
687  GEOSPreparedGeom_destroy(prepared_edge);
688  GEOSGeom_destroy(edgegg);
689  _lwt_release_edges(edges, num_edges);
690  lwerror("Could not convert edge geometry to GEOS: %s", lwgeom_geos_errmsg);
691  return -1;
692  }
693 
694  LWDEBUGF(2, "Edge %d converted to GEOS", edge_id);
695 
696  /* check if the edge crosses our edge (not boundary-boundary) */
697 
698  relate = GEOSRelateBoundaryNodeRule(eegg, edgegg, 2);
699  if ( ! relate ) {
700  GEOSGeom_destroy(eegg);
701  GEOSPreparedGeom_destroy(prepared_edge);
702  GEOSGeom_destroy(edgegg);
703  _lwt_release_edges(edges, num_edges);
704  lwerror("GEOSRelateBoundaryNodeRule error: %s", lwgeom_geos_errmsg);
705  return -1;
706  }
707 
708  LWDEBUGF(2, "Edge %d relate pattern is %s", edge_id, relate);
709 
710  match = GEOSRelatePatternMatch(relate, "F********");
711  if ( match ) {
712  /* error or no interior intersection */
713  GEOSGeom_destroy(eegg);
714  GEOSFree(relate);
715  if ( match == 2 ) {
716  _lwt_release_edges(edges, num_edges);
717  GEOSPreparedGeom_destroy(prepared_edge);
718  GEOSGeom_destroy(edgegg);
719  lwerror("GEOSRelatePatternMatch error: %s", lwgeom_geos_errmsg);
720  return -1;
721  }
722  else continue; /* no interior intersection */
723  }
724 
725  match = GEOSRelatePatternMatch(relate, "1FFF*FFF2");
726  if ( match ) {
727  _lwt_release_edges(edges, num_edges);
728  GEOSPreparedGeom_destroy(prepared_edge);
729  GEOSGeom_destroy(edgegg);
730  GEOSGeom_destroy(eegg);
731  GEOSFree(relate);
732  if ( match == 2 ) {
733  lwerror("GEOSRelatePatternMatch error: %s", lwgeom_geos_errmsg);
734  } else {
735  lwerror("SQL/MM Spatial exception - coincident edge %" LWTFMT_ELEMID,
736  edge_id);
737  }
738  return -1;
739  }
740 
741  match = GEOSRelatePatternMatch(relate, "1********");
742  if ( match ) {
743  _lwt_release_edges(edges, num_edges);
744  GEOSPreparedGeom_destroy(prepared_edge);
745  GEOSGeom_destroy(edgegg);
746  GEOSGeom_destroy(eegg);
747  GEOSFree(relate);
748  if ( match == 2 ) {
749  lwerror("GEOSRelatePatternMatch error: %s", lwgeom_geos_errmsg);
750  } else {
751  lwerror("Spatial exception - geometry intersects edge %"
752  LWTFMT_ELEMID, edge_id);
753  }
754  return -1;
755  }
756 
757  match = GEOSRelatePatternMatch(relate, "T********");
758  if ( match ) {
759  _lwt_release_edges(edges, num_edges);
760  GEOSPreparedGeom_destroy(prepared_edge);
761  GEOSGeom_destroy(edgegg);
762  GEOSGeom_destroy(eegg);
763  GEOSFree(relate);
764  if ( match == 2 ) {
765  lwerror("GEOSRelatePatternMatch error: %s", lwgeom_geos_errmsg);
766  } else {
767  lwerror("SQL/MM Spatial exception - geometry crosses edge %"
768  LWTFMT_ELEMID, edge_id);
769  }
770  return -1;
771  }
772 
773  LWDEBUGF(2, "Edge %d analisys completed, it does no harm", edge_id);
774 
775  GEOSFree(relate);
776  GEOSGeom_destroy(eegg);
777  }
778  if ( edges ) _lwt_release_edges(edges, num_edges);
779  /* would be NULL if num_edges was 0 */
780 
781  GEOSPreparedGeom_destroy(prepared_edge);
782  GEOSGeom_destroy(edgegg);
783 
784  return 0;
785 }
786 
787 
790  LWT_ELEMID endNode, const LWLINE* geom )
791 {
792  int num_nodes;
793  int i;
794  LWT_ISO_EDGE newedge;
795  LWT_ISO_NODE *endpoints;
796  LWT_ELEMID containing_face = -1;
797  LWT_ELEMID node_ids[2];
798  LWT_ISO_NODE updated_nodes[2];
799  int skipISOChecks = 0;
800  POINT2D p1, p2;
801 
802  /* NOT IN THE SPECS:
803  * A closed edge is never isolated (as it forms a face)
804  */
805  if ( startNode == endNode )
806  {
807  lwerror("Closed edges would not be isolated, try lwt_AddEdgeNewFaces");
808  return -1;
809  }
810 
811  if ( ! skipISOChecks )
812  {
813  /* Acurve must be simple */
814  if ( ! lwgeom_is_simple(lwline_as_lwgeom(geom)) )
815  {
816  lwerror("SQL/MM Spatial exception - curve not simple");
817  return -1;
818  }
819  }
820 
821  /*
822  * Check for:
823  * existence of nodes
824  * nodes faces match
825  * Extract:
826  * nodes face id
827  * nodes geoms
828  */
829  num_nodes = 2;
830  node_ids[0] = startNode;
831  node_ids[1] = endNode;
832  endpoints = lwt_be_getNodeById( topo, node_ids, &num_nodes,
834  if ( num_nodes < 0 )
835  {
836  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
837  return -1;
838  }
839  else if ( num_nodes < 2 )
840  {
841  if ( num_nodes ) _lwt_release_nodes(endpoints, num_nodes);
842  lwerror("SQL/MM Spatial exception - non-existent node");
843  return -1;
844  }
845  for ( i=0; i<num_nodes; ++i )
846  {
847  const LWT_ISO_NODE *n = &(endpoints[i]);
848  if ( n->containing_face == -1 )
849  {
850  _lwt_release_nodes(endpoints, num_nodes);
851  lwerror("SQL/MM Spatial exception - not isolated node");
852  return -1;
853  }
854  if ( containing_face == -1 ) containing_face = n->containing_face;
855  else if ( containing_face != n->containing_face )
856  {
857  _lwt_release_nodes(endpoints, num_nodes);
858  lwerror("SQL/MM Spatial exception - nodes in different faces");
859  return -1;
860  }
861 
862  if ( ! skipISOChecks )
863  {
864  if ( n->node_id == startNode )
865  {
866  /* l) Check that start point of acurve match start node geoms. */
867  getPoint2d_p(geom->points, 0, &p1);
868  getPoint2d_p(n->geom->point, 0, &p2);
869  if ( ! p2d_same(&p1, &p2) )
870  {
871  _lwt_release_nodes(endpoints, num_nodes);
872  lwerror("SQL/MM Spatial exception - "
873  "start node not geometry start point.");
874  return -1;
875  }
876  }
877  else
878  {
879  /* m) Check that end point of acurve match end node geoms. */
880  getPoint2d_p(geom->points, geom->points->npoints-1, &p1);
881  getPoint2d_p(n->geom->point, 0, &p2);
882  if ( ! p2d_same(&p1, &p2) )
883  {
884  _lwt_release_nodes(endpoints, num_nodes);
885  lwerror("SQL/MM Spatial exception - "
886  "end node not geometry end point.");
887  return -1;
888  }
889  }
890  }
891  }
892 
893  if ( num_nodes ) _lwt_release_nodes(endpoints, num_nodes);
894 
895  if ( ! skipISOChecks )
896  {
897  if ( _lwt_CheckEdgeCrossing( topo, startNode, endNode, geom, 0 ) )
898  {
899  /* would have called lwerror already, leaking :( */
900  return -1;
901  }
902  }
903 
904  /*
905  * All checks passed, time to prepare the new edge
906  */
907 
908  newedge.edge_id = lwt_be_getNextEdgeId( topo );
909  if ( newedge.edge_id == -1 ) {
910  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
911  return -1;
912  }
913 
914  /* TODO: this should likely be an exception instead ! */
915  if ( containing_face == -1 ) containing_face = 0;
916 
917  newedge.start_node = startNode;
918  newedge.end_node = endNode;
919  newedge.face_left = newedge.face_right = containing_face;
920  newedge.next_left = -newedge.edge_id;
921  newedge.next_right = newedge.edge_id;
922  newedge.geom = (LWLINE *)geom; /* const cast.. */
923 
924  int ret = lwt_be_insertEdges(topo, &newedge, 1);
925  if ( ret == -1 ) {
926  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
927  return -1;
928  } else if ( ret == 0 ) {
929  lwerror("Insertion of split edge failed (no reason)");
930  return -1;
931  }
932 
933  /*
934  * Update Node containing_face values
935  *
936  * the nodes anode and anothernode are no more isolated
937  * because now there is an edge connecting them
938  */
939  updated_nodes[0].node_id = startNode;
940  updated_nodes[0].containing_face = -1;
941  updated_nodes[1].node_id = endNode;
942  updated_nodes[1].containing_face = -1;
943  ret = lwt_be_updateNodesById(topo, updated_nodes, 2,
945  if ( ret == -1 ) {
946  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
947  return -1;
948  }
949 
950  return newedge.edge_id;
951 }
952 
953 static LWCOLLECTION *
954 _lwt_EdgeSplit( LWT_TOPOLOGY* topo, LWT_ELEMID edge, LWPOINT* pt, int skipISOChecks, LWT_ISO_EDGE** oldedge )
955 {
956  LWGEOM *split;
957  LWCOLLECTION *split_col;
958  int i;
959 
960  /* Get edge */
961  i = 1;
962  LWDEBUG(1, "calling lwt_be_getEdgeById");
963  *oldedge = lwt_be_getEdgeById(topo, &edge, &i, LWT_COL_EDGE_ALL);
964  LWDEBUGF(1, "lwt_be_getEdgeById returned %p", *oldedge);
965  if ( ! *oldedge )
966  {
967  LWDEBUGF(1, "lwt_be_getEdgeById returned NULL and set i=%d", i);
968  if ( i == -1 )
969  {
970  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
971  return NULL;
972  }
973  else if ( i == 0 )
974  {
975  lwerror("SQL/MM Spatial exception - non-existent edge");
976  return NULL;
977  }
978  else
979  {
980  lwerror("Backend coding error: getEdgeById callback returned NULL "
981  "but numelements output parameter has value %d "
982  "(expected 0 or 1)", i);
983  return NULL;
984  }
985  }
986 
987 
988  /*
989  * - check if a coincident node already exists
990  */
991  if ( ! skipISOChecks )
992  {
993  LWDEBUG(1, "calling lwt_be_ExistsCoincidentNode");
994  if ( lwt_be_ExistsCoincidentNode(topo, pt) ) /*x*/
995  {
996  LWDEBUG(1, "lwt_be_ExistsCoincidentNode returned");
997  _lwt_release_edges(*oldedge, 1);
998  lwerror("SQL/MM Spatial exception - coincident node");
999  return NULL;
1000  }
1001  LWDEBUG(1, "lwt_be_ExistsCoincidentNode returned");
1002  }
1003 
1004  /* Split edge */
1005  split = lwgeom_split((LWGEOM*)(*oldedge)->geom, (LWGEOM*)pt);
1006  if ( ! split )
1007  {
1008  _lwt_release_edges(*oldedge, 1);
1009  lwerror("could not split edge by point ?");
1010  return NULL;
1011  }
1012  split_col = lwgeom_as_lwcollection(split);
1013  if ( ! split_col ) {
1014  _lwt_release_edges(*oldedge, 1);
1015  lwgeom_free(split);
1016  lwerror("lwgeom_as_lwcollection returned NULL");
1017  return NULL;
1018  }
1019  if (split_col->ngeoms < 2) {
1020  _lwt_release_edges(*oldedge, 1);
1021  lwgeom_free(split);
1022  lwerror("SQL/MM Spatial exception - point not on edge");
1023  return NULL;
1024  }
1025 
1026 #if 0
1027  {
1028  size_t sz;
1029  char *wkt = lwgeom_to_wkt((LWGEOM*)split_col, WKT_EXTENDED, 2, &sz);
1030  LWDEBUGF(1, "returning split col: %s", wkt);
1031  lwfree(wkt);
1032  }
1033 #endif
1034  return split_col;
1035 }
1036 
1037 LWT_ELEMID
1039  LWPOINT* pt, int skipISOChecks )
1040 {
1041  LWT_ISO_NODE node;
1042  LWT_ISO_EDGE* oldedge = NULL;
1043  LWCOLLECTION *split_col;
1044  const LWGEOM *oldedge_geom;
1045  const LWGEOM *newedge_geom;
1046  LWT_ISO_EDGE newedge1;
1047  LWT_ISO_EDGE seledge, updedge, excedge;
1048  int ret;
1049 
1050  split_col = _lwt_EdgeSplit( topo, edge, pt, skipISOChecks, &oldedge );
1051  if ( ! split_col ) return -1; /* should have raised an exception */
1052  oldedge_geom = split_col->geoms[0];
1053  newedge_geom = split_col->geoms[1];
1054  /* Make sure the SRID is set on the subgeom */
1055  ((LWGEOM*)oldedge_geom)->srid = split_col->srid;
1056  ((LWGEOM*)newedge_geom)->srid = split_col->srid;
1057 
1058  /* Add new node, getting new id back */
1059  node.node_id = -1;
1060  node.containing_face = -1; /* means not-isolated */
1061  node.geom = pt;
1062  if ( ! lwt_be_insertNodes(topo, &node, 1) )
1063  {
1064  _lwt_release_edges(oldedge, 1);
1065  lwcollection_free(split_col);
1066  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1067  return -1;
1068  }
1069  if (node.node_id == -1) {
1070  /* should have been set by backend */
1071  _lwt_release_edges(oldedge, 1);
1072  lwcollection_free(split_col);
1073  lwerror("Backend coding error: "
1074  "insertNodes callback did not return node_id");
1075  return -1;
1076  }
1077 
1078  /* Insert the new edge */
1079  newedge1.edge_id = lwt_be_getNextEdgeId(topo);
1080  if ( newedge1.edge_id == -1 ) {
1081  _lwt_release_edges(oldedge, 1);
1082  lwcollection_free(split_col);
1083  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1084  return -1;
1085  }
1086  newedge1.start_node = node.node_id;
1087  newedge1.end_node = oldedge->end_node;
1088  newedge1.face_left = oldedge->face_left;
1089  newedge1.face_right = oldedge->face_right;
1090  newedge1.next_left = oldedge->next_left == -oldedge->edge_id ?
1091  -newedge1.edge_id : oldedge->next_left;
1092  newedge1.next_right = -oldedge->edge_id;
1093  newedge1.geom = lwgeom_as_lwline(newedge_geom);
1094  /* lwgeom_split of a line should only return lines ... */
1095  if ( ! newedge1.geom ) {
1096  _lwt_release_edges(oldedge, 1);
1097  lwcollection_free(split_col);
1098  lwerror("first geometry in lwgeom_split output is not a line");
1099  return -1;
1100  }
1101  ret = lwt_be_insertEdges(topo, &newedge1, 1);
1102  if ( ret == -1 ) {
1103  _lwt_release_edges(oldedge, 1);
1104  lwcollection_free(split_col);
1105  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1106  return -1;
1107  } else if ( ret == 0 ) {
1108  _lwt_release_edges(oldedge, 1);
1109  lwcollection_free(split_col);
1110  lwerror("Insertion of split edge failed (no reason)");
1111  return -1;
1112  }
1113 
1114  /* Update the old edge */
1115  updedge.geom = lwgeom_as_lwline(oldedge_geom);
1116  /* lwgeom_split of a line should only return lines ... */
1117  if ( ! updedge.geom ) {
1118  _lwt_release_edges(oldedge, 1);
1119  lwcollection_free(split_col);
1120  lwerror("second geometry in lwgeom_split output is not a line");
1121  return -1;
1122  }
1123  updedge.next_left = newedge1.edge_id;
1124  updedge.end_node = node.node_id;
1125  ret = lwt_be_updateEdges(topo,
1126  oldedge, LWT_COL_EDGE_EDGE_ID,
1128  NULL, 0);
1129  if ( ret == -1 ) {
1130  _lwt_release_edges(oldedge, 1);
1131  lwcollection_free(split_col);
1132  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1133  return -1;
1134  } else if ( ret == 0 ) {
1135  _lwt_release_edges(oldedge, 1);
1136  lwcollection_free(split_col);
1137  lwerror("Edge being split (%d) disappeared during operations?", oldedge->edge_id);
1138  return -1;
1139  } else if ( ret > 1 ) {
1140  _lwt_release_edges(oldedge, 1);
1141  lwcollection_free(split_col);
1142  lwerror("More than a single edge found with id %d !", oldedge->edge_id);
1143  return -1;
1144  }
1145 
1146  /* Update all next edge references to match new layout (ST_ModEdgeSplit) */
1147 
1148  updedge.next_right = -newedge1.edge_id;
1149  excedge.edge_id = newedge1.edge_id;
1150  seledge.next_right = -oldedge->edge_id;
1151  seledge.start_node = oldedge->end_node;
1152  ret = lwt_be_updateEdges(topo,
1154  &updedge, LWT_COL_EDGE_NEXT_RIGHT,
1155  &excedge, LWT_COL_EDGE_EDGE_ID);
1156  if ( ret == -1 ) {
1157  _lwt_release_edges(oldedge, 1);
1158  lwcollection_free(split_col);
1159  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1160  return -1;
1161  }
1162 
1163  updedge.next_left = -newedge1.edge_id;
1164  excedge.edge_id = newedge1.edge_id;
1165  seledge.next_left = -oldedge->edge_id;
1166  seledge.end_node = oldedge->end_node;
1167  ret = lwt_be_updateEdges(topo,
1169  &updedge, LWT_COL_EDGE_NEXT_LEFT,
1170  &excedge, LWT_COL_EDGE_EDGE_ID);
1171  if ( ret == -1 ) {
1172  _lwt_release_edges(oldedge, 1);
1173  lwcollection_free(split_col);
1174  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1175  return -1;
1176  }
1177 
1178  /* Update TopoGeometries composition */
1179  ret = lwt_be_updateTopoGeomEdgeSplit(topo, oldedge->edge_id, newedge1.edge_id, -1);
1180  if ( ! ret ) {
1181  _lwt_release_edges(oldedge, 1);
1182  lwcollection_free(split_col);
1183  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1184  return -1;
1185  }
1186 
1187  _lwt_release_edges(oldedge, 1);
1188  lwcollection_free(split_col);
1189 
1190  /* return new node id */
1191  return node.node_id;
1192 }
1193 
1194 LWT_ELEMID
1196  LWPOINT* pt, int skipISOChecks )
1197 {
1198  LWT_ISO_NODE node;
1199  LWT_ISO_EDGE* oldedge = NULL;
1200  LWCOLLECTION *split_col;
1201  const LWGEOM *oldedge_geom;
1202  const LWGEOM *newedge_geom;
1203  LWT_ISO_EDGE newedges[2];
1204  LWT_ISO_EDGE seledge, updedge;
1205  int ret;
1206 
1207  split_col = _lwt_EdgeSplit( topo, edge, pt, skipISOChecks, &oldedge );
1208  if ( ! split_col ) return -1; /* should have raised an exception */
1209  oldedge_geom = split_col->geoms[0];
1210  newedge_geom = split_col->geoms[1];
1211  /* Make sure the SRID is set on the subgeom */
1212  ((LWGEOM*)oldedge_geom)->srid = split_col->srid;
1213  ((LWGEOM*)newedge_geom)->srid = split_col->srid;
1214 
1215  /* Add new node, getting new id back */
1216  node.node_id = -1;
1217  node.containing_face = -1; /* means not-isolated */
1218  node.geom = pt;
1219  if ( ! lwt_be_insertNodes(topo, &node, 1) )
1220  {
1221  _lwt_release_edges(oldedge, 1);
1222  lwcollection_free(split_col);
1223  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1224  return -1;
1225  }
1226  if (node.node_id == -1) {
1227  _lwt_release_edges(oldedge, 1);
1228  lwcollection_free(split_col);
1229  /* should have been set by backend */
1230  lwerror("Backend coding error: "
1231  "insertNodes callback did not return node_id");
1232  return -1;
1233  }
1234 
1235  /* Delete the old edge */
1236  seledge.edge_id = edge;
1237  ret = lwt_be_deleteEdges(topo, &seledge, LWT_COL_EDGE_EDGE_ID);
1238  if ( ret == -1 ) {
1239  _lwt_release_edges(oldedge, 1);
1240  lwcollection_free(split_col);
1241  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1242  return -1;
1243  }
1244 
1245  /* Get new edges identifiers */
1246  newedges[0].edge_id = lwt_be_getNextEdgeId(topo);
1247  if ( newedges[0].edge_id == -1 ) {
1248  _lwt_release_edges(oldedge, 1);
1249  lwcollection_free(split_col);
1250  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1251  return -1;
1252  }
1253  newedges[1].edge_id = lwt_be_getNextEdgeId(topo);
1254  if ( newedges[1].edge_id == -1 ) {
1255  _lwt_release_edges(oldedge, 1);
1256  lwcollection_free(split_col);
1257  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1258  return -1;
1259  }
1260 
1261  /* Define the first new edge (to new node) */
1262  newedges[0].start_node = oldedge->start_node;
1263  newedges[0].end_node = node.node_id;
1264  newedges[0].face_left = oldedge->face_left;
1265  newedges[0].face_right = oldedge->face_right;
1266  newedges[0].next_left = newedges[1].edge_id;
1267  if ( oldedge->next_right == edge )
1268  newedges[0].next_right = newedges[0].edge_id;
1269  else if ( oldedge->next_right == -edge )
1270  newedges[0].next_right = -newedges[1].edge_id;
1271  else
1272  newedges[0].next_right = oldedge->next_right;
1273  newedges[0].geom = lwgeom_as_lwline(oldedge_geom);
1274  /* lwgeom_split of a line should only return lines ... */
1275  if ( ! newedges[0].geom ) {
1276  _lwt_release_edges(oldedge, 1);
1277  lwcollection_free(split_col);
1278  lwerror("first geometry in lwgeom_split output is not a line");
1279  return -1;
1280  }
1281 
1282  /* Define the second new edge (from new node) */
1283  newedges[1].start_node = node.node_id;
1284  newedges[1].end_node = oldedge->end_node;
1285  newedges[1].face_left = oldedge->face_left;
1286  newedges[1].face_right = oldedge->face_right;
1287  newedges[1].next_right = -newedges[0].edge_id;
1288  if ( oldedge->next_left == -edge )
1289  newedges[1].next_left = -newedges[1].edge_id;
1290  else if ( oldedge->next_left == edge )
1291  newedges[1].next_left = newedges[0].edge_id;
1292  else
1293  newedges[1].next_left = oldedge->next_left;
1294  newedges[1].geom = lwgeom_as_lwline(newedge_geom);
1295  /* lwgeom_split of a line should only return lines ... */
1296  if ( ! newedges[1].geom ) {
1297  _lwt_release_edges(oldedge, 1);
1298  lwcollection_free(split_col);
1299  lwerror("second geometry in lwgeom_split output is not a line");
1300  return -1;
1301  }
1302 
1303  /* Insert both new edges */
1304  ret = lwt_be_insertEdges(topo, newedges, 2);
1305  if ( ret == -1 ) {
1306  _lwt_release_edges(oldedge, 1);
1307  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1308  return -1;
1309  } else if ( ret == 0 ) {
1310  _lwt_release_edges(oldedge, 1);
1311  lwcollection_free(split_col);
1312  lwerror("Insertion of split edge failed (no reason)");
1313  return -1;
1314  }
1315 
1316  /* Update all next edge references pointing to old edge id */
1317 
1318  updedge.next_right = newedges[1].edge_id;
1319  seledge.next_right = edge;
1320  seledge.start_node = oldedge->start_node;
1321  ret = lwt_be_updateEdges(topo,
1323  &updedge, LWT_COL_EDGE_NEXT_RIGHT,
1324  NULL, 0);
1325  if ( ret == -1 ) {
1326  _lwt_release_edges(oldedge, 1);
1327  lwcollection_free(split_col);
1328  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1329  return -1;
1330  }
1331 
1332  updedge.next_right = -newedges[0].edge_id;
1333  seledge.next_right = -edge;
1334  seledge.start_node = oldedge->end_node;
1335  ret = lwt_be_updateEdges(topo,
1337  &updedge, LWT_COL_EDGE_NEXT_RIGHT,
1338  NULL, 0);
1339  if ( ret == -1 ) {
1340  _lwt_release_edges(oldedge, 1);
1341  lwcollection_free(split_col);
1342  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1343  return -1;
1344  }
1345 
1346  updedge.next_left = newedges[0].edge_id;
1347  seledge.next_left = edge;
1348  seledge.end_node = oldedge->start_node;
1349  ret = lwt_be_updateEdges(topo,
1351  &updedge, LWT_COL_EDGE_NEXT_LEFT,
1352  NULL, 0);
1353  if ( ret == -1 ) {
1354  _lwt_release_edges(oldedge, 1);
1355  lwcollection_free(split_col);
1356  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1357  return -1;
1358  }
1359 
1360  updedge.next_left = -newedges[1].edge_id;
1361  seledge.next_left = -edge;
1362  seledge.end_node = oldedge->end_node;
1363  ret = lwt_be_updateEdges(topo,
1365  &updedge, LWT_COL_EDGE_NEXT_LEFT,
1366  NULL, 0);
1367  if ( ret == -1 ) {
1368  _lwt_release_edges(oldedge, 1);
1369  lwcollection_release(split_col);
1370  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1371  return -1;
1372  }
1373 
1374  /* Update TopoGeometries composition */
1375  ret = lwt_be_updateTopoGeomEdgeSplit(topo, oldedge->edge_id, newedges[0].edge_id, newedges[1].edge_id);
1376  if ( ! ret ) {
1377  _lwt_release_edges(oldedge, 1);
1378  lwcollection_free(split_col);
1379  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1380  return -1;
1381  }
1382 
1383  _lwt_release_edges(oldedge, 1);
1384  lwcollection_free(split_col);
1385 
1386  /* return new node id */
1387  return node.node_id;
1388 }
1389 
1390 /* Data structure used by AddEdgeX functions */
1391 typedef struct edgeend_t {
1392  /* Signed identifier of next clockwise edge (+outgoing,-incoming) */
1394  /* Identifier of face between myaz and next CW edge */
1396  /* Signed identifier of next counterclockwise edge (+outgoing,-incoming) */
1398  /* Identifier of face between myaz and next CCW edge */
1401  double myaz; /* azimuth of edgeend geometry */
1402 } edgeend;
1403 
1404 /*
1405  * Get first distinct vertex from endpoint
1406  * @param pa the pointarray to seek points in
1407  * @param ref the point we want to search a distinct one
1408  * @param from vertex index to start from
1409  * @param dir 1 to go forward
1410  * -1 to go backward
1411  * @return 0 if edge is collapsed (no distinct points)
1412  */
1413 static int
1414 _lwt_FirstDistinctVertex2D(const POINTARRAY* pa, POINT2D *ref, int from, int dir, POINT2D *op)
1415 {
1416  int i, toofar, inc;
1417  POINT2D fp;
1418 
1419  if ( dir > 0 )
1420  {
1421  toofar = pa->npoints;
1422  inc = 1;
1423  }
1424  else
1425  {
1426  toofar = -1;
1427  inc = -1;
1428  }
1429 
1430  LWDEBUGF(1, "first point is index %d", from);
1431  fp = *ref; /* getPoint2d_p(pa, from, &fp); */
1432  for ( i = from+inc; i != toofar; i += inc )
1433  {
1434  LWDEBUGF(1, "testing point %d", i);
1435  getPoint2d_p(pa, i, op); /* pick next point */
1436  if ( p2d_same(op, &fp) ) continue; /* equal to startpoint */
1437  /* this is a good one, neither same of start nor of end point */
1438  return 1; /* found */
1439  }
1440 
1441  /* no distinct vertices found */
1442  return 0;
1443 }
1444 
1445 
1446 /*
1447  * Return non-zero on failure (lwerror is invoked in that case)
1448  * Possible failures:
1449  * -1 no two distinct vertices exist
1450  * -2 azimuth computation failed for first edge end
1451  */
1452 static int
1454  POINT2D *fp, POINT2D *lp)
1455 {
1456  POINTARRAY *pa = edge->points;
1457  POINT2D pt;
1458 
1459  fee->nextCW = fee->nextCCW =
1460  lee->nextCW = lee->nextCCW = 0;
1461  fee->cwFace = fee->ccwFace =
1462  lee->cwFace = lee->ccwFace = -1;
1463 
1464  /* Compute azimuth of first edge end */
1465  LWDEBUG(1, "computing azimuth of first edge end");
1466  if ( ! _lwt_FirstDistinctVertex2D(pa, fp, 0, 1, &pt) )
1467  {
1468  lwerror("Invalid edge (no two distinct vertices exist)");
1469  return -1;
1470  }
1471  if ( ! azimuth_pt_pt(fp, &pt, &(fee->myaz)) ) {
1472  lwerror("error computing azimuth of first edgeend [%g %g,%g %g]",
1473  fp->x, fp->y, pt.x, pt.y);
1474  return -2;
1475  }
1476  LWDEBUGF(1, "azimuth of first edge end [%g %g,%g %g] is %g",
1477  fp->x, fp->y, pt.x, pt.y, fee->myaz);
1478 
1479  /* Compute azimuth of second edge end */
1480  LWDEBUG(1, "computing azimuth of second edge end");
1481  if ( ! _lwt_FirstDistinctVertex2D(pa, lp, pa->npoints-1, -1, &pt) )
1482  {
1483  lwerror("Invalid edge (no two distinct vertices exist)");
1484  return -1;
1485  }
1486  if ( ! azimuth_pt_pt(lp, &pt, &(lee->myaz)) ) {
1487  lwerror("error computing azimuth of last edgeend [%g %g,%g %g]",
1488  lp->x, lp->y, pt.x, pt.y);
1489  return -2;
1490  }
1491  LWDEBUGF(1, "azimuth of last edge end [%g %g,%g %g] is %g",
1492  lp->x, lp->y, pt.x, pt.y, lee->myaz);
1493 
1494  return 0;
1495 }
1496 
1497 /*
1498  * Find the first edges encountered going clockwise and counterclockwise
1499  * around a node, starting from the given azimuth, and take
1500  * note of the face on the both sides.
1501  *
1502  * @param topo the topology to act upon
1503  * @param node the identifier of the node to analyze
1504  * @param data input (myaz) / output (nextCW, nextCCW) parameter
1505  * @param other edgeend, if also incident to given node (closed edge).
1506  * @param myedge_id identifier of the edge id that data->myaz belongs to
1507  * @return number of incident edges found
1508  *
1509  */
1510 static int
1512  edgeend *other, int myedge_id )
1513 {
1514  LWT_ISO_EDGE *edges;
1515  int numedges = 1;
1516  int i;
1517  double minaz, maxaz;
1518  double az, azdif;
1519 
1520  data->nextCW = data->nextCCW = 0;
1521  data->cwFace = data->ccwFace = -1;
1522 
1523  if ( other ) {
1524  azdif = other->myaz - data->myaz;
1525  if ( azdif < 0 ) azdif += 2 * M_PI;
1526  minaz = maxaz = azdif;
1527  /* TODO: set nextCW/nextCCW/cwFace/ccwFace to other->something ? */
1528  LWDEBUGF(1, "Other edge end has cwFace=%d and ccwFace=%d",
1529  other->cwFace, other->ccwFace);
1530  } else {
1531  minaz = maxaz = -1;
1532  }
1533 
1534  LWDEBUGF(1, "Looking for edges incident to node %" LWTFMT_ELEMID
1535  " and adjacent to azimuth %g", node, data->myaz);
1536 
1537  /* Get incident edges */
1538  edges = lwt_be_getEdgeByNode( topo, &node, &numedges, LWT_COL_EDGE_ALL );
1539  if ( numedges == -1 ) {
1540  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1541  return 0;
1542  }
1543 
1544  LWDEBUGF(1, "getEdgeByNode returned %d edges, minaz=%g, maxaz=%g",
1545  numedges, minaz, maxaz);
1546 
1547  /* For each incident edge-end (1 or 2): */
1548  for ( i = 0; i < numedges; ++i )
1549  {
1550  LWT_ISO_EDGE *edge;
1551  LWGEOM *g;
1552  LWGEOM *cleangeom;
1553  POINT2D p1, p2;
1554  POINTARRAY *pa;
1555 
1556  edge = &(edges[i]);
1557 
1558  if ( edge->edge_id == myedge_id ) continue;
1559 
1560  g = lwline_as_lwgeom(edge->geom);
1561  /* NOTE: remove_repeated_points call could be replaced by
1562  * some other mean to pick two distinct points for endpoints */
1563  cleangeom = lwgeom_remove_repeated_points( g, 0 );
1564  pa = lwgeom_as_lwline(cleangeom)->points;
1565 
1566  if ( pa->npoints < 2 ) {{
1567  LWT_ELEMID id = edge->edge_id;
1568  _lwt_release_edges(edges, numedges);
1569  lwgeom_free(cleangeom);
1570  lwerror("corrupted topology: edge %" LWTFMT_ELEMID
1571  " does not have two distinct points", id);
1572  return -1;
1573  }}
1574 
1575  if ( edge->start_node == node ) {
1576  getPoint2d_p(pa, 0, &p1);
1577  getPoint2d_p(pa, 1, &p2);
1578  LWDEBUGF(1, "edge %" LWTFMT_ELEMID
1579  " starts on node %" LWTFMT_ELEMID
1580  ", edgeend is %g,%g-%g,%g",
1581  edge->edge_id, node, p1.x, p1.y, p2.x, p2.y);
1582  if ( ! azimuth_pt_pt(&p1, &p2, &az) ) {{
1583  LWT_ELEMID id = edge->edge_id;
1584  _lwt_release_edges(edges, numedges);
1585  lwgeom_free(cleangeom);
1586  lwerror("error computing azimuth of edge %d first edgeend [%g,%g-%g,%g]",
1587  id, p1.x, p1.y, p2.x, p2.y);
1588  return -1;
1589  }}
1590  azdif = az - data->myaz;
1591  LWDEBUGF(1, "azimuth of edge %" LWTFMT_ELEMID
1592  ": %g (diff: %g)", edge->edge_id, az, azdif);
1593 
1594  if ( azdif < 0 ) azdif += 2 * M_PI;
1595  if ( minaz == -1 ) {
1596  minaz = maxaz = azdif;
1597  data->nextCW = data->nextCCW = edge->edge_id; /* outgoing */
1598  data->cwFace = edge->face_left;
1599  data->ccwFace = edge->face_right;
1600  LWDEBUGF(1, "new nextCW and nextCCW edge is %" LWTFMT_ELEMID
1601  ", outgoing, "
1602  "with face_left %" LWTFMT_ELEMID " and face_right %" LWTFMT_ELEMID
1603  " (face_right is new ccwFace, face_left is new cwFace)",
1604  edge->edge_id, edge->face_left,
1605  edge->face_right);
1606  } else {
1607  if ( azdif < minaz ) {
1608  data->nextCW = edge->edge_id; /* outgoing */
1609  data->cwFace = edge->face_left;
1610  LWDEBUGF(1, "new nextCW edge is %" LWTFMT_ELEMID
1611  ", outgoing, "
1612  "with face_left %" LWTFMT_ELEMID " and face_right %" LWTFMT_ELEMID
1613  " (previous had minaz=%g, face_left is new cwFace)",
1614  edge->edge_id, edge->face_left,
1615  edge->face_right, minaz);
1616  minaz = azdif;
1617  }
1618  else if ( azdif > maxaz ) {
1619  data->nextCCW = edge->edge_id; /* outgoing */
1620  data->ccwFace = edge->face_right;
1621  LWDEBUGF(1, "new nextCCW edge is %" LWTFMT_ELEMID
1622  ", outgoing, "
1623  "with face_left %" LWTFMT_ELEMID " and face_right %" LWTFMT_ELEMID
1624  " (previous had maxaz=%g, face_right is new ccwFace)",
1625  edge->edge_id, edge->face_left,
1626  edge->face_right, maxaz);
1627  maxaz = azdif;
1628  }
1629  }
1630  }
1631 
1632  if ( edge->end_node == node ) {
1633  getPoint2d_p(pa, pa->npoints-1, &p1);
1634  getPoint2d_p(pa, pa->npoints-2, &p2);
1635  LWDEBUGF(1, "edge %" LWTFMT_ELEMID " ends on node %" LWTFMT_ELEMID
1636  ", edgeend is %g,%g-%g,%g",
1637  edge->edge_id, node, p1.x, p1.y, p2.x, p2.y);
1638  if ( ! azimuth_pt_pt(&p1, &p2, &az) ) {{
1639  LWT_ELEMID id = edge->edge_id;
1640  _lwt_release_edges(edges, numedges);
1641  lwgeom_free(cleangeom);
1642  lwerror("error computing azimuth of edge %d last edgeend [%g,%g-%g,%g]",
1643  id, p1.x, p1.y, p2.x, p2.y);
1644  return -1;
1645  }}
1646  azdif = az - data->myaz;
1647  LWDEBUGF(1, "azimuth of edge %" LWTFMT_ELEMID
1648  ": %g (diff: %g)", edge->edge_id, az, azdif);
1649  if ( azdif < 0 ) azdif += 2 * M_PI;
1650  if ( minaz == -1 ) {
1651  minaz = maxaz = azdif;
1652  data->nextCW = data->nextCCW = -edge->edge_id; /* incoming */
1653  data->cwFace = edge->face_right;
1654  data->ccwFace = edge->face_left;
1655  LWDEBUGF(1, "new nextCW and nextCCW edge is %" LWTFMT_ELEMID
1656  ", incoming, "
1657  "with face_left %" LWTFMT_ELEMID " and face_right %" LWTFMT_ELEMID
1658  " (face_right is new cwFace, face_left is new ccwFace)",
1659  edge->edge_id, edge->face_left,
1660  edge->face_right);
1661  } else {
1662  if ( azdif < minaz ) {
1663  data->nextCW = -edge->edge_id; /* incoming */
1664  data->cwFace = edge->face_right;
1665  LWDEBUGF(1, "new nextCW edge is %" LWTFMT_ELEMID
1666  ", incoming, "
1667  "with face_left %" LWTFMT_ELEMID " and face_right %" LWTFMT_ELEMID
1668  " (previous had minaz=%g, face_right is new cwFace)",
1669  edge->edge_id, edge->face_left,
1670  edge->face_right, minaz);
1671  minaz = azdif;
1672  }
1673  else if ( azdif > maxaz ) {
1674  data->nextCCW = -edge->edge_id; /* incoming */
1675  data->ccwFace = edge->face_left;
1676  LWDEBUGF(1, "new nextCCW edge is %" LWTFMT_ELEMID
1677  ", outgoing, from start point, "
1678  "with face_left %" LWTFMT_ELEMID " and face_right %" LWTFMT_ELEMID
1679  " (previous had maxaz=%g, face_left is new ccwFace)",
1680  edge->edge_id, edge->face_left,
1681  edge->face_right, maxaz);
1682  maxaz = azdif;
1683  }
1684  }
1685  }
1686 
1687  lwgeom_free(cleangeom);
1688  }
1689  if ( numedges ) _lwt_release_edges(edges, numedges);
1690 
1691  LWDEBUGF(1, "edges adjacent to azimuth %g"
1692  " (incident to node %" LWTFMT_ELEMID ")"
1693  ": CW:%" LWTFMT_ELEMID "(%g) CCW:%" LWTFMT_ELEMID "(%g)",
1694  data->myaz, node, data->nextCW, minaz,
1695  data->nextCCW, maxaz);
1696 
1697  if ( myedge_id < 1 && numedges && data->cwFace != data->ccwFace )
1698  {
1699  if ( data->cwFace != -1 && data->ccwFace != -1 ) {
1700  lwerror("Corrupted topology: adjacent edges %" LWTFMT_ELEMID " and %" LWTFMT_ELEMID
1701  " bind different face (%" LWTFMT_ELEMID " and %" LWTFMT_ELEMID ")",
1702  data->nextCW, data->nextCCW,
1703  data->cwFace, data->ccwFace);
1704  return -1;
1705  }
1706  }
1707 
1708  /* Return number of incident edges found */
1709  return numedges;
1710 }
1711 
1712 /*
1713  * Get a point internal to the line and write it into the "ip"
1714  * parameter
1715  *
1716  * return 0 on failure (line is empty or collapsed), 1 otherwise
1717  */
1718 static int
1720 {
1721  int i;
1722  POINT2D fp, lp, tp;
1723  POINTARRAY *pa = edge->points;
1724 
1725  if ( pa->npoints < 2 ) return 0; /* empty or structurally collapsed */
1726 
1727  getPoint2d_p(pa, 0, &fp); /* save first point */
1728  getPoint2d_p(pa, pa->npoints-1, &lp); /* save last point */
1729  for (i=1; i<pa->npoints-1; ++i)
1730  {
1731  getPoint2d_p(pa, i, &tp); /* pick next point */
1732  if ( p2d_same(&tp, &fp) ) continue; /* equal to startpoint */
1733  if ( p2d_same(&tp, &lp) ) continue; /* equal to endpoint */
1734  /* this is a good one, neither same of start nor of end point */
1735  *ip = tp;
1736  return 1; /* found */
1737  }
1738 
1739  /* no distinct vertex found */
1740 
1741  /* interpolate if start point != end point */
1742 
1743  if ( p2d_same(&fp, &lp) ) return 0; /* no distinct points in edge */
1744 
1745  ip->x = fp.x + ( (lp.x - fp.x) * 0.5 );
1746  ip->y = fp.y + ( (lp.y - fp.y) * 0.5 );
1747 
1748  return 1;
1749 }
1750 
1751 /*
1752  * Add a split face by walking on the edge side.
1753  *
1754  * @param topo the topology to act upon
1755  * @param sedge edge id and walking side and direction
1756  * (forward,left:positive backward,right:negative)
1757  * @param face the face in which the edge identifier is known to be
1758  * @param mbr_only do not create a new face but update MBR of the current
1759  *
1760  * @return:
1761  * -1: if mbr_only was requested
1762  * 0: if the edge does not form a ring
1763  * -1: if it is impossible to create a face on the requested side
1764  * ( new face on the side is the universe )
1765  * -2: error
1766  * >0 : id of newly added face
1767  */
1768 static LWT_ELEMID
1770  LWT_ELEMID sedge, LWT_ELEMID face,
1771  int mbr_only )
1772 {
1773  int numedges, numfaceedges, i, j;
1774  int newface_outside;
1775  int num_signed_edge_ids;
1776  LWT_ELEMID *signed_edge_ids;
1777  LWT_ELEMID *edge_ids;
1778  LWT_ISO_EDGE *edges;
1779  LWT_ISO_EDGE *ring_edges;
1780  LWT_ISO_EDGE *forward_edges = NULL;
1781  int forward_edges_count = 0;
1782  LWT_ISO_EDGE *backward_edges = NULL;
1783  int backward_edges_count = 0;
1784 
1785  signed_edge_ids = lwt_be_getRingEdges(topo, sedge,
1786  &num_signed_edge_ids, 0);
1787  if ( ! signed_edge_ids ) {
1788  lwerror("Backend error (no ring edges for edge %" LWTFMT_ELEMID "): %s",
1789  sedge, lwt_be_lastErrorMessage(topo->be_iface));
1790  return -2;
1791  }
1792  LWDEBUGF(1, "getRingEdges returned %d edges", num_signed_edge_ids);
1793 
1794  /* You can't get to the other side of an edge forming a ring */
1795  for (i=0; i<num_signed_edge_ids; ++i) {
1796  if ( signed_edge_ids[i] == -sedge ) {
1797  /* No split here */
1798  LWDEBUG(1, "not a ring");
1799  lwfree( signed_edge_ids );
1800  return 0;
1801  }
1802  }
1803 
1804  LWDEBUGF(1, "Edge %" LWTFMT_ELEMID " split face %" LWTFMT_ELEMID " (mbr_only:%d)",
1805  sedge, face, mbr_only);
1806 
1807  /* Construct a polygon using edges of the ring */
1808  numedges = 0;
1809  edge_ids = lwalloc(sizeof(LWT_ELEMID)*num_signed_edge_ids);
1810  for (i=0; i<num_signed_edge_ids; ++i) {
1811  int absid = llabs(signed_edge_ids[i]);
1812  int found = 0;
1813  /* Do not add the same edge twice */
1814  for (j=0; j<numedges; ++j) {
1815  if ( edge_ids[j] == absid ) {
1816  found = 1;
1817  break;
1818  }
1819  }
1820  if ( ! found ) edge_ids[numedges++] = absid;
1821  }
1822  i = numedges;
1823  ring_edges = lwt_be_getEdgeById(topo, edge_ids, &i,
1825  lwfree( edge_ids );
1826  if ( i == -1 )
1827  {
1828  lwfree( signed_edge_ids );
1829  /* ring_edges should be NULL */
1830  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1831  return -2;
1832  }
1833  else if ( i != numedges )
1834  {
1835  lwfree( signed_edge_ids );
1836  _lwt_release_edges(ring_edges, numedges);
1837  lwerror("Unexpected error: %d edges found when expecting %d", i, numedges);
1838  return -2;
1839  }
1840 
1841  /* Should now build a polygon with those edges, in the order
1842  * given by GetRingEdges.
1843  */
1844  POINTARRAY *pa = NULL;
1845  for ( i=0; i<num_signed_edge_ids; ++i )
1846  {
1847  LWT_ELEMID eid = signed_edge_ids[i];
1848  LWDEBUGF(1, "Edge %d in ring of edge %" LWTFMT_ELEMID " is edge %" LWTFMT_ELEMID,
1849  i, sedge, eid);
1850  LWT_ISO_EDGE *edge = NULL;
1851  POINTARRAY *epa;
1852  for ( j=0; j<numedges; ++j )
1853  {
1854  if ( ring_edges[j].edge_id == llabs(eid) )
1855  {
1856  edge = &(ring_edges[j]);
1857  break;
1858  }
1859  }
1860  if ( edge == NULL )
1861  {
1862  lwfree( signed_edge_ids );
1863  _lwt_release_edges(ring_edges, numedges);
1864  lwerror("missing edge that was found in ring edges loop");
1865  return -2;
1866  }
1867 
1868  if ( pa == NULL )
1869  {
1870  pa = ptarray_clone_deep(edge->geom->points);
1871  if ( eid < 0 ) ptarray_reverse(pa);
1872  }
1873  else
1874  {
1875  if ( eid < 0 )
1876  {
1877  epa = ptarray_clone_deep(edge->geom->points);
1878  ptarray_reverse(epa);
1879  ptarray_append_ptarray(pa, epa, 0);
1880  ptarray_free(epa);
1881  }
1882  else
1883  {
1884  /* avoid a clone here */
1885  ptarray_append_ptarray(pa, edge->geom->points, 0);
1886  }
1887  }
1888  }
1889  POINTARRAY **points = lwalloc(sizeof(POINTARRAY*));
1890  points[0] = pa;
1891  /* NOTE: the ring may very well have collapsed components,
1892  * which would make it topologically invalid
1893  */
1894  LWPOLY* shell = lwpoly_construct(0, 0, 1, points);
1895 
1896  int isccw = ptarray_isccw(pa);
1897  LWDEBUGF(1, "Ring of edge %" LWTFMT_ELEMID " is %sclockwise",
1898  sedge, isccw ? "counter" : "");
1899  const GBOX* shellbox = lwgeom_get_bbox(lwpoly_as_lwgeom(shell));
1900 
1901  if ( face == 0 )
1902  {
1903  /* Edge split the universe face */
1904  if ( ! isccw )
1905  {
1906  lwpoly_free(shell);
1907  lwfree( signed_edge_ids );
1908  _lwt_release_edges(ring_edges, numedges);
1909  /* Face on the left side of this ring is the universe face.
1910  * Next call (for the other side) should create the split face
1911  */
1912  LWDEBUG(1, "The left face of this clockwise ring is the universe, "
1913  "won't create a new face there");
1914  return -1;
1915  }
1916  }
1917 
1918  if ( mbr_only && face != 0 )
1919  {
1920  if ( isccw )
1921  {{
1922  LWT_ISO_FACE updface;
1923  updface.face_id = face;
1924  updface.mbr = (GBOX *)shellbox; /* const cast, we won't free it, later */
1925  int ret = lwt_be_updateFacesById( topo, &updface, 1 );
1926  if ( ret == -1 )
1927  {
1928  lwfree( signed_edge_ids );
1929  _lwt_release_edges(ring_edges, numedges);
1930  lwpoly_free(shell); /* NOTE: owns shellbox above */
1931  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1932  return -2;
1933  }
1934  if ( ret != 1 )
1935  {
1936  lwfree( signed_edge_ids );
1937  _lwt_release_edges(ring_edges, numedges);
1938  lwpoly_free(shell); /* NOTE: owns shellbox above */
1939  lwerror("Unexpected error: %d faces found when expecting 1", ret);
1940  return -2;
1941  }
1942  }}
1943  lwfree( signed_edge_ids );
1944  _lwt_release_edges(ring_edges, numedges);
1945  lwpoly_free(shell); /* NOTE: owns shellbox above */
1946  return -1; /* mbr only was requested */
1947  }
1948 
1949  LWT_ISO_FACE *oldface = NULL;
1950  LWT_ISO_FACE newface;
1951  newface.face_id = -1;
1952  if ( face != 0 && ! isccw)
1953  {{
1954  /* Face created an hole in an outer face */
1955  int nfaces = 1;
1956  oldface = lwt_be_getFaceById(topo, &face, &nfaces, LWT_COL_FACE_ALL);
1957  if ( nfaces == -1 )
1958  {
1959  lwfree( signed_edge_ids );
1960  lwpoly_free(shell); /* NOTE: owns shellbox */
1961  _lwt_release_edges(ring_edges, numedges);
1962  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1963  return -2;
1964  }
1965  if ( nfaces != 1 )
1966  {
1967  lwfree( signed_edge_ids );
1968  lwpoly_free(shell); /* NOTE: owns shellbox */
1969  _lwt_release_edges(ring_edges, numedges);
1970  lwerror("Unexpected error: %d faces found when expecting 1", nfaces);
1971  return -2;
1972  }
1973  newface.mbr = oldface->mbr;
1974  }}
1975  else
1976  {
1977  newface.mbr = (GBOX *)shellbox; /* const cast, we won't free it, later */
1978  }
1979 
1980  /* Insert the new face */
1981  int ret = lwt_be_insertFaces( topo, &newface, 1 );
1982  if ( ret == -1 )
1983  {
1984  lwfree( signed_edge_ids );
1985  lwpoly_free(shell); /* NOTE: owns shellbox */
1986  _lwt_release_edges(ring_edges, numedges);
1987  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
1988  return -2;
1989  }
1990  if ( ret != 1 )
1991  {
1992  lwfree( signed_edge_ids );
1993  lwpoly_free(shell); /* NOTE: owns shellbox */
1994  _lwt_release_edges(ring_edges, numedges);
1995  lwerror("Unexpected error: %d faces inserted when expecting 1", ret);
1996  return -2;
1997  }
1998  if ( oldface ) {
1999  newface.mbr = NULL; /* it is a reference to oldface mbr... */
2000  _lwt_release_faces(oldface, 1);
2001  }
2002 
2003  /* Update side location of new face edges */
2004 
2005  /* We want the new face to be on the left, if possible */
2006  if ( face != 0 && ! isccw ) { /* ring is clockwise in a real face */
2007  /* face shrinked, must update all non-contained edges and nodes */
2008  LWDEBUG(1, "New face is on the outside of the ring, updating rings in former shell");
2009  newface_outside = 1;
2010  /* newface is outside */
2011  } else {
2012  LWDEBUG(1, "New face is on the inside of the ring, updating forward edges in new ring");
2013  newface_outside = 0;
2014  /* newface is inside */
2015  }
2016 
2017  /* Update edges bounding the old face */
2018  /* (1) fetch all edges where left_face or right_face is = oldface */
2019  int fields = LWT_COL_EDGE_EDGE_ID |
2023  ;
2024  numfaceedges = 1;
2025  edges = lwt_be_getEdgeByFace( topo, &face, &numfaceedges, fields, newface.mbr );
2026  if ( numfaceedges == -1 ) {
2027  lwfree( signed_edge_ids );
2028  _lwt_release_edges(ring_edges, numedges);
2029  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2030  return -2;
2031  }
2032  LWDEBUGF(1, "lwt_be_getEdgeByFace returned %d edges", numfaceedges);
2033  GEOSGeometry *shellgg = 0;
2034  const GEOSPreparedGeometry* prepshell = 0;
2035  shellgg = LWGEOM2GEOS( lwpoly_as_lwgeom(shell), 0);
2036  if ( ! shellgg ) {
2037  lwpoly_free(shell);
2038  lwfree(signed_edge_ids);
2039  _lwt_release_edges(ring_edges, numedges);
2040  _lwt_release_edges(edges, numfaceedges);
2041  lwerror("Could not convert shell geometry to GEOS: %s", lwgeom_geos_errmsg);
2042  return -2;
2043  }
2044  prepshell = GEOSPrepare( shellgg );
2045  if ( ! prepshell ) {
2046  GEOSGeom_destroy(shellgg);
2047  lwpoly_free(shell);
2048  lwfree(signed_edge_ids);
2049  _lwt_release_edges(ring_edges, numedges);
2050  _lwt_release_edges(edges, numfaceedges);
2051  lwerror("Could not prepare shell geometry: %s", lwgeom_geos_errmsg);
2052  return -2;
2053  }
2054 
2055  if ( numfaceedges )
2056  {
2057  forward_edges = lwalloc(sizeof(LWT_ISO_EDGE)*numfaceedges);
2058  forward_edges_count = 0;
2059  backward_edges = lwalloc(sizeof(LWT_ISO_EDGE)*numfaceedges);
2060  backward_edges_count = 0;
2061 
2062  /* (2) loop over the results and: */
2063  for ( i=0; i<numfaceedges; ++i )
2064  {
2065  LWT_ISO_EDGE *e = &(edges[i]);
2066  int found = 0;
2067  int contains;
2068  GEOSGeometry *egg;
2069  LWPOINT *epgeom;
2070  POINT2D ep;
2071 
2072  /* (2.1) skip edges whose ID is in the list of boundary edges ? */
2073  for ( j=0; j<num_signed_edge_ids; ++j )
2074  {
2075  int seid = signed_edge_ids[j];
2076  if ( seid == e->edge_id )
2077  {
2078  /* IDEA: remove entry from signed_edge_ids ? */
2079  LWDEBUGF(1, "Edge %d is a forward edge of the new ring", e->edge_id);
2080  forward_edges[forward_edges_count].edge_id = e->edge_id;
2081  forward_edges[forward_edges_count++].face_left = newface.face_id;
2082  found++;
2083  if ( found == 2 ) break;
2084  }
2085  else if ( -seid == e->edge_id )
2086  {
2087  /* IDEA: remove entry from signed_edge_ids ? */
2088  LWDEBUGF(1, "Edge %d is a backward edge of the new ring", e->edge_id);
2089  backward_edges[backward_edges_count].edge_id = e->edge_id;
2090  backward_edges[backward_edges_count++].face_right = newface.face_id;
2091  found++;
2092  if ( found == 2 ) break;
2093  }
2094  }
2095  if ( found ) continue;
2096 
2097  /* We need to check only a single point
2098  * (to avoid collapsed elements of the shell polygon
2099  * giving false positive).
2100  * The point but must not be an endpoint.
2101  */
2102  if ( ! _lwt_GetInteriorEdgePoint(e->geom, &ep) )
2103  {
2104  GEOSPreparedGeom_destroy(prepshell);
2105  GEOSGeom_destroy(shellgg);
2106  lwfree(signed_edge_ids);
2107  lwpoly_free(shell);
2108  lwfree(forward_edges); /* contents owned by ring_edges */
2109  lwfree(backward_edges); /* contents owned by ring_edges */
2110  _lwt_release_edges(ring_edges, numedges);
2111  _lwt_release_edges(edges, numfaceedges);
2112  lwerror("Could not find interior point for edge %d: %s",
2114  return -2;
2115  }
2116 
2117  epgeom = lwpoint_make2d(0, ep.x, ep.y);
2118  egg = LWGEOM2GEOS( lwpoint_as_lwgeom(epgeom) , 0);
2119  lwpoint_free(epgeom);
2120  if ( ! egg ) {
2121  GEOSPreparedGeom_destroy(prepshell);
2122  GEOSGeom_destroy(shellgg);
2123  lwfree(signed_edge_ids);
2124  lwpoly_free(shell);
2125  lwfree(forward_edges); /* contents owned by ring_edges */
2126  lwfree(backward_edges); /* contents owned by ring_edges */
2127  _lwt_release_edges(ring_edges, numedges);
2128  _lwt_release_edges(edges, numfaceedges);
2129  lwerror("Could not convert edge geometry to GEOS: %s",
2131  return -2;
2132  }
2133  /* IDEA: can be optimized by computing this on our side rather
2134  * than on GEOS (saves conversion of big edges) */
2135  /* IDEA: check that bounding box shortcut is taken, or use
2136  * shellbox to do it here */
2137  contains = GEOSPreparedContains( prepshell, egg );
2138  GEOSGeom_destroy(egg);
2139  if ( contains == 2 )
2140  {
2141  GEOSPreparedGeom_destroy(prepshell);
2142  GEOSGeom_destroy(shellgg);
2143  lwfree(signed_edge_ids);
2144  lwpoly_free(shell);
2145  lwfree(forward_edges); /* contents owned by ring_edges */
2146  lwfree(backward_edges); /* contents owned by ring_edges */
2147  _lwt_release_edges(ring_edges, numedges);
2148  _lwt_release_edges(edges, numfaceedges);
2149  lwerror("GEOS exception on PreparedContains: %s", lwgeom_geos_errmsg);
2150  return -2;
2151  }
2152  LWDEBUGF(1, "Edge %d %scontained in new ring", e->edge_id,
2153  (contains?"":"not "));
2154 
2155  /* (2.2) skip edges (NOT, if newface_outside) contained in shell */
2156  if ( newface_outside )
2157  {
2158  if ( contains )
2159  {
2160  LWDEBUGF(1, "Edge %d contained in an hole of the new face",
2161  e->edge_id);
2162  continue;
2163  }
2164  }
2165  else
2166  {
2167  if ( ! contains )
2168  {
2169  LWDEBUGF(1, "Edge %d not contained in the face shell",
2170  e->edge_id);
2171  continue;
2172  }
2173  }
2174 
2175  /* (2.3) push to forward_edges if left_face = oface */
2176  if ( e->face_left == face )
2177  {
2178  LWDEBUGF(1, "Edge %d has new face on the left side", e->edge_id);
2179  forward_edges[forward_edges_count].edge_id = e->edge_id;
2180  forward_edges[forward_edges_count++].face_left = newface.face_id;
2181  }
2182 
2183  /* (2.4) push to backward_edges if right_face = oface */
2184  if ( e->face_right == face )
2185  {
2186  LWDEBUGF(1, "Edge %d has new face on the right side", e->edge_id);
2187  backward_edges[backward_edges_count].edge_id = e->edge_id;
2188  backward_edges[backward_edges_count++].face_right = newface.face_id;
2189  }
2190  }
2191 
2192  /* Update forward edges */
2193  if ( forward_edges_count )
2194  {
2195  ret = lwt_be_updateEdgesById(topo, forward_edges,
2196  forward_edges_count,
2198  if ( ret == -1 )
2199  {
2200  lwfree( signed_edge_ids );
2201  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2202  return -2;
2203  }
2204  if ( ret != forward_edges_count )
2205  {
2206  lwfree( signed_edge_ids );
2207  lwerror("Unexpected error: %d edges updated when expecting %d",
2208  ret, forward_edges_count);
2209  return -2;
2210  }
2211  }
2212 
2213  /* Update backward edges */
2214  if ( backward_edges_count )
2215  {
2216  ret = lwt_be_updateEdgesById(topo, backward_edges,
2217  backward_edges_count,
2219  if ( ret == -1 )
2220  {
2221  lwfree( signed_edge_ids );
2222  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2223  return -2;
2224  }
2225  if ( ret != backward_edges_count )
2226  {
2227  lwfree( signed_edge_ids );
2228  lwerror("Unexpected error: %d edges updated when expecting %d",
2229  ret, backward_edges_count);
2230  return -2;
2231  }
2232  }
2233 
2234  lwfree(forward_edges);
2235  lwfree(backward_edges);
2236 
2237  }
2238 
2239  _lwt_release_edges(ring_edges, numedges);
2240  _lwt_release_edges(edges, numfaceedges);
2241 
2242  /* Update isolated nodes which are now in new face */
2243  int numisonodes = 1;
2245  LWT_ISO_NODE *nodes = lwt_be_getNodeByFace(topo, &face,
2246  &numisonodes, fields, newface.mbr);
2247  if ( numisonodes == -1 ) {
2248  lwfree( signed_edge_ids );
2249  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2250  return -2;
2251  }
2252  if ( numisonodes ) {
2253  LWT_ISO_NODE *updated_nodes = lwalloc(sizeof(LWT_ISO_NODE)*numisonodes);
2254  int nodes_to_update = 0;
2255  for (i=0; i<numisonodes; ++i)
2256  {
2257  LWT_ISO_NODE *n = &(nodes[i]);
2258  GEOSGeometry *ngg;
2259  ngg = LWGEOM2GEOS( lwpoint_as_lwgeom(n->geom), 0 );
2260  int contains;
2261  if ( ! ngg ) {
2262  _lwt_release_nodes(nodes, numisonodes);
2263  if ( prepshell ) GEOSPreparedGeom_destroy(prepshell);
2264  if ( shellgg ) GEOSGeom_destroy(shellgg);
2265  lwfree(signed_edge_ids);
2266  lwpoly_free(shell);
2267  lwerror("Could not convert node geometry to GEOS: %s",
2269  return -2;
2270  }
2271  contains = GEOSPreparedContains( prepshell, ngg );
2272  GEOSGeom_destroy(ngg);
2273  if ( contains == 2 )
2274  {
2275  _lwt_release_nodes(nodes, numisonodes);
2276  if ( prepshell ) GEOSPreparedGeom_destroy(prepshell);
2277  if ( shellgg ) GEOSGeom_destroy(shellgg);
2278  lwfree(signed_edge_ids);
2279  lwpoly_free(shell);
2280  lwerror("GEOS exception on PreparedContains: %s", lwgeom_geos_errmsg);
2281  return -2;
2282  }
2283  LWDEBUGF(1, "Node %d is %scontained in new ring, newface is %s",
2284  n->node_id, contains ? "" : "not ",
2285  newface_outside ? "outside" : "inside" );
2286  if ( newface_outside )
2287  {
2288  if ( contains )
2289  {
2290  LWDEBUGF(1, "Node %d contained in an hole of the new face",
2291  n->node_id);
2292  continue;
2293  }
2294  }
2295  else
2296  {
2297  if ( ! contains )
2298  {
2299  LWDEBUGF(1, "Node %d not contained in the face shell",
2300  n->node_id);
2301  continue;
2302  }
2303  }
2304  updated_nodes[nodes_to_update].node_id = n->node_id;
2305  updated_nodes[nodes_to_update++].containing_face =
2306  newface.face_id;
2307  LWDEBUGF(1, "Node %d will be updated", n->node_id);
2308  }
2309  _lwt_release_nodes(nodes, numisonodes);
2310  if ( nodes_to_update )
2311  {
2312  int ret = lwt_be_updateNodesById(topo, updated_nodes,
2313  nodes_to_update,
2315  if ( ret == -1 ) {
2316  lwfree( signed_edge_ids );
2317  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2318  return -2;
2319  }
2320  }
2321  lwfree(updated_nodes);
2322  }
2323 
2324  GEOSPreparedGeom_destroy(prepshell);
2325  GEOSGeom_destroy(shellgg);
2326  lwfree(signed_edge_ids);
2327  lwpoly_free(shell);
2328 
2329  return newface.face_id;
2330 }
2331 
2332 static LWT_ELEMID
2334  LWT_ELEMID start_node, LWT_ELEMID end_node,
2335  LWLINE *geom, int skipChecks, int modFace )
2336 {
2337  LWT_ISO_EDGE newedge;
2338  LWGEOM *cleangeom;
2339  edgeend span; /* start point analisys */
2340  edgeend epan; /* end point analisys */
2341  POINT2D p1, pn, p2;
2342  POINTARRAY *pa;
2343  LWT_ELEMID node_ids[2];
2344  const LWPOINT *start_node_geom = NULL;
2345  const LWPOINT *end_node_geom = NULL;
2346  int num_nodes;
2347  LWT_ISO_NODE *endpoints;
2348  int i;
2349  int prev_left;
2350  int prev_right;
2351  LWT_ISO_EDGE seledge;
2352  LWT_ISO_EDGE updedge;
2353 
2354  if ( ! skipChecks )
2355  {
2356  /* curve must be simple */
2357  if ( ! lwgeom_is_simple(lwline_as_lwgeom(geom)) )
2358  {
2359  lwerror("SQL/MM Spatial exception - curve not simple");
2360  return -1;
2361  }
2362  }
2363 
2364  newedge.start_node = start_node;
2365  newedge.end_node = end_node;
2366  newedge.geom = geom;
2367  newedge.face_left = -1;
2368  newedge.face_right = -1;
2369  cleangeom = lwgeom_remove_repeated_points( lwline_as_lwgeom(geom), 0 );
2370 
2371  pa = lwgeom_as_lwline(cleangeom)->points;
2372  if ( pa->npoints < 2 ) {
2373  lwgeom_free(cleangeom);
2374  lwerror("Invalid edge (no two distinct vertices exist)");
2375  return -1;
2376  }
2377 
2378  /* Initialize endpoint info (some of that ) */
2379  span.cwFace = span.ccwFace =
2380  epan.cwFace = epan.ccwFace = -1;
2381 
2382  /* Compute azimut of first edge end on start node */
2383  getPoint2d_p(pa, 0, &p1);
2384  getPoint2d_p(pa, 1, &pn);
2385  if ( p2d_same(&p1, &pn) ) {
2386  lwgeom_free(cleangeom);
2387  /* Can still happen, for 2-point lines */
2388  lwerror("Invalid edge (no two distinct vertices exist)");
2389  return -1;
2390  }
2391  if ( ! azimuth_pt_pt(&p1, &pn, &span.myaz) ) {
2392  lwgeom_free(cleangeom);
2393  lwerror("error computing azimuth of first edgeend [%g,%g-%g,%g]",
2394  p1.x, p1.y, pn.x, pn.y);
2395  return -1;
2396  }
2397  LWDEBUGF(1, "edge's start node is %g,%g", p1.x, p1.y);
2398 
2399  /* Compute azimuth of last edge end on end node */
2400  getPoint2d_p(pa, pa->npoints-1, &p2);
2401  getPoint2d_p(pa, pa->npoints-2, &pn);
2402  lwgeom_free(cleangeom);
2403  if ( ! azimuth_pt_pt(&p2, &pn, &epan.myaz) ) {
2404  lwerror("error computing azimuth of last edgeend [%g,%g-%g,%g]",
2405  p2.x, p2.y, pn.x, pn.y);
2406  return -1;
2407  }
2408  LWDEBUGF(1, "edge's end node is %g,%g", p2.x, p2.y);
2409 
2410  /*
2411  * Check endpoints existance, match with Curve geometry
2412  * and get face information (if any)
2413  */
2414 
2415  if ( start_node != end_node ) {
2416  num_nodes = 2;
2417  node_ids[0] = start_node;
2418  node_ids[1] = end_node;
2419  } else {
2420  num_nodes = 1;
2421  node_ids[0] = start_node;
2422  }
2423 
2424  endpoints = lwt_be_getNodeById( topo, node_ids, &num_nodes, LWT_COL_NODE_ALL );
2425  if ( num_nodes < 0 ) {
2426  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2427  return -1;
2428  }
2429  for ( i=0; i<num_nodes; ++i )
2430  {
2431  LWT_ISO_NODE* node = &(endpoints[i]);
2432  if ( node->containing_face != -1 )
2433  {
2434  if ( newedge.face_left == -1 )
2435  {
2436  newedge.face_left = newedge.face_right = node->containing_face;
2437  }
2438  else if ( newedge.face_left != node->containing_face )
2439  {
2440  _lwt_release_nodes(endpoints, num_nodes);
2441  lwerror("SQL/MM Spatial exception - geometry crosses an edge"
2442  " (endnodes in faces %" LWTFMT_ELEMID " and %" LWTFMT_ELEMID ")",
2443  newedge.face_left, node->containing_face);
2444  }
2445  }
2446 
2447  LWDEBUGF(1, "Node %d, with geom %p (looking for %d and %d)",
2448  node->node_id, node->geom, start_node, end_node);
2449  if ( node->node_id == start_node ) {
2450  start_node_geom = node->geom;
2451  }
2452  if ( node->node_id == end_node ) {
2453  end_node_geom = node->geom;
2454  }
2455  }
2456 
2457  if ( ! skipChecks )
2458  {
2459  if ( ! start_node_geom )
2460  {
2461  if ( num_nodes ) _lwt_release_nodes(endpoints, num_nodes);
2462  lwerror("SQL/MM Spatial exception - non-existent node");
2463  return -1;
2464  }
2465  else
2466  {
2467  pa = start_node_geom->point;
2468  getPoint2d_p(pa, 0, &pn);
2469  if ( ! p2d_same(&pn, &p1) )
2470  {
2471  if ( num_nodes ) _lwt_release_nodes(endpoints, num_nodes);
2472  lwerror("SQL/MM Spatial exception"
2473  " - start node not geometry start point."
2474  //" - start node not geometry start point (%g,%g != %g,%g).", pn.x, pn.y, p1.x, p1.y
2475  );
2476  return -1;
2477  }
2478  }
2479 
2480  if ( ! end_node_geom )
2481  {
2482  if ( num_nodes ) _lwt_release_nodes(endpoints, num_nodes);
2483  lwerror("SQL/MM Spatial exception - non-existent node");
2484  return -1;
2485  }
2486  else
2487  {
2488  pa = end_node_geom->point;
2489  getPoint2d_p(pa, 0, &pn);
2490  if ( ! p2d_same(&pn, &p2) )
2491  {
2492  if ( num_nodes ) _lwt_release_nodes(endpoints, num_nodes);
2493  lwerror("SQL/MM Spatial exception"
2494  " - end node not geometry end point."
2495  //" - end node not geometry end point (%g,%g != %g,%g).", pn.x, pn.y, p2.x, p2.y
2496  );
2497  return -1;
2498  }
2499  }
2500 
2501  if ( num_nodes ) _lwt_release_nodes(endpoints, num_nodes);
2502 
2503  if ( _lwt_CheckEdgeCrossing( topo, start_node, end_node, geom, 0 ) )
2504  return -1;
2505 
2506  } /* ! skipChecks */
2507 
2508  /*
2509  * All checks passed, time to prepare the new edge
2510  */
2511 
2512  newedge.edge_id = lwt_be_getNextEdgeId( topo );
2513  if ( newedge.edge_id == -1 ) {
2514  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2515  return -1;
2516  }
2517 
2518  /* Find adjacent edges to each endpoint */
2519  int isclosed = start_node == end_node;
2520  int found;
2521  found = _lwt_FindAdjacentEdges( topo, start_node, &span,
2522  isclosed ? &epan : NULL, -1 );
2523  if ( found ) {
2524  span.was_isolated = 0;
2525  newedge.next_right = span.nextCW ? span.nextCW : -newedge.edge_id;
2526  prev_left = span.nextCCW ? -span.nextCCW : newedge.edge_id;
2527  LWDEBUGF(1, "New edge %d is connected on start node, "
2528  "next_right is %d, prev_left is %d",
2529  newedge.edge_id, newedge.next_right, prev_left);
2530  if ( newedge.face_right == -1 ) {
2531  newedge.face_right = span.cwFace;
2532  }
2533  if ( newedge.face_left == -1 ) {
2534  newedge.face_left = span.ccwFace;
2535  }
2536  } else {
2537  span.was_isolated = 1;
2538  newedge.next_right = isclosed ? -newedge.edge_id : newedge.edge_id;
2539  prev_left = isclosed ? newedge.edge_id : -newedge.edge_id;
2540  LWDEBUGF(1, "New edge %d is isolated on start node, "
2541  "next_right is %d, prev_left is %d",
2542  newedge.edge_id, newedge.next_right, prev_left);
2543  }
2544 
2545  found = _lwt_FindAdjacentEdges( topo, end_node, &epan,
2546  isclosed ? &span : NULL, -1 );
2547  if ( found ) {
2548  epan.was_isolated = 0;
2549  newedge.next_left = epan.nextCW ? epan.nextCW : newedge.edge_id;
2550  prev_right = epan.nextCCW ? -epan.nextCCW : -newedge.edge_id;
2551  LWDEBUGF(1, "New edge %d is connected on end node, "
2552  "next_left is %d, prev_right is %d",
2553  newedge.edge_id, newedge.next_left, prev_right);
2554  if ( newedge.face_right == -1 ) {
2555  newedge.face_right = span.ccwFace;
2556  }
2557  if ( newedge.face_left == -1 ) {
2558  newedge.face_left = span.cwFace;
2559  }
2560  } else {
2561  epan.was_isolated = 1;
2562  newedge.next_left = isclosed ? newedge.edge_id : -newedge.edge_id;
2563  prev_right = isclosed ? -newedge.edge_id : newedge.edge_id;
2564  LWDEBUGF(1, "New edge %d is isolated on end node, "
2565  "next_left is %d, prev_right is %d",
2566  newedge.edge_id, newedge.next_left, prev_right);
2567  }
2568 
2569  /*
2570  * If we don't have faces setup by now we must have encountered
2571  * a malformed topology (no containing_face on isolated nodes, no
2572  * left/right faces on adjacent edges or mismatching values)
2573  */
2574  if ( newedge.face_left != newedge.face_right )
2575  {
2576  lwerror("Left(%" LWTFMT_ELEMID ")/right(%" LWTFMT_ELEMID ")"
2577  "faces mismatch: invalid topology ?",
2578  newedge.face_left, newedge.face_right);
2579  return -1;
2580  }
2581  else if ( newedge.face_left == -1 )
2582  {
2583  lwerror("Could not derive edge face from linked primitives:"
2584  " invalid topology ?");
2585  return -1;
2586  }
2587 
2588  /*
2589  * Insert the new edge, and update all linking
2590  */
2591 
2592  int ret = lwt_be_insertEdges(topo, &newedge, 1);
2593  if ( ret == -1 ) {
2594  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2595  return -1;
2596  } else if ( ret == 0 ) {
2597  lwerror("Insertion of split edge failed (no reason)");
2598  return -1;
2599  }
2600 
2601  int updfields;
2602 
2603  /* Link prev_left to us
2604  * (if it's not us already) */
2605  if ( llabs(prev_left) != newedge.edge_id )
2606  {
2607  if ( prev_left > 0 )
2608  {
2609  /* its next_left_edge is us */
2610  updfields = LWT_COL_EDGE_NEXT_LEFT;
2611  updedge.next_left = newedge.edge_id;
2612  seledge.edge_id = prev_left;
2613  }
2614  else
2615  {
2616  /* its next_right_edge is us */
2617  updfields = LWT_COL_EDGE_NEXT_RIGHT;
2618  updedge.next_right = newedge.edge_id;
2619  seledge.edge_id = -prev_left;
2620  }
2621 
2622  ret = lwt_be_updateEdges(topo,
2623  &seledge, LWT_COL_EDGE_EDGE_ID,
2624  &updedge, updfields,
2625  NULL, 0);
2626  if ( ret == -1 ) {
2627  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2628  return -1;
2629  }
2630  }
2631 
2632  /* Link prev_right to us
2633  * (if it's not us already) */
2634  if ( llabs(prev_right) != newedge.edge_id )
2635  {
2636  if ( prev_right > 0 )
2637  {
2638  /* its next_left_edge is -us */
2639  updfields = LWT_COL_EDGE_NEXT_LEFT;
2640  updedge.next_left = -newedge.edge_id;
2641  seledge.edge_id = prev_right;
2642  }
2643  else
2644  {
2645  /* its next_right_edge is -us */
2646  updfields = LWT_COL_EDGE_NEXT_RIGHT;
2647  updedge.next_right = -newedge.edge_id;
2648  seledge.edge_id = -prev_right;
2649  }
2650 
2651  ret = lwt_be_updateEdges(topo,
2652  &seledge, LWT_COL_EDGE_EDGE_ID,
2653  &updedge, updfields,
2654  NULL, 0);
2655  if ( ret == -1 ) {
2656  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2657  return -1;
2658  }
2659  }
2660 
2661  /* NOT IN THE SPECS...
2662  * set containing_face = null for start_node and end_node
2663  * if they where isolated
2664  *
2665  */
2666  LWT_ISO_NODE updnode, selnode;
2667  updnode.containing_face = -1;
2668  if ( span.was_isolated )
2669  {
2670  selnode.node_id = start_node;
2671  ret = lwt_be_updateNodes(topo,
2672  &selnode, LWT_COL_NODE_NODE_ID,
2673  &updnode, LWT_COL_NODE_CONTAINING_FACE,
2674  NULL, 0);
2675  if ( ret == -1 ) {
2676  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2677  return -1;
2678  }
2679  }
2680  if ( epan.was_isolated )
2681  {
2682  selnode.node_id = end_node;
2683  ret = lwt_be_updateNodes(topo,
2684  &selnode, LWT_COL_NODE_NODE_ID,
2685  &updnode, LWT_COL_NODE_CONTAINING_FACE,
2686  NULL, 0);
2687  if ( ret == -1 ) {
2688  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2689  return -1;
2690  }
2691  }
2692 
2693  /* Check face splitting */
2694 
2695  if ( ! isclosed && ( epan.was_isolated || span.was_isolated ) )
2696  {
2697  LWDEBUG(1, "New edge is dangling, so it cannot split any face");
2698  return newedge.edge_id; /* no split */
2699  }
2700 
2701  int newface1 = -1;
2702 
2703  /* IDEA: avoid building edge ring if input is closed, which means we
2704  * know in advance it splits a face */
2705 
2706  if ( ! modFace )
2707  {
2708  newface1 = _lwt_AddFaceSplit( topo, -newedge.edge_id, newedge.face_left, 0 );
2709  if ( newface1 == 0 ) {
2710  LWDEBUG(1, "New edge does not split any face");
2711  return newedge.edge_id; /* no split */
2712  }
2713  }
2714 
2715  int newface = _lwt_AddFaceSplit( topo, newedge.edge_id,
2716  newedge.face_left, 0 );
2717  if ( modFace )
2718  {
2719  if ( newface == 0 ) {
2720  LWDEBUG(1, "New edge does not split any face");
2721  return newedge.edge_id; /* no split */
2722  }
2723 
2724  if ( newface < 0 )
2725  {
2726  /* face on the left is the universe face */
2727  /* must be forming a maximal ring in universal face */
2728  newface = _lwt_AddFaceSplit( topo, -newedge.edge_id,
2729  newedge.face_left, 0 );
2730  if ( newface < 0 ) return newedge.edge_id; /* no split */
2731  }
2732  else
2733  {
2734  _lwt_AddFaceSplit( topo, -newedge.edge_id, newedge.face_left, 1 );
2735  }
2736  }
2737 
2738  /*
2739  * Update topogeometries, if needed
2740  */
2741  if ( newedge.face_left != 0 )
2742  {
2743  ret = lwt_be_updateTopoGeomFaceSplit(topo, newedge.face_left,
2744  newface, newface1);
2745  if ( ret == 0 ) {
2746  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2747  return -1;
2748  }
2749 
2750  if ( ! modFace )
2751  {
2752  /* drop old face from the face table */
2753  ret = lwt_be_deleteFacesById(topo, &(newedge.face_left), 1);
2754  if ( ret == -1 ) {
2755  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2756  return -1;
2757  }
2758  }
2759  }
2760 
2761  return newedge.edge_id;
2762 }
2763 
2764 LWT_ELEMID
2766  LWT_ELEMID start_node, LWT_ELEMID end_node,
2767  LWLINE *geom, int skipChecks )
2768 {
2769  return _lwt_AddEdge( topo, start_node, end_node, geom, skipChecks, 1 );
2770 }
2771 
2772 LWT_ELEMID
2774  LWT_ELEMID start_node, LWT_ELEMID end_node,
2775  LWLINE *geom, int skipChecks )
2776 {
2777  return _lwt_AddEdge( topo, start_node, end_node, geom, skipChecks, 0 );
2778 }
2779 
2780 static LWGEOM *
2781 _lwt_FaceByEdges(LWT_TOPOLOGY *topo, LWT_ISO_EDGE *edges, int numfaceedges)
2782 {
2783  LWGEOM *outg;
2784  LWCOLLECTION *bounds;
2785  LWGEOM **geoms = lwalloc( sizeof(LWGEOM*) * numfaceedges );
2786  int i, validedges = 0;
2787 
2788  for ( i=0; i<numfaceedges; ++i )
2789  {
2790  /* NOTE: skipping edges with same face on both sides, although
2791  * correct, results in a failure to build faces from
2792  * invalid topologies as expected by legacy tests.
2793  * TODO: update legacy tests expectances/unleash this skipping ?
2794  */
2795  /* if ( edges[i].face_left == edges[i].face_right ) continue; */
2796  geoms[validedges++] = lwline_as_lwgeom(edges[i].geom);
2797  }
2798  if ( ! validedges )
2799  {
2800  /* Face has no valid boundary edges, we'll return EMPTY, see
2801  * https://trac.osgeo.org/postgis/ticket/3221 */
2802  if ( numfaceedges ) lwfree(geoms);
2803  LWDEBUG(1, "_lwt_FaceByEdges returning empty polygon");
2804  return lwpoly_as_lwgeom(
2805  lwpoly_construct_empty(topo->srid, topo->hasZ, 0)
2806  );
2807  }
2809  topo->srid,
2810  NULL, /* gbox */
2811  validedges,
2812  geoms);
2813  outg = lwgeom_buildarea( lwcollection_as_lwgeom(bounds) );
2814  lwcollection_release(bounds);
2815  lwfree(geoms);
2816 #if 0
2817  {
2818  size_t sz;
2819  char *wkt = lwgeom_to_wkt(outg, WKT_EXTENDED, 2, &sz);
2820  LWDEBUGF(1, "_lwt_FaceByEdges returning area: %s", wkt);
2821  lwfree(wkt);
2822  }
2823 #endif
2824  return outg;
2825 }
2826 
2827 LWGEOM*
2829 {
2830  int numfaceedges;
2831  LWT_ISO_EDGE *edges;
2832  LWT_ISO_FACE *face;
2833  LWPOLY *out;
2834  LWGEOM *outg;
2835  int i;
2836  int fields;
2837 
2838  if ( faceid == 0 )
2839  {
2840  lwerror("SQL/MM Spatial exception - universal face has no geometry");
2841  return NULL;
2842  }
2843 
2844  /* Construct the face geometry */
2845  numfaceedges = 1;
2846  fields = LWT_COL_EDGE_GEOM |
2849  ;
2850  edges = lwt_be_getEdgeByFace( topo, &faceid, &numfaceedges, fields, NULL );
2851  if ( numfaceedges == -1 ) {
2852  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2853  return NULL;
2854  }
2855 
2856  if ( numfaceedges == 0 )
2857  {
2858  i = 1;
2859  face = lwt_be_getFaceById(topo, &faceid, &i, LWT_COL_FACE_FACE_ID);
2860  if ( i == -1 ) {
2861  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
2862  return NULL;
2863  }
2864  if ( i == 0 ) {
2865  lwerror("SQL/MM Spatial exception - non-existent face.");
2866  return NULL;
2867  }
2868  lwfree( face );
2869  if ( i > 1 ) {
2870  lwerror("Corrupted topology: multiple face records have face_id=%"
2871  PRId64, faceid);
2872  return NULL;
2873  }
2874  /* Face has no boundary edges, we'll return EMPTY, see
2875  * https://trac.osgeo.org/postgis/ticket/3221 */
2876  out = lwpoly_construct_empty(topo->srid, topo->hasZ, 0);
2877  return lwpoly_as_lwgeom(out);
2878  }
2879 
2880  outg = _lwt_FaceByEdges( topo, edges, numfaceedges );
2881  _lwt_release_edges(edges, numfaceedges);
2882 
2883  return outg;
2884 }
2885 
2886 /* Find which edge from the "edges" set defines the next
2887  * portion of the given "ring".
2888  *
2889  * The edge might be either forward or backward.
2890  *
2891  * @param ring The ring to find definition of.
2892  * It is assumed it does not contain duplicated vertices.
2893  * @param from offset of the ring point to start looking from
2894  * @param edges array of edges to search into
2895  * @param numedges number of edges in the edges array
2896  *
2897  * @return index of the edge defining the next ring portion or
2898  * -1 if no edge was found to be part of the ring
2899  */
2900 static int
2901 _lwt_FindNextRingEdge(const POINTARRAY *ring, int from,
2902  const LWT_ISO_EDGE *edges, int numedges)
2903 {
2904  int i;
2905  POINT2D p1;
2906 
2907  /* Get starting ring point */
2908  getPoint2d_p(ring, from, &p1);
2909 
2910  LWDEBUGF(1, "Ring's 'from' point (%d) is %g,%g", from, p1.x, p1.y);
2911 
2912  /* find the edges defining the next portion of ring starting from
2913  * vertex "from" */
2914  for ( i=0; i<numedges; ++i )
2915  {
2916  const LWT_ISO_EDGE *isoe = &(edges[i]);
2917  LWLINE *edge = isoe->geom;
2918  POINTARRAY *epa = edge->points;
2919  POINT2D p2, pt;
2920  int match = 0;
2921  int j;
2922 
2923  /* Skip if the edge is a dangling one */
2924  if ( isoe->face_left == isoe->face_right )
2925  {
2926  LWDEBUGF(3, "_lwt_FindNextRingEdge: edge %" LWTFMT_ELEMID
2927  " has same face (%" LWTFMT_ELEMID
2928  ") on both sides, skipping",
2929  isoe->edge_id, isoe->face_left);
2930  continue;
2931  }
2932 
2933 #if 0
2934  size_t sz;
2935  LWDEBUGF(1, "Edge %" LWTFMT_ELEMID " is %s",
2936  isoe->edge_id,
2937  lwgeom_to_wkt(lwline_as_lwgeom(edge), WKT_EXTENDED, 2, &sz));
2938 #endif
2939 
2940  /* ptarray_remove_repeated_points ? */
2941 
2942  getPoint2d_p(epa, 0, &p2);
2943  LWDEBUGF(1, "Edge %" LWTFMT_ELEMID " 'first' point is %g,%g",
2944  isoe->edge_id, p2.x, p2.y);
2945  LWDEBUGF(1, "Rings's 'from' point is still %g,%g", p1.x, p1.y);
2946  if ( p2d_same(&p1, &p2) )
2947  {
2948  LWDEBUG(1, "p2d_same(p1,p2) returned true");
2949  LWDEBUGF(1, "First point of edge %" LWTFMT_ELEMID
2950  " matches ring vertex %d", isoe->edge_id, from);
2951  /* first point matches, let's check next non-equal one */
2952  for ( j=1; j<epa->npoints; ++j )
2953  {
2954  getPoint2d_p(epa, j, &p2);
2955  LWDEBUGF(1, "Edge %" LWTFMT_ELEMID " 'next' point %d is %g,%g",
2956  isoe->edge_id, j, p2.x, p2.y);
2957  /* we won't check duplicated edge points */
2958  if ( p2d_same(&p1, &p2) ) continue;
2959  /* we assume there are no duplicated points in ring */
2960  getPoint2d_p(ring, from+1, &pt);
2961  LWDEBUGF(1, "Ring's point %d is %g,%g",
2962  from+1, pt.x, pt.y);
2963  match = p2d_same(&pt, &p2);
2964  break; /* we want to check a single non-equal next vertex */
2965  }
2966 #if POSTGIS_DEBUG_LEVEL > 0
2967  if ( match ) {
2968  LWDEBUGF(1, "Prev point of edge %" LWTFMT_ELEMID
2969  " matches ring vertex %d", isoe->edge_id, from+1);
2970  } else {
2971  LWDEBUGF(1, "Prev point of edge %" LWTFMT_ELEMID
2972  " does not match ring vertex %d", isoe->edge_id, from+1);
2973  }
2974 #endif
2975  }
2976 
2977  if ( ! match )
2978  {
2979  LWDEBUGF(1, "Edge %" LWTFMT_ELEMID " did not match as forward",
2980  isoe->edge_id);
2981  getPoint2d_p(epa, epa->npoints-1, &p2);
2982  LWDEBUGF(1, "Edge %" LWTFMT_ELEMID " 'last' point is %g,%g",
2983  isoe->edge_id, p2.x, p2.y);
2984  if ( p2d_same(&p1, &p2) )
2985  {
2986  LWDEBUGF(1, "Last point of edge %" LWTFMT_ELEMID
2987  " matches ring vertex %d", isoe->edge_id, from);
2988  /* last point matches, let's check next non-equal one */
2989  for ( j=epa->npoints-2; j>=0; --j )
2990  {
2991  getPoint2d_p(epa, j, &p2);
2992  LWDEBUGF(1, "Edge %" LWTFMT_ELEMID " 'prev' point %d is %g,%g",
2993  isoe->edge_id, j, p2.x, p2.y);
2994  /* we won't check duplicated edge points */
2995  if ( p2d_same(&p1, &p2) ) continue;
2996  /* we assume there are no duplicated points in ring */
2997  getPoint2d_p(ring, from+1, &pt);
2998  LWDEBUGF(1, "Ring's point %d is %g,%g",
2999  from+1, pt.x, pt.y);
3000  match = p2d_same(&pt, &p2);
3001  break; /* we want to check a single non-equal next vertex */
3002  }
3003  }
3004 #if POSTGIS_DEBUG_LEVEL > 0
3005  if ( match ) {
3006  LWDEBUGF(1, "Prev point of edge %" LWTFMT_ELEMID
3007  " matches ring vertex %d", isoe->edge_id, from+1);
3008  } else {
3009  LWDEBUGF(1, "Prev point of edge %" LWTFMT_ELEMID
3010  " does not match ring vertex %d", isoe->edge_id, from+1);
3011  }
3012 #endif
3013  }
3014 
3015  if ( match ) return i;
3016 
3017  }
3018 
3019  return -1;
3020 }
3021 
3022 /* Reverse values in array between "from" (inclusive)
3023  * and "to" (exclusive) indexes */
3024 static void
3025 _lwt_ReverseElemidArray(LWT_ELEMID *ary, int from, int to)
3026 {
3027  LWT_ELEMID t;
3028  while (from < to)
3029  {
3030  t = ary[from];
3031  ary[from++] = ary[to];
3032  ary[to--] = t;
3033  }
3034 }
3035 
3036 /* Rotate values in array between "from" (inclusive)
3037  * and "to" (exclusive) indexes, so that "rotidx" is
3038  * the new value at "from" */
3039 static void
3040 _lwt_RotateElemidArray(LWT_ELEMID *ary, int from, int to, int rotidx)
3041 {
3042  _lwt_ReverseElemidArray(ary, from, rotidx-1);
3043  _lwt_ReverseElemidArray(ary, rotidx, to-1);
3044  _lwt_ReverseElemidArray(ary, from, to-1);
3045 }
3046 
3047 
3048 int
3050 {
3051  LWGEOM *face;
3052  LWPOLY *facepoly;
3053  LWT_ISO_EDGE *edges;
3054  int numfaceedges;
3055  int fields, i;
3056  int nseid = 0; /* number of signed edge ids */
3057  int prevseid;
3058  LWT_ELEMID *seid; /* signed edge ids */
3059 
3060  /* Get list of face edges */
3061  numfaceedges = 1;
3062  fields = LWT_COL_EDGE_EDGE_ID |
3066  ;
3067  edges = lwt_be_getEdgeByFace( topo, &face_id, &numfaceedges, fields, NULL );
3068  if ( numfaceedges == -1 ) {
3069  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3070  return -1;
3071  }
3072  if ( ! numfaceedges ) return 0; /* no edges in output */
3073 
3074  /* order edges by occurrence in face */
3075 
3076  face = _lwt_FaceByEdges(topo, edges, numfaceedges);
3077  if ( ! face )
3078  {
3079  /* _lwt_FaceByEdges should have already invoked lwerror in this case */
3080  _lwt_release_edges(edges, numfaceedges);
3081  return -1;
3082  }
3083 
3084  if ( lwgeom_is_empty(face) )
3085  {
3086  /* no edges in output */
3087  _lwt_release_edges(edges, numfaceedges);
3088  lwgeom_free(face);
3089  return 0;
3090  }
3091 
3092  /* force_lhr, if the face is not the universe */
3093  /* _lwt_FaceByEdges seems to guaranteed RHR */
3094  /* lwgeom_force_clockwise(face); */
3095  if ( face_id ) lwgeom_reverse(face);
3096 
3097 #if 0
3098  {
3099  size_t sz;
3100  char *wkt = lwgeom_to_wkt(face, WKT_EXTENDED, 6, &sz);
3101  LWDEBUGF(1, "Geometry of face %" LWTFMT_ELEMID " is: %s",
3102  face_id, wkt);
3103  lwfree(wkt);
3104  }
3105 #endif
3106 
3107  facepoly = lwgeom_as_lwpoly(face);
3108  if ( ! facepoly )
3109  {
3110  _lwt_release_edges(edges, numfaceedges);
3111  lwgeom_free(face);
3112  lwerror("Geometry of face %" LWTFMT_ELEMID " is not a polygon", face_id);
3113  return -1;
3114  }
3115 
3116  nseid = prevseid = 0;
3117  seid = lwalloc( sizeof(LWT_ELEMID) * numfaceedges );
3118 
3119  /* for each ring of the face polygon... */
3120  for ( i=0; i<facepoly->nrings; ++i )
3121  {
3122  const POINTARRAY *ring = facepoly->rings[i];
3123  int j = 0;
3124  LWT_ISO_EDGE *nextedge;
3125  LWLINE *nextline;
3126 
3127  LWDEBUGF(1, "Ring %d has %d points", i, ring->npoints);
3128 
3129  while ( j < ring->npoints-1 )
3130  {
3131  LWDEBUGF(1, "Looking for edge covering ring %d from vertex %d",
3132  i, j);
3133 
3134  int edgeno = _lwt_FindNextRingEdge(ring, j, edges, numfaceedges);
3135  if ( edgeno == -1 )
3136  {
3137  /* should never happen */
3138  _lwt_release_edges(edges, numfaceedges);
3139  lwgeom_free(face);
3140  lwfree(seid);
3141  lwerror("No edge (among %d) found to be defining geometry of face %"
3142  LWTFMT_ELEMID, numfaceedges, face_id);
3143  return -1;
3144  }
3145 
3146  nextedge = &(edges[edgeno]);
3147  nextline = nextedge->geom;
3148 
3149  LWDEBUGF(1, "Edge %" LWTFMT_ELEMID
3150  " covers ring %d from vertex %d to %d",
3151  nextedge->edge_id, i, j, j + nextline->points->npoints - 1);
3152 
3153 #if 0
3154  {
3155  size_t sz;
3156  char *wkt = lwgeom_to_wkt(lwline_as_lwgeom(nextline), WKT_EXTENDED, 6, &sz);
3157  LWDEBUGF(1, "Edge %" LWTFMT_ELEMID " is %s",
3158  nextedge->edge_id, wkt);
3159  lwfree(wkt);
3160  }
3161 #endif
3162 
3163  j += nextline->points->npoints - 1;
3164 
3165  /* Add next edge to the output array */
3166  seid[nseid++] = nextedge->face_left == face_id ?
3167  nextedge->edge_id :
3168  -nextedge->edge_id;
3169 
3170  /* avoid checking again on next time turn */
3171  nextedge->face_left = nextedge->face_right = -1;
3172  }
3173 
3174  /* now "scroll" the list of edges so that the one
3175  * with smaller absolute edge_id is first */
3176  /* Range is: [prevseid, nseid) -- [inclusive, exclusive) */
3177  if ( (nseid - prevseid) > 1 )
3178  {{
3179  LWT_ELEMID minid = 0;
3180  int minidx = 0;
3181  LWDEBUGF(1, "Looking for smallest id among the %d edges "
3182  "composing ring %d", (nseid-prevseid), i);
3183  for ( j=prevseid; j<nseid; ++j )
3184  {
3185  LWT_ELEMID id = llabs(seid[j]);
3186  LWDEBUGF(1, "Abs id of edge in pos %d is %" LWTFMT_ELEMID, j, id);
3187  if ( ! minid || id < minid )
3188  {
3189  minid = id;
3190  minidx = j;
3191  }
3192  }
3193  LWDEBUGF(1, "Smallest id is %" LWTFMT_ELEMID
3194  " at position %d", minid, minidx);
3195  if ( minidx != prevseid )
3196  {
3197  _lwt_RotateElemidArray(seid, prevseid, nseid, minidx);
3198  }
3199  }}
3200 
3201  prevseid = nseid;
3202  }
3203 
3204  lwgeom_free(face);
3205  _lwt_release_edges(edges, numfaceedges);
3206 
3207  *out = seid;
3208  return nseid;
3209 }
3210 
3211 static GEOSGeometry *
3212 _lwt_EdgeMotionArea(LWLINE *geom, int isclosed)
3213 {
3214  GEOSGeometry *gg;
3215  POINT4D p4d;
3216  POINTARRAY *pa;
3217  POINTARRAY **pas;
3218  LWPOLY *poly;
3219  LWGEOM *g;
3220 
3221  pas = lwalloc(sizeof(POINTARRAY*));
3222 
3223  initGEOS(lwnotice, lwgeom_geos_error);
3224 
3225  if ( isclosed )
3226  {
3227  pas[0] = ptarray_clone_deep( geom->points );
3228  poly = lwpoly_construct(0, 0, 1, pas);
3229  gg = LWGEOM2GEOS( lwpoly_as_lwgeom(poly), 0 );
3230  lwpoly_free(poly); /* should also delete the pointarrays */
3231  }
3232  else
3233  {
3234  pa = geom->points;
3235  getPoint4d_p(pa, 0, &p4d);
3236  pas[0] = ptarray_clone_deep( pa );
3237  /* don't bother dup check */
3238  if ( LW_FAILURE == ptarray_append_point(pas[0], &p4d, LW_TRUE) )
3239  {
3240  ptarray_free(pas[0]);
3241  lwfree(pas);
3242  lwerror("Could not append point to pointarray");
3243  return NULL;
3244  }
3245  poly = lwpoly_construct(0, NULL, 1, pas);
3246  /* make valid, in case the edge self-intersects on its first-last
3247  * vertex segment */
3249  lwpoly_free(poly); /* should also delete the pointarrays */
3250  if ( ! g )
3251  {
3252  lwerror("Could not make edge motion area valid");
3253  return NULL;
3254  }
3255  gg = LWGEOM2GEOS(g, 0);
3256  lwgeom_free(g);
3257  }
3258  if ( ! gg )
3259  {
3260  lwerror("Could not convert old edge area geometry to GEOS: %s",
3262  return NULL;
3263  }
3264  return gg;
3265 }
3266 
3267 int
3269 {
3270  LWT_ISO_EDGE *oldedge;
3271  LWT_ISO_EDGE newedge;
3272  POINT2D p1, p2, pt;
3273  int i;
3274  int isclosed = 0;
3275 
3276  /* curve must be simple */
3277  if ( ! lwgeom_is_simple(lwline_as_lwgeom(geom)) )
3278  {
3279  lwerror("SQL/MM Spatial exception - curve not simple");
3280  return -1;
3281  }
3282 
3283  i = 1;
3284  oldedge = lwt_be_getEdgeById(topo, &edge_id, &i, LWT_COL_EDGE_ALL);
3285  if ( ! oldedge )
3286  {
3287  LWDEBUGF(1, "lwt_ChangeEdgeGeom: "
3288  "lwt_be_getEdgeById returned NULL and set i=%d", i);
3289  if ( i == -1 )
3290  {
3291  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3292  return -1;
3293  }
3294  else if ( i == 0 )
3295  {
3296  lwerror("SQL/MM Spatial exception - non-existent edge %"
3297  LWTFMT_ELEMID, edge_id);
3298  return -1;
3299  }
3300  else
3301  {
3302  lwerror("Backend coding error: getEdgeById callback returned NULL "
3303  "but numelements output parameter has value %d "
3304  "(expected 0 or 1)", i);
3305  return -1;
3306  }
3307  }
3308 
3309  LWDEBUGF(1, "lwt_ChangeEdgeGeom: "
3310  "old edge has %d points, new edge has %d points",
3311  oldedge->geom->points->npoints, geom->points->npoints);
3312 
3313  /*
3314  * e) Check StartPoint consistency
3315  */
3316  getPoint2d_p(oldedge->geom->points, 0, &p1);
3317  getPoint2d_p(geom->points, 0, &pt);
3318  if ( ! p2d_same(&p1, &pt) )
3319  {
3320  _lwt_release_edges(oldedge, 1);
3321  lwerror("SQL/MM Spatial exception - "
3322  "start node not geometry start point.");
3323  return -1;
3324  }
3325 
3326  /*
3327  * f) Check EndPoint consistency
3328  */
3329  if ( oldedge->geom->points->npoints < 2 )
3330  {
3331  _lwt_release_edges(oldedge, 1);
3332  lwerror("Corrupted topology: edge %" LWTFMT_ELEMID
3333  " has less than 2 vertices", oldedge->edge_id);
3334  return -1;
3335  }
3336  getPoint2d_p(oldedge->geom->points, oldedge->geom->points->npoints-1, &p2);
3337  if ( geom->points->npoints < 2 )
3338  {
3339  _lwt_release_edges(oldedge, 1);
3340  lwerror("Invalid edge: less than 2 vertices");
3341  return -1;
3342  }
3343  getPoint2d_p(geom->points, geom->points->npoints-1, &pt);
3344  if ( ! p2d_same(&pt, &p2) )
3345  {
3346  _lwt_release_edges(oldedge, 1);
3347  lwerror("SQL/MM Spatial exception - "
3348  "end node not geometry end point.");
3349  return -1;
3350  }
3351 
3352  /* Not in the specs:
3353  * if the edge is closed, check we didn't change winding !
3354  * (should be part of isomorphism checking)
3355  */
3356  if ( oldedge->start_node == oldedge->end_node )
3357  {
3358  isclosed = 1;
3359 #if 1 /* TODO: this is actually bogus as a test */
3360  /* check for valid edge (distinct vertices must exist) */
3361  if ( ! _lwt_GetInteriorEdgePoint(geom, &pt) )
3362  {
3363  _lwt_release_edges(oldedge, 1);
3364  lwerror("Invalid edge (no two distinct vertices exist)");
3365  return -1;
3366  }
3367 #endif
3368 
3369  if ( ptarray_isccw(oldedge->geom->points) !=
3370  ptarray_isccw(geom->points) )
3371  {
3372  _lwt_release_edges(oldedge, 1);
3373  lwerror("Edge twist at node POINT(%g %g)", p1.x, p1.y);
3374  return -1;
3375  }
3376  }
3377 
3378  if ( _lwt_CheckEdgeCrossing(topo, oldedge->start_node,
3379  oldedge->end_node, geom, edge_id ) )
3380  {
3381  /* would have called lwerror already, leaking :( */
3382  _lwt_release_edges(oldedge, 1);
3383  return -1;
3384  }
3385 
3386  LWDEBUG(1, "lwt_ChangeEdgeGeom: "
3387  "edge crossing check passed ");
3388 
3389  /*
3390  * Not in the specs:
3391  * Check topological isomorphism
3392  */
3393 
3394  /* Check that the "motion range" doesn't include any node */
3395  // 1. compute combined bbox of old and new edge
3396  GBOX mbox; /* motion box */
3397  lwgeom_add_bbox((LWGEOM*)oldedge->geom); /* just in case */
3398  lwgeom_add_bbox((LWGEOM*)geom); /* just in case */
3399  gbox_union(oldedge->geom->bbox, geom->bbox, &mbox);
3400  // 2. fetch all nodes in the combined box
3401  LWT_ISO_NODE *nodes;
3402  int numnodes;
3403  nodes = lwt_be_getNodeWithinBox2D(topo, &mbox, &numnodes,
3404  LWT_COL_NODE_ALL, 0);
3405  LWDEBUGF(1, "lwt_be_getNodeWithinBox2D returned %d nodes", numnodes);
3406  if ( numnodes == -1 ) {
3407  _lwt_release_edges(oldedge, 1);
3408  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3409  return -1;
3410  }
3411  // 3. if any node beside endnodes are found:
3412  if ( numnodes > ( 1 + isclosed ? 0 : 1 ) )
3413  {{
3414  GEOSGeometry *oarea, *narea;
3415  const GEOSPreparedGeometry *oareap, *nareap;
3416 
3417  initGEOS(lwnotice, lwgeom_geos_error);
3418 
3419  oarea = _lwt_EdgeMotionArea(oldedge->geom, isclosed);
3420  if ( ! oarea )
3421  {
3422  _lwt_release_edges(oldedge, 1);
3423  lwerror("Could not compute edge motion area for old edge");
3424  return -1;
3425  }
3426 
3427  narea = _lwt_EdgeMotionArea(geom, isclosed);
3428  if ( ! narea )
3429  {
3430  GEOSGeom_destroy(oarea);
3431  _lwt_release_edges(oldedge, 1);
3432  lwerror("Could not compute edge motion area for new edge");
3433  return -1;
3434  }
3435 
3436  // 3.2. bail out if any node is in one and not the other
3437  oareap = GEOSPrepare( oarea );
3438  nareap = GEOSPrepare( narea );
3439  for (i=0; i<numnodes; ++i)
3440  {
3441  LWT_ISO_NODE *n = &(nodes[i]);
3442  GEOSGeometry *ngg;
3443  int ocont, ncont;
3444  size_t sz;
3445  char *wkt;
3446  if ( n->node_id == oldedge->start_node ) continue;
3447  if ( n->node_id == oldedge->end_node ) continue;
3448  ngg = LWGEOM2GEOS( lwpoint_as_lwgeom(n->geom) , 0);
3449  ocont = GEOSPreparedContains( oareap, ngg );
3450  ncont = GEOSPreparedContains( nareap, ngg );
3451  GEOSGeom_destroy(ngg);
3452  if (ocont == 2 || ncont == 2)
3453  {
3454  _lwt_release_nodes(nodes, numnodes);
3455  GEOSPreparedGeom_destroy(oareap);
3456  GEOSGeom_destroy(oarea);
3457  GEOSPreparedGeom_destroy(nareap);
3458  GEOSGeom_destroy(narea);
3459  lwerror("GEOS exception on PreparedContains: %s", lwgeom_geos_errmsg);
3460  return -1;
3461  }
3462  if (ocont != ncont)
3463  {
3464  GEOSPreparedGeom_destroy(oareap);
3465  GEOSGeom_destroy(oarea);
3466  GEOSPreparedGeom_destroy(nareap);
3467  GEOSGeom_destroy(narea);
3468  wkt = lwgeom_to_wkt(lwpoint_as_lwgeom(n->geom), WKT_ISO, 15, &sz);
3469  _lwt_release_nodes(nodes, numnodes);
3470  lwerror("Edge motion collision at %s", wkt);
3471  lwfree(wkt); /* would not necessarely reach this point */
3472  return -1;
3473  }
3474  }
3475  GEOSPreparedGeom_destroy(oareap);
3476  GEOSGeom_destroy(oarea);
3477  GEOSPreparedGeom_destroy(nareap);
3478  GEOSGeom_destroy(narea);
3479  }}
3480  if ( numnodes ) _lwt_release_nodes(nodes, numnodes);
3481 
3482  LWDEBUG(1, "nodes containment check passed");
3483 
3484  /*
3485  * Check edge adjacency before
3486  * TODO: can be optimized to gather azimuths of all edge ends once
3487  */
3488 
3489  edgeend span_pre, epan_pre;
3490  /* initialize span_pre.myaz and epan_pre.myaz with existing edge */
3491  i = _lwt_InitEdgeEndByLine(&span_pre, &epan_pre,
3492  oldedge->geom, &p1, &p2);
3493  if ( i ) return -1; /* lwerror should have been raised */
3494  _lwt_FindAdjacentEdges( topo, oldedge->start_node, &span_pre,
3495  isclosed ? &epan_pre : NULL, edge_id );
3496  _lwt_FindAdjacentEdges( topo, oldedge->end_node, &epan_pre,
3497  isclosed ? &span_pre : NULL, edge_id );
3498 
3499  LWDEBUGF(1, "edges adjacent to old edge are %" LWTFMT_ELEMID
3500  " and %" LWTFMT_ELEMID " (first point), %" LWTFMT_ELEMID
3501  " and %" LWTFMT_ELEMID " (last point)",
3502  span_pre.nextCW, span_pre.nextCCW,
3503  epan_pre.nextCW, epan_pre.nextCCW);
3504 
3505  /* update edge geometry */
3506  newedge.edge_id = edge_id;
3507  newedge.geom = geom;
3508  i = lwt_be_updateEdgesById(topo, &newedge, 1, LWT_COL_EDGE_GEOM);
3509  if ( i == -1 )
3510  {
3511  _lwt_release_edges(oldedge, 1);
3512  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3513  return -1;
3514  }
3515  if ( ! i )
3516  {
3517  _lwt_release_edges(oldedge, 1);
3518  lwerror("Unexpected error: %d edges updated when expecting 1", i);
3519  return -1;
3520  }
3521 
3522  /*
3523  * Check edge adjacency after
3524  */
3525  edgeend span_post, epan_post;
3526  i = _lwt_InitEdgeEndByLine(&span_post, &epan_post, geom, &p1, &p2);
3527  if ( i ) return -1; /* lwerror should have been raised */
3528  /* initialize epan_post.myaz and epan_post.myaz */
3529  i = _lwt_InitEdgeEndByLine(&span_post, &epan_post,
3530  geom, &p1, &p2);
3531  if ( i ) return -1; /* lwerror should have been raised */
3532  _lwt_FindAdjacentEdges( topo, oldedge->start_node, &span_post,
3533  isclosed ? &epan_post : NULL, edge_id );
3534  _lwt_FindAdjacentEdges( topo, oldedge->end_node, &epan_post,
3535  isclosed ? &span_post : NULL, edge_id );
3536 
3537  LWDEBUGF(1, "edges adjacent to new edge are %" LWTFMT_ELEMID
3538  " and %" LWTFMT_ELEMID " (first point), %" LWTFMT_ELEMID
3539  " and %" LWTFMT_ELEMID " (last point)",
3540  span_pre.nextCW, span_pre.nextCCW,
3541  epan_pre.nextCW, epan_pre.nextCCW);
3542 
3543 
3544  /* Bail out if next CW or CCW edge on start node changed */
3545  if ( span_pre.nextCW != span_post.nextCW ||
3546  span_pre.nextCCW != span_post.nextCCW )
3547  {{
3548  LWT_ELEMID nid = oldedge->start_node;
3549  _lwt_release_edges(oldedge, 1);
3550  lwerror("Edge changed disposition around start node %"
3551  LWTFMT_ELEMID, nid);
3552  return -1;
3553  }}
3554 
3555  /* Bail out if next CW or CCW edge on end node changed */
3556  if ( epan_pre.nextCW != epan_post.nextCW ||
3557  epan_pre.nextCCW != epan_post.nextCCW )
3558  {{
3559  LWT_ELEMID nid = oldedge->end_node;
3560  _lwt_release_edges(oldedge, 1);
3561  lwerror("Edge changed disposition around end node %"
3562  LWTFMT_ELEMID, nid);
3563  return -1;
3564  }}
3565 
3566  /*
3567  -- Update faces MBR of left and right faces
3568  -- TODO: think about ways to optimize this part, like see if
3569  -- the old edge geometry partecipated in the definition
3570  -- of the current MBR (for shrinking) or the new edge MBR
3571  -- would be larger than the old face MBR...
3572  --
3573  */
3574  int facestoupdate = 0;
3575  LWT_ISO_FACE faces[2];
3576  LWGEOM *nface1 = NULL;
3577  LWGEOM *nface2 = NULL;
3578  if ( oldedge->face_left != 0 )
3579  {
3580  nface1 = lwt_GetFaceGeometry(topo, oldedge->face_left);
3581  if ( ! nface1 )
3582  {
3583  lwerror("lwt_ChangeEdgeGeom could not construct face %"
3584  PRId64 ", on the left of edge %" PRId64,
3585  oldedge->face_left, edge_id);
3586  return -1;
3587  }
3588 #if 0
3589  {
3590  size_t sz;
3591  char *wkt = lwgeom_to_wkt(nface1, WKT_EXTENDED, 2, &sz);
3592  LWDEBUGF(1, "new geometry of face left (%d): %s", (int)oldedge->face_left, wkt);
3593  lwfree(wkt);
3594  }
3595 #endif
3596  lwgeom_add_bbox(nface1);
3597  faces[facestoupdate].face_id = oldedge->face_left;
3598  /* ownership left to nface */
3599  faces[facestoupdate++].mbr = nface1->bbox;
3600  }
3601  if ( oldedge->face_right != 0
3602  /* no need to update twice the same face.. */
3603  && oldedge->face_right != oldedge->face_left )
3604  {
3605  nface2 = lwt_GetFaceGeometry(topo, oldedge->face_right);
3606  if ( ! nface2 )
3607  {
3608  lwerror("lwt_ChangeEdgeGeom could not construct face %"
3609  PRId64 ", on the right of edge %" PRId64,
3610  oldedge->face_right, edge_id);
3611  return -1;
3612  }
3613 #if 0
3614  {
3615  size_t sz;
3616  char *wkt = lwgeom_to_wkt(nface2, WKT_EXTENDED, 2, &sz);
3617  LWDEBUGF(1, "new geometry of face right (%d): %s", (int)oldedge->face_right, wkt);
3618  lwfree(wkt);
3619  }
3620 #endif
3621  lwgeom_add_bbox(nface2);
3622  faces[facestoupdate].face_id = oldedge->face_right;
3623  faces[facestoupdate++].mbr = nface2->bbox; /* ownership left to nface */
3624  }
3625  LWDEBUGF(1, "%d faces to update", facestoupdate);
3626  if ( facestoupdate )
3627  {
3628  i = lwt_be_updateFacesById( topo, &(faces[0]), facestoupdate );
3629  if ( i != facestoupdate )
3630  {
3631  if ( nface1 ) lwgeom_free(nface1);
3632  if ( nface2 ) lwgeom_free(nface2);
3633  _lwt_release_edges(oldedge, 1);
3634  if ( i == -1 )
3635  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3636  else
3637  lwerror("Unexpected error: %d faces found when expecting 1", i);
3638  return -1;
3639  }
3640  }
3641  if ( nface1 ) lwgeom_free(nface1);
3642  if ( nface2 ) lwgeom_free(nface2);
3643 
3644  LWDEBUG(1, "all done, cleaning up edges");
3645 
3646  _lwt_release_edges(oldedge, 1);
3647  return 0; /* success */
3648 }
3649 
3650 /* Only return CONTAINING_FACE in the node object */
3651 static LWT_ISO_NODE *
3653 {
3654  LWT_ISO_NODE *node;
3655  int n = 1;
3656 
3657  node = lwt_be_getNodeById( topo, &nid, &n, LWT_COL_NODE_CONTAINING_FACE );
3658  if ( n < 0 ) {
3659  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3660  return 0;
3661  }
3662  if ( n < 1 ) {
3663  lwerror("SQL/MM Spatial exception - non-existent node");
3664  return 0;
3665  }
3666  if ( node->containing_face == -1 )
3667  {
3668  lwfree(node);
3669  lwerror("SQL/MM Spatial exception - not isolated node");
3670  return 0;
3671  }
3672 
3673  return node;
3674 }
3675 
3676 int
3678 {
3679  LWT_ISO_NODE *node;
3680  int ret;
3681 
3682  node = _lwt_GetIsoNode( topo, nid );
3683  if ( ! node ) return -1;
3684 
3685  if ( lwt_be_ExistsCoincidentNode(topo, pt) )
3686  {
3687  lwfree(node);
3688  lwerror("SQL/MM Spatial exception - coincident node");
3689  return -1;
3690  }
3691 
3692  if ( lwt_be_ExistsEdgeIntersectingPoint(topo, pt) )
3693  {
3694  lwfree(node);
3695  lwerror("SQL/MM Spatial exception - edge crosses node.");
3696  return -1;
3697  }
3698 
3699  /* TODO: check that the new point is in the same containing face !
3700  * See https://trac.osgeo.org/postgis/ticket/3232
3701  */
3702 
3703  node->node_id = nid;
3704  node->geom = pt;
3705  ret = lwt_be_updateNodesById(topo, node, 1,
3707  if ( ret == -1 ) {
3708  lwfree(node);
3709  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3710  return -1;
3711  }
3712 
3713  lwfree(node);
3714  return 0;
3715 }
3716 
3717 int
3719 {
3720  LWT_ISO_NODE *node;
3721  int n = 1;
3722 
3723  node = _lwt_GetIsoNode( topo, nid );
3724  if ( ! node ) return -1;
3725 
3726  n = lwt_be_deleteNodesById( topo, &nid, n );
3727  if ( n == -1 )
3728  {
3729  lwfree(node);
3730  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3731  return -1;
3732  }
3733  if ( n != 1 )
3734  {
3735  lwfree(node);
3736  lwerror("Unexpected error: %d nodes deleted when expecting 1", n);
3737  return -1;
3738  }
3739 
3740  /* TODO: notify to caller about node being removed ?
3741  * See https://trac.osgeo.org/postgis/ticket/3231
3742  */
3743 
3744  lwfree(node);
3745  return 0; /* success */
3746 }
3747 
3748 int
3750 {
3751  LWT_ISO_EDGE deledge;
3752  LWT_ISO_EDGE *edge;
3753  LWT_ELEMID nid[2];
3754  LWT_ISO_NODE upd_node[2];
3755  LWT_ELEMID containing_face;
3756  int n = 1;
3757  int i;
3758 
3759  edge = lwt_be_getEdgeById( topo, &id, &n, LWT_COL_EDGE_START_NODE|
3763  if ( ! edge )
3764  {
3765  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3766  return -1;
3767  }
3768  if ( ! n )
3769  {
3770  lwerror("SQL/MM Spatial exception - non-existent edge");
3771  return -1;
3772  }
3773  if ( n > 1 )
3774  {
3775  lwfree(edge);
3776  lwerror("Corrupted topology: more than a single edge have id %"
3777  LWTFMT_ELEMID, id);
3778  return -1;
3779  }
3780 
3781  if ( edge[0].face_left != edge[0].face_right )
3782  {
3783  lwfree(edge);
3784  lwerror("SQL/MM Spatial exception - not isolated edge");
3785  return -1;
3786  }
3787  containing_face = edge[0].face_left;
3788 
3789  nid[0] = edge[0].start_node;
3790  nid[1] = edge[0].end_node;
3791  lwfree(edge);
3792 
3793  n = 2;
3794  edge = lwt_be_getEdgeByNode( topo, nid, &n, LWT_COL_EDGE_EDGE_ID );
3795  if ( n == -1 )
3796  {
3797  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3798  return -1;
3799  }
3800  for ( i=0; i<n; ++i )
3801  {
3802  if ( edge[i].edge_id == id ) continue;
3803  lwfree(edge);
3804  lwerror("SQL/MM Spatial exception - not isolated edge");
3805  return -1;
3806  }
3807  if ( edge ) lwfree(edge);
3808 
3809  deledge.edge_id = id;
3810  n = lwt_be_deleteEdges( topo, &deledge, LWT_COL_EDGE_EDGE_ID );
3811  if ( n == -1 )
3812  {
3813  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3814  return -1;
3815  }
3816  if ( n != 1 )
3817  {
3818  lwerror("Unexpected error: %d edges deleted when expecting 1", n);
3819  return -1;
3820  }
3821 
3822  upd_node[0].node_id = nid[0];
3823  upd_node[0].containing_face = containing_face;
3824  n = 1;
3825  if ( nid[1] != nid[0] ) {
3826  upd_node[1].node_id = nid[1];
3827  upd_node[1].containing_face = containing_face;
3828  ++n;
3829  }
3830  n = lwt_be_updateNodesById(topo, upd_node, n,
3832  if ( n == -1 )
3833  {
3834  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3835  return -1;
3836  }
3837 
3838  /* TODO: notify to caller about edge being removed ?
3839  * See https://trac.osgeo.org/postgis/ticket/3248
3840  */
3841 
3842  return 0; /* success */
3843 }
3844 
3845 /* Used by _lwt_RemEdge to update edge face ref on healing
3846  *
3847  * @param of old face id (never 0 as you cannot remove face 0)
3848  * @param nf new face id
3849  * @return 0 on success, -1 on backend error
3850  */
3851 static int
3853 {
3854  LWT_ISO_EDGE sel_edge, upd_edge;
3855  int ret;
3856 
3857  assert( of != 0 );
3858 
3859  /* Update face_left for all edges still referencing old face */
3860  sel_edge.face_left = of;
3861  upd_edge.face_left = nf;
3862  ret = lwt_be_updateEdges(topo, &sel_edge, LWT_COL_EDGE_FACE_LEFT,
3863  &upd_edge, LWT_COL_EDGE_FACE_LEFT,
3864  NULL, 0);
3865  if ( ret == -1 ) return -1;
3866 
3867  /* Update face_right for all edges still referencing old face */
3868  sel_edge.face_right = of;
3869  upd_edge.face_right = nf;
3870  ret = lwt_be_updateEdges(topo, &sel_edge, LWT_COL_EDGE_FACE_RIGHT,
3871  &upd_edge, LWT_COL_EDGE_FACE_RIGHT,
3872  NULL, 0);
3873  if ( ret == -1 ) return -1;
3874 
3875  return 0;
3876 }
3877 
3878 /* Used by _lwt_RemEdge to update node face ref on healing
3879  *
3880  * @param of old face id (never 0 as you cannot remove face 0)
3881  * @param nf new face id
3882  * @return 0 on success, -1 on backend error
3883  */
3884 static int
3886 {
3887  LWT_ISO_NODE sel, upd;
3888  int ret;
3889 
3890  assert( of != 0 );
3891 
3892  /* Update face_left for all edges still referencing old face */
3893  sel.containing_face = of;
3894  upd.containing_face = nf;
3897  NULL, 0);
3898  if ( ret == -1 ) return -1;
3899 
3900  return 0;
3901 }
3902 
3903 /* Used by lwt_RemEdgeModFace and lwt_RemEdgeNewFaces
3904  *
3905  * Returns -1 on error, identifier of the face that takes up the space
3906  * previously occupied by the removed edge if modFace is 1, identifier of
3907  * the created face (0 if none) if modFace is 0.
3908  */
3909 static LWT_ELEMID
3910 _lwt_RemEdge( LWT_TOPOLOGY* topo, LWT_ELEMID edge_id, int modFace )
3911 {
3912  int i, nedges, nfaces, fields;
3913  LWT_ISO_EDGE *edge = NULL;
3914  LWT_ISO_EDGE *upd_edge = NULL;
3915  LWT_ISO_EDGE upd_edge_left[2];
3916  int nedge_left = 0;
3917  LWT_ISO_EDGE upd_edge_right[2];
3918  int nedge_right = 0;
3919  LWT_ISO_NODE upd_node[2];
3920  int nnode = 0;
3921  LWT_ISO_FACE *faces = NULL;
3922  LWT_ISO_FACE newface;
3923  LWT_ELEMID node_ids[2];
3924  LWT_ELEMID face_ids[2];
3925  int fnode_edges = 0; /* number of edges on the first node (excluded
3926  * the one being removed ) */
3927  int lnode_edges = 0; /* number of edges on the last node (excluded
3928  * the one being removed ) */
3929 
3930  newface.face_id = 0;
3931 
3932  i = 1;
3933  edge = lwt_be_getEdgeById(topo, &edge_id, &i, LWT_COL_EDGE_ALL);
3934  if ( ! edge )
3935  {
3936  LWDEBUGF(1, "lwt_be_getEdgeById returned NULL and set i=%d", i);
3937  if ( i == -1 )
3938  {
3939  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3940  return -1;
3941  }
3942  else if ( i == 0 )
3943  {
3944  lwerror("SQL/MM Spatial exception - non-existent edge %"
3945  LWTFMT_ELEMID, edge_id);
3946  return -1;
3947  }
3948  else
3949  {
3950  lwerror("Backend coding error: getEdgeById callback returned NULL "
3951  "but numelements output parameter has value %d "
3952  "(expected 0 or 1)", i);
3953  return -1;
3954  }
3955  }
3956 
3957  if ( ! lwt_be_checkTopoGeomRemEdge(topo, edge_id,
3958  edge->face_left, edge->face_right) )
3959  {
3961  return -1;
3962  }
3963 
3964  LWDEBUG(1, "Updating next_{right,left}_face of ring edges...");
3965 
3966  /* Update edge linking */
3967 
3968  nedges = 0;
3969  node_ids[nedges++] = edge->start_node;
3970  if ( edge->end_node != edge->start_node )
3971  {
3972  node_ids[nedges++] = edge->end_node;
3973  }
3977  upd_edge = lwt_be_getEdgeByNode( topo, &(node_ids[0]), &nedges, fields );
3978  if ( nedges == -1 ) {
3979  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
3980  return -1;
3981  }
3982  nedge_left = nedge_right = 0;
3983  for ( i=0; i<nedges; ++i )
3984  {
3985  LWT_ISO_EDGE *e = &(upd_edge[i]);
3986  if ( e->edge_id == edge_id ) continue;
3987  if ( e->start_node == edge->start_node || e->end_node == edge->start_node )
3988  {
3989  ++fnode_edges;
3990  }
3991  if ( e->start_node == edge->end_node || e->end_node == edge->end_node )
3992  {
3993  ++lnode_edges;
3994  }
3995  if ( e->next_left == -edge_id )
3996  {
3997  upd_edge_left[nedge_left].edge_id = e->edge_id;
3998  upd_edge_left[nedge_left++].next_left =
3999  edge->next_left != edge_id ? edge->next_left : edge->next_right;
4000  }
4001  else if ( e->next_left == edge_id )
4002  {
4003  upd_edge_left[nedge_left].edge_id = e->edge_id;
4004  upd_edge_left[nedge_left++].next_left =
4005  edge->next_right != -edge_id ? edge->next_right : edge->next_left;
4006  }
4007 
4008  if ( e->next_right == -edge_id )
4009  {
4010  upd_edge_right[nedge_right].edge_id = e->edge_id;
4011  upd_edge_right[nedge_right++].next_right =
4012  edge->next_left != edge_id ? edge->next_left : edge->next_right;
4013  }
4014  else if ( e->next_right == edge_id )
4015  {
4016  upd_edge_right[nedge_right].edge_id = e->edge_id;
4017  upd_edge_right[nedge_right++].next_right =
4018  edge->next_right != -edge_id ? edge->next_right : edge->next_left;
4019  }
4020  }
4021 
4022  if ( nedge_left )
4023  {
4024  LWDEBUGF(1, "updating %d 'next_left' edges", nedge_left);
4025  /* update edges in upd_edge_left set next_left */
4026  i = lwt_be_updateEdgesById(topo, &(upd_edge_left[0]), nedge_left,
4028  if ( i == -1 )
4029  {
4030  _lwt_release_edges(edge, 1);
4031  lwfree(upd_edge);
4032  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4033  return -1;
4034  }
4035  }
4036  if ( nedge_right )
4037  {
4038  LWDEBUGF(1, "updating %d 'next_right' edges", nedge_right);
4039  /* update edges in upd_edge_right set next_right */
4040  i = lwt_be_updateEdgesById(topo, &(upd_edge_right[0]), nedge_right,
4042  if ( i == -1 )
4043  {
4044  _lwt_release_edges(edge, 1);
4045  lwfree(upd_edge);
4046  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4047  return -1;
4048  }
4049  }
4050  LWDEBUGF(1, "releasing %d updateable edges in %p", nedges, upd_edge);
4051  lwfree(upd_edge);
4052 
4053  /* Id of face that will take up all the space previously
4054  * taken by left and right faces of the edge */
4055  LWT_ELEMID floodface;
4056 
4057  /* Find floodface, and update its mbr if != 0 */
4058  if ( edge->face_left == edge->face_right )
4059  {
4060  floodface = edge->face_right;
4061  }
4062  else
4063  {
4064  /* Two faces healed */
4065  if ( edge->face_left == 0 || edge->face_right == 0 )
4066  {
4067  floodface = 0;
4068  LWDEBUG(1, "floodface is universe");
4069  }
4070  else
4071  {
4072  /* we choose right face as the face that will remain
4073  * to be symmetric with ST_AddEdgeModFace */
4074  floodface = edge->face_right;
4075  LWDEBUGF(1, "floodface is %" LWTFMT_ELEMID, floodface);
4076  /* update mbr of floodface as union of mbr of both faces */
4077  face_ids[0] = edge->face_left;
4078  face_ids[1] = edge->face_right;
4079  nfaces = 2;
4080  fields = LWT_COL_FACE_ALL;
4081  faces = lwt_be_getFaceById(topo, face_ids, &nfaces, fields);
4082  if ( nfaces == -1 ) {
4083  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4084  return -1;
4085  }
4086  GBOX *box1=NULL;
4087  GBOX *box2=NULL;
4088  for ( i=0; i<nfaces; ++i )
4089  {
4090  if ( faces[i].face_id == edge->face_left )
4091  {
4092  if ( ! box1 ) box1 = faces[i].mbr;
4093  else
4094  {
4095  i = edge->face_left;
4096  _lwt_release_edges(edge, 1);
4097  _lwt_release_faces(faces, nfaces);
4098  lwerror("corrupted topology: more than 1 face have face_id=%"
4099  LWTFMT_ELEMID, i);
4100  return -1;
4101  }
4102  }
4103  else if ( faces[i].face_id == edge->face_right )
4104  {
4105  if ( ! box2 ) box2 = faces[i].mbr;
4106  else
4107  {
4108  i = edge->face_right;
4109  _lwt_release_edges(edge, 1);
4110  _lwt_release_faces(faces, nfaces);
4111  lwerror("corrupted topology: more than 1 face have face_id=%"
4112  LWTFMT_ELEMID, i);
4113  return -1;
4114  }
4115  }
4116  else
4117  {
4118  i = faces[i].face_id;
4119  _lwt_release_edges(edge, 1);
4120  _lwt_release_faces(faces, nfaces);
4121  lwerror("Backend coding error: getFaceById returned face "
4122  "with non-requested id %" LWTFMT_ELEMID, i);
4123  return -1;
4124  }
4125  }
4126  if ( ! box1 ) {
4127  i = edge->face_left;
4128  _lwt_release_edges(edge, 1);
4129  _lwt_release_faces(faces, nfaces);
4130  lwerror("corrupted topology: no face have face_id=%"
4131  LWTFMT_ELEMID " (left face for edge %"
4132  LWTFMT_ELEMID ")", i, edge_id);
4133  return -1;
4134  }
4135  if ( ! box2 ) {
4136  i = edge->face_right;
4137  _lwt_release_edges(edge, 1);
4138  _lwt_release_faces(faces, nfaces);
4139  lwerror("corrupted topology: no face have face_id=%"
4140  LWTFMT_ELEMID " (right face for edge %"
4141  LWTFMT_ELEMID ")", i, edge_id);
4142  return -1;
4143  }
4144  gbox_merge(box2, box1); /* box1 is now the union of the two */
4145  newface.mbr = box1;
4146  if ( modFace )
4147  {
4148  newface.face_id = floodface;
4149  i = lwt_be_updateFacesById( topo, &newface, 1 );
4150  _lwt_release_faces(faces, 2);
4151  if ( i == -1 )
4152  {
4153  _lwt_release_edges(edge, 1);
4154  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4155  return -1;
4156  }
4157  if ( i != 1 )
4158  {
4159  _lwt_release_edges(edge, 1);
4160  lwerror("Unexpected error: %d faces updated when expecting 1", i);
4161  return -1;
4162  }
4163  }
4164  else
4165  {
4166  /* New face replaces the old two faces */
4167  newface.face_id = -1;
4168  i = lwt_be_insertFaces( topo, &newface, 1 );
4169  _lwt_release_faces(faces, 2);
4170  if ( i == -1 )
4171  {
4172  _lwt_release_edges(edge, 1);
4173  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4174  return -1;
4175  }
4176  if ( i != 1 )
4177  {
4178  _lwt_release_edges(edge, 1);
4179  lwerror("Unexpected error: %d faces inserted when expecting 1", i);
4180  return -1;
4181  }
4182  floodface = newface.face_id;
4183  }
4184  }
4185 
4186  /* Update face references for edges and nodes still referencing
4187  * the removed face(s) */
4188 
4189  if ( edge->face_left != floodface )
4190  {
4191  if ( -1 == _lwt_UpdateEdgeFaceRef(topo, edge->face_left, floodface) )
4192  {
4193  _lwt_release_edges(edge, 1);
4194  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4195  return -1;
4196  }
4197  if ( -1 == _lwt_UpdateNodeFaceRef(topo, edge->face_left, floodface) )
4198  {
4199  _lwt_release_edges(edge, 1);
4200  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4201  return -1;
4202  }
4203  }
4204 
4205  if ( edge->face_right != floodface )
4206  {
4207  if ( -1 == _lwt_UpdateEdgeFaceRef(topo, edge->face_right, floodface) )
4208  {
4209  _lwt_release_edges(edge, 1);
4210  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4211  return -1;
4212  }
4213  if ( -1 == _lwt_UpdateNodeFaceRef(topo, edge->face_right, floodface) )
4214  {
4215  _lwt_release_edges(edge, 1);
4216  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4217  return -1;
4218  }
4219  }
4220 
4221  /* Update topogeoms on heal */
4222  if ( ! lwt_be_updateTopoGeomFaceHeal(topo,
4223  edge->face_right, edge->face_left,
4224  floodface) )
4225  {
4226  _lwt_release_edges(edge, 1);
4228  return -1;
4229  }
4230  } /* two faces healed */
4231 
4232  /* Delete the edge */
4233  i = lwt_be_deleteEdges(topo, edge, LWT_COL_EDGE_EDGE_ID);
4234  if ( i == -1 ) {
4235  _lwt_release_edges(edge, 1);
4236  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4237  return -1;
4238  }
4239 
4240  /* If any of the edge nodes remained isolated, set
4241  * containing_face = floodface
4242  */
4243  if ( ! fnode_edges )
4244  {
4245  upd_node[nnode].node_id = edge->start_node;
4246  upd_node[nnode].containing_face = floodface;
4247  ++nnode;
4248  }
4249  if ( edge->end_node != edge->start_node && ! lnode_edges )
4250  {
4251  upd_node[nnode].node_id = edge->end_node;
4252  upd_node[nnode].containing_face = floodface;
4253  ++nnode;
4254  }
4255  if ( nnode )
4256  {
4257  i = lwt_be_updateNodesById(topo, upd_node, nnode,
4259  if ( i == -1 ) {
4260  _lwt_release_edges(edge, 1);
4261  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4262  return -1;
4263  }
4264  }
4265 
4266  if ( edge->face_left != edge->face_right )
4267  /* or there'd be no face to remove */
4268  {
4269  LWT_ELEMID ids[2];
4270  int nids = 0;
4271  if ( edge->face_right != floodface )
4272  ids[nids++] = edge->face_right;
4273  if ( edge->face_left != floodface )
4274  ids[nids++] = edge->face_left;
4275  i = lwt_be_deleteFacesById(topo, ids, nids);
4276  if ( i == -1 ) {
4277  _lwt_release_edges(edge, 1);
4278  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4279  return -1;
4280  }
4281  }
4282 
4283  _lwt_release_edges(edge, 1);
4284  return modFace ? floodface : newface.face_id;
4285 }
4286 
4287 LWT_ELEMID
4289 {
4290  return _lwt_RemEdge( topo, edge_id, 1 );
4291 }
4292 
4293 LWT_ELEMID
4295 {
4296  return _lwt_RemEdge( topo, edge_id, 0 );
4297 }
4298 
4299 static LWT_ELEMID
4301  int modEdge )
4302 {
4303  LWT_ELEMID ids[2];
4304  LWT_ELEMID commonnode = -1;
4305  int caseno = 0;
4306  LWT_ISO_EDGE *node_edges;
4307  int num_node_edges;
4308  LWT_ISO_EDGE *edges;
4309  LWT_ISO_EDGE *e1 = NULL;
4310  LWT_ISO_EDGE *e2 = NULL;
4311  LWT_ISO_EDGE newedge, updedge, seledge;
4312  int nedges, i;
4313  int e1freenode;
4314  int e2sign, e2freenode;
4315  POINTARRAY *pa;
4316  char buf[256];
4317  char *ptr;
4318  size_t bufleft = 256;
4319 
4320  ptr = buf;
4321 
4322  /* NOT IN THE SPECS: see if the same edge is given twice.. */
4323  if ( eid1 == eid2 )
4324  {
4325  lwerror("Cannot heal edge %" LWTFMT_ELEMID
4326  " with itself, try with another", eid1);
4327  return -1;
4328  }
4329  ids[0] = eid1;
4330  ids[1] = eid2;
4331  nedges = 2;
4332  edges = lwt_be_getEdgeById(topo, ids, &nedges, LWT_COL_EDGE_ALL);
4333  if ( nedges == -1 )
4334  {
4335  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4336  return -1;
4337  }
4338  for ( i=0; i<nedges; ++i )
4339  {
4340  if ( edges[i].edge_id == eid1 ) {
4341  if ( e1 ) {
4342  _lwt_release_edges(edges, nedges);
4343  lwerror("Corrupted topology: multiple edges have id %"
4344  LWTFMT_ELEMID, eid1);
4345  return -1;
4346  }
4347  e1 = &(edges[i]);
4348  }
4349  else if ( edges[i].edge_id == eid2 ) {
4350  if ( e2 ) {
4351  _lwt_release_edges(edges, nedges);
4352  lwerror("Corrupted topology: multiple edges have id %"
4353  LWTFMT_ELEMID, eid2);
4354  return -1;
4355  }
4356  e2 = &(edges[i]);
4357  }
4358  }
4359  if ( ! e1 )
4360  {
4361  if ( edges ) _lwt_release_edges(edges, nedges);
4362  lwerror("SQL/MM Spatial exception - non-existent edge %"
4363  LWTFMT_ELEMID, eid1);
4364  return -1;
4365  }
4366  if ( ! e2 )
4367  {
4368  if ( edges ) _lwt_release_edges(edges, nedges);
4369  lwerror("SQL/MM Spatial exception - non-existent edge %"
4370  LWTFMT_ELEMID, eid2);
4371  return -1;
4372  }
4373 
4374  /* NOT IN THE SPECS: See if any of the two edges are closed. */
4375  if ( e1->start_node == e1->end_node )
4376  {
4377  _lwt_release_edges(edges, nedges);
4378  lwerror("Edge %" LWTFMT_ELEMID " is closed, cannot heal to edge %"
4379  LWTFMT_ELEMID, eid1, eid2);
4380  return -1;
4381  }
4382  if ( e2->start_node == e2->end_node )
4383  {
4384  _lwt_release_edges(edges, nedges);
4385  lwerror("Edge %" LWTFMT_ELEMID " is closed, cannot heal to edge %"
4386  LWTFMT_ELEMID, eid2, eid1);
4387  return -1;
4388  }
4389 
4390  /* Find common node */
4391 
4392  if ( e1->end_node == e2->start_node )
4393  {
4394  commonnode = e1->end_node;
4395  caseno = 1;
4396  }
4397  else if ( e1->end_node == e2->end_node )
4398  {
4399  commonnode = e1->end_node;
4400  caseno = 2;
4401  }
4402  /* Check if any other edge is connected to the common node, if found */
4403  if ( commonnode != -1 )
4404  {
4405  num_node_edges = 1;
4406  node_edges = lwt_be_getEdgeByNode( topo, &commonnode,
4407  &num_node_edges, LWT_COL_EDGE_EDGE_ID );
4408  if ( num_node_edges == -1 ) {
4409  _lwt_release_edges(edges, nedges);
4410  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4411  return -1;
4412  }
4413  for (i=0; i<num_node_edges; ++i)
4414  {
4415  int r;
4416  if ( node_edges[i].edge_id == eid1 ) continue;
4417  if ( node_edges[i].edge_id == eid2 ) continue;
4418  commonnode = -1;
4419  /* append to string, for error message */
4420  if ( bufleft ) {
4421  r = snprintf(ptr, bufleft, "%s%" LWTFMT_ELEMID,
4422  ( ptr==buf ? "" : "," ), node_edges[i].edge_id);
4423  if ( r >= bufleft )
4424  {
4425  bufleft = 0;
4426  buf[252] = '.';
4427  buf[253] = '.';
4428  buf[254] = '.';
4429  buf[255] = '\0';
4430  }
4431  else
4432  {
4433  bufleft -= r;
4434  ptr += r;
4435  }
4436  }
4437  }
4438  lwfree(node_edges);
4439  }
4440 
4441  if ( commonnode == -1 )
4442  {
4443  if ( e1->start_node == e2->start_node )
4444  {
4445  commonnode = e1->start_node;
4446  caseno = 3;
4447  }
4448  else if ( e1->start_node == e2->end_node )
4449  {
4450  commonnode = e1->start_node;
4451  caseno = 4;
4452  }
4453  /* Check if any other edge is connected to the common node, if found */
4454  if ( commonnode != -1 )
4455  {
4456  num_node_edges = 1;
4457  node_edges = lwt_be_getEdgeByNode( topo, &commonnode,
4458  &num_node_edges, LWT_COL_EDGE_EDGE_ID );
4459  if ( num_node_edges == -1 ) {
4460  _lwt_release_edges(edges, nedges);
4461  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4462  return -1;
4463  }
4464  for (i=0; i<num_node_edges; ++i)
4465  {
4466  int r;
4467  if ( node_edges[i].edge_id == eid1 ) continue;
4468  if ( node_edges[i].edge_id == eid2 ) continue;
4469  commonnode = -1;
4470  /* append to string, for error message */
4471  if ( bufleft ) {
4472  r = snprintf(ptr, bufleft, "%s%" LWTFMT_ELEMID,
4473  ( ptr==buf ? "" : "," ), node_edges[i].edge_id);
4474  if ( r >= bufleft )
4475  {
4476  bufleft = 0;
4477  buf[252] = '.';
4478  buf[253] = '.';
4479  buf[254] = '.';
4480  buf[255] = '\0';
4481  }
4482  else
4483  {
4484  bufleft -= r;
4485  ptr += r;
4486  }
4487  }
4488  }
4489  if ( num_node_edges ) lwfree(node_edges);
4490  }
4491  }
4492 
4493  if ( commonnode == -1 )
4494  {
4495  _lwt_release_edges(edges, nedges);
4496  if ( ptr != buf )
4497  {
4498  lwerror("SQL/MM Spatial exception - other edges connected (%s)",
4499  buf);
4500  }
4501  else
4502  {
4503  lwerror("SQL/MM Spatial exception - non-connected edges");
4504  }
4505  return -1;
4506  }
4507 
4508  if ( ! lwt_be_checkTopoGeomRemNode(topo, commonnode,
4509  eid1, eid2 ) )
4510  {
4511  _lwt_release_edges(edges, nedges);
4513  return -1;
4514  }
4515 
4516  /* Construct the geometry of the new edge */
4517  switch (caseno)
4518  {
4519  case 1: /* e1.end = e2.start */
4520  pa = ptarray_clone_deep(e1->geom->points);
4521  //pa = ptarray_merge(pa, e2->geom->points);
4522  ptarray_append_ptarray(pa, e2->geom->points, 0);
4523  newedge.start_node = e1->start_node;
4524  newedge.end_node = e2->end_node;
4525  newedge.next_left = e2->next_left;
4526  newedge.next_right = e1->next_right;
4527  e1freenode = 1;
4528  e2freenode = -1;
4529  e2sign = 1;
4530  break;
4531  case 2: /* e1.end = e2.end */
4532  {
4533  POINTARRAY *pa2;
4534  pa2 = ptarray_clone_deep(e2->geom->points);
4535  ptarray_reverse(pa2);
4536  pa = ptarray_clone_deep(e1->geom->points);
4537  //pa = ptarray_merge(e1->geom->points, pa);
4538  ptarray_append_ptarray(pa, pa2, 0);
4539  ptarray_free(pa2);
4540  newedge.start_node = e1->start_node;
4541  newedge.end_node = e2->start_node;
4542  newedge.next_left = e2->next_right;
4543  newedge.next_right = e1->next_right;
4544  e1freenode = 1;
4545  e2freenode = 1;
4546  e2sign = -1;
4547  break;
4548  }
4549  case 3: /* e1.start = e2.start */
4550  pa = ptarray_clone_deep(e2->geom->points);
4551  ptarray_reverse(pa);
4552  //pa = ptarray_merge(pa, e1->geom->points);
4553  ptarray_append_ptarray(pa, e1->geom->points, 0);
4554  newedge.end_node = e1->end_node;
4555  newedge.start_node = e2->end_node;
4556  newedge.next_left = e1->next_left;
4557  newedge.next_right = e2->next_left;
4558  e1freenode = -1;
4559  e2freenode = -1;
4560  e2sign = -1;
4561  break;
4562  case 4: /* e1.start = e2.end */
4563  pa = ptarray_clone_deep(e2->geom->points);
4564  //pa = ptarray_merge(pa, e1->geom->points);
4565  ptarray_append_ptarray(pa, e1->geom->points, 0);
4566  newedge.end_node = e1->end_node;
4567  newedge.start_node = e2->start_node;
4568  newedge.next_left = e1->next_left;
4569  newedge.next_right = e2->next_right;
4570  e1freenode = -1;
4571  e2freenode = 1;
4572  e2sign = 1;
4573  break;
4574  default:
4575  pa = NULL;
4576  e1freenode = 0;
4577  e2freenode = 0;
4578  e2sign = 0;
4579  _lwt_release_edges(edges, nedges);
4580  lwerror("Coding error: caseno=%d should never happen", caseno);
4581  break;
4582  }
4583  newedge.geom = lwline_construct(topo->srid, NULL, pa);
4584 
4585  if ( modEdge )
4586  {
4587  /* Update data of the first edge */
4588  newedge.edge_id = eid1;
4589  i = lwt_be_updateEdgesById(topo, &newedge, 1,
4595  if ( i == -1 )
4596  {
4597  lwline_free(newedge.geom);
4598  _lwt_release_edges(edges, nedges);
4599  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4600  return -1;
4601  }
4602  else if ( i != 1 )
4603  {
4604  lwline_free(newedge.geom);
4605  if ( edges ) _lwt_release_edges(edges, nedges);
4606  lwerror("Unexpected error: %d edges updated when expecting 1", i);
4607  return -1;
4608  }
4609  }
4610  else
4611  {
4612  /* Add new edge */
4613  newedge.edge_id = -1;
4614  newedge.face_left = e1->face_left;
4615  newedge.face_right = e1->face_right;
4616  i = lwt_be_insertEdges(topo, &newedge, 1);
4617  if ( i == -1 ) {
4618  lwline_free(newedge.geom);
4619  _lwt_release_edges(edges, nedges);
4620  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4621  return -1;
4622  } else if ( i == 0 ) {
4623  lwline_free(newedge.geom);
4624  _lwt_release_edges(edges, nedges);
4625  lwerror("Insertion of split edge failed (no reason)");
4626  return -1;
4627  }
4628  }
4629  lwline_free(newedge.geom);
4630 
4631  /*
4632  -- Update next_left_edge/next_right_edge for
4633  -- any edge having them still pointing at the edge being removed
4634  -- (eid2 only when modEdge, or both otherwise)
4635  --
4636  -- NOTE:
4637  -- e#freenode is 1 when edge# end node was the common node
4638  -- and -1 otherwise. This gives the sign of possibly found references
4639  -- to its "free" (non connected to other edge) endnode.
4640  -- e2sign is -1 if edge1 direction is opposite to edge2 direction,
4641  -- or 1 otherwise.
4642  --
4643  */
4644 
4645  /* update edges connected to e2's boundary from their end node */
4646  seledge.next_left = e2freenode * eid2;
4647  updedge.next_left = e2freenode * newedge.edge_id * e2sign;
4648  i = lwt_be_updateEdges(topo, &seledge, LWT_COL_EDGE_NEXT_LEFT,
4649  &updedge, LWT_COL_EDGE_NEXT_LEFT,
4650  NULL, 0);
4651  if ( i == -1 )
4652  {
4653  _lwt_release_edges(edges, nedges);
4654  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4655  return -1;
4656  }
4657 
4658  /* update edges connected to e2's boundary from their start node */
4659  seledge.next_right = e2freenode * eid2;
4660  updedge.next_right = e2freenode * newedge.edge_id * e2sign;
4661  i = lwt_be_updateEdges(topo, &seledge, LWT_COL_EDGE_NEXT_RIGHT,
4662  &updedge, LWT_COL_EDGE_NEXT_RIGHT,
4663  NULL, 0);
4664  if ( i == -1 )
4665  {
4666  _lwt_release_edges(edges, nedges);
4667  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4668  return -1;
4669  }
4670 
4671  if ( ! modEdge )
4672  {
4673  /* update edges connected to e1's boundary from their end node */
4674  seledge.next_left = e1freenode * eid1;
4675  updedge.next_left = e1freenode * newedge.edge_id;
4676  i = lwt_be_updateEdges(topo, &seledge, LWT_COL_EDGE_NEXT_LEFT,
4677  &updedge, LWT_COL_EDGE_NEXT_LEFT,
4678  NULL, 0);
4679  if ( i == -1 )
4680  {
4681  _lwt_release_edges(edges, nedges);
4682  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4683  return -1;
4684  }
4685 
4686  /* update edges connected to e1's boundary from their start node */
4687  seledge.next_right = e1freenode * eid1;
4688  updedge.next_right = e1freenode * newedge.edge_id;
4689  i = lwt_be_updateEdges(topo, &seledge, LWT_COL_EDGE_NEXT_RIGHT,
4690  &updedge, LWT_COL_EDGE_NEXT_RIGHT,
4691  NULL, 0);
4692  if ( i == -1 )
4693  {
4694  _lwt_release_edges(edges, nedges);
4695  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4696  return -1;
4697  }
4698  }
4699 
4700  /* delete the edges (only second on modEdge or both) */
4702  if ( i == -1 )
4703  {
4704  _lwt_release_edges(edges, nedges);
4705  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4706  return -1;
4707  }
4708  if ( ! modEdge ) {
4710  if ( i == -1 )
4711  {
4712  _lwt_release_edges(edges, nedges);
4713  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4714  return -1;
4715  }
4716  }
4717 
4718  _lwt_release_edges(edges, nedges);
4719 
4720  /* delete the common node */
4721  i = lwt_be_deleteNodesById( topo, &commonnode, 1 );
4722  if ( i == -1 )
4723  {
4724  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4725  return -1;
4726  }
4727 
4728  /*
4729  --
4730  -- NOT IN THE SPECS:
4731  -- Drop composition rows involving second
4732  -- edge, as the first edge took its space,
4733  -- and all affected TopoGeom have been previously checked
4734  -- for being composed by both edges.
4735  */
4736  if ( ! lwt_be_updateTopoGeomEdgeHeal(topo,
4737  eid1, eid2, newedge.edge_id) )
4738  {
4740  return -1;
4741  }
4742 
4743  return modEdge ? commonnode : newedge.edge_id;
4744 }
4745 
4746 LWT_ELEMID
4748 {
4749  return _lwt_HealEdges( topo, e1, e2, 1 );
4750 }
4751 
4752 LWT_ELEMID
4754 {
4755  return _lwt_HealEdges( topo, e1, e2, 0 );
4756 }
4757 
4758 LWT_ELEMID
4759 lwt_GetNodeByPoint(LWT_TOPOLOGY *topo, LWPOINT *pt, double tol)
4760 {
4761  LWT_ISO_NODE *elem;
4762  int num;
4763  int flds = LWT_COL_NODE_NODE_ID|LWT_COL_NODE_GEOM; /* geom not needed */
4764  LWT_ELEMID id = 0;
4765  POINT2D qp; /* query point */
4766 
4767  if ( ! getPoint2d_p(pt->point, 0, &qp) )
4768  {
4769  lwerror("Empty query point");
4770  return -1;
4771  }
4772  elem = lwt_be_getNodeWithinDistance2D(topo, pt, tol, &num, flds, 0);
4773  if ( num == -1 )
4774  {
4775  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4776  return -1;
4777  }
4778  else if ( num )
4779  {
4780  if ( num > 1 )
4781  {
4782  _lwt_release_nodes(elem, num);
4783  lwerror("Two or more nodes found");
4784  return -1;
4785  }
4786  id = elem[0].node_id;
4787  _lwt_release_nodes(elem, num);
4788  }
4789 
4790  return id;
4791 }
4792 
4793 LWT_ELEMID
4794 lwt_GetEdgeByPoint(LWT_TOPOLOGY *topo, LWPOINT *pt, double tol)
4795 {
4796  LWT_ISO_EDGE *elem;
4797  int num, i;
4798  int flds = LWT_COL_EDGE_EDGE_ID|LWT_COL_EDGE_GEOM; /* GEOM is not needed */
4799  LWT_ELEMID id = 0;
4800  LWGEOM *qp = lwpoint_as_lwgeom(pt); /* query point */
4801 
4802  if ( lwgeom_is_empty(qp) )
4803  {
4804  lwerror("Empty query point");
4805  return -1;
4806  }
4807  elem = lwt_be_getEdgeWithinDistance2D(topo, pt, tol, &num, flds, 0);
4808  if ( num == -1 )
4809  {
4810  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4811  return -1;
4812  }
4813  for (i=0; i<num;++i)
4814  {
4815  LWT_ISO_EDGE *e = &(elem[i]);
4816 #if 0
4817  LWGEOM* geom;
4818  double dist;
4819 
4820  if ( ! e->geom )
4821  {
4822  _lwt_release_edges(elem, num);
4823  lwnotice("Corrupted topology: edge %" LWTFMT_ELEMID
4824  " has null geometry", e->edge_id);
4825  continue;
4826  }
4827 
4828  /* Should we check for intersection not being on an endpoint
4829  * as documented ? */
4830  geom = lwline_as_lwgeom(e->geom);
4831  dist = lwgeom_mindistance2d_tolerance(geom, qp, tol);
4832  if ( dist > tol ) continue;
4833 #endif
4834 
4835  if ( id )
4836  {
4837  _lwt_release_edges(elem, num);
4838  lwerror("Two or more edges found");
4839  return -1;
4840  }
4841  else id = e->edge_id;
4842  }
4843 
4844  if ( num ) _lwt_release_edges(elem, num);
4845 
4846  return id;
4847 }
4848 
4849 LWT_ELEMID
4850 lwt_GetFaceByPoint(LWT_TOPOLOGY *topo, LWPOINT *pt, double tol)
4851 {
4852  LWT_ELEMID id = 0;
4853  LWT_ISO_EDGE *elem;
4854  int num, i;
4855  int flds = LWT_COL_EDGE_EDGE_ID |
4859  LWGEOM *qp = lwpoint_as_lwgeom(pt);
4860 
4861  id = lwt_be_getFaceContainingPoint(topo, pt);
4862  if ( id == -2 ) {
4863  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4864  return -1;
4865  }
4866 
4867  if ( id > 0 )
4868  {
4869  return id;
4870  }
4871  id = 0; /* or it'll be -1 for not found */
4872 
4873  LWDEBUG(1, "No face properly contains query point,"
4874  " looking for edges");
4875 
4876  /* Not in a face, may be in universe or on edge, let's check
4877  * for distance */
4878  /* NOTE: we never pass a tolerance of 0 to avoid ever using
4879  * ST_Within, which doesn't include endpoints matches */
4880  elem = lwt_be_getEdgeWithinDistance2D(topo, pt, tol?tol:1e-5, &num, flds, 0);
4881  if ( num == -1 )
4882  {
4883  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
4884  return -1;
4885  }
4886  for (i=0; i<num; ++i)
4887  {
4888  LWT_ISO_EDGE *e = &(elem[i]);
4889  LWT_ELEMID eface = 0;
4890  LWGEOM* geom;
4891  double dist;
4892 
4893  if ( ! e->geom )
4894  {
4895  _lwt_release_edges(elem, num);
4896  lwnotice("Corrupted topology: edge %" LWTFMT_ELEMID
4897  " has null geometry", e->edge_id);
4898  continue;
4899  }
4900 
4901  /* don't consider dangling edges */
4902  if ( e->face_left == e->face_right )
4903  {
4904  LWDEBUGF(1, "Edge %" LWTFMT_ELEMID
4905  " is dangling, won't consider it", e->edge_id);
4906  continue;
4907  }
4908 
4909  geom = lwline_as_lwgeom(e->geom);
4910  dist = lwgeom_mindistance2d_tolerance(geom, qp, tol);
4911 
4912  LWDEBUGF(1, "Distance from edge %" LWTFMT_ELEMID
4913  " is %g (tol=%g)", e->edge_id, dist, tol);
4914 
4915  /* we won't consider edges too far */
4916  if ( dist > tol ) continue;
4917  if ( e->face_left == 0 ) {
4918  eface = e->face_right;
4919  }
4920  else if ( e->face_right == 0 ) {
4921  eface = e->face_left;
4922  }
4923  else {
4924  _lwt_release_edges(elem, num);
4925  lwerror("Two or more faces found");
4926  return -1;
4927  }
4928 
4929  if ( id && id != eface )
4930  {
4931  _lwt_release_edges(elem, num);
4932  lwerror("Two or more faces found"
4933 #if 0 /* debugging */
4934  " (%" LWTFMT_ELEMID
4935  " and %" LWTFMT_ELEMID ")", id, eface
4936 #endif
4937  );
4938  return -1;
4939  }
4940  else id = eface;
4941  }
4942  if ( num ) _lwt_release_edges(elem, num);
4943 
4944  return id;
4945 }
4946 
4947 /* Return the smallest delta that can perturbate
4948  * the maximum absolute value of a geometry ordinate
4949  */
4950 static double
4952 {
4953  const GBOX* gbox;
4954  double max;
4955  double ret;
4956 
4957  gbox = lwgeom_get_bbox(g);
4958  if ( ! gbox ) return 0; /* empty */
4959  max = FP_ABS(gbox->xmin);
4960  if ( max < FP_ABS(gbox->xmax) ) max = FP_ABS(gbox->xmax);
4961  if ( max < FP_ABS(gbox->ymin) ) max = FP_ABS(gbox->ymin);
4962  if ( max < FP_ABS(gbox->ymax) ) max = FP_ABS(gbox->ymax);
4963 
4964  ret = 3.6 * pow(10, - ( 15 - log10(max?max:1.0) ) );
4965 
4966  return ret;
4967 }
4968 
4969 #define _LWT_MINTOLERANCE( topo, geom ) ( \
4970  topo->precision ? topo->precision : _lwt_minTolerance(geom) )
4971 
4972 typedef struct scored_pointer_t {
4973  void *ptr;
4974  double score;
4975 } scored_pointer;
4976 
4977 static int
4978 compare_scored_pointer(const void *si1, const void *si2)
4979 {
4980  double a = ((scored_pointer *)si1)->score;
4981  double b = ((scored_pointer *)si2)->score;
4982  if ( a < b )
4983  return -1;
4984  else if ( a > b )
4985  return 1;
4986  else
4987  return 0;
4988 }
4989 
4990 LWT_ELEMID
4991 lwt_AddPoint(LWT_TOPOLOGY* topo, LWPOINT* point, double tol)
4992 {
4993  int num, i;
4994  double mindist = FLT_MAX;
4995  LWT_ISO_NODE *nodes, *nodes2;
4996  LWT_ISO_EDGE *edges, *edges2;
4997  LWGEOM *pt = lwpoint_as_lwgeom(point);
4998  int flds;
4999  LWT_ELEMID id = 0;
5000  scored_pointer *sorted;
5001 
5002  /* Get tolerance, if 0 was given */
5003  if ( ! tol ) tol = _LWT_MINTOLERANCE( topo, pt );
5004 
5005  LWDEBUGG(1, pt, "Adding point");
5006 
5007  /*
5008  -- 1. Check if any existing node is closer than the given precision
5009  -- and if so pick the closest
5010  TODO: use WithinBox2D
5011  */
5013  nodes = lwt_be_getNodeWithinDistance2D(topo, point, tol, &num, flds, 0);
5014  if ( num == -1 )
5015  {
5016  lwerror("Backend error: %s", lwt_be_lastErrorMessage(topo->be_iface));
5017  return -1;
5018  }
5019  if ( num )
5020  {
5021  LWDEBUGF(1, "New point is within %.15g units of %d nodes", tol, num);
5022  /* Order by distance if there are more than a single return */
5023  if ( num > 1 )
5024  {{
5025  sorted= lwalloc(sizeof(scored_pointer)*num);
5026  for (i=0; i<num; ++i)
5027  {
5028  sorted[i].ptr = nodes+i;
5029  sorted[i].score = lwgeom_mindistance2d(lwpoint_as_lwgeom(nodes[i].geom), pt);
5030  LWDEBUGF(1, "Node %" LWTFMT_ELEMID " distance: %.15g",
5031  ((LWT_ISO_NODE*)(sorted[i].ptr))->node_id, sorted[i].score);
5032  }
5033  qsort(sorted, num, sizeof(scored_pointer), compare_scored_pointer);
5034  nodes2 = lwalloc(sizeof(LWT_ISO_NODE)*num);
5035  for (i=0; i<num; ++i)
5036  {
5037  nodes2[i] = *((LWT_ISO_NODE*)sorted[i].ptr);
5038  }
5039  lwfree(sorted);
5040  lwfree(nodes);
5041  nodes = nodes2;
5042  }}
5043 
5044  for ( i=0; i<num; ++i )
5045  {
5046  LWT_ISO_NODE *n = &(nodes[i]);
5047  LWGEOM *g = lwpoint_as_lwgeom(n->geom);
5048  double dist = lwgeom_mindistance2d(g, pt);
5049  /* TODO: move this check in the previous sort scan ... */
5050  if ( dist >= tol ) continue; /* must be closer than tolerated */
5051  if ( ! id || dist < mindist )
5052  {
5053  id = n->node_id;
5054  mindist = dist;
5055  }
5056  }
5057  if ( id )
5058  {
5059  /* found an existing node */
5060  if ( nodes ) _lwt_release_nodes(nodes, num);
5061  return id;
5062  }
5063  }
5064 
5065  initGEOS(lwnotice, lwgeom_geos_error);
5066 
5067  /*
5068  -- 2. Check if any existing edge falls within tolerance
5069  -- and if so split it by a point projected on it
5070  TODO: use WithinBox2D
5071  */