Visualizing orbital swath data¶

Visualizing L2 and L3 data products that are generated from sensors that follow orbital tracks around the earth requires some special care to ensure that the underlying data are combined in reasonable ways. The TiTiler-CMR API surfaces several query parameters from the CMR API that make it possible to filter granules down by attributes such as TRACK_NUMBER or other properties that define logical groups of granules.

This notebook contains visualization examples that use granule filters to render coherent views of L2 collections like NISAR Beta Geocoded Polarimetric Covariance (GCOV) and Harmonize Landsat Sentinel (HLS), and L3 collections like OPERA Surface Displacement from Sentinel-1.

In [1]:

import json
from datetime import datetime, timedelta, UTC

import httpx
import shapely
from folium import GeoJson, GeoJsonPopup, GeoJsonTooltip, LayerControl, Map, TileLayer
from shapely.geometry import box, mapping

titiler_endpoint = "https://staging.openveda.cloud/api/titiler-cmr"  # staging endpoint

import json from datetime import datetime, timedelta, UTC import httpx import shapely from folium import GeoJson, GeoJsonPopup, GeoJsonTooltip, LayerControl, Map, TileLayer from shapely.geometry import box, mapping titiler_endpoint = "https://staging.openveda.cloud/api/titiler-cmr" # staging endpoint

NISAR GCOV in India¶

Shortly before this notebook was written (March 2026) ASF published ~22k granules from the NISAR Beta GCOV collection. It is early days for NISAR but TiTiler-CMR can be used to visualize a mosaic of granules.

The /bbox/.../granules endpoint can be used to retrieve the granule metadata that matches a TiTiler-CMR query as a GeoJSON feature collection. This is useful because it allows you to look directly at the metadata from the granules that TiTiler-CMR will use for a set of CMR search terms.

In [2]:

india_bbox = (66.191, 6.930, 91.297, 34.657)

nisar_assets_req = httpx.get(
    f"{titiler_endpoint}/bbox/{','.join(str(b) for b in india_bbox)}/granules",
    params={
        "collection_concept_id": "C3622214170-ASF",
        "temporal": "2026-01-01T00:00:00Z/2026-02-01T00:00:00Z",  # January 2026
        "f": "geojson",  # return granules as a GeoJSON FeatureCollection
        "sortkey": "-start_date",  # sort granules in descending order by time
        "skipcovered": True,  # do not return multiple granules with identical geometries
        "items_limit": 400,  # global granule query limit
    },
    timeout=None,
)

nisar_assets_req.raise_for_status()
nisar_assets_geojson = nisar_assets_req.json()

print("found", len(nisar_assets_geojson["features"]), "granules")

india_bbox = (66.191, 6.930, 91.297, 34.657) nisar_assets_req = httpx.get( f"{titiler_endpoint}/bbox/{','.join(str(b) for b in india_bbox)}/granules", params={ "collection_concept_id": "C3622214170-ASF", "temporal": "2026-01-01T00:00:00Z/2026-02-01T00:00:00Z", # January 2026 "f": "geojson", # return granules as a GeoJSON FeatureCollection "sortkey": "-start_date", # sort granules in descending order by time "skipcovered": True, # do not return multiple granules with identical geometries "items_limit": 400, # global granule query limit }, timeout=None, ) nisar_assets_req.raise_for_status() nisar_assets_geojson = nisar_assets_req.json() print("found", len(nisar_assets_geojson["features"]), "granules")

found 267 granules

In [3]:

for feature in nisar_assets_geojson["features"]:
    # move additional_attributes into properties so we can use
    # them more easily in maps
    feature["properties"].update(
        {
            attribute["name"]: (
                attribute["values"][0]
                if len(attribute["values"]) == 1
                else attribute["values"]
            )
            for attribute in feature["properties"]["additional_attributes"]
        }
    )

    # drop additional_attributes after moving values out
    feature["properties"].pop("additional_attributes")

    feature["properties"]["datetime"] = datetime.fromisoformat(
        feature["properties"]["temporal_extent"]["range_date_time"][
            "beginning_date_time"
        ]
    ).isoformat()

for feature in nisar_assets_geojson["features"]: # move additional_attributes into properties so we can use # them more easily in maps feature["properties"].update( { attribute["name"]: ( attribute["values"][0] if len(attribute["values"]) == 1 else attribute["values"] ) for attribute in feature["properties"]["additional_attributes"] } ) # drop additional_attributes after moving values out feature["properties"].pop("additional_attributes") feature["properties"]["datetime"] = datetime.fromisoformat( feature["properties"]["temporal_extent"]["range_date_time"][ "beginning_date_time" ] ).isoformat()

In [4]:

m = Map(
    location=((india_bbox[1] + india_bbox[3]) / 2, (india_bbox[0] + india_bbox[2]) / 2),
    zoom_start=5,
)

GeoJson(
    nisar_assets_geojson,
    name="granule footprints",
    tooltip=GeoJsonTooltip(fields=["granule_ur"], localize=True),
    popup=GeoJsonPopup(fields=["datetime", "ASCENDING_DESCENDING", "TRACK_NUMBER"]),
    style_function=lambda x: {
        "color": "magenta"
        if x["properties"]["ASCENDING_DESCENDING"] == "ASCENDING"
        else "blue",
        "weight": 1,
        "fillOpacity": 0,
    },
).add_to(m)

LayerControl().add_to(m)

m

m = Map( location=((india_bbox[1] + india_bbox[3]) / 2, (india_bbox[0] + india_bbox[2]) / 2), zoom_start=5, ) GeoJson( nisar_assets_geojson, name="granule footprints", tooltip=GeoJsonTooltip(fields=["granule_ur"], localize=True), popup=GeoJsonPopup(fields=["datetime", "ASCENDING_DESCENDING", "TRACK_NUMBER"]), style_function=lambda x: { "color": "magenta" if x["properties"]["ASCENDING_DESCENDING"] == "ASCENDING" else "blue", "weight": 1, "fillOpacity": 0, }, ).add_to(m) LayerControl().add_to(m) m

Out[4]:

Make this Notebook Trusted to load map: File -> Trust Notebook

In [5]:

nisar_tilejson_params = {
    "collection_concept_id": "C3622214170-ASF",
    "temporal": "2026-01-01T00:00:00Z/2026-02-01T00:00:00Z",
    "bounding_box": ",".join(str(b) for b in india_bbox),  # optional
    "sortkey": "-start_date",  # prioritize most recent granules during rendering
    "skipcovered": True,  # skip granules with geometry identical to collected granules
    "exitwhenfull": True,  # exit CMR search once tile geometry is covered by assets
    "variables": ["HHHH", "HVHV"],  # two variables to be combined via expression
    "expression": "10 * log10(b1); 10 * log10(b2); 10 * log10(b1/b2)",
    "rescale": [[-20, 0], [-30, 5], [2, 18]],
}

nisar_frequencyA_params = {
    "group": "/science/LSAR/GCOV/grids/frequencyA",
    "minzoom": 10,
}

nisar_frequencyB_params = {
    "group": "/science/LSAR/GCOV/grids/frequencyB",
    "minzoom": 6,
    "maxzoom": 10,
}

nisar_tilejson_params = { "collection_concept_id": "C3622214170-ASF", "temporal": "2026-01-01T00:00:00Z/2026-02-01T00:00:00Z", "bounding_box": ",".join(str(b) for b in india_bbox), # optional "sortkey": "-start_date", # prioritize most recent granules during rendering "skipcovered": True, # skip granules with geometry identical to collected granules "exitwhenfull": True, # exit CMR search once tile geometry is covered by assets "variables": ["HHHH", "HVHV"], # two variables to be combined via expression "expression": "10 * log10(b1); 10 * log10(b2); 10 * log10(b1/b2)", "rescale": [[-20, 0], [-30, 5], [2, 18]], } nisar_frequencyA_params = { "group": "/science/LSAR/GCOV/grids/frequencyA", "minzoom": 10, } nisar_frequencyB_params = { "group": "/science/LSAR/GCOV/grids/frequencyB", "minzoom": 6, "maxzoom": 10, }

Visualize a single track¶

The NISAR granules can be overwhelming for the TiTiler-CMR application. One way to limit the number of granules that any single tile request attempts to open is to apply an orbit/track filter to the CMR granule queries. This is done by specifying an attribute query parameter (docs).

The NISAR granules contain HHHH and VVVV arrays at two resolutions, represented by frequencyA (high res) and frequencyB (low res). We can use this feature to compose a set of tile layers that will create a visually seamless view of the data by rendering tiles for frequencyB at lower zoom levels (zoomed out) and frequencyA at high zoom levels (zoomed in).

In [6]:

# pick a track number then visualize the track
track_number = 149

track_params = {"attribute": f"int,TRACK_NUMBER,{track_number}"}

nisar_track_frequencyA_tilejson_req = httpx.get(
    f"{titiler_endpoint}/xarray/WebMercatorQuad/tilejson.json",
    params={
        **nisar_tilejson_params,
        **nisar_frequencyA_params,
        **track_params,
    },
    timeout=None,
)

nisar_track_frequencyA_tilejson_req.raise_for_status()
nisar_track_frequencyA_tilejson = nisar_track_frequencyA_tilejson_req.json()

nisar_track_frequencyB_tilejson_req = httpx.get(
    f"{titiler_endpoint}/xarray/WebMercatorQuad/tilejson.json",
    params={
        **nisar_tilejson_params,
        **nisar_frequencyB_params,
        **track_params,
    },
    timeout=None,
)

nisar_track_frequencyB_tilejson_req.raise_for_status()
nisar_track_frequencyB_tilejson = nisar_track_frequencyB_tilejson_req.json()

# pick a track number then visualize the track track_number = 149 track_params = {"attribute": f"int,TRACK_NUMBER,{track_number}"} nisar_track_frequencyA_tilejson_req = httpx.get( f"{titiler_endpoint}/xarray/WebMercatorQuad/tilejson.json", params={ **nisar_tilejson_params, **nisar_frequencyA_params, **track_params, }, timeout=None, ) nisar_track_frequencyA_tilejson_req.raise_for_status() nisar_track_frequencyA_tilejson = nisar_track_frequencyA_tilejson_req.json() nisar_track_frequencyB_tilejson_req = httpx.get( f"{titiler_endpoint}/xarray/WebMercatorQuad/tilejson.json", params={ **nisar_tilejson_params, **nisar_frequencyB_params, **track_params, }, timeout=None, ) nisar_track_frequencyB_tilejson_req.raise_for_status() nisar_track_frequencyB_tilejson = nisar_track_frequencyB_tilejson_req.json()

In [7]:

nisar_track_features = {
    "type": "FeatureCollection",
    "features": [
        feature
        for feature in nisar_assets_geojson["features"]
        if int(feature["properties"]["TRACK_NUMBER"]) == track_number
    ],
}

aoi = shapely.union_all(
    [box(*feature["bbox"]) for feature in nisar_track_features["features"]]
)

m = Map(
    location=((aoi.bounds[1] + aoi.bounds[3]) / 2, (aoi.bounds[0] + aoi.bounds[2]) / 2),
    zoom_start=7,
)

GeoJson(
    nisar_track_features,
    name="granule footprints",
    popup=GeoJsonPopup(
        fields=["granule_ur", "datetime", "ASCENDING_DESCENDING", "TRACK_NUMBER"]
    ),
    style_function=lambda x: {
        "color": "orange",
        "weight": 2,
        "fillOpacity": 0,
    },
).add_to(m)

TileLayer(
    tiles=nisar_track_frequencyB_tilejson["tiles"][0],
    name="NISAR GCOV frequencyB (lower res)",
    opacity=1,
    attr="NASA/ISRO",
    overlay=True,
    min_zoom=6,
    max_zoom=9,
).add_to(m)

TileLayer(
    tiles=nisar_track_frequencyA_tilejson["tiles"][0],
    name="NISAR GCOV frequencyA (higher res)",
    opacity=1,
    attr="NASA/ISRO",
    overlay=True,
    min_zoom=10,
).add_to(m)

LayerControl().add_to(m)

m

nisar_track_features = { "type": "FeatureCollection", "features": [ feature for feature in nisar_assets_geojson["features"] if int(feature["properties"]["TRACK_NUMBER"]) == track_number ], } aoi = shapely.union_all( [box(*feature["bbox"]) for feature in nisar_track_features["features"]] ) m = Map( location=((aoi.bounds[1] + aoi.bounds[3]) / 2, (aoi.bounds[0] + aoi.bounds[2]) / 2), zoom_start=7, ) GeoJson( nisar_track_features, name="granule footprints", popup=GeoJsonPopup( fields=["granule_ur", "datetime", "ASCENDING_DESCENDING", "TRACK_NUMBER"] ), style_function=lambda x: { "color": "orange", "weight": 2, "fillOpacity": 0, }, ).add_to(m) TileLayer( tiles=nisar_track_frequencyB_tilejson["tiles"][0], name="NISAR GCOV frequencyB (lower res)", opacity=1, attr="NASA/ISRO", overlay=True, min_zoom=6, max_zoom=9, ).add_to(m) TileLayer( tiles=nisar_track_frequencyA_tilejson["tiles"][0], name="NISAR GCOV frequencyA (higher res)", opacity=1, attr="NASA/ISRO", overlay=True, min_zoom=10, ).add_to(m) LayerControl().add_to(m) m

Out[7]:

Make this Notebook Trusted to load map: File -> Trust Notebook

Visualizing multiple orbits/tracks¶

It is possible to visualize a mosaic across multiple orbits/tracks, but this will result in TiTiler-CMR opening more granules per tile request and due to the structure of the granule files this is a resource intensive process! These visualizations should be used with a higher minimum zoom level than the orbit-wise tile layers.

In [8]:

# write an attribute filter to select only granules on ascending orbits
ascending_params = {"attribute": "string,ASCENDING_DESCENDING,ASCENDING"}

nisar_ascending_frequencyA_tilejson_req = httpx.get(
    f"{titiler_endpoint}/xarray/WebMercatorQuad/tilejson.json",
    params={
        **nisar_tilejson_params,
        **nisar_frequencyA_params,
        **ascending_params,
    },
    timeout=None,
)

nisar_ascending_frequencyA_tilejson_req.raise_for_status()
nisar_ascending_frequencyA_tilejson = nisar_ascending_frequencyA_tilejson_req.json()

nisar_ascending_frequencyB_tilejson_req = httpx.get(
    f"{titiler_endpoint}/xarray/WebMercatorQuad/tilejson.json",
    params={
        **nisar_tilejson_params,
        **nisar_frequencyB_params,
        **ascending_params,
        # minzoom must be higher since TiTiler-CMR will struggle with too many of these granules!
        "minzoom": 8,
    },
    timeout=None,
)

nisar_ascending_frequencyB_tilejson_req.raise_for_status()
nisar_ascemdomg_frequencyB_tilejson = nisar_ascending_frequencyB_tilejson_req.json()

# write an attribute filter to select only granules on ascending orbits ascending_params = {"attribute": "string,ASCENDING_DESCENDING,ASCENDING"} nisar_ascending_frequencyA_tilejson_req = httpx.get( f"{titiler_endpoint}/xarray/WebMercatorQuad/tilejson.json", params={ **nisar_tilejson_params, **nisar_frequencyA_params, **ascending_params, }, timeout=None, ) nisar_ascending_frequencyA_tilejson_req.raise_for_status() nisar_ascending_frequencyA_tilejson = nisar_ascending_frequencyA_tilejson_req.json() nisar_ascending_frequencyB_tilejson_req = httpx.get( f"{titiler_endpoint}/xarray/WebMercatorQuad/tilejson.json", params={ **nisar_tilejson_params, **nisar_frequencyB_params, **ascending_params, # minzoom must be higher since TiTiler-CMR will struggle with too many of these granules! "minzoom": 8, }, timeout=None, ) nisar_ascending_frequencyB_tilejson_req.raise_for_status() nisar_ascemdomg_frequencyB_tilejson = nisar_ascending_frequencyB_tilejson_req.json()

In [9]:

m = Map(
    location=((india_bbox[1] + india_bbox[3]) / 2, (india_bbox[0] + india_bbox[2]) / 2),
    zoom_start=8,
)

nisar_ascending_features = {
    "type": "FeatureCollection",
    "features": [
        feature
        for feature in nisar_assets_geojson["features"]
        if feature["properties"]["ASCENDING_DESCENDING"] == "ASCENDING"
    ],
}

GeoJson(
    nisar_ascending_features,
    name="granule footprints",
    popup=GeoJsonPopup(
        fields=["granule_ur", "datetime", "ASCENDING_DESCENDING", "TRACK_NUMBER"]
    ),
    style_function=lambda x: {
        "color": "orange",
        "weight": 2,
        "fillOpacity": 0,
    },
).add_to(m)

TileLayer(
    tiles=nisar_ascemdomg_frequencyB_tilejson["tiles"][0],
    name="NISAR GCOV frequencyB (lower res)",
    opacity=1,
    attr="NASA/ISRO",
    overlay=True,
    min_zoom=8,
    max_zoom=10,
).add_to(m)

TileLayer(
    tiles=nisar_ascending_frequencyA_tilejson["tiles"][0],
    name="NISAR GCOV frequencyA (higher res)",
    opacity=1,
    attr="NASA/ISRO",
    overlay=True,
    min_zoom=11,
).add_to(m)

LayerControl().add_to(m)

m

m = Map( location=((india_bbox[1] + india_bbox[3]) / 2, (india_bbox[0] + india_bbox[2]) / 2), zoom_start=8, ) nisar_ascending_features = { "type": "FeatureCollection", "features": [ feature for feature in nisar_assets_geojson["features"] if feature["properties"]["ASCENDING_DESCENDING"] == "ASCENDING" ], } GeoJson( nisar_ascending_features, name="granule footprints", popup=GeoJsonPopup( fields=["granule_ur", "datetime", "ASCENDING_DESCENDING", "TRACK_NUMBER"] ), style_function=lambda x: { "color": "orange", "weight": 2, "fillOpacity": 0, }, ).add_to(m) TileLayer( tiles=nisar_ascemdomg_frequencyB_tilejson["tiles"][0], name="NISAR GCOV frequencyB (lower res)", opacity=1, attr="NASA/ISRO", overlay=True, min_zoom=8, max_zoom=10, ).add_to(m) TileLayer( tiles=nisar_ascending_frequencyA_tilejson["tiles"][0], name="NISAR GCOV frequencyA (higher res)", opacity=1, attr="NASA/ISRO", overlay=True, min_zoom=11, ).add_to(m) LayerControl().add_to(m) m

Out[9]:

Make this Notebook Trusted to load map: File -> Trust Notebook

OPERA Surface Displacement from Sentinel-1¶

Some collections require very careful filtering in order to produce valid mosaics. The OPERA Surface Displacement from Sentinel-1 collection contains granules that represent a comparison between two specific points in time. For this reason it does not make much sense to create mosaics of granules where those points in time do not match very closely.

The granules describe the reference and secondary timestamps in the additional attributes REFERENCE_ZERO_DOPPLER_START_TIME and SECONDARY_ZERO_DOPPLER_END_TIME. These values are represented as strings in the metadata but the CMR API allows you to filter down to granules where the additional attribute value falls in a range. This works in this case because the string datetimes are ordered alphabetically/chronologically.

In [10]:

opera_disp_collection_concept_id = "C3294057315-ASF"

ca_bbox = (-127, 31, -113, 42)

# we want granules that share the same start/end date in the interval attributes
start_datetime = datetime(2016, 7, 5, tzinfo=UTC)
end_datetime = datetime(2017, 1, 25, tzinfo=UTC)

reference_interval = (
    start_datetime,
    start_datetime + timedelta(days=1) - timedelta(seconds=1),
)

secondary_interval = (
    end_datetime,
    end_datetime + timedelta(days=1) - timedelta(seconds=1),
)

dt_fmt = "%Y-%m-%dT%H:%M:%SZ"

opera_disp_params = {
    "collection_concept_id": opera_disp_collection_concept_id,
    "attribute": [
        f"string,REFERENCE_ZERO_DOPPLER_START_TIME,{reference_interval[0].strftime(dt_fmt)},{reference_interval[1].strftime(dt_fmt)}",
        f"string,SECONDARY_ZERO_DOPPLER_END_TIME,{secondary_interval[0].strftime(dt_fmt)},{secondary_interval[1].strftime(dt_fmt)}",
    ],
    "sortkey": "-start_date",
}

opera_disp_assets_req = httpx.get(
    f"{titiler_endpoint}/bbox/{','.join(str(b) for b in ca_bbox)}/granules",
    params={
        "skipcovered": False,
        "limit": 500,
        "f": "geojson",
        **opera_disp_params,
    },
    timeout=None,
)

opera_disp_assets_req.raise_for_status()
opera_disp_assets_geojson = opera_disp_assets_req.json()
print("found", len(opera_disp_assets_geojson["features"]), "granules")

opera_disp_collection_concept_id = "C3294057315-ASF" ca_bbox = (-127, 31, -113, 42) # we want granules that share the same start/end date in the interval attributes start_datetime = datetime(2016, 7, 5, tzinfo=UTC) end_datetime = datetime(2017, 1, 25, tzinfo=UTC) reference_interval = ( start_datetime, start_datetime + timedelta(days=1) - timedelta(seconds=1), ) secondary_interval = ( end_datetime, end_datetime + timedelta(days=1) - timedelta(seconds=1), ) dt_fmt = "%Y-%m-%dT%H:%M:%SZ" opera_disp_params = { "collection_concept_id": opera_disp_collection_concept_id, "attribute": [ f"string,REFERENCE_ZERO_DOPPLER_START_TIME,{reference_interval[0].strftime(dt_fmt)},{reference_interval[1].strftime(dt_fmt)}", f"string,SECONDARY_ZERO_DOPPLER_END_TIME,{secondary_interval[0].strftime(dt_fmt)},{secondary_interval[1].strftime(dt_fmt)}", ], "sortkey": "-start_date", } opera_disp_assets_req = httpx.get( f"{titiler_endpoint}/bbox/{','.join(str(b) for b in ca_bbox)}/granules", params={ "skipcovered": False, "limit": 500, "f": "geojson", **opera_disp_params, }, timeout=None, ) opera_disp_assets_req.raise_for_status() opera_disp_assets_geojson = opera_disp_assets_req.json() print("found", len(opera_disp_assets_geojson["features"]), "granules")

found 2 granules

You can check the granule_ur values for the returned granules:

In [11]:

for granule in opera_disp_assets_geojson["features"]:
    print(granule["properties"]["granule_ur"])

for granule in opera_disp_assets_geojson["features"]: print(granule["properties"]["granule_ur"])

OPERA_L3_DISP-S1_IW_F09154_VV_20160705T020626Z_20170125T020626Z_v1.0_20250408T162732Z
OPERA_L3_DISP-S1_IW_F09155_VV_20160705T020645Z_20170125T020646Z_v1.0_20250408T163031Z

Each granule comes with a list of additional_attributes that are useful for filtering.

In [12]:

print(
    json.dumps(
        opera_disp_assets_geojson["features"][0]["properties"]["additional_attributes"],
        indent=2,
    )
)

print( json.dumps( opera_disp_assets_geojson["features"][0]["properties"]["additional_attributes"], indent=2, ) )

[
  {
    "name": "ASCENDING_DESCENDING",
    "values": [
      "ASCENDING"
    ]
  },
  {
    "name": "FRAME_NUMBER",
    "values": [
      "9154"
    ]
  },
  {
    "name": "PATH_NUMBER",
    "values": [
      "35"
    ]
  },
  {
    "name": "POLARIZATION",
    "values": [
      "VV"
    ]
  },
  {
    "name": "PROCESSING_TYPE",
    "values": [
      "DISP-S1"
    ]
  },
  {
    "name": "PRODUCT_VERSION",
    "values": [
      "1.0"
    ]
  },
  {
    "name": "REFERENCE_ZERO_DOPPLER_END_TIME",
    "values": [
      "2016-07-05T02:06:48Z"
    ]
  },
  {
    "name": "REFERENCE_ZERO_DOPPLER_START_TIME",
    "values": [
      "2016-07-05T02:06:26Z"
    ]
  },
  {
    "name": "SECONDARY_ZERO_DOPPLER_END_TIME",
    "values": [
      "2017-01-25T02:06:49Z"
    ]
  },
  {
    "name": "SECONDARY_ZERO_DOPPLER_START_TIME",
    "values": [
      "2017-01-25T02:06:26Z"
    ]
  },
  {
    "name": "STACK_ID",
    "values": [
      "IW_F09154_VV"
    ]
  }
]

We can request a tilejson that can be used to render tiles for this mosaic of two granules like this:

In [13]:

opera_disp_tilejson_req = httpx.get(
    f"{titiler_endpoint}/xarray/WebMercatorQuad/tilejson.json",
    params={
        **opera_disp_params,
        "bounding_box": ",".join(str(b) for b in ca_bbox),
        "exitwhenfull": True,
        "variables": ["displacement"],
        "rescale": [[-0.05, 0.05]],
        "colormap_name": "rainbow",
        "minzoom": 6,
    },
    timeout=None,
)

opera_disp_tilejson_req.raise_for_status()
opera_disp_tilejson = opera_disp_tilejson_req.json()
print(opera_disp_tilejson)

opera_disp_tilejson_req = httpx.get( f"{titiler_endpoint}/xarray/WebMercatorQuad/tilejson.json", params={ **opera_disp_params, "bounding_box": ",".join(str(b) for b in ca_bbox), "exitwhenfull": True, "variables": ["displacement"], "rescale": [[-0.05, 0.05]], "colormap_name": "rainbow", "minzoom": 6, }, timeout=None, ) opera_disp_tilejson_req.raise_for_status() opera_disp_tilejson = opera_disp_tilejson_req.json() print(opera_disp_tilejson)

{'tilejson': '3.0.0', 'version': '1.0.0', 'scheme': 'xyz', 'tiles': ['https://staging.openveda.cloud/api/titiler-cmr/xarray/tiles/WebMercatorQuad/{z}/{x}/{y}?collection_concept_id=C3294057315-ASF&attribute=string%2CREFERENCE_ZERO_DOPPLER_START_TIME%2C2016-07-05T00%3A00%3A00Z%2C2016-07-05T23%3A59%3A59Z&attribute=string%2CSECONDARY_ZERO_DOPPLER_END_TIME%2C2017-01-25T00%3A00%3A00Z%2C2017-01-25T23%3A59%3A59Z&sortkey=-start_date&bounding_box=-127%2C31%2C-113%2C42&exitwhenfull=true&variables=displacement&rescale=%5B-0.05%2C+0.05%5D&colormap_name=rainbow&tilesize=512'], 'minzoom': 6, 'maxzoom': 18, 'bounds': [-127.0, 31.0, -113.0, 42.0], 'center': [-120.0, 36.5, 6]}

In [14]:

aoi = shapely.union_all(
    [box(*feature["bbox"]) for feature in opera_disp_assets_geojson["features"]]
)

m = Map(
    location=((aoi.bounds[1] + aoi.bounds[3]) / 2, (aoi.bounds[0] + aoi.bounds[2]) / 2),
    zoom_start=8,
)

GeoJson(
    opera_disp_assets_geojson,
    name="granule footprints",
    popup=GeoJsonPopup(fields=["granule_ur"]),
    style_function=lambda x: {
        "color": "magenta",
        "weight": 2,
        "fillOpacity": 0,
    },
).add_to(m)

TileLayer(
    tiles=opera_disp_tilejson["tiles"][0],
    name="OPERA-DISP",
    opacity=1,
    attr="NASA",
    min_zoom=6,
).add_to(m)

LayerControl().add_to(m)


m

aoi = shapely.union_all( [box(*feature["bbox"]) for feature in opera_disp_assets_geojson["features"]] ) m = Map( location=((aoi.bounds[1] + aoi.bounds[3]) / 2, (aoi.bounds[0] + aoi.bounds[2]) / 2), zoom_start=8, ) GeoJson( opera_disp_assets_geojson, name="granule footprints", popup=GeoJsonPopup(fields=["granule_ur"]), style_function=lambda x: { "color": "magenta", "weight": 2, "fillOpacity": 0, }, ).add_to(m) TileLayer( tiles=opera_disp_tilejson["tiles"][0], name="OPERA-DISP", opacity=1, attr="NASA", min_zoom=6, ).add_to(m) LayerControl().add_to(m) m

Out[14]:

Make this Notebook Trusted to load map: File -> Trust Notebook

Harmonized Landsat Sentinel-2¶

Visualizing L2 satellite imagery collections comes with a different set of considerations. Users will want to be able to browse seamless mosaics of imagery probably with minimal cloud cover. As of v1.0 this can be accomplished with TiTiler-CMR.

The key to efficient visualization of satellite imagery datasets from CMR is to use the sort_key field to return the granules in ascending order by cloud cover so the least cloudy granules come back first, then add exitwhenfull to stop returning granules in the search when the tile geometry is fully covered by the returned granules.

See also: HLS configuration docs

In [15]:

hls_bbox = (-93, 45, -89, 49)
minzoom = 7

hls_geojson = GeoJson(
    mapping(box(*hls_bbox)),
    style_function=lambda x: {
        "color": "orange",
        "weight": 1,
        "fillOpacity": 0,
    },
    name="AOI",
)

m = Map(
    location=((hls_bbox[1] + hls_bbox[3]) / 2, (hls_bbox[0] + hls_bbox[2]) / 2),
    zoom_start=minzoom,
)

hls_geojson.add_to(m)

hls_s30_tilejson = httpx.get(
    f"{titiler_endpoint}/rasterio/WebMercatorQuad/tilejson.json",
    params={
        # CMR Search Parameters
        "collection_concept_id": "C2021957295-LPCLOUD",  # HLSS30 2.0
        "bounding_box": ",".join(str(b) for b in hls_bbox),
        "temporal": "2025-07-01T00:00:00Z/2025-07-30T23:59:59Z",  # a whole month!
        "sort_key": "cloud_cover",  # sort by cloud cover in ascending order
        "exitwhenfull": True,
        "skipcovered": True,
        # asset selection logic
        "assets_regex": "B[0-9][0-9A-Za-z]",
        "assets": ["B04", "B03", "B02"],  # true color
        # The data is in type of Uint16 so we need to apply some
        # rescaling/color_formula in order to create PNGs
        "color_formula": "Gamma RGB 3.5 Saturation 1.2 Sigmoidal RGB 15 0.35",
    },
    timeout=None,
).json()


TileLayer(
    tiles=hls_s30_tilejson["tiles"][0],
    opacity=1,
    name="HLSS30 true color",
    overlay=True,
    show=True,
    control=True,
    attr="NASA",
).add_to(m)

LayerControl(position="topright", collapsed=False).add_to(m)

m

hls_bbox = (-93, 45, -89, 49) minzoom = 7 hls_geojson = GeoJson( mapping(box(*hls_bbox)), style_function=lambda x: { "color": "orange", "weight": 1, "fillOpacity": 0, }, name="AOI", ) m = Map( location=((hls_bbox[1] + hls_bbox[3]) / 2, (hls_bbox[0] + hls_bbox[2]) / 2), zoom_start=minzoom, ) hls_geojson.add_to(m) hls_s30_tilejson = httpx.get( f"{titiler_endpoint}/rasterio/WebMercatorQuad/tilejson.json", params={ # CMR Search Parameters "collection_concept_id": "C2021957295-LPCLOUD", # HLSS30 2.0 "bounding_box": ",".join(str(b) for b in hls_bbox), "temporal": "2025-07-01T00:00:00Z/2025-07-30T23:59:59Z", # a whole month! "sort_key": "cloud_cover", # sort by cloud cover in ascending order "exitwhenfull": True, "skipcovered": True, # asset selection logic "assets_regex": "B[0-9][0-9A-Za-z]", "assets": ["B04", "B03", "B02"], # true color # The data is in type of Uint16 so we need to apply some # rescaling/color_formula in order to create PNGs "color_formula": "Gamma RGB 3.5 Saturation 1.2 Sigmoidal RGB 15 0.35", }, timeout=None, ).json() TileLayer( tiles=hls_s30_tilejson["tiles"][0], opacity=1, name="HLSS30 true color", overlay=True, show=True, control=True, attr="NASA", ).add_to(m) LayerControl(position="topright", collapsed=False).add_to(m) m

Out[15]:

Make this Notebook Trusted to load map: File -> Trust Notebook

To demonstrate why the sort_key=cloud_cover query parameter is important, look at what happens if we skip it and just let the query return granules regardless of cloud cover.

In [16]:

m = Map(
    location=((hls_bbox[1] + hls_bbox[3]) / 2, (hls_bbox[0] + hls_bbox[2]) / 2),
    zoom_start=minzoom,
)

hls_geojson.add_to(m)

hls_s30_tilejson = httpx.get(
    f"{titiler_endpoint}/rasterio/WebMercatorQuad/tilejson.json",
    params={
        # CMR Search Parameters
        "collection_concept_id": "C2021957295-LPCLOUD",  # HLSS30 2.0
        "bounding_box": ",".join(str(b) for b in hls_bbox),
        "temporal": "2025-07-01T00:00:00Z/2025-07-30T23:59:59Z",  # a whole month!
        "exitwhenfull": True,
        "skipcovered": True,
        # asset selection logic
        "assets_regex": "B[0-9][0-9A-Za-z]",
        "assets": ["B04", "B03", "B02"],  # true color
        # The data is in type of Uint16 so we need to apply some
        # rescaling/color_formula in order to create PNGs
        "color_formula": "Gamma RGB 3.5 Saturation 1.2 Sigmoidal RGB 15 0.35",
    },
    timeout=None,
).json()


TileLayer(
    tiles=hls_s30_tilejson["tiles"][0],
    opacity=1,
    name="HLSS30 true color",
    overlay=True,
    show=True,
    control=True,
    attr="NASA",
).add_to(m)

LayerControl(position="topright", collapsed=False).add_to(m)

m

m = Map( location=((hls_bbox[1] + hls_bbox[3]) / 2, (hls_bbox[0] + hls_bbox[2]) / 2), zoom_start=minzoom, ) hls_geojson.add_to(m) hls_s30_tilejson = httpx.get( f"{titiler_endpoint}/rasterio/WebMercatorQuad/tilejson.json", params={ # CMR Search Parameters "collection_concept_id": "C2021957295-LPCLOUD", # HLSS30 2.0 "bounding_box": ",".join(str(b) for b in hls_bbox), "temporal": "2025-07-01T00:00:00Z/2025-07-30T23:59:59Z", # a whole month! "exitwhenfull": True, "skipcovered": True, # asset selection logic "assets_regex": "B[0-9][0-9A-Za-z]", "assets": ["B04", "B03", "B02"], # true color # The data is in type of Uint16 so we need to apply some # rescaling/color_formula in order to create PNGs "color_formula": "Gamma RGB 3.5 Saturation 1.2 Sigmoidal RGB 15 0.35", }, timeout=None, ).json() TileLayer( tiles=hls_s30_tilejson["tiles"][0], opacity=1, name="HLSS30 true color", overlay=True, show=True, control=True, attr="NASA", ).add_to(m) LayerControl(position="topright", collapsed=False).add_to(m) m

Out[16]:

Make this Notebook Trusted to load map: File -> Trust Notebook

It is also possible to apply a cloud cover filter to limit the returned granules to a defined range of cloud cover values. This will often result in areas being completely uncovered by rendered data, so use with caution.

In [17]:

m = Map(
    location=((hls_bbox[1] + hls_bbox[3]) / 2, (hls_bbox[0] + hls_bbox[2]) / 2),
    zoom_start=minzoom,
)

hls_geojson.add_to(m)

hls_s30_tilejson = httpx.get(
    f"{titiler_endpoint}/rasterio/WebMercatorQuad/tilejson.json",
    params={
        # CMR Search Parameters
        "collection_concept_id": "C2021957295-LPCLOUD",  # HLSS30 2.0
        "bounding_box": ",".join(str(b) for b in hls_bbox),
        "temporal": "2025-07-01T00:00:00Z/2025-07-30T23:59:59Z",  # a whole month!
        "cloud_cover": "0,20",
        "exitwhenfull": True,
        "skipcovered": True,
        # asset selection logic
        "assets_regex": "B[0-9][0-9A-Za-z]",
        "assets": ["B04", "B03", "B02"],  # true color
        # The data is in type of Uint16 so we need to apply some
        # rescaling/color_formula in order to create PNGs
        "color_formula": "Gamma RGB 3.5 Saturation 1.2 Sigmoidal RGB 15 0.35",
    },
    timeout=None,
).json()


TileLayer(
    tiles=hls_s30_tilejson["tiles"][0],
    opacity=1,
    name="HLSS30 true color",
    overlay=True,
    show=True,
    control=True,
    attr="NASA",
).add_to(m)

LayerControl(position="topright", collapsed=False).add_to(m)

m

m = Map( location=((hls_bbox[1] + hls_bbox[3]) / 2, (hls_bbox[0] + hls_bbox[2]) / 2), zoom_start=minzoom, ) hls_geojson.add_to(m) hls_s30_tilejson = httpx.get( f"{titiler_endpoint}/rasterio/WebMercatorQuad/tilejson.json", params={ # CMR Search Parameters "collection_concept_id": "C2021957295-LPCLOUD", # HLSS30 2.0 "bounding_box": ",".join(str(b) for b in hls_bbox), "temporal": "2025-07-01T00:00:00Z/2025-07-30T23:59:59Z", # a whole month! "cloud_cover": "0,20", "exitwhenfull": True, "skipcovered": True, # asset selection logic "assets_regex": "B[0-9][0-9A-Za-z]", "assets": ["B04", "B03", "B02"], # true color # The data is in type of Uint16 so we need to apply some # rescaling/color_formula in order to create PNGs "color_formula": "Gamma RGB 3.5 Saturation 1.2 Sigmoidal RGB 15 0.35", }, timeout=None, ).json() TileLayer( tiles=hls_s30_tilejson["tiles"][0], opacity=1, name="HLSS30 true color", overlay=True, show=True, control=True, attr="NASA", ).add_to(m) LayerControl(position="topright", collapsed=False).add_to(m) m

Out[17]:

Make this Notebook Trusted to load map: File -> Trust Notebook

In [ ]: