Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

Is there a PostGIS implementation of the Brabyn and Skelly road network travel time model for the New Zealand road network? This model is outlined in:

Brabyn L, Skelly C. Modeling population access to New Zealand public hospitals. Int J Health Geogr. 2002;1(1):3

and

Lauder C, Skelly C, Brabyn L. Developing and validating a road travel time network for cost path analysis. The 13th Annual Colloquium of the Spatial Information Research Centre. 2001.

share|improve this question

1 Answer 1

up vote 1 down vote accepted

I have created this function as a PostGIS implementation of the travel time model. Improvements and suggestions welcomed.

create or replace function public.brabynskellytraveltime(
    /*
    this is an implementation of Brabyn and Skelly's travel time model for New
    Zealand road networks. 1. Brabyn L, Skelly C. Modeling population access to
    New Zealand public hospitals. Int J Health Geogr [Internet]. 2002;1(1):3.
    Available from: http://www.pubmedcentral.nih.gov/articlerender.fcgi
    ?artid=149398&tool=pmcentrez&rendertype=abstract.

    Aspects of this model are also described in 1. Lauder C, Skelly C, Brabyn L.
    Developing and validating a road travel time network for cost path analysis.
    The 13th Annual Colloquium of the Spatial Information Research Centre. 2001.

        This function expects a SQL statement that contains the following fields
    for each road segment:
        sealed::bool - indicates if the road is sealed
      urban::bool - indicates if the road lies within an urban area
      motorway::bool - indicates if the road is a motorway
      lanes::int - indicates the number of lanes the road has
      the_geom::geometry - the line representing each road segment

    Brabyn and Skelly use a fixed minimum segment length of 500m and a threshold
    for bendy and straight roads at 1.02 so these is values are not adjustable
    in this function.

    Note that Brabyn and Skelly use a concept of removing nodes until the
    segment length is 500m or more, rather than re-creating the segments
    to obtain exactly 500m segments.

    Brabyn and Skelly in their literature were not clear about what they did in
    the event that a super segment contained multiple classes. This function
    assumes the lowest common class for the super segment where mixed class
    exists i.e. if the super segment contains any unsealed roads, the whole
    super segment is calculated to be unsealed.

    authored: Jayden MacRae : 7 Oct 2012
  */

  sql text,                     -- sql containing the required structure
  id OUT int,                   -- returned record id
  speed OUT int,            -- returned speed for segment
  time_seconds OUT double precision,    -- returned time to travel segment
  segment OUT geometry
) returns setof record as
$BODY$
declare
  counter int;                                          -- internal counter returned from function as record ID
    current_length double precision;    -- current length of the super segment
  is_motorway bool;                                 -- super segment motorway flag
  is_sealed bool;                                       -- super segment sealed flag
  is_urban bool;                                        -- super segment urban flag
  lane_count int;                                       -- super segment lane count
  linear_length double precision;       -- length of start point to end point of super segment
  point_counter int;                                -- number of points in super segment
  points geometry[];                                -- points that will make up the super segment
  r record;                                                 -- records returned from sql statement
  sinuosity double precision;               -- sinuosity of the super segment
  speed int;                                                -- derived super segment speed
  super_segment geometry;
  time_seconds double precision;        -- time to travel super segment
  temp_point geometry;                          -- stores temporary point for maninuplation
begin
    counter := 0;
    for r in execute sql                        -- iterate through each record in passed sql
  loop


    if array_length(points,1) is null then      -- initialise variables
        current_length := 0;
      point_counter := -1;
      is_motorway := true;
      is_sealed := true;
      is_urban := false;
      lane_count := 2;
    end if;

    point_counter := point_counter + 2; -- increment to increase array size
        points[point_counter - 1] := ST_STARTPOINT(r.the_geom);
    points[point_counter] := ST_ENDPOINT(r.the_geom);
    is_motorway := is_motorway and r.motorway;
    is_sealed := is_sealed and r.sealed;
    is_urban := is_urban or r.urban;
    if r.lanes < lane_count then lane_count = 1; end if;

    /*
        Check the ordering of points to ensure that they are being added to the
      array in the correct order. Checking the distance to a reference point
      will determine whether the order makes sense.
    */

    if point_counter = 2 then
      if ST_DISTANCE(points[0], points[2]) < ST_DISTANCE(points[1], points[2]) then
        temp_point := points[0];        -- switches 1st and 2nd points around
        points[0] := points[1];
        points[1] := temp_point;
      end if;
    end if;

        -- reorder the last two points based on the previous point
    if ST_DISTANCE(points[point_counter],points[point_counter - 2]) < ST_DISTANCE(points[point_counter -1], points[point_counter - 2]) then
      temp_point := points[point_counter];
      points[point_counter] := points[point_counter -1];
      points[point_counter - 1] := temp_point;
    end if;


    current_length:= ST_LENGTH(ST_MAKELINE(points));

    if current_length >= 500 then
        counter := counter + 1;
        super_segment:= ST_MAKELINE(points);
      linear_length := ST_DISTANCE(points[0],points[array_length(points,1)-1]);
      sinuosity := current_length / linear_length;

      if is_sealed and is_urban then speed := 30;
      elseif is_motorway then speed := 80;
      elseif is_sealed and not is_urban and lane_count = 2 and sinuosity < 1.02 then speed := 80;
      elseif is_sealed and not is_urban and lane_count = 2 and sinuosity >= 1.02 then speed := 60;
      elseif is_sealed and not is_urban and lane_count = 1 and sinuosity < 1.02 then speed := 80;
      elseif is_sealed and not is_urban and lane_count = 1 and sinuosity >= 1.02 then speed := 40;
      elseif not is_sealed and sinuosity < 1.02 then speed := 50;
      elseif not is_sealed and sinuosity >= 1.02 then speed := 30;
      else speed := 30;             -- catch all, not in Brabyn Skelly model, but added just in case
      end if;

      points := '{}';                   -- clear the points array
      time_seconds := current_length / ((speed * 1000) / 3600);

            return query select counter, speed, time_seconds, super_segment;

    end if;
  end loop;

    -- deal with any remaining segement that hasn't fitted into the prior 500 m
  -- limit , don't know a better more elegent way to deal with this
    if array_length(points,1) is not null then
      counter := counter + 1;
        super_segment:= ST_MAKELINE(points);
      linear_length := ST_DISTANCE(points[0],points[array_length(points,1)-1]);
      sinuosity := current_length / linear_length;

      if is_sealed and is_urban then speed := 30;
      elseif is_motorway then speed := 80;
      elseif is_sealed and not is_urban and lane_count = 2 and sinuosity < 1.02 then speed := 80;
      elseif is_sealed and not is_urban and lane_count = 2 and sinuosity >= 1.02 then speed := 60;
      elseif is_sealed and not is_urban and lane_count = 1 and sinuosity < 1.02 then speed := 80;
      elseif is_sealed and not is_urban and lane_count = 1 and sinuosity >= 1.02 then speed := 40;
      elseif not is_sealed and sinuosity < 1.02 then speed := 50;
      elseif not is_sealed and sinuosity >= 1.02 then speed := 30;
      else speed := 30;             -- catch all, not in Brabyn Skelly model, but added just in case
      end if;

      time_seconds := current_length / ((speed * 1000) / 3600);

            return query select counter, speed, time_seconds, super_segment;
  end if;

end

$BODY$
language plpgsql;
share|improve this answer

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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