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