10

I was checking the best way to cut linestrings by points.

The scenario is: lots of streets, need the segments cut by intersect points, kind of this:

example

I got

  • linestrings (full uncuted by points) table

  • st_intersection points table

I need to get the independent linestring segments cut by the intersection points table.

I'm using PostGIS functions, and found several approaches, but every one of them has giving me some kind of issue.

Here's what I already tested:

1

Line table: 1 row , st_memunion of 1200 lines Point table: 1700 rows (points)

What's bad: really takes really lot of time and memory flush. Can't create multiple tables at same time cause memory just can't handle it. And the result is dirty and messy. Instead of giving me the correct row number, and I need to clean it up later (well explained here Split lines at intersection points)

CREATE TABLE lines_with_messy_result AS (
SELECT
((ST_DUMP(ST_SPLIT(a.geom,b.ix))).geom) as geom
FROM st_union_lines a
INNER JOIN lots_of_points b
ON ST_INTERSECTS(a.geom, b.ix)
);

--then need to clean this up
create table lines_segments_cleaned as (
SELECT DISTINCT ON (ST_AsBinary(geom)) 
geom 
FROM 
lines_with_messy_result
);

source of this way/approach: https://stackoverflow.com/questions/25753348/how-do-i-divide-city-streets-by-intersection-using-postgis


2

Same lines/points table. Still messy results and need to clean this up. Still lots of time to finish query.

--messy table    
Create table messy_lines_segments as 
Select 
row_number() over() as id,
(st_dump(st_split(input.geom, blade.ix))).geom as geom
from st_union_lines input
inner join lots_of_points blade on st_intersects(input.geom, blade.ix);

--then cleaning the messy table
delete from messy_lines_segments a
where exists 
(
select 1 from messy_lines_segments b where a.id != b.id and st_coveredby(b.geom,a.geom)
);

source of this way / approach : Split lines at intersection points


3

I also found this function: https://github.com/Remi-C/PPPP_utilities/blob/master/postgis/rc_Split_Line_By_Points.sql

which has the good thing that it does not leave a messy_result that then I need to clean it up. But you do need st_memunion from both sides (lines table and points table)

It's kind of:

create table osm.calles_cortadas_03segmentos_sanluis as (
SELECT result.geom
FROM 
osm.calles_cortadas_01uniones_sanluis AS line, 
osm.calles_cortadas_00intersecciones_sanluis AS point,
rc_split_line_by_points(
input_line:=line.geom
,input_points:=point.ix
,tolerance:=4
) AS result
);

But its also super-long hours for getting the results. And I also tried with longer tables (10k lines, 14k points) and I just get memory out issues.

I also tried Esri's ArcGIS with also bad results...

So, what's the best way of getting this done with PostGIS geom functions?

I mean, without stepping into topology.

Or what's your best recommendation?

3
  • 1
    buffer points with a very small tolerance, erase the polylines inside the buffers, snap the new endpoints together. . Jul 18, 2015 at 1:42
  • 4
    You should put the final update as an answer, since it answers the question
    – raphael
    Aug 18, 2015 at 19:22
  • It is fine to self-answer questions and can earn you reputation to unlock additional functionality on this site.
    – PolyGeo
    Oct 23, 2015 at 7:05

2 Answers 2

5

Best non-topology solution ever.. it REALLY works fast and easy (believe me i tested lots of ways to do this):

--please add this function:
https://github.com/Remi-C/PPPP_utilities/blob/master/postgis/rc_split_multi.sql

-- please create a universal sequence for unique id multipurpose
CREATE SEQUENCE select_id
  INCREMENT 1
  MINVALUE 1
  MAXVALUE 9223372036854775807
  START 1
  CACHE 1
  CYCLE;



-- lines (non st_unioned linestrings, with gist index) : linetable
-- points (non st_unioned points, with gist index) : pointtable

-- first, lines with cuts (pointtable)
create table segment_lines_00 as (
WITH points_intersecting_lines AS( 
   SELECT lines.id AS lines_id, lines.geom AS line_geom,  ST_Collect(points.ix) AS blade
   FROM 
   linetable as lines, 
   pointtable as points
   WHERE st_dwithin(lines.geom, points.ix, 4) = true
   GROUP BY lines.id, lines.geom
)
SELECT lines_id, rc_Split_multi(line_geom, blade, 4)
FROM points_intersecting_lines
);
CREATE INDEX ON segment_lines_00 USING gist(rc_split_multi);


-- then, segments cutted by points
create table segment_lines_01 as (
select 
(ST_Dump(rc_split_multi)).geom as geom,  nextval('SELECT_id'::regclass) AS gid  from segment_lines_00
);
CREATE INDEX ON segment_lines_01 USING gist(geom);
1

I also add, in addition, my own method, without using st_split():

For each line I check if there are some intersecting points.

If yes my temptable will look like this:

[line_id, fraction]
[1      , 0       ]
[1      , 0.3     ]
[1      , 1       ]
[2      , 0       ]
[2      , 0.88    ]
[2      , 1       ]
...

A table that contains the line ID and the fraction of the line length where a point intersect the line.

Then I combine the fraction by pair to create the new cutted lines

Recquirement:

  • The table "points" contains the endpoints of each line.
  • The points perfectly intersect the lines.

The query:

DROP TABLE IF EXISTS temptable;

CREATE TABLE temptable as 
(
    SELECT ST_LineLocatePoint(a.geom,b.geom) as fraction,a.id, a.geom as shape
    FROM lines a, points b
    WHERE st_intersects(a.geom,b.geom)
    UNION
    SELECT 1 as fraction ,a.id, a.geom as shape -- It avoid to get an incomplete line when the original line is circular
    FROM lines a
    ORDER BY a.id, fraction -- important to keep the same order
);

ALTER TABLE temptable  ADD COLUMN gid SERIAL PRIMARY KEY;

DROP TABLE IF EXISTS LINE_SPLIT;
CREATE TABLE LINE_SPLIT AS
(
SELECT b.fend,a.*, ST_LineSubstring(a.shape,a.fstart,b.fend) as geom
FROM
    (
    SELECT *, fraction as fstart FROM temptable
    LIMIT (SELECT COUNT(*) FROM temptable)-1
    ) a,
    (
    SELECT fraction as fend,id-1 as id FROM temptable
    OFFSET 1
    ) b
WHERE a.id = b.id
AND a.fstart < b.fend
);

ALTER TABLE LINE_SPLIT DROP COLUMN IF EXISTS shape;
DROP TABLE IF EXISTS temptable;
3
  • nice! txs! ill check this one too May 9, 2019 at 16:17
  • @obchardon The second part of the query doesn't work. Table comes up blank.
    – MJM
    Aug 23 at 16:41
  • @ MJM I've mixed id and gid. Of course gid is identical as id.
    – obchardon
    Aug 24 at 12:52

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.