PostGIS 3.7.0dev-r@@SVN_REVISION@@
Loading...
Searching...
No Matches

◆ _lwt_AddLine()

static LWT_ELEMID * _lwt_AddLine ( LWT_TOPOLOGY topo,
LWLINE line,
double  tol,
int *  nedges,
int  handleFaceSplit,
int  maxNewEdges 
)
static

Definition at line 7253 of file lwgeom_topo.c.

7255{
7256 LWGEOM *geomsbuf[1];
7257 LWGEOM **geoms;
7258 uint32_t ngeoms;
7259 LWGEOM *noded, *tmp;
7260 LWCOLLECTION *col;
7261 LWT_ELEMID *ids;
7262 LWT_ISO_EDGE *edges;
7263 LWT_ISO_NODE *nodes;
7264 uint64_t num, numedges = 0, numnodes = 0;
7265 int num_new_edges = 0;
7266 uint64_t i;
7267 GBOX qbox;
7268 int forward;
7269 int input_was_closed = 0;
7270 POINT4D originalStartPoint;
7271
7272 if ( lwline_is_empty(line) )
7273 {
7274 *nedges = 0;
7275 return NULL;
7276 }
7277
7278 if ( lwline_is_closed(line) )
7279 {
7280 input_was_closed = 1;
7281 getPoint4d_p( line->points, 0, &originalStartPoint);
7282 LWDEBUGF(1, "Input line is closed, original point is %g,%g", originalStartPoint.x, originalStartPoint.y);
7283 }
7284
7285 *nedges = -1; /* error condition, by default */
7286
7287 /* Get tolerance, if 0 was given */
7288 if ( tol == -1 ) tol = _LWT_MINTOLERANCE( topo, (LWGEOM*)line );
7289 LWDEBUGF(1, "Working tolerance:%.15g", tol);
7290 LWDEBUGF(1, "Input line has srid=%d", line->srid);
7291
7292 /* Remove consecutive vertices below given tolerance upfront */
7293 if ( tol )
7294 {{
7296 tmp = lwline_as_lwgeom(clean); /* NOTE: might collapse to non-simple */
7297 LWDEBUGG(1, tmp, "Repeated-point removed");
7298 }} else tmp=(LWGEOM*)line;
7299
7300 /* 1. Self-node */
7301 noded = lwgeom_node((LWGEOM*)tmp);
7302 if ( tmp != (LWGEOM*)line ) lwgeom_free(tmp);
7303 if ( ! noded ) return NULL; /* should have called lwerror already */
7304 LWDEBUGG(1, noded, "Noded");
7305
7306 qbox = *lwgeom_get_bbox( lwline_as_lwgeom(line) );
7307 LWDEBUGF(1, "Line BOX is %.15g %.15g, %.15g %.15g", qbox.xmin, qbox.ymin,
7308 qbox.xmax, qbox.ymax);
7309 gbox_expand(&qbox, tol);
7310 LWDEBUGF(1, "BOX expanded by %g is %.15g %.15g, %.15g %.15g",
7311 tol, qbox.xmin, qbox.ymin, qbox.xmax, qbox.ymax);
7312
7313 LWGEOM **nearby = 0;
7314 int nearbyindex = 0;
7315 int nearbycount = 0;
7316
7317 /* 2.0. Find edges falling within tol distance */
7318 edges = lwt_be_getEdgeWithinBox2D( topo, &qbox, &numedges, LWT_COL_EDGE_ALL, 0 );
7319 if (numedges == UINT64_MAX)
7320 {
7321 lwgeom_free(noded);
7323 return NULL;
7324 }
7325 LWDEBUGF(1, "Line has %u points, its bbox intersects %llu edges bboxes",
7326 line->points->npoints, numedges);
7327 if ( numedges )
7328 {{
7329 /* TODO: compute earlier for reuse later ? */
7330 GEOSGeometry *noded_g = LWGEOM2GEOS(noded, 0);
7331 /* collect those whose distance from us is < tol */
7332 nearbycount += numedges;
7333 nearby = lwalloc(numedges * sizeof(LWGEOM *));
7334 for (i=0; i<numedges; ++i)
7335 {{
7336 LW_ON_INTERRUPT(return NULL);
7337 LWT_ISO_EDGE *e = &(edges[i]);
7338 LWGEOM *g = lwline_as_lwgeom(e->geom);
7339 GEOSGeometry *edge_g = LWGEOM2GEOS(g, 0);
7340 LWDEBUGF(2, "Computing distance from edge %" LWTFMT_ELEMID " with %u points", i, e->geom->points->npoints);
7341 double dist;
7342 if ( 0 == GEOSDistanceIndexed(edge_g, noded_g, &dist) ) {
7343 GEOSGeom_destroy(edge_g);
7344 GEOSGeom_destroy(noded_g);
7345 lwgeom_free(noded);
7346 lwerror("GEOSDistanceIndexed error: %s", lwgeom_geos_errmsg);
7347 return NULL;
7348 }
7349 GEOSGeom_destroy(edge_g);
7350 if ( dist && dist >= tol ) continue;
7351 nearby[nearbyindex++] = g;
7352 }}
7353 LWDEBUGF(1, "Found %d edges closer than tolerance (%g)", nearbyindex, tol);
7354 GEOSGeom_destroy(noded_g);
7355 }}
7356 int nearbyedgecount = nearbyindex;
7357
7358 /* 2.1. Find isolated nodes falling within tol distance
7359 *
7360 * TODO: add backend-interface support for only getting isolated nodes
7361 */
7362 nodes = lwt_be_getNodeWithinBox2D( topo, &qbox, &numnodes, LWT_COL_NODE_ALL, 0 );
7363 if (numnodes == UINT64_MAX)
7364 {
7365 lwgeom_free(noded);
7367 return NULL;
7368 }
7369 LWDEBUGF(1, "Line bbox intersects %llu nodes bboxes", numnodes);
7370 if ( numnodes )
7371 {{
7372 /* collect those whose distance from us is < tol */
7373 nearbycount = nearbyedgecount + numnodes;
7374 nearby = nearby ?
7375 lwrealloc(nearby, nearbycount * sizeof(LWGEOM *))
7376 :
7377 lwalloc(nearbycount * sizeof(LWGEOM *))
7378 ;
7379 int nn = 0;
7380 for (i=0; i<numnodes; ++i)
7381 {
7382 LWT_ISO_NODE *n = &(nodes[i]);
7383 if ( n->containing_face == -1 ) continue; /* skip not-isolated nodes */
7384 LWGEOM *g = lwpoint_as_lwgeom(n->geom);
7385 double dist = lwgeom_mindistance2d(g, noded);
7386 /* must be closer than tolerated, unless distance is zero */
7387 if ( dist && dist >= tol )
7388 {
7389 LWDEBUGF(1, "Node %" LWTFMT_ELEMID " is %g units away, we tolerate only %g", n->node_id, dist, tol);
7390 continue;
7391 }
7392 nearby[nearbyindex++] = g;
7393 nn = nn + 1;
7394 }
7395 LWDEBUGF(1, "Found %d isolated nodes closer than tolerance (%g)", nn, tol);
7396 }}
7397 int nearbynodecount = nearbyindex - nearbyedgecount;
7398 nearbycount = nearbyindex;
7399
7400 LWDEBUGF(1, "Number of nearby elements is %d", nearbycount);
7401
7402 /* 2.2. Snap to nearby elements */
7403 if ( nearbycount )
7404 {{
7405 LWCOLLECTION *col;
7406 LWGEOM *elems;
7407
7409 NULL, nearbycount, nearby);
7410 elems = lwcollection_as_lwgeom(col);
7411
7412 LWDEBUGG(1, elems, "Collected nearby elements");
7413
7414 tmp = _lwt_toposnap(noded, elems, tol);
7415 lwgeom_free(noded);
7416 noded = tmp;
7417 LWDEBUGG(1, noded, "Elements-snapped");
7418 if ( input_was_closed )
7419 {{
7420 /* Recompute start point in case it moved */
7421 LWLINE *scrolled = lwgeom_as_lwline(noded);
7422 if (scrolled)
7423 {
7424 getPoint4d_p( scrolled->points, 0, &originalStartPoint);
7425 LWDEBUGF(1, "Closed input line start point after snap %g,%g", originalStartPoint.x, originalStartPoint.y);
7426 }
7427 }}
7428
7429 /* will not release the geoms array */
7431
7432 /*
7433 -- re-node to account for ST_Snap introduced self-intersections
7434 -- See http://trac.osgeo.org/postgis/ticket/1714
7435 -- TODO: consider running UnaryUnion once after all noding
7436 */
7437 tmp = lwgeom_unaryunion(noded);
7438 lwgeom_free(noded);
7439 noded = tmp;
7440 LWDEBUGG(1, noded, "Unary-unioned");
7441 }}
7442
7443 /* 2.3. Node with nearby edges */
7444 if ( nearbyedgecount )
7445 {{
7446 LWCOLLECTION *col;
7447 LWGEOM *iedges; /* just an alias for col */
7448 LWGEOM *diff, *xset;
7449
7450 LWDEBUGF(1, "Line intersects %d edges", nearbyedgecount);
7451
7453 NULL, nearbyedgecount, nearby);
7454 iedges = lwcollection_as_lwgeom(col);
7455 LWDEBUGG(1, iedges, "Collected edges");
7456
7457 LWDEBUGF(1, "Diffing noded, with srid=%d "
7458 "and intersecting edges, with srid=%d",
7459 noded->srid, iedges->srid);
7460 diff = lwgeom_difference(noded, iedges);
7461 LWDEBUGG(1, diff, "Differenced");
7462
7463 LWDEBUGF(1, "Intersecting noded, with srid=%d "
7464 "and intersecting edges, with srid=%d",
7465 noded->srid, iedges->srid);
7466 xset = lwgeom_intersection(noded, iedges);
7467 LWDEBUGG(1, xset, "Intersected");
7468 lwgeom_free(noded);
7469
7470 /* We linemerge here because INTERSECTION, as of GEOS 3.8,
7471 * will result in shared segments being output as multiple
7472 * lines rather than a single line. Example:
7473
7474 INTERSECTION(
7475 'LINESTRING(0 0, 5 0, 8 0, 10 0,12 0)',
7476 'LINESTRING(5 0, 8 0, 10 0)'
7477 )
7478 ==
7479 MULTILINESTRING((5 0,8 0),(8 0,10 0))
7480
7481 * We will re-split in a subsequent step, by splitting
7482 * the final line with pre-existing nodes
7483 */
7484 LWDEBUG(1, "Linemerging intersection");
7485 tmp = lwgeom_linemerge(xset);
7486 LWDEBUGG(1, tmp, "Linemerged");
7487 lwgeom_free(xset);
7488 xset = tmp;
7489
7490 /*
7491 * Here we union the (linemerged) intersection with
7492 * the difference (new lines)
7493 */
7494 LWDEBUG(1, "Unioning difference and (linemerged) intersection");
7495 noded = lwgeom_union(diff, xset);
7496 LWDEBUGG(1, noded, "Diff-Xset Unioned");
7497 lwgeom_free(xset);
7498 lwgeom_free(diff);
7499
7500 /* will not release the geoms array */
7502
7503 if ( input_was_closed )
7504 {{
7505 LWLINE *scrolled = lwgeom_as_lwline(noded);
7506 if (scrolled) {
7507 if ( lwline_is_closed(scrolled) ) {
7508 ptarray_scroll_in_place(scrolled->points, &originalStartPoint);
7509 }
7510 else {
7511 LWDEBUGG(1, lwline_as_lwgeom(scrolled), "Closed input became non closed");
7512 }
7513 }
7514 else {
7515 LWDEBUGG(1, noded, "Diff-Xset Unioned cannot be scrolled");
7516 }
7517 }}
7518
7519
7520 }}
7521
7522
7523 /* 2.4. Split by pre-existing nodes
7524 *
7525 * Pre-existing nodes are isolated nodes AND endpoints
7526 * of intersecting edges
7527 */
7528 if ( nearbyedgecount )
7529 {
7530 nearbycount += nearbyedgecount * 2; /* make space for endpoints */
7531 nearby = lwrealloc(nearby, nearbycount * sizeof(LWGEOM *));
7532 for (int i=0; i<nearbyedgecount; i++)
7533 {
7534 LWLINE *edge = lwgeom_as_lwline(nearby[i]);
7535 LWPOINT *startNode = lwline_get_lwpoint(edge, 0);
7536 LWPOINT *endNode = lwline_get_lwpoint(edge, edge->points->npoints-1);
7537 /* TODO: only add if within distance to noded AND if not duplicated */
7538 nearby[nearbyindex++] = lwpoint_as_lwgeom(startNode);
7539 nearbynodecount++;
7540 nearby[nearbyindex++] = lwpoint_as_lwgeom(endNode);
7541 nearbynodecount++;
7542 }
7543 }
7544 if ( nearbynodecount )
7545 {
7547 NULL, nearbynodecount,
7548 nearby + nearbyedgecount);
7549 LWGEOM *inodes = lwcollection_as_lwgeom(col);
7550 /* TODO: use lwgeom_split of lwgeom_union ... */
7551 tmp = _lwt_split_by_nodes(noded, inodes);
7552 lwgeom_free(noded);
7553 noded = tmp;
7554 LWDEBUGG(1, noded, "Node-split");
7555 /* will not release the geoms array */
7557 }
7558
7559
7560 LWDEBUG(1, "Freeing up nearby elements");
7561
7562 /* TODO: free up endpoints of nearbyedges */
7563 if ( nearby ) lwfree(nearby);
7564 if ( nodes ) _lwt_release_nodes(nodes, numnodes);
7565 if ( edges ) _lwt_release_edges(edges, numedges);
7566
7567 LWDEBUGG(2, noded, "Finally-noded");
7568
7569 /* 3. For each (now-noded) segment, insert an edge */
7570 col = lwgeom_as_lwcollection(noded);
7571 if ( col )
7572 {
7573 LWDEBUG(1, "Noded line was a collection");
7574 geoms = col->geoms;
7575 ngeoms = col->ngeoms;
7576 }
7577 else
7578 {
7579 LWDEBUG(1, "Noded line was a single geom");
7580 geomsbuf[0] = noded;
7581 geoms = geomsbuf;
7582 ngeoms = 1;
7583 }
7584
7585 LWDEBUGF(1, "Line was split into %d edges", ngeoms);
7586
7587 /* TODO: refactor to first add all nodes (re-snapping edges if
7588 * needed) and then check all edges for existing already
7589 * ( so to save a DB scan for each edge to be added )
7590 */
7591 ids = lwalloc(sizeof(LWT_ELEMID)*ngeoms);
7592 num = 0;
7593 for ( i=0; i<ngeoms; ++i )
7594 {
7595 int edgeNewEdges;
7596 LWT_ELEMID id;
7597 LWGEOM *g = geoms[i];
7598 g->srid = noded->srid;
7599
7600#if POSTGIS_DEBUG_LEVEL > 0
7601 {
7602 size_t sz;
7603 char *wkt1 = lwgeom_to_wkt(g, WKT_EXTENDED, 15, &sz);
7604 LWDEBUGF(1, "Component %llu of split line is: %s", i, wkt1);
7605 lwfree(wkt1);
7606 }
7607#endif
7608
7609 forward = -1; /* will be set to either 0 or 1 if the edge already existed */
7610 id = _lwt_AddLineEdge( topo, lwgeom_as_lwline(g), tol, handleFaceSplit, &forward, &edgeNewEdges );
7611 num_new_edges += edgeNewEdges;
7612 /* if forward is still == -1 this was NOT an existing edge ? */
7613 if ( forward == -1 )
7614 {
7615 ++num_new_edges;
7616 }
7617
7618 LWDEBUGF(1, "_lwt_AddLineEdge returned %" LWTFMT_ELEMID
7619 " (forward ? %d), reported to create %d new edges (total new edges: %d)",
7620 id, forward, edgeNewEdges, num_new_edges);
7621 if ( id < 0 )
7622 {
7623 lwgeom_free(noded);
7624 lwfree(ids);
7625 return NULL;
7626 }
7627
7628 if ( maxNewEdges >= 0 && num_new_edges > maxNewEdges )
7629 {
7630 lwgeom_free(noded);
7631 lwfree(ids);
7632 lwerror("Adding line to topology requires creating more edges than the requested limit of %d", maxNewEdges);
7633 return NULL;
7634 }
7635
7636 if ( ! id )
7637 {
7638 LWDEBUGF(1, "Component %llu of split line collapsed", i);
7639 continue;
7640 }
7641
7642 LWDEBUGF(1, "Component %llu of split line is %s edge %" LWTFMT_ELEMID,
7643 i, forward ? "forward" : "backward", id);
7644 ids[num++] = forward ? id : -id; /* TODO: skip duplicates */
7645 }
7646
7647 LWDEBUGG(2, noded, "Noded before free");
7648 lwgeom_free(noded);
7649
7650 /* TODO: XXX remove duplicated ids if not done before */
7651
7652 *nedges = num;
7653 return ids;
7654}
void gbox_expand(GBOX *g, double d)
Move the box minimums down and the maximums up by the distance provided.
Definition gbox.c:97
char lwgeom_geos_errmsg[LWGEOM_GEOS_ERRMSG_MAXSIZE]
GEOSGeometry * LWGEOM2GEOS(const LWGEOM *lwgeom, uint8_t autofix)
LWCOLLECTION * lwcollection_construct(uint8_t type, int32_t srid, GBOX *bbox, uint32_t ngeoms, LWGEOM **geoms)
LWGEOM * lwpoint_as_lwgeom(const LWPOINT *obj)
Definition lwgeom.c:372
LWGEOM * lwgeom_unaryunion(const LWGEOM *geom1)
#define COLLECTIONTYPE
Definition liblwgeom.h:108
void * lwrealloc(void *mem, size_t size)
Definition lwutil.c:242
LWGEOM * lwgeom_node(const LWGEOM *lwgeom_in)
LWGEOM * lwgeom_linemerge(const LWGEOM *geom1)
void lwgeom_free(LWGEOM *geom)
Definition lwgeom.c:1246
LWPOINT * lwline_get_lwpoint(const LWLINE *line, uint32_t where)
Returns freshly allocated LWPOINT that corresponds to the index where.
Definition lwline.c:319
LWCOLLECTION * lwgeom_as_lwcollection(const LWGEOM *lwgeom)
Definition lwgeom.c:261
#define WKT_EXTENDED
Definition liblwgeom.h:2221
char * lwgeom_to_wkt(const LWGEOM *geom, uint8_t variant, int precision, size_t *size_out)
WKT emitter function.
Definition lwout_wkt.c:708
#define MULTIPOINTTYPE
Definition liblwgeom.h:105
void * lwalloc(size_t size)
Definition lwutil.c:227
LWGEOM * lwgeom_intersection(const LWGEOM *geom1, const LWGEOM *geom2)
void lwfree(void *mem)
Definition lwutil.c:248
LWGEOM * lwline_as_lwgeom(const LWLINE *obj)
Definition lwgeom.c:367
LWGEOM * lwgeom_difference(const LWGEOM *geom1, const LWGEOM *geom2)
void lwcollection_release(LWCOLLECTION *lwcollection)
double lwgeom_mindistance2d(const LWGEOM *lw1, const LWGEOM *lw2)
Function initializing min distance calculation.
Definition measures.c:212
int getPoint4d_p(const POINTARRAY *pa, uint32_t n, POINT4D *point)
Definition lwgeom_api.c:125
LWLINE * lwgeom_as_lwline(const LWGEOM *lwgeom)
Definition lwgeom.c:207
LWGEOM * lwgeom_union(const LWGEOM *geom1, const LWGEOM *geom2)
const GBOX * lwgeom_get_bbox(const LWGEOM *lwgeom)
Get a non-empty geometry bounding box, computing and caching it if not already there.
Definition lwgeom.c:771
LWGEOM * lwcollection_as_lwgeom(const LWCOLLECTION *obj)
Definition lwgeom.c:337
LWGEOM * lwline_remove_repeated_points(const LWLINE *in, double tolerance)
Definition lwline.c:449
int lwline_is_empty(const LWLINE *line)
#define LW_ON_INTERRUPT(x)
int lwline_is_closed(const LWLINE *line)
Definition lwline.c:455
int ptarray_scroll_in_place(POINTARRAY *pa, const POINT4D *newbase)
Definition ptarray.c:2337
LWT_INT64 LWT_ELEMID
Identifier of topology element.
#define LWT_COL_EDGE_ALL
#define LWT_COL_NODE_ALL
#define PGTOPO_BE_ERROR()
#define LWTFMT_ELEMID
#define LWDEBUG(level, msg)
Definition lwgeom_log.h:101
#define LWDEBUGF(level, msg,...)
Definition lwgeom_log.h:106
void void lwerror(const char *fmt,...) __attribute__((format(printf
Write a notice out to the error handler.
#define LWDEBUGG(level, geom, msg)
Definition lwgeom_log.h:111
static LWGEOM * _lwt_toposnap(LWGEOM *src, LWGEOM *tgt, double tol)
static void _lwt_release_nodes(LWT_ISO_NODE *nodes, int num_nodes)
static LWGEOM * _lwt_split_by_nodes(const LWGEOM *g, const LWGEOM *nodes)
void _lwt_release_edges(LWT_ISO_EDGE *edges, int num_edges)
#define _LWT_MINTOLERANCE(topo, geom)
static LWT_ISO_NODE * lwt_be_getNodeWithinBox2D(const LWT_TOPOLOGY *topo, const GBOX *box, uint64_t *numelems, int fields, uint64_t limit)
static LWT_ELEMID _lwt_AddLineEdge(LWT_TOPOLOGY *topo, LWLINE *edge, double tol, int handleFaceSplit, int *forward, int *numNewEdges)
LWT_ISO_EDGE * lwt_be_getEdgeWithinBox2D(const LWT_TOPOLOGY *topo, const GBOX *box, uint64_t *numelems, int fields, uint64_t limit)
double ymax
Definition liblwgeom.h:357
double xmax
Definition liblwgeom.h:355
double ymin
Definition liblwgeom.h:356
double xmin
Definition liblwgeom.h:354
uint32_t ngeoms
Definition liblwgeom.h:580
LWGEOM ** geoms
Definition liblwgeom.h:575
int32_t srid
Definition liblwgeom.h:460
POINTARRAY * points
Definition liblwgeom.h:483
int32_t srid
Definition liblwgeom.h:484
LWT_ELEMID node_id
LWT_ELEMID containing_face
LWPOINT * geom
double x
Definition liblwgeom.h:414
double y
Definition liblwgeom.h:414
uint32_t npoints
Definition liblwgeom.h:427

References _lwt_AddLineEdge(), _LWT_MINTOLERANCE, _lwt_release_edges(), _lwt_release_nodes(), _lwt_split_by_nodes(), _lwt_toposnap(), COLLECTIONTYPE, LWT_ISO_NODE::containing_face, gbox_expand(), LWT_ISO_NODE::geom, LWT_ISO_EDGE::geom, LWCOLLECTION::geoms, getPoint4d_p(), LW_ON_INTERRUPT, lwalloc(), lwcollection_as_lwgeom(), lwcollection_construct(), lwcollection_release(), LWDEBUG, LWDEBUGF, LWDEBUGG, lwerror(), lwfree(), LWGEOM2GEOS(), lwgeom_as_lwcollection(), lwgeom_as_lwline(), lwgeom_difference(), lwgeom_free(), lwgeom_geos_errmsg, lwgeom_get_bbox(), lwgeom_intersection(), lwgeom_linemerge(), lwgeom_mindistance2d(), lwgeom_node(), lwgeom_to_wkt(), lwgeom_unaryunion(), lwgeom_union(), lwline_as_lwgeom(), lwline_get_lwpoint(), lwline_is_closed(), lwline_is_empty(), lwline_remove_repeated_points(), lwpoint_as_lwgeom(), lwrealloc(), lwt_be_getEdgeWithinBox2D(), lwt_be_getNodeWithinBox2D(), LWT_COL_EDGE_ALL, LWT_COL_NODE_ALL, LWTFMT_ELEMID, MULTIPOINTTYPE, LWCOLLECTION::ngeoms, LWT_ISO_NODE::node_id, POINTARRAY::npoints, PGTOPO_BE_ERROR, LWLINE::points, ptarray_scroll_in_place(), LWGEOM::srid, LWLINE::srid, LWT_TOPOLOGY_T::srid, WKT_EXTENDED, POINT4D::x, GBOX::xmax, GBOX::xmin, POINT4D::y, GBOX::ymax, and GBOX::ymin.

Referenced by lwt_AddLine(), and lwt_AddLineNoFace().

Here is the call graph for this function:
Here is the caller graph for this function: