Integrating third party raster data

As of OSRM version 4.8.0, OSRM supports incorporating third-party data into a Lua profile with a native (non-PostGIS) API.

Data Format

As of the initial release, this feature only supports plaintext ASCII grid sources — looking something like

3 4 56 66 6  7
3 4 4  56 20 3
5 6 5  3  14 1

Data available immediately in this format includes, for example, ArcInfo ASCII grids of SRTM data. These files also include six lines of metadata, which are necessary for loading a source file; for example:

ncols         6001
nrows         6001
xllcorner     -125.00041606132
yllcorner     34.999583357538
cellsize      0.00083333333333333
NODATA_value  -9999

These lines should be stripped before loading this source file.

Profile API

To use a raster source in a profile, include a source_function in your profile. When present, this function is called once by a specific lua state that can later be used for per-segment updates. A source_function takes no arguments in its signature and should be used to call sources:load as such:

function source_function ()
  mysource = sources:load(
    "../path/to/raster_source.asc",
    0,    -- longitude min
    0.1,  -- longitude max
    0,    -- latitude min
    0.1,  -- latitude max
    5,    -- number of rows
    4     -- number of cols
  )
end

• sources is a single instance of a SourceContainer class that is bound to this lua state. load is a binding to SourceContainer::loadRasterSource.
• mysource is an int ID that is saved to the global scope of this lua state, meaning it will be available for use in other functions later. If you are loading multiple sources, you can assign them to multiple variables.

mysource is now available for use in a segment_function, a separate function to write into the profile. Its function signature is segment_function(segment) where segment has fields:

• source and target both have lat,lon properties (premultiplied by COORDINATE_PRECISION)
  • It is advisable, then, to either hardcode or read from environment the lat and lon min/max in the source_function — the latter can be done like

    LON_MIN=tonumber(os.getenv("LON_MIN"))
    

    — and then multiplying all of them by the precision constant, bound to the lua state as constants.precision, before exiting the scope of source_function
• distance is a const double in meters
• weight is a segment weight in properties.weight_name units
• duration is a segment duration in seconds

In a segment_function you can query loaded raster sources using a nearest-neighbor query or a bilinear interpolation query, which have signatures like

-- nearest neighbor:
  local sourceData = sources:query(mysource, segment.source.lon, segment.source.lat)
-- bilinear interpolation:
  local sourceData = sources:interpolate(mysource, segment.target.lon, segment.target.lat)

where the signatures both look like (unsigned int source_id, int lon, int lat). They both return an instance of RasterDatum with a member std::int32_t datum. Out-of-bounds queries return a specific int std::numeric_limits<std::int32_t>::max(), which is bound as a member function to the lua state as invalid_data(). So you could use this to find whether a query was out of bounds:

  if sourceData.datum ~= sourceData.invalid_data() then
    -- this is in bounds
  end

…though it is faster to compare source and target to the raster bounds before querying the source.

The complete additions of both functions then might look like so (since [OSRM version 5.6.0] the parameters of segment_function () changed):

api_version = 1
properties.force_split_edges = true

local LON_MIN=tonumber(os.getenv("LON_MIN"))
local LON_MAX=tonumber(os.getenv("LON_MAX"))
local LAT_MIN=tonumber(os.getenv("LAT_MIN"))
local LAT_MAX=tonumber(os.getenv("LAT_MAX"))

function source_function()
    raster_source = sources:load(
        os.getenv("RASTER_SOURCE_PATH"),
        LON_MIN,
        LON_MAX,
        LAT_MIN,
        LAT_MAX,
        tonumber(os.getenv("NROWS")),
        tonumber(os.getenv("NCOLS")))
    LON_MIN = LON_MIN * constants.precision
    LON_MAX = LON_MAX * constants.precision
    LAT_MIN = LAT_MIN * constants.precision
    LAT_MAX = LAT_MAX * constants.precision
end

function segment_function (segment)
  local out_of_bounds = false
  if segment.source.lon < LON_MIN or segment.source.lon > LON_MAX or
     segment.source.lat < LAT_MIN or segment.source.lat > LAT_MAX or
     segment.target.lon < LON_MIN or segment.target.lon > LON_MAX or
     segment.target.lat < LAT_MIN or segment.target.lat > LAT_MAX then
        out_of_bounds = true
  end

  if out_of_bounds == false then
    local sourceData = sources:interpolate(raster_source, segment.source.lon, segment.source.lat)
    local targetData = sources:interpolate(raster_source, segment.target.lon, segment.target.lat)
    local elev_delta = targetData.datum - sourceData.datum

    local slope = 0
    local penalize = 0

    if distance ~= 0 and targetData.datum > 0 and sourceData.datum > 0 then
      slope = elev_delta / segment.distance
    end
    
    -- these types of heuristics are fairly arbitrary and take some trial and error
    if slope > 0.08 then
        penalize = 0.6
    end

    if slope ~= 0 then
      segment.weight = segment.weight * (1 - penalize)
      -- a very rought estimate of duration of the edge so that the time estimate is more accurate
      segment.duration = segment.duration * (1+slope)
    end
  end
end

⚠️ Default behavior during the extraction stage is to merge forward and backward edges into one edge if travel modes, speeds and rates are equal. If the profile modifies in segment_function forward and backward weights or durations independently then the global properties flag force_split_edges must be set to true. Otherwise segment_function will be called only once for a segment that corresponds to the forward edge.