Разбить линии на непересекающиеся подмножества на основе точек

Учитывая таблицу с геометрией линии и одну или несколько точек, которые привязаны к этой линии в отдельной таблице, я хотел бы разбить каждую линию с одной или несколькими пересекающимися точками в каждом из мест, где линия пересекает точку.

Например, есть линия L с тремя пересекающимися точками A, B и C по порядку вдоль геометрии линии. Я хотел бы вернуть L в виде четырех различных геометрий: от начала L до A, от A до B вдоль L, от B до C вдоль L и от C до конца L.

В прошлом я использовал стройные для этой задачи, которая является проблемой линейных ссылок ( http://sgillies.net/blog/1040/shapely-recipes/ ). Тем не менее, это не может быть осуществимо в этом случае, который имеет много миллионов линий и точек. Вместо этого я ищу решение, использующее PostgreSQL / PostGIS.

Обратите внимание, что точки должны находиться на одной линии. Кроме того, точка может действительно находиться в начале или конце линии, и в этом случае линия не должна быть разделена (если нет других точек, которые не совпадают с начальной или конечной точками той же линии). Линии подмножества должны сохранять свое направление и свои атрибуты, но атрибуты точечных объектов не имеют значения.

Функция PostGIS ST_Split , вероятно, то, что вы хотите.

PostGIS 2.2+ теперь поддерживает геометрию Multi * в ST_Split.

Для более старых версий PostGIS читайте дальше:

Чтобы получить одну строку, разделенную на несколько точек, вы можете использовать что-то вроде этой функции многоточечной обертки plpgsql. Я упростил это до случая «разделить (несколько) линий с (несколькими) точками» ниже:

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;

Затем, чтобы создать многоточечную геометрию для резки, используйте ST_Collect и создайте ее вручную из входных данных:

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

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

Или соберите его из подзапроса:

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

Обновление до PostGIS 2.2 , где ST_Split был расширен для поддержки разделения по многолинейной, многоточечной или (много) полигональной границе.

postgis=# SELECT postgis_version(),
                  ST_AsText(ST_Split('LINESTRING(0 0, 2 0)', 'MULTIPOINT(0 0, 1 0)'));
-[ RECORD 1 ]---+------------------------------------------------------------
postgis_version | 2.2 USE_GEOS=1 USE_PROJ=1 USE_STATS=1
st_astext       | GEOMETRYCOLLECTION(LINESTRING(1 0,2 0),LINESTRING(0 0,1 0))

У меня нет полного ответа для вас, но ST_Line_Locate_Point принимает линию и точку в качестве аргументов и возвращает число от 0 до 1, представляющее расстояние вдоль линии до позиции, ближайшей к этой точке.

ST_Line_Substring принимает в качестве аргументов строку и два числа, каждое от 0 до 1. Числа представляют позиции на линии в виде дробных расстояний. Функция возвращает отрезок, который проходит между этими двумя позициями.

Работая с этими двумя функциями, вы сможете достичь того, чего хотите.

Меня уже дважды просили об этом, извините за задержку. Это вряд ли можно считать кратким решением; Я написал это, когда чуть дальше по кривой обучения, чем сейчас. Любые советы приветствуются, даже стилистические.

--Inputs:
--walkingNetwork = Line features representing edges pedestrians can walk on
--stops = Bus stops
--NOTE: stops.geom is already constrained to be coincident with line features
--from walkingNetwork. They may be on a vertex or between two vertices.

--This series of queries returns a version of walkingNetwork, with edges split
--into separate features where they intersect stops.

CREATE TABLE tmp_lineswithstops AS (
    WITH subq AS (
        SELECT
        ST_Line_Locate_Point(
            roads.geom,
            ST_ClosestPoint(roads.geom, stops.geom)
        ) AS LR,
        rank() OVER (
            PARTITION BY roads.gid
            ORDER BY ST_Line_Locate_Point(
                roads.geom,
                ST_ClosestPoint(roads.geom, stops.geom)
            )
        ) AS LRRank,
        ST_ClosestPoint(roads.geom, stops.geom),
        roads.*
        FROM walkingNetwork AS roads
        LEFT OUTER JOIN stops
        ON ST_Distance(roads.geom, stops.geom) < 0.01
        WHERE ST_Equals(ST_StartPoint(roads.geom), stops.geom) IS false
        AND ST_Equals(ST_EndPoint(roads.geom), stops.geom) IS false
        ORDER BY gid, LRRank
    )
    SELECT * FROM subq
);

-- Calculate the interior edges with a join
--If the match is null, calculate the line to the end
CREATE TABLE tmp_testsplit AS (
    SELECT
    l1.gid,
    l1.geom,
    l1.lr AS LR1,
    l1.st_closestpoint AS LR1geom,
    l1.lrrank AS lr1rank,
    l2.lr AS LR2,
    l2.st_closestpoint AS LR2geom,
    l2.lrrank AS lr2rank,
    CASE WHEN l2.lrrank IS NULL -- When the point is the last along the line
        THEN ST_Line_Substring(l1.geom, l1.lr, 1) --get the substring line to the end
        ELSE ST_Line_Substring(l1.geom, l1.lr, l2.lr) --get the substring between the two points
    END AS sublinegeom
    FROM tmp_lineswithstops AS l1
    LEFT OUTER JOIN tmp_lineswithstops AS l2
    ON l1.gid = l2.gid
    AND l2.lrrank = (l1.lrrank + 1)
);

--Calculate the start to first stop edge
INSERT INTO tmp_testsplit (gid, geom, lr1, lr1geom, lr1rank, lr2, lr2geom, lr2rank, sublinegeom)
SELECT gid, geom,
0 as lr1,
ST_StartPoint(geom) as lr1geom,
0 as lr1rank,
lr AS lr2,
st_closestpoint AS lr2geom,
lrrank AS lr2rank,
ST_Line_Substring(l1.geom, 0, lr) AS sublinegeom --Start to point
FROM tmp_lineswithstops AS l1
WHERE l1.lrrank = 1;

--Now match back to the original road features, both modified and unmodified
CREATE TABLE walkingNetwork_split AS (
    SELECT
    roadssplit.sublinegeom,
    roadssplit.gid AS sgid, --split-gid
    roads.*
    FROM tmp_testsplit AS roadssplit
    JOIN walkingNetwork AS r
    ON r.gid = roadssplit.gid
    RIGHT OUTER JOIN walkingNetwork AS roads --Original edges with null if unchanged, original edges with split geom otherwise
    ON roads.gid = roadssplit.gid
);

--Now update the necessary columns, and drop the temporary columns
--You'll probably need to work on your own length and cost functions
--Here I assume it's valid to just multiply the old cost by the fraction of
--the length the now-split line represents of the non-split line
UPDATE walkingNetwork_split
SET geom = sublinegeom,
lengthz = lengthz*(ST_Length(sublinegeom)/ST_Length(geom)),
walk_seconds_ft = walk_seconds_ft*(ST_Length(sublinegeom)/ST_Length(geom)),
walk_seconds_tf = walk_seconds_tf*(ST_Length(sublinegeom)/ST_Length(geom))
WHERE sublinegeom IS NOT NULL
AND ST_Length(sublinegeom) > 0;
ALTER TABLE walkingNetwork_split
DROP COLUMN sublinegeom,
DROP COLUMN sgid;

--Drop intermediate tables
--You probably could use actual temporary tables;
--I prefer to have a sanity check at each stage
DROP TABLE IF EXISTS tmp_testsplit;
DROP TABLE IF EXISTS tmp_lineswithstops;

--Assign the edges a new unique id, so we can use this as source/target columns in pgRouting
ALTER TABLE walkingNetwork_split
DROP COLUMN IF EXISTS fid;
ALTER TABLE walkingNetwork_split
ADD COLUMN fid INTEGER;
CREATE SEQUENCE roads_seq;
UPDATE walkingNetwork_split
SET fid = nextval('roads_seq');
ALTER TABLE walkingNetwork_split
ADD PRIMARY KEY ("fid");

Я хочу расширить ответы выше с точки зрения начинающего. В этом сценарии у вас есть ряд точек, и вы смотрите, чтобы использовать их как «лезвие», чтобы разрезать линии на сегменты. Весь этот пример предполагает, что вы сначала привязали свои точки к линии и что точки имеют атрибут уникального идентификатора из их привязанной линии. Я использую «column_id» для представления уникального идентификатора строки.

Во-первых , вы хотите сгруппировать свои точки в несколько точек, когда на линию выпадает более одного лезвия. В противном случае функция split_line_multipoint действует как функция ST_Split, что не является желаемым результатом.

CREATE TABLE multple_terminal_lines AS
SELECT ST_Multi(ST_Union(the_geom)) as the_geom, a.matched_alid
FROM    point_table a
        INNER JOIN
        (
            SELECT  column_id
            FROM    point_table
            GROUP   BY column_id
            HAVING  COUNT(*) > 1
        ) b ON a.column_id = b.column_id
GROUP BY a.column_id;

Затем вы хотите разделить свою сеть на основе этих многоточечных.

CREATE TABLE split_multi AS
SELECT (ST_Dump(split_line_multipoint(ST_Snap(a.the_geometry, b.the_geom, 0.00001),b.the_geom))).geom as the_geom
FROM line_table a
JOIN multple_terminal_lines b 
ON a.column_id = b.column_id;

Повторите шаги 1 и 2 со своими линиями, которые имеют только одну точку пересечения. Для этого необходимо обновить код с шага 1 до «HAVING COUNT (*) = 1». Переименуйте таблицы соответственно.

Затем создайте дубликат таблицы строк и удалите записи с точками на них.

CREATE TABLE line_dup AS
SELECT * FROM line_table;
-- Delete shared entries
DELETE FROM line_dup
WHERE column_id in (SELECT DISTINCT column_id FROM split_single) OR column_id in (SELECT DISTINCT column_id FROM split_multi) ;

Наконец , объедините три таблицы вместе, используя UNION ALL:

CREATE TABLE line_filtered AS 
SELECT the_geom
FROM split_single
UNION ALL 
SELECT the_geom
FROM split_multi
UNION ALL 
SELECT the_geom
FROM line_dup;

BAM!