[GIS] Splitting lines at intersections in PostGIS

linear-referencingNetworkpostgissplitting

I have oriented polylines (Montreal's train network) and stops (train stations).

My network looks like this:
Train Network

(closer look):
Train stations

I'd like to split my tracks at each stations and ideally add fields FromStation, ToStation to my new Train tracks table.

I've used ST_Snap to snap my points to my tracks so there's an intersection. There's no guaranty though that there's a node where it intersects.

I tried some variation using ST_Split. No error, but no results either. Just empty geometry in the end.

My data comes from the gtfs, so for now, I have a shape_id, trip_headsign and route_short_name on each line.

If I could have a resulting table like

shape_id fromstation tostation
2_1      VRM1B       SLR3B
2_1      SLR3B       SLR2B
2_1      SLR2B       SLR1B
etc.

(shape_id 2_1 is the blue line west bound), that would be ideal.

I need somehow to use that to create a personnel allocation model for our inspectors and agents on the network. Each segments will be considered a road in Network Analyst with usage restriction to match the train schedule at each train station… Never found anything like that done, but I've been given the task to at least try.

Best Answer

You might need to use the linear referencing functions to get what you want. For each point, use ST_LineLocatePoint to find out what proportion along the line the point is, then for each pairing of points, use ST_LocateBetween to extract the portion of the line between the points. So rather than splitting, you're using linear referencing to extract.

http://postgis.net/workshops/postgis-intro/linear_referencing.html

Related Solutions

PostGIS Splitting – Splitting Lines into Non-Overlapping Subsets Based on Points Using PostGIS

The ST_Split PostGIS function is probably what you want.

PostGIS 2.2+ now supports Multi* geometries in ST_Split.

For older versions of PostGIS, read on:

To get a single line split by multiple points, you could use something like this multipoint wrapper plpgsql function. I've simplified it just to the "split (multi)lines with (multi)points" case below:

DROP FUNCTION IF EXISTS split_line_multipoint(input_geom geometry, blade geometry);
CREATE FUNCTION split_line_multipoint(input_geom geometry, blade geometry)
  RETURNS geometry AS
$BODY$
    -- this function is a wrapper around the function ST_Split 
    -- to allow splitting multilines with multipoints
    --
    DECLARE
        result geometry;
        simple_blade geometry;
        blade_geometry_type text := GeometryType(blade);
        geom_geometry_type text := GeometryType(input_geom);
    BEGIN
        IF blade_geometry_type NOT ILIKE 'MULTI%' THEN
            RETURN ST_Split(input_geom, blade);
        ELSIF blade_geometry_type NOT ILIKE '%POINT' THEN
            RAISE NOTICE 'Need a Point/MultiPoint blade';
            RETURN NULL;
        END IF;

        IF geom_geometry_type NOT ILIKE '%LINESTRING' THEN
            RAISE NOTICE 'Need a LineString/MultiLineString input_geom';
            RETURN NULL;
        END IF;

        result := input_geom;           
        -- Loop on all the points in the blade
        FOR simple_blade IN SELECT (ST_Dump(ST_CollectionExtract(blade, 1))).geom
        LOOP
            -- keep splitting the previous result
            result := ST_CollectionExtract(ST_Split(result, simple_blade), 2);
        END LOOP;
        RETURN result;
    END;
$BODY$
LANGUAGE plpgsql IMMUTABLE;

-- testing
SELECT ST_AsText(split_line_multipoint(geom, blade))
    FROM (
        SELECT ST_GeomFromText('Multilinestring((-3 0, 3 0),(-1 0, 1 0))') AS geom,
        ST_GeomFromText('MULTIPOINT((-0.5 0),(0.5 0))') AS blade
        --ST_GeomFromText('POINT(-0.5 0)') AS blade
    ) AS T;

Then to create a multipoint geometry to cut by, use ST_Collect and either create it manually from inputs:

SELECT ST_AsText(ST_Collect(
  ST_GeomFromText('POINT(1 2)'),
  ST_GeomFromText('POINT(-2 3)')
));

st_astext
----------
MULTIPOINT(1 2,-2 3)

Or aggregate it from a subquery:

SELECT stusps,
  ST_Multi(ST_Collect(f.the_geom)) as singlegeom
FROM (SELECT stusps, (ST_Dump(the_geom)).geom As the_geom
      FROM somestatetable ) As f
GROUP BY stusps

[GIS] Splitting lines with PostGIS and view result in QGIS

There is a plpgsql function is this post that splits a Line using multiple points.

However, it can also be done in a query, using generate_series, ST_Line_Locate_Point and ST_GeometryN to get knife points and ST_Line_Substring to split the line. Whichever you prefer, I leave up to you, as this query is not exactly terse.

 CREATE TABLE sometable AS
 WITH 
 lines as (SELECT ST_GeomFromText('LINESTRING(0 0, 2 2, 10 10, 15 15, 20 20)') as geom),
 points as (SELECT ST_GeomFromText('MultiPoint((0 0), (2 2), (10 10), (15 15), (20 20))') as geom),
 numgeoms (pt_count) as (SELECT st_numgeometries(geom) + 1 as pt_count FROM points),
 knife_points as (
   SELECT x as segment, CASE WHEN x = 0 then 0 
     WHEN x = (select pt_count from numgeoms) THEN 1 ELSE  
     ST_Line_Locate_Point(l.geom, ST_GeometryN(p.geom,x)) END as line_fraction_end, 
     CASE WHEN x = 1 THEN 0  else 
     ST_Line_Locate_Point(l.geom, ST_GeometryN(p.geom,x-1)) END as line_fraction_start 
   FROM points p, lines l, (SELECT generate_series(0, (SELECT pt_count from numgeoms)) as x ) g),
   segments as 
    (SELECT ST_Line_Substring(geom, line_fraction_start, line_fraction_end) as geom, segment 
     FROM knife_points, lines WHERE segment >0 )
   SELECT geom, segment from segments 
   WHERE ST_GeometryType(geom) = 'ST_LineString';

where you would replace the initial lines and points tables with your own. I suspect this query could be simplified quite a bit by using lead or lag, to get calculate each successive pair of knife points, but it does work.

The output of the above is,

 LINESTRING(0 0,2 2)     |       2
 LINESTRING(2 2,10 10)   |       3
 LINESTRING(10 10,15 15) |       4
 LINESTRING(15 15,20 20) |       5

Best Answer

Related Solutions

PostGIS Splitting – Splitting Lines into Non-Overlapping Subsets Based on Points Using PostGIS

[GIS] Splitting lines with PostGIS and view result in QGIS

Related Question