Item Search with Intersects

This notebook shows the use of pystac-client to perform item search with the intersects parameter, to restrict the results to an Area of Interest (AOI).

Client

We first connect to an API by retrieving the root catalog, or landing page, of the API with the Client.open function.

from pystac_client import Client
from typing import Any, Dict
import json

# STAC API root URL
URL = 'https://planetarycomputer.microsoft.com/api/stac/v1'

client = Client.open(URL)

Item Search

When the Catalog is a STAC API, we have the ability to search for items based on spatio-temporal properties.

The STAC API endpoint /search accepts a parameter intersects that is a GeoJSON Geometry representing the AOI of the search.

The search method of the pystac_client Client class accepts several different types of objects:

1. a string representing a GeoJSON geometry

2. a dictionary representing a GeoJSON geometry

3. any object that implements a __geo_interface__ property, an informal specification supported by several libraries for generating a GeoJSON representation from an object. Several prominent libraries support this protocol for their objects that represent geometries, including Shapely, ArcPy, PySAL, geojson, pyshp, descartes, and pygeoif

Item Search with AOI as a Dictionary

# AOI around Delfzijl, in northern Netherlands
aoi_as_dict: Dict[str, Any] = {
    "type": "Polygon",
    "coordinates": [
      [
        [
          6,
          53
        ],
        [
          7,
          53
        ],
        [
          7,
          54
        ],
        [
          6,
          54
        ],
        [
          6,
          53
        ]
      ]
    ]
}

search = client.search(
    max_items = 25,
    collections = "aster-l1t",
    intersects = aoi_as_dict,
)

print(f"AOI as dictionary, found {len(list(search.items()))} items")

AOI as dictionary, found 25 items

Item Search with AOI as a String

aoi_as_str: str  = json.dumps(aoi_as_dict)

search = client.search(
    max_items = 25,
    collections = "aster-l1t",
    intersects = aoi_as_str,
)

print(f"AOI as string, found {len(list(search.items()))} items")

AOI as string, found 25 items

Item Search with AOI as a Shapely Geometry Object

import shapely.geometry

aoi_as_shapely_shape = shapely.geometry.shape(aoi_as_dict)

search = client.search(
    max_items = 25,
    collections = "aster-l1t",
    intersects = aoi_as_shapely_shape,
)

print(f"AOI as Shapely Geometry object from shape(), found {len(list(search.items()))} items")

aoi_as_shapely_polygon = shapely.geometry.Polygon(aoi_as_dict["coordinates"][0])

search = client.search(
    max_items = 25,
    collections = "aster-l1t",
    intersects = aoi_as_shapely_polygon,
)

print(f"AOI as Shapely Geometry object with Polygon, found {len(list(search.items()))} items")

AOI as Shapely Geometry object from shape(), found 25 items

AOI as Shapely Geometry object with Polygon, found 25 items

Item Search with AOI as a “geojson” library object

import geojson

aoi_as_geojson_polygon = geojson.Polygon(aoi_as_dict["coordinates"])

search = client.search(
    max_items = 25,
    collections = "aster-l1t",
    intersects = aoi_as_geojson_polygon,
)

print(f"AOI as geojson Polygon, found {len(list(search.items()))} items")

AOI as geojson Polygon, found 25 items

Item Search with AOI as a pygeoif Object

import pygeoif
aoi_as_pygeoif_polygon = pygeoif.geometry.Polygon(aoi_as_dict["coordinates"][0])

search = client.search(
    max_items = 25,
    collections = "aster-l1t",
    intersects = aoi_as_pygeoif_polygon,
)

print(f"AOI as pygeoif Polygon, found {len(list(search.items()))} items")

AOI as pygeoif Polygon, found 25 items