Spatial Operations

This page provides an overview of spatial operations available in openEO for manipulating and analysing geospatial data cubes. These processes answer the question of where in space the analysis should be applied.

NoteBackend support overview

The buttons on this page are a documentation aid, not a live capability registry. Use the openEO Hub for a cross-backend overview: open Filters, select a process under Processes, and the Hub shows only matching services. The Hub data is crawled and cached.

Filter an area of interest

openEO provides several spatial operations to manipulate the geospatial extent of the datacube. It allows users to load only the relevant spatial subset of the data, reducing unnecessary data transfer and processing. There are two main approaches: filtering by a bounding box or by an irregular geometry, using either of the following processes:

• filter_bbox keeps only pixels within a bounding box.
• filter_spatial keeps only pixels within an irregular GeoJSON geometry.

Filter by bounding box

filter_bbox is typically the first spatial operation applied because it reduces data transfer and processing by focusing only on the area of interest. It is suitable for rectangular regions and is computationally efficient.

Python

import openeo

connection = openeo.connect("openeofed.dataspace.copernicus.eu").authenticate_oidc()
cube = connection.load_collection(
    "SENTINEL2_L2A",
    temporal_extent=["2024-06-01", "2024-06-30"],
    bands=["B04", "B08"],
)

study_area = cube.filter_bbox(
    west=4.30, east=4.55, south=50.80, north=50.98, crs="EPSG:4326"
)

R

library(openeo)

connection <- connect("openeofed.dataspace.copernicus.eu") |> authenticate_oidc()
cube <- connection$load_collection(
    "SENTINEL2_L2A",
    temporal_extent=c("2024-06-01", "2024-06-30"),
    bands=c("B04", "B08")
)

study_area <- cube$filter_bbox(
    west=4.30, east=4.55, south=50.80, north=50.98, crs="EPSG:4326"
)

JavaScript

import openeo from "openeo";

const connection = await openeo.connect("openeofed.dataspace.copernicus.eu").authenticate_oidc();
const cube = connection.load_collection(
    "SENTINEL2_L2A",
    {temporal_extent: ["2024-06-01", "2024-06-30"], bands: ["B04", "B08"]}
);

const study_area = cube.filter_bbox({
 west: 4.30, east: 4.55, south: 50.80, north: 50.98, crs: "EPSG:4326"
});

Filter by irregular geometry

Alternatively, when the area of interest is not rectangular, use filter_spatial, which lets the user pass a GeoJSON geometry to define it. This preserves the original grid but marks pixels outside the geometry as no data. For a quick rectangular crop, prefer filter_bbox.

Python

boundary = {
    "type": "Polygon",
    "coordinates": [[[4.31, 50.81], [4.54, 50.81], [4.54, 50.97],
                     [4.31, 50.97], [4.31, 50.81]]],
}
study_area = cube.filter_spatial(boundary)

R

boundary <- list(
    type = "Polygon",
    coordinates = list(list(c(4.31, 50.81), c(4.54, 50.81), c(4.54, 50.97),
                            c(4.31, 50.97), c(4.31, 50.81)))
)
study_area <- cube$filter_spatial(boundary)

JavaScript

const boundary = {
 type: "Polygon",
 coordinates: [[[4.31, 50.81], [4.54, 50.81], [4.54, 50.97],
 [4.31, 50.97], [4.31, 50.81]]],
};
const study_area = cube.filter_spatial(boundary);

Additionally, if the geometry is stored and hosted externally, you can fetch it using the standard load_url process and then pass it to filter_spatial. This enables dynamic, flexible spatial filtering using external GeoJSON resources.

Python

import openeo

connection = openeo.connect("openeofed.dataspace.copernicus.eu").authenticate_oidc()
cube = connection.load_collection(
    "SENTINEL2_L2A",
    {temporal_extent: ["2024-06-01", "2024-06-30"], bands: ["B04", "B08"]}
)

geojson_url = "https://example.com/geometry.geojson"
geometry = cube.load_url(geojson_url)
study_area = cube.filter_spatial(geometry)

R

library(openeo)

connection <- openeo::connect("openeofed.dataspace.copernicus.eu")$authenticate_oidc()
cube <- connection$load_collection(
    "SENTINEL2_L2A",
    list(temporal_extent = c("2024-06-01", "2024-06-30"), bands = c("B04", "B08"))
)

geojson_url <- "https://example.com/geometry.geojson"
geometry <- cube$load_url(geojson_url)
study_area <- cube$filter_spatial(geometry)

JavaScript

import openeo from "openeo";

const connection = await openeo.connect("openeofed.dataspace.copernicus.eu").authenticate_oidc();
const cube = connection.load_collection(
    "SENTINEL2_L2A",
    {temporal_extent: ["2024-06-01", "2024-06-30"], bands: ["B04", "B08"]}
);

const geojson_url = "https://example.com/geometry.geojson";
const geometry = await cube.load_url(geojson_url);
const study_area = cube.filter_spatial(geometry);

While the filter_bbox and filter_spatial processes can be used once the datacube is loaded, it is recommended to apply spatial filtering when loading the datacube to reduce the amount of data being processed and improve performance. This is particularly important for large datasets or when working with limited computational resources.

Resample a cube spatially

With the resample_spatial process, you can adjust the spatial resolution and projection of a datacube to match your analysis requirements. In other words, resample_spatial changes the spatial resolution and, optionally, the coordinate reference system. This is useful when combining datasets with different grids or when a coarser output is sufficient. The interpolation method matters: use near for classes and masks, and a continuous-data method such as bilinear for reflectance or temperature.

Python

aligned = study_area.resample_spatial(
    resolution=20,
    projection="EPSG:3857",
    method="bilinear",
)

R

aligned <- study_area$resample_spatial(
    resolution = 20,
    projection = "EPSG:3857",
    method = "bilinear"
)

JavaScript

const aligned = study_area.resample_spatial({
 resolution: 20,
 projection: "EPSG:3857",
 method: "bilinear"
});

NoteAlign to another cube

When a second cube is the target grid, use the resample_cube_spatial process rather than guessing a projection and resolution. It preserves the target cube’s spatial layout.

resample_cube_spatial aligns a cube with the spatial grid of a target cube. This process can be used as an alternative to specifying the resolution and projection, ensuring perfect alignment with the target cube.

Python

aligned = source.resample_cube_spatial(target)

R

aligned <- source$resample_cube_spatial(target)

JavaScript

const aligned = source.resample_cube_spatial(target);

Zonal statistics

To calculate numerical summaries (e.g., mean, sum) for each polygon in a vector dataset, you can use the aggregate_spatial process. This is commonly referred to as performing zonal statistics in geospatial analysis.

aggregate_spatial summarises raster values for each feature in a GeoJSON or vector collection. A common example is the mean NDVI per agricultural parcel. The result is a vector cube: each input feature receives its own value or time series.

Python

parcels = {
    "type": "FeatureCollection",
    "features": [{
        "type": "Feature",
        "properties": {"parcel_id": "A01"},
        "geometry": {
            "type": "Polygon",
            "coordinates": [[[4.35, 50.85], [4.38, 50.85], [4.38, 50.87],
                             [4.35, 50.87], [4.35, 50.85]]],
        },
    }],
}

parcel_means = aligned.aggregate_spatial(geometries=parcels, reducer="mean")
parcel_means.download("parcel_means.geojson", format="GeoJSON")

R

parcel_means <- aligned$aggregate_spatial(geometries = parcels, reducer = "mean")
parcel_means$download("parcel_means.geojson", format = "GeoJSON")

JavaScript

const parcelMeans = aligned.aggregate_spatial({ geometries: parcels, reducer: "mean" });
parcelMeans.download("parcel_means.geojson", { format: "GeoJSON" });

TipExample notebook

Parcel Delineation using openEO API uses aggregate_spatial to compute zonal statistics per parcel.

Filter vector features

The filter_vector process limits the vector data cube to the specified geometries. The process operates on geometries as defined in the OGC Simple Features standard. All geometries that were empty or have become empty will be removed from the data cube. Alternatively, use filter_bbox to filter with a bounding box or filter_spatial to filter a raster data cube based on geometries.

Python

selected = vector_cube.filter_vector(geometries=study_area)

R

selected <- vector_cube$filter_vector(geometries = study_area)

JavaScript

const selected = vector_cube.filter_vector({ geometries: study_area });

Reduce pixels within one geometry

The user can use reduce_spatial to summarise raster values for a single geometry, for example, to calculate the mean value for a single study area. However, use aggregate_spatial instead when you have multiple features and want a separate result for each one.

Python

mean_value = cube.reduce_spatial(geometries=study_area, reducer="mean")

R

mean_value <- cube$reduce_spatial(geometries = study_area, reducer = "mean")

JavaScript

const meanValue = cube.reduce_spatial({ geometries: study_area, reducer: "mean" });

Aggregate spatial windows

If the intention is to aggregate values over regularly sized spatial windows rather than irregular vector geometries, aggregate_spatial_window is the appropriate choice. It is particularly useful for grid-based analyses where each window has the same dimensions.

Python

windowed = cube.aggregate_spatial_window(size=5, reducer="mean")

R

windowed <- cube$aggregate_spatial_window(size = 5, reducer = "mean")

JavaScript

const windowed = cube.aggregate_spatial_window({ size: 5, reducer: "mean" });

Apply a spatial kernel

Usecases that require a fixed convolution over the spatial neighbourhood, such as smoothing, sharpening, or edge detection, can benefit from the apply_kernel process. The kernel defines the weights applied around each pixel and is preferable to a custom neighbourhood callback when the weights are known in advance.

Python

smoothed = cube.apply_kernel(kernel=[[1, 1, 1], [1, 1, 1], [1, 1, 1]])

R

smoothed <- cube$apply_kernel(kernel = matrix(1, nrow = 3, ncol = 3))

JavaScript

const smoothed = cube.apply_kernel({ kernel: [[1, 1, 1], [1, 1, 1], [1, 1, 1]] });

Process spatial neighborhoods

With the apply_neighborhood process, you can define a custom operation for each moving window over the spatial neighbourhood, for example, a local percentile or a conditional texture measure.

In other words, this process applies a process to a neighbourhood of pixels in a sliding-window fashion with (optional) overlap. In this case, data chunking is explicitly controlled by the user. Dimensions and number of labels are fully preserved. This is the most versatile and widely used function to work with UDF’s.

Python

udf = openeo.UDF(
"""
import xarray

def apply_datacube(cube: xarray.DataArray, context: dict) -> xarray.DataArray:
 #### PYTHON UDF EXAMPLE with an operation on the data cube
 return cube
"""
)
neighborhood_result = cube.apply_neighborhood(
    process=udf,
    size=[
        {"dimension": "x", "value": 384, "unit": "px"},
        {"dimension": "y", "value": 384, "unit": "px"},
    ],
    overlap=[
        {"dimension": "x", "value": 64, "unit": "px"},
        {"dimension": "y", "value": 64, "unit": "px"},
    ]
)

R

udf <- openeo::UDF(
"
library(xarray)

apply_datacube <- function(cube, context) {
 #### R UDF EXAMPLE with an operation on the data cube
 return(cube)
}
"
)
neighborhood_result <- cube$apply_neighborhood(
  process = udf,
  size = list(
    list(dimension = "x", value = 384, unit = "px"),
    list(dimension = "y", value = 384, unit = "px")
  ),
  overlap = list(
    list(dimension = "x", value = 64, unit = "px"),
    list(dimension = "y", value = 64, unit = "px")
  )
)

TipExample notebook

Hillshade from Copernicus 30 m DEM derives slope, aspect, and hillshade with apply_neighborhood-style band math over a DEM.

Mask outside a polygon

While mask mentioned in the preprocessing section is used for applying a mask on a raster cube, mask_polygon is specifically used to restrict the cube to an irregular study area defined by vector geometries. All pixels for which the point at the pixel centre does not intersect with any polygon (as defined in the Simple Features standard by the OGC) are replaced. This behaviour can be inverted by setting the parameter inside to true.

Python

masked = cube.mask_polygon(geometries=study_area)

R

masked <- cube$mask_polygon(geometries = study_area)

JavaScript

let masked = cube.mask_polygon({geometries: study_area});

Buffer vector geometries

A well-known task when working with vector geometries is to create buffers around features to account for spatial uncertainty or to define areas of influence. Users can use vector_buffer to expand or contract geometries before filtering or aggregating. A positive distance creates an area around a feature; a negative distance can create an interior buffer where supported.

Python

buffered = openeo.processes.vector_buffer(data=geometry, distance=100)

R

buffered <- openeo.processes.vector_buffer(data = geometry, distance = 100)

JavaScript

let buffered = openeo.processes.vector_buffer({data: geometry, distance: 100});

Reproject vector geometries

When working with vector geometry-related tasks, such as vector buffering using vector_buffer, it is important to ensure that the geometries are in the correct coordinate reference system (CRS). If the spatial reference system unit is not meters, a UnitMismatch error is thrown. Use vector_reproject() to convert the geometries to a suitable spatial reference system.

Python

reprojected = openeo.processes.vector_reproject(
    data=geometry, target_crs="EPSG:3035"
)

R

reprojected <- openeo.processes.vector_reproject(data = geometry, target_crs = "EPSG:3035")

JavaScript

let reprojected = openeo.processes.vector_reproject({data: geometry, target_crs: "EPSG:3035"});

Apply a process per polygon

apply_polygon applies a process to the sub data cube covered by each polygon in a set of geometries, rather than masking the cube once and running a separate process afterwards. This is useful when a per-parcel operation (e.g. a custom smoothing or classification step) needs to run independently for every polygon, while pixels outside of each polygon are replaced with mask_value.

Python

per_parcel = cube.apply_polygon(
    geometries=parcels,
    process=lambda data: data.linear_scale_range(0, 1, 0, 255),
)

R

per_parcel <- cube$apply_polygon(
    geometries = parcels,
    process = function(data) linear_scale_range(data, inputMin = 0, inputMax = 1, outputMin = 0, outputMax = 255)
)

JavaScript

const perParcel = cube.apply_polygon({
  geometries: parcels,
  process: data => data.linearScaleRange(0, 1, 0, 255),
});

Sample random points from geometries

vector_to_random_points generates a vector data cube of points by sampling random locations from within input geometries. It is useful for building training or validation point sets from polygons, for example to extract pixel values for a classifier.

Python

sample_points = openeo.processes.vector_to_random_points(
    data=parcels, geometry_count=10, seed=42
)

R

sample_points <- openeo.processes.vector_to_random_points(data = parcels, geometry_count = 10, seed = 42)

JavaScript

let samplePoints = openeo.processes.vector_to_random_points({data: parcels, geometry_count: 10, seed: 42});

Sample regularly spaced points from geometries

vector_to_regular_points generates a vector data cube of points spaced at a fixed minimum distance within input geometries, which is useful for building a regular sampling grid over a study area instead of random points.

Python

grid_points = openeo.processes.vector_to_regular_points(data=parcels, distance=100)

R

grid_points <- openeo.processes.vector_to_regular_points(data = parcels, distance = 100)

JavaScript

let gridPoints = openeo.processes.vector_to_regular_points({data: parcels, distance: 100});