Chapter 19. Raster Sampling

Table of Contents

Sampling Operator
Raquet / QUADBIN Sampling
Restriction and Predicate Operators
Build and Install
Production Roadmap

The PostGIS raster extension stores gridded coverage data (elevation, temperature, satellite imagery, …) in the raster type. Each raster has one or more bands; a band holds a 2D array of pixel values, a spatial extent, a pixel size, and an optional nodata sentinel value.

MobilityDB's raster sampling operator evaluates a raster band at the instants of a moving-point trajectory and returns a tfloat. Instants whose position falls outside the raster extent or on a nodata pixel are silently dropped. This enables queries such as what elevation profile did each vehicle follow? or which trips entered a temperature threshold band?

To use the operators described here, MobilityDB must be built with -DRASTER=ON. On such a build the generated mobilitydb.control declares requires = 'postgis, postgis_raster', so a single CASCADE creates the full stack:

CREATE EXTENSION mobilitydb CASCADE;
-- NOTICE:  installing required extension "postgis"
-- NOTICE:  installing required extension "postgis_raster"

Sampling Operator

  • Sample a raster band at each instant of a trajectory

    raster_value(raster,tgeompoint,band integer DEFAULT 1) → tfloat

    Evaluates pixel band at every instant of traj using bilinear-nearest sampling (ST_Value with resample=false). Returns an instant-set tfloat whose values are rounded to four decimal places. Returns NULL when no instant intersects the raster.

    -- Build a 3×3 synthetic elevation raster (SRID 4326, 1° pixels)
    WITH rast AS (
      SELECT ST_SetValues(
        ST_AddBand(
          ST_MakeEmptyRaster(3, 3, 0.0, 3.0, 1.0, -1.0, 0.0, 0.0, 4326),
          '32BF'::text, 0.0::float8, NULL::float8),
        1, 1, 1,
        ARRAY[[10.0::float4, 20.0::float4, 30.0::float4],
              [40.0::float4, 50.0::float4, 60.0::float4],
              [70.0::float4, 80.0::float4, 90.0::float4]]) AS r
    )
    SELECT raster_value(r,
      tgeompoint 'SRID=4326;[POINT(0.5 2.5)@2000-01-01 00:00:00+00,
                              POINT(2.5 2.5)@2000-01-02 00:00:00+00,
                              POINT(0.5 0.5)@2000-01-03 00:00:00+00]')::text
    FROM rast;
    -- {10@2000-01-01 00:00:00+00, 30@2000-01-02 00:00:00+00, 70@2000-01-03 00:00:00+00}
    

    Applied to the BerlinMOD trajectory dataset with a terrain elevation raster:

    -- Elevation profile per trip (requires a loaded elevation raster)
    SELECT tripid, raster_value(elev, trip4326) AS elevation
    FROM trips, elevation_raster;
    
    -- Average elevation per trip
    SELECT tripid, avgValue(raster_value(elev, trip4326)) AS mean_elev
    FROM trips, elevation_raster;
    
    -- Trips that crossed terrain above 200 m
    SELECT tripid
    FROM trips, elevation_raster
    WHERE atValues(raster_value(elev, trip4326),
                   floatspan '[200, 9999]') IS NOT NULL;