Calculate all points in network that are x distance from a starting point in PostGIS with pgrouting

Question

I am trying to write an SQL that calculates all points at a given distance from a point within a network (most likely the wrong term) where the network is stored in a PostgreSQL/PostGIS database as geography objects.

To try to simplify the problem I have created a network of MultilineStrings (as I get multilinestring from the database):

"MULTILINESTRING((-3395658.95218449 9999999.63159823,-3395629.87976969 9999999.12378312))"
"MULTILINESTRING((-3395640.58204213 9999999.31072238,-3395605.60319134 10000047.4596158))"
"MULTILINESTRING((-3395629.87976969 9999999.12378312,-3395601.369858 9999960.94086569))"
"MULTILINESTRING((-3395654.02194139 9999999.54548045,-3395654.02194139 9999977.87419904))"
"MULTILINESTRING((-3395654.02194139 9999977.87419904,-3395568.2969413 9999977.87419904))"
"MULTILINESTRING((-3395591.32900369 9999977.87419904,-3395591.58027466 10000017.8262824))"
"MULTILINESTRING((-3395591.46379143 9999999.30544906,-3395591.46379143 9999999.30544906,-3395618.83235802 9999999.30544906))"
"MULTILINESTRING((-3395608.51010251 9999999.30544906,-3395608.51360801 10000016.5033657))"
"MULTILINESTRING((-3395591.46379143 9999999.30544906,-3395567.93313922 9999999.33852198))"
"MULTILINESTRING((-3395582.44960782 9999999.31830406,-3395582.22063923 9999985.58018863))"
"MULTILINESTRING((-3395658.95218449 9999999.63159823,-3395684.61438934 9999974.46768862))"
"MULTILINESTRING((-3395637.01606056 9999977.94918962,-3395637.25397262 9999992.45935531))"
"MULTILINESTRING((-3395625.67040625 10000019.9754834,-3395672.70813933 10000020.5051887))"

with a starting point @ -3395649.31,9999999.42 trying to find distances 50m away.

This should give a layout like below, where the green point is the starting point and the red points are the ones I would like to calculate:

Note: if while following a path and I hit a dead-end I would do a 180 and keep going till I reach the distance and at intersection split and follow each new path.

i.e:

Some errors in the data I have found so far, Not all lines intersect where they should be an intersection, so I also need to work out a way to extend lines to make them intersect when following a path (normally the end/start of the line in within 500mm).

I have converted my data set to topology and noded with pgrouting (@geozelot pointed me to pgrouting), what I'm still missing is how to add a starting point to the existing noded data to do the calculations from?

SQL for noded data set if it helps:

CREATE TABLE public.gis_test_noded (
    id bigint NOT NULL,
    old_id integer,
    sub_id integer,
    source bigint,
    target bigint,
    cost double precision DEFAULT 0,
    geom public.geometry(LineString,32755)
);
COPY public.gis_test_noded (id, old_id, sub_id, source, target, cost, geom) FROM stdin;
1   1   1   1   2   4.930995161664557   0102000020F37F0000020000006E2EE17925E849C17E0D36F4CF126341BBF9CE0223E849C1699374F1CF126341
2   1   2   2   3   13.441949393534069  0102000020F37F000002000000BBF9CE0223E849C1699374F1CF126341455B804A1CE849C11070F1E9CF126341
3   1   3   3   4   10.703904974912364  0102000020F37F000002000000455B804A1CE849C11070F1E9CF1263410B4B9CF016E849C10708F6E3CF126341
4   2   1   3   22  25.48298626673298   0102000020F37F000002000000465B804A1CE849C11070F1E9CF12634111075FCD14E849C123B6AE7DD2126341
5   2   2   22  17  34.03033945258891   0102000020F37F00000200000011075FCD14E849C123B6AE7DD2126341B45F35CD0AE849C1272CB5EED5126341
6   3   1   4   5   26.519541477588916  0102000020F37F0000020000000B4B9CF016E849C10708F6E3CF1263413189B7010FE849C14470F93BCD126341
7   3   2   5   6   21.132848264454513  0102000020F37F0000020000003189B7010FE849C14470F93BCD126341CD8157AF08E849C15D921B1ECB126341
8   5   1   7   20  17.005880832206458  0102000020F37F000002000000BBF9CE0223E849C14470F93BCD126341BC450E821AE849C14470F93BCD126341
9   5   2   20  5   23.002646987792104  0102000020F37F000002000000BC450E821AE849C14470F93BCD1263413189B7010FE849C14470F93BCD126341
10  5   3   5   8   22.6844098768197    0102000020F37F0000020000003189B7010FE849C14470F93BCD12634106CB1CAA03E849C14470F93BCD126341
11  5   4   8   9   23.0320623931475    0102000020F37F00000200000006CB1CAA03E849C14470F93BCD1263412D2C0226F8E749C14470F93BCD126341
12  6   1   8   10  21.431673879337882  0102000020F37F00000200000006CB1CAA03E849C14470F93BCD12634181845DBB03E849C11D3DC6E9CF126341
13  6   3   10  11  18.521199648420286  0102000020F37F00000200000081845DBB03E849C11D3DC6E9CF126341A67046CA03E849C1D2E7703AD2126341
14  7   1   10  12  17.046311081852764  0102000020F37F00000200000081845DBB03E849C11D3DC6E9CF126341040A4B410CE849C11D3DC6E9CF126341
15  7   2   12  13  10.322255508508533  0102000020F37F000002000000040A4B410CE849C11D3DC6E9CF12634127B58A6A11E849C11D3DC6E9CF126341
16  9   1   10  14  9.014195599909053   0102000020F37F00000200000081845DBB03E849C11D3DC6E9CF126341DFA58C39FFE749C153072EEACF126341
17  9   2   14  15  14.516479855393287  0102000020F37F000002000000DFA58C39FFE749C153072EEACF1263411B1B71F7F7E749C10C2CD5EACF126341
19  10  2   14  16  13.739838008188064  0102000020F37F000002000000DFA58C39FFE749C153072EEACF12634108E83D1CFFE749C1C0E79032CE126341
20  4   1   2   7   21.671281406655908  0102000020F37F000002000000BBF9CE0223E849C1699374F1CF126341BBF9CE0223E849C14470F93BCD126341
21  8   1   12  18  17.197917041939753  0102000020F37F000002000000040A4B410CE849C11D3DC6E9CF12634143E8BD410CE849C17C921B10D2126341
22  11  1   1   19  35.94121735091404   0102000020F37F0000020000006E2EE17925E849C17E0D36F4CF126341564FA44E32E849C1204EF7CECC126341
23  12  1   20  21  14.51211598932906   0102000020F37F000002000000BC450E821AE849C1E7C25F3ECD126341C92C82A01AE849C1E609B30ECF126341
24  13  1   22  23  47.04071556131955   0102000020F37F00000200000041DFCFD514E849C1F928377FD21263413B4FA45A2CE849C16D812A90D2126341
.

Graph theory algorithms for routing applications are implemented in pgRouting, a PostgreSQL extension utilizing the PostGIS spatial types and functions. — geozelot, May 04 '21 at 06:05
Thanks @geozelot, I have read up a little on prgrouting and have converted the test data to the formats it needs, I believe I should be using pgr_drivingDistance but I'm struggling in starting from my startpoint and node the closest node, is the a way to add the startnode in to the noded dataset? — Yendor, May 06 '21 at 04:38
Look into pgr_withPointsDD of the still proposed, but supported withPoints - Family of functions - note that, in order to refer to the pids of the points_sql as start/end vertices, you'd have to prefix the actual values with - (minus) in the function call! See https://gis.stackexchange.com/questions/377630/pgrouting-shortest-path-starting-from-the-closest-vertex-not-the-closest-node/377641#377641 — geozelot, May 06 '21 at 10:04
...and https://gis.stackexchange.com/questions/379120/how-to-identify-pgrouting-nodes-i-am-interested-after-splitting-road-network-in/379124#379124 where you find a link to a post with yet some more helpful links to make use of an efficient (K)NN query to get the required edge_id and frac parameters. And after you made it work, self-answer you question! — geozelot, May 06 '21 at 10:11

score 2 · Accepted Answer · answered May 10 '21 at 05:02

With using pgr_withPointDD and from: https://stackoverflow.com/questions/51572251/how-to-get-a-road-path-which-can-be-travelled-in-10-minutes-from-a-location I came up with:

WITH
  dd AS (
    SELECT pg.node,
           pg.edge,
           pg.cost
    FROM pgr_withPointsDD(
        'SELECT id, source, target, cost FROM gis_test_noded ORDER BY id',
        'WITH 
                start_vertex AS (
                        SELECT ST_SetSRID(ST_MakePoint(-3395649.31, 9999999.42), 32755) AS pt
                ),
                nearest_line AS (
                        SELECT id, ST_Distance(geom, start_vertex.pt)
                        FROM gis_test_noded, start_vertex
                        ORDER BY 2 ASC LIMIT 1
                ),
                new_vertex AS (
                        SELECT nearest_line.id AS edge_id, ST_LineLocatePoint(gis_test_noded.geom, start_vertex.pt) as fraction
                        FROM gis_test_noded, start_vertex, nearest_line
                        WHERE gis_test_noded.id = nearest_line.id
                )
        SELECT edge_id, fraction FROM new_vertex',
        -1, 50,
        directed := false) AS pg
  ),
  dd_edgs AS (
    SELECT edg.id,
           edg.geom
    FROM gis_test_noded AS edg
    JOIN dd AS d1
      ON edg.source = d1.node
    JOIN dd AS d2
      ON edg.target = d2.node
  ),
  dd_ext AS (
    SELECT edg.id,
             CASE
               WHEN dd.node = edg.source
                  THEN ST_LineInterpolatePoint(edg.geom, (50 - dd.cost) / 50)
                  ELSE ST_LineInterpolatePoint(edg.geom, 1 - ((50 - dd.cost) / 50))
             END AS geom
    FROM dd
    JOIN gis_test_noded AS edg
      ON dd.node IN (edg.source, edg.target) AND edg.id NOT IN (SELECT id FROM dd_edgs)
  )
SELECT id,
       geom
FROM dd_ext;

which gives me:

"id"	"geom"
5	"0101000020F37F0000603305E60FE849C18543AF2DD4126341"
24	"0101000020F37F0000283E115E20E849C17AEA8687D2126341"
23	"0101000020F37F0000B2A626961AE849C1012FC670CE126341"
7	"0101000020F37F00008990B6090CE849C19541823DCC126341"
10	"0101000020F37F0000F95B28AE09E849C14470F93BCD126341"

Where -1 green dot is the starting point, yellow dots are the calculated endpoints and red dots are where I would have expected to see a yellow. But for now this answers my original question of what type of SQL I would need with the added help of using pgrouting.

Calculate all points in network that are x distance from a starting point in PostGIS with pgrouting

1 Answers1