1

Let me start by saying I'm assuming I'm missing something very silly but I can't see it.

I have a table of points and I'm trying to generate a polygon around each point that estimates a maximum of 1000m driving distance from the point.

I'm using the current Ordnance Survey Vector Map District road layer loaded into postgis as multilinestrings, with Postgis 2.2 and Postgres 9.3. All layers are using the OSGB 27700 projection. The point data has a column which contains the nearest road node calculated from the source road data.

When I run the SQL for a single point to generate the list of individual roads using pgrouting driving distance it works perfectly:

SELECT seq, 
       id1 AS node, 
       id2 AS edge, 
       cost, 
       geom
  FROM pgr_drivingdistance(
       'SELECT gid AS id, 
               source, 
               target, 
               st_length(geom) as cost 
          FROM public.os_vd_2014_road',
       1529980, -- The node
       1000,    -- The driving distance
       false, false
                          ) as di
  JOIN os_vd_2014_road_nodes AS pt
  ON di.id1 = pt.id;

(Since this was essentially Anita Graser's code from her examples I expected nothing else).

However, I want to create a general case, ideally returning a polygon displaying the limits of each route but at the least a single geometry for each point.

To do this I created a function that would hopefully do that for a given point, and then call it in a simple SELECT statement:

CREATE OR REPLACE FUNCTION vdsl_hull (integer, integer)
RETURNS geometry AS $$

SELECT St_Collect(geom) AS the_geom
  FROM pgr_drivingdistance(
       'SELECT gid AS id, 
               source, 
               target, 
               st_length(geom) as cost 
          FROM public.os_vd_2014_road',
       $1,    -- The node
       $2,    -- The driving distance
       false, false
                          ) as di
  JOIN os_vd_2014_road_nodes AS pt
  ON di.id1 = pt.id;
$$ LANGUAGE SQL;

SELECT bt_pcp_list.pcp_id, 
       bt_pcp_list.id AS gid,
       vdsl_hull(bt_pcp_list.nearest_node, 1000) AS the_geom
FROM   public.bt_pcp_list
WHERE  district_code = 'SMBRD';

This "appears" to work in that it generates output for each point, a valid "gid" value, and a column labelled "the_geom" but the geometry column is empty. The next stage was to surround the ST_Collect with ST_ConcaveHull() but clearly that's a stage too far for now.

I'm assuming the error in my thinking is something to do with the resulting geometry from the ST_Collect but could be completely wrong.

Any help would be very welcome!

1

I think the trouble is that function pgr_drivingdistance is returning:

set of pgr_costResult[]:
seq: row sequence
id1: node ID
id2: edge ID
cost: cost to get to this node ID

You don't have any geometry here...

You have to join this result again with os_vd_2014_road table to reach geometry and then st_collect it and st_concavehull

  • Thank you - it was indeed a silly mistake on my part. I have now linked the result back to the vertices table generated when I created the topology, and created the hull from the points. I'll post the result shortly just for completeness. – Adrian Jun 24 '14 at 7:48
  • You're welcome :) The only people who do nothing do not make mistakes :) – Jendrusk Jun 24 '14 at 8:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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