Get map from COGs - NO2

Demonstrates generating a map for a given area.
Author

Leo Thomas, Julia Signell

Published

February 7, 2023

Run this notebook

You can launch this notbook using mybinder, by clicking the button below.

Binder

Approach

  1. Fetch STAC item for a particular date and collection - NO2
  2. Pass STAC item in to the raster API /stac/tilejson.json endpoint
  3. Visualize tiles using folium
import requests
import folium

Declare your collection of interest

You can discover available collections the following ways:

  • Programmatically: see example in the list-collections.ipynb notebook
  • JSON API: https://staging-stac.delta-backend.com/collections
  • STAC Browser: http://veda-staging-stac-browser.s3-website-us-west-2.amazonaws.com
STAC_API_URL = "https://openveda.cloud/api/stac"
RASTER_API_URL = "https://openveda.cloud/api/raster"

collection_id = "no2-monthly"

Fetch STAC collection

We will use requests to fetch all the metadata about the collection of interest from STAC.

collection = requests.get(f"{STAC_API_URL}/collections/{collection_id}").json()
collection
{'id': 'no2-monthly',
 'type': 'Collection',
 'links': [{'rel': 'items',
   'type': 'application/geo+json',
   'href': 'https://openveda.cloud/api/stac/collections/no2-monthly/items'},
  {'rel': 'parent',
   'type': 'application/json',
   'href': 'https://openveda.cloud/api/stac/'},
  {'rel': 'root',
   'type': 'application/json',
   'href': 'https://openveda.cloud/api/stac/'},
  {'rel': 'self',
   'type': 'application/json',
   'href': 'https://openveda.cloud/api/stac/collections/no2-monthly'}],
 'title': 'NO₂',
 'assets': {'thumbnail': {'href': 'https://thumbnails.openveda.cloud/no2--dataset-cover.jpg',
   'type': 'image/jpeg',
   'roles': ['thumbnail'],
   'title': 'Thumbnail',
   'description': 'Photo by [Mick Truyts](https://unsplash.com/photos/x6WQeNYJC1w) (Power plant shooting steam at the sky)'}},
 'extent': {'spatial': {'bbox': [[-180.0, -90.0, 180.0, 90.0]]},
  'temporal': {'interval': [['2016-01-01T00:00:00+00:00',
     '2022-12-31T00:00:00+00:00']]}},
 'license': 'MIT',
 'renders': {'dashboard': {'bidx': [1],
   'title': 'VEDA Dashboard Render Parameters',
   'assets': ['cog_default'],
   'rescale': [[0, 15000000000000000]],
   'resampling': 'bilinear',
   'color_formula': 'gamma r 1.05',
   'colormap_name': 'rdbu_r'}},
 'providers': [{'url': 'https://disc.gsfc.nasa.gov/',
   'name': 'NASA Goddard Earth Sciences Data and Information Services Center',
   'roles': ['producer', 'processor']},
  {'url': 'https://www.earthdata.nasa.gov/dashboard/',
   'name': 'NASA VEDA',
   'roles': ['host']}],
 'summaries': {'datetime': ['2016-01-01T00:00:00Z', '2023-09-30T00:00:00Z']},
 'description': 'Darker colors indicate higher nitrogen dioxide (NO₂) levels and more activity. Lighter colors indicate lower levels of NO₂ and less activity. Missing pixels indicate areas of no data most likely associated with cloud cover or snow.',
 'item_assets': {'cog_default': {'type': 'image/tiff; application=geotiff; profile=cloud-optimized',
   'roles': ['data', 'layer'],
   'title': 'Default COG Layer',
   'description': 'Cloud optimized default layer to display on map'}},
 'stac_version': '1.0.0',
 'stac_extensions': ['https://stac-extensions.github.io/item-assets/v1.0.0/schema.json',
  'https://stac-extensions.github.io/render/v1.0.0/schema.json'],
 'dashboard:is_periodic': True,
 'dashboard:time_density': 'month'}

Fetch STAC item for a particular time

We can use the search API to find the item that matches exactly our time of interest.

response = requests.post(
    f"{STAC_API_URL}/search",
    json={
        "collections": [collection_id],
        "query": {"datetime": {"eq": "2021-01-01T00:00:00"}},
        "limit": 100,
    },
).json()
items = response["features"]
len(items)
1

Let’s take a look at that one item.

item = items[0]
item
{'id': 'OMI_trno2_0.10x0.10_202101_Col3_V4.nc',
 'bbox': [-180.0, -90.0, 180.0, 90.0],
 'type': 'Feature',
 'links': [{'rel': 'collection',
   'type': 'application/json',
   'href': 'https://openveda.cloud/api/stac/collections/no2-monthly'},
  {'rel': 'parent',
   'type': 'application/json',
   'href': 'https://openveda.cloud/api/stac/collections/no2-monthly'},
  {'rel': 'root',
   'type': 'application/json',
   'href': 'https://openveda.cloud/api/stac/'},
  {'rel': 'self',
   'type': 'application/geo+json',
   'href': 'https://openveda.cloud/api/stac/collections/no2-monthly/items/OMI_trno2_0.10x0.10_202101_Col3_V4.nc'},
  {'title': 'Map of Item',
   'href': 'https://openveda.cloud/api/raster/collections/no2-monthly/items/OMI_trno2_0.10x0.10_202101_Col3_V4.nc/map?bidx=1&assets=cog_default&rescale=0%2C15000000000000000&resampling=bilinear&color_formula=gamma+r+1.05&colormap_name=rdbu_r',
   'rel': 'preview',
   'type': 'text/html'}],
 'assets': {'cog_default': {'href': 's3://veda-data-store/no2-monthly/OMI_trno2_0.10x0.10_202101_Col3_V4.nc.tif',
   'type': 'image/tiff; application=geotiff',
   'roles': ['data', 'layer'],
   'title': 'Default COG Layer',
   'proj:bbox': [-180.0, -90.0, 180.0, 90.0],
   'proj:epsg': 4326,
   'proj:wkt2': 'GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]',
   'proj:shape': [1800, 3600],
   'description': 'Cloud optimized default layer to display on map',
   'raster:bands': [{'scale': 1.0,
     'nodata': -1.2676506002282294e+30,
     'offset': 0.0,
     'sampling': 'area',
     'data_type': 'float32',
     'histogram': {'max': 35781585143857150,
      'min': -4107596126486528.0,
      'count': 11,
      'buckets': [7437, 432387, 2866, 699, 356, 207, 76, 27, 7, 1]},
     'statistics': {'mean': 367152773066762.6,
      'stddev': 961254458662885.4,
      'maximum': 35781585143857150,
      'minimum': -4107596126486528.0,
      'valid_percent': 84.69829559326172}}],
   'proj:geometry': {'type': 'Polygon',
    'coordinates': [[[-180.0, -90.0],
      [180.0, -90.0],
      [180.0, 90.0],
      [-180.0, 90.0],
      [-180.0, -90.0]]]},
   'proj:projjson': {'id': {'code': 4326, 'authority': 'EPSG'},
    'name': 'WGS 84',
    'type': 'GeographicCRS',
    'datum': {'name': 'World Geodetic System 1984',
     'type': 'GeodeticReferenceFrame',
     'ellipsoid': {'name': 'WGS 84',
      'semi_major_axis': 6378137,
      'inverse_flattening': 298.257223563}},
    '$schema': 'https://proj.org/schemas/v0.4/projjson.schema.json',
    'coordinate_system': {'axis': [{'name': 'Geodetic latitude',
       'unit': 'degree',
       'direction': 'north',
       'abbreviation': 'Lat'},
      {'name': 'Geodetic longitude',
       'unit': 'degree',
       'direction': 'east',
       'abbreviation': 'Lon'}],
     'subtype': 'ellipsoidal'}},
   'proj:transform': [0.1, 0.0, -180.0, 0.0, -0.1, 90.0, 0.0, 0.0, 1.0]},
  'rendered_preview': {'title': 'Rendered preview',
   'href': 'https://openveda.cloud/api/raster/collections/no2-monthly/items/OMI_trno2_0.10x0.10_202101_Col3_V4.nc/preview.png?bidx=1&assets=cog_default&rescale=0%2C15000000000000000&resampling=bilinear&color_formula=gamma+r+1.05&colormap_name=rdbu_r',
   'rel': 'preview',
   'roles': ['overview'],
   'type': 'image/png'}},
 'geometry': {'type': 'Polygon',
  'coordinates': [[[-180, -90],
    [180, -90],
    [180, 90],
    [-180, 90],
    [-180, -90]]]},
 'collection': 'no2-monthly',
 'properties': {'end_datetime': '2021-01-31T00:00:00+00:00',
  'start_datetime': '2021-01-01T00:00:00+00:00'},
 'stac_version': '1.0.0',
 'stac_extensions': ['https://stac-extensions.github.io/raster/v1.1.0/schema.json',
  'https://stac-extensions.github.io/projection/v1.1.0/schema.json']}
item_stats = item['assets']['cog_default']['raster:bands'][0]['statistics']
rescale_values = item_stats['minimum'], item_stats['maximum']

Use /stac/tilejson.json to get tiles

We pass the, item id, collection name, and the rescale_values in to the RASTER API endpoint and get back a tile.

tiles = requests.get(
    f"{RASTER_API_URL}/collections/{collection_id}/items/{item['id']}/tilejson.json?"
    "&assets=cog_default"
    "&color_formula=gamma+r+1.05&colormap_name=rdbu_r"
    f"&rescale={rescale_values[0]},{rescale_values[1]}",
).json()
tiles
{'tilejson': '2.2.0',
 'version': '1.0.0',
 'scheme': 'xyz',
 'tiles': ['https://openveda.cloud/api/raster/collections/no2-monthly/items/OMI_trno2_0.10x0.10_202101_Col3_V4.nc/tiles/WebMercatorQuad/{z}/{x}/{y}@1x?assets=cog_default&color_formula=gamma+r+1.05&colormap_name=rdbu_r&rescale=-4107596126486528.0%2C35781585143857150'],
 'minzoom': 0,
 'maxzoom': 24,
 'bounds': [-180.0, -90.0, 180.0, 90.0],
 'center': [0.0, 0.0, 0]}

With that tile url in hand we can create a simple visualization using folium.

folium.Map(
    tiles=tiles["tiles"][0],
    min_zoom=3,
    attr="VEDA",
)
Make this Notebook Trusted to load map: File -> Trust Notebook