STAC Item Creation

Once you have created a collection (see stac-collection-creation.ipynb) use this notebook to create items for that collection.

Run this notebook

You can launch this notebook in VEDA JupyterHub by clicking the link below.

Launch in VEDA JupyterHub (requires access)

Learn more

Inside the Hub

This notebook was written on a VEDA JupyterHub instance

See (VEDA Analytics JupyterHub Access)[https://nasa-impact.github.io/veda-docs/veda-jh-access.html] for information about how to gain access.

Outside the Hub

You are welcome to run this anywhere you like (Note: alternatively you can run this on https://daskhub.veda.smce.nasa.gov/, MAAP, locally, …), just make sure that the data is accessible, or get in contact with the VEDA team to enable access.

Install extra packages

!pip install -U rio_stac parse xpystac pystac nbss-upload --quiet

from datetime import datetime

import rio_stac
import xarray as xr
from parse import parse

Create pystac.Item

In this section we will be creating a pystac.Item object. This is the part of that notebook that you should update.

Declare constants

Start by declaring some string fields.

COLLECTION_ID = "no2-monthly-diff"
ITEM_ID = "OMI_trno2_0.10x0.10_202212_Col3_V4.nc"
SOURCE = "s3://veda-data-store-staging/no2-monthly-diff/OMI_trno2_0.10x0.10_202212_Col3_V4.nc.tif"

Calculate datetime

Create a function that calculates datetime when given an item_id. You can change this to depend on the source instead if that works better.

def datetime_func(item_id: str) -> datetime:
    """Given the item_id, figure out the datetime"""

    fields = parse("OMI_trno2_0.10x0.10_{year:4}{month:2}_Col3_V4.nc", item_id)
    year = int(fields["year"])
    month = int(fields["month"])
    day = 1
    return datetime(year, month, day)

Test out the datetime function:

datetime_func(ITEM_ID)

datetime.datetime(2022, 12, 1, 0, 0)

Put it together

Now take your constants and datetime function and create the STAC Item using rio_stac.

item = rio_stac.stac.create_stac_item(
    id=ITEM_ID,
    source=SOURCE,
    collection=COLLECTION_ID,
    input_datetime=datetime_func(ITEM_ID),
    with_proj=True,
    with_raster=True,
    asset_name="cog_default",
    asset_roles=["data", "layer"],
    asset_media_type="image/tiff; application=geotiff; profile=cloud-optimized",
)

Try it out!

Now that you have an item you can try it out and make sure it looks good and passes validation checks.

item.validate()

['https://schemas.stacspec.org/v1.0.0/item-spec/json-schema/item.json',
 'https://stac-extensions.github.io/projection/v1.1.0/schema.json',
 'https://stac-extensions.github.io/raster/v1.1.0/schema.json']

item.to_dict()

{'type': 'Feature',
 'stac_version': '1.0.0',
 'id': 'OMI_trno2_0.10x0.10_202212_Col3_V4.nc',
 'properties': {'proj:epsg': 4326,
  '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:bbox': [-180.0, -90.0, 180.0, 90.0],
  'proj:shape': [1800, 3600],
  'proj:transform': [0.1, 0.0, -180.0, 0.0, -0.1, 90.0, 0.0, 0.0, 1.0],
  'proj:projjson': {'$schema': 'https://proj.org/schemas/v0.7/projjson.schema.json',
   'type': 'GeographicCRS',
   'name': 'WGS 84',
   'datum': {'type': 'GeodeticReferenceFrame',
    'name': 'World Geodetic System 1984',
    'ellipsoid': {'name': 'WGS 84',
     'semi_major_axis': 6378137,
     'inverse_flattening': 298.257223563}},
   'coordinate_system': {'subtype': 'ellipsoidal',
    'axis': [{'name': 'Geodetic latitude',
      'abbreviation': 'Lat',
      'direction': 'north',
      'unit': 'degree'},
     {'name': 'Geodetic longitude',
      'abbreviation': 'Lon',
      'direction': 'east',
      'unit': 'degree'}]},
   'id': {'authority': 'EPSG', 'code': 4326}},
  'proj:wkt2': 'GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]',
  'datetime': '2022-12-01T00:00:00Z'},
 'geometry': {'type': 'Polygon',
  'coordinates': [[(-180.0, -90.0),
    (180.0, -90.0),
    (180.0, 90.0),
    (-180.0, 90.0),
    (-180.0, -90.0)]]},
 'links': [{'rel': 'collection',
   'href': 'no2-monthly-diff',
   'type': 'application/json'}],
 'assets': {'cog_default': {'href': 's3://veda-data-store-staging/no2-monthly-diff/OMI_trno2_0.10x0.10_202212_Col3_V4.nc.tif',
   'type': 'image/tiff; application=geotiff; profile=cloud-optimized',
   'raster:bands': [{'data_type': 'float32',
     'scale': 1.0,
     'offset': 0.0,
     'sampling': 'area',
     'nodata': -1.2676506002282294e+30,
     'statistics': {'mean': 12233282717799.0,
      'minimum': -1.30282195779584e+16,
      'maximum': 2.082349180465971e+16,
      'stddev': 416857512760678.5,
      'valid_percent': 82.7056884765625},
     'histogram': {'count': 11,
      'min': -1.30282195779584e+16,
      'max': 2.082349180465971e+16,
      'buckets': [20, 138, 881, 421049, 11300, 203, 20, 3, 1, 1]}}],
   'roles': ['data', 'layer']}},
 'bbox': [-180.0, -90.0, 180.0, 90.0],
 'stac_extensions': ['https://stac-extensions.github.io/projection/v1.1.0/schema.json',
  'https://stac-extensions.github.io/raster/v1.1.0/schema.json'],
 'collection': 'no2-monthly-diff'}

Plot it (optional)

Create a quick visual to make sure that data loads and visualizes properly.

data = xr.open_dataset(item).cog_default.isel(time=0)
data

<xarray.DataArray 'cog_default' (latitude: 1800, longitude: 3600)> Size: 26MB
[6480000 values with dtype=float32]
Coordinates:
  * latitude     (latitude) float64 14kB 89.95 89.85 89.75 ... -89.85 -89.95
  * longitude    (longitude) float64 29kB -179.9 -179.8 -179.8 ... 179.9 180.0
    spatial_ref  int32 4B ...
    time         datetime64[ns] 8B 2022-12-01
Attributes:
    nodata:   -1.2676506002282294e+30

xarray.DataArray

'cog_default'

• latitude: 1800
• longitude: 3600

• ...

  [6480000 values with dtype=float32]

• Coordinates: (4)
  • latitude
    (latitude)
    float64
    89.95 89.85 89.75 ... -89.85 -89.95

    units :

    degrees_north

    resolution :

    -0.1

    crs :

    EPSG:4326

    array([ 89.95,  89.85,  89.75, ..., -89.75, -89.85, -89.95])

  • longitude
    (longitude)
    float64
    -179.9 -179.8 ... 179.9 180.0

    units :

    degrees_east

    resolution :

    0.1

    crs :

    EPSG:4326

    array([-179.95, -179.85, -179.75, ...,  179.75,  179.85,  179.95])

  • spatial_ref
    ()
    int32
    ...

    spatial_ref :

    GEOGCRS["WGS 84",ENSEMBLE["World Geodetic System 1984 ensemble",MEMBER["World Geodetic System 1984 (Transit)"],MEMBER["World Geodetic System 1984 (G730)"],MEMBER["World Geodetic System 1984 (G873)"],MEMBER["World Geodetic System 1984 (G1150)"],MEMBER["World Geodetic System 1984 (G1674)"],MEMBER["World Geodetic System 1984 (G1762)"],MEMBER["World Geodetic System 1984 (G2139)"],ELLIPSOID["WGS 84",6378137,298.257223563,LENGTHUNIT["metre",1]],ENSEMBLEACCURACY[2.0]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433]],CS[ellipsoidal,2],AXIS["geodetic latitude (Lat)",north,ORDER[1],ANGLEUNIT["degree",0.0174532925199433]],AXIS["geodetic longitude (Lon)",east,ORDER[2],ANGLEUNIT["degree",0.0174532925199433]],USAGE[SCOPE["Horizontal component of 3D system."],AREA["World."],BBOX[-90,-180,90,180]],ID["EPSG",4326]]

    crs_wkt :

    GEOGCRS["WGS 84",ENSEMBLE["World Geodetic System 1984 ensemble",MEMBER["World Geodetic System 1984 (Transit)"],MEMBER["World Geodetic System 1984 (G730)"],MEMBER["World Geodetic System 1984 (G873)"],MEMBER["World Geodetic System 1984 (G1150)"],MEMBER["World Geodetic System 1984 (G1674)"],MEMBER["World Geodetic System 1984 (G1762)"],MEMBER["World Geodetic System 1984 (G2139)"],ELLIPSOID["WGS 84",6378137,298.257223563,LENGTHUNIT["metre",1]],ENSEMBLEACCURACY[2.0]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433]],CS[ellipsoidal,2],AXIS["geodetic latitude (Lat)",north,ORDER[1],ANGLEUNIT["degree",0.0174532925199433]],AXIS["geodetic longitude (Lon)",east,ORDER[2],ANGLEUNIT["degree",0.0174532925199433]],USAGE[SCOPE["Horizontal component of 3D system."],AREA["World."],BBOX[-90,-180,90,180]],ID["EPSG",4326]]

    semi_major_axis :

    6378137.0

    semi_minor_axis :

    6356752.314245179

    inverse_flattening :

    298.257223563

    reference_ellipsoid_name :

    WGS 84

    longitude_of_prime_meridian :

    0.0

    prime_meridian_name :

    Greenwich

    geographic_crs_name :

    WGS 84

    horizontal_datum_name :

    World Geodetic System 1984 ensemble

    grid_mapping_name :

    latitude_longitude

    GeoTransform :

    -180 0.100000000000000005551115 0 90 0 -0.100000000000000005551115

    [1 values with dtype=int32]

  • time
    ()
    datetime64[ns]
    2022-12-01

    array('2022-12-01T00:00:00.000000000', dtype='datetime64[ns]')

• Indexes: (2)
  • latitude
    PandasIndex

    PandasIndex(Index([             89.95,  89.85000000000001,              89.75,
                        89.65,              89.55,              89.45,
            89.35000000000001,              89.25,              89.15,
                        89.05,
           ...
                       -89.05, -89.15000000000002, -89.25000000000001,
           -89.35000000000001,             -89.45,             -89.55,
           -89.65000000000002, -89.75000000000001, -89.85000000000001,
                       -89.95],
          dtype='float64', name='latitude', length=1800))

  • longitude
    PandasIndex

    PandasIndex(Index([            -179.95,             -179.85,             -179.75,
           -179.64999999999998, -179.54999999999998,             -179.45,
                       -179.35,             -179.25, -179.14999999999998,
           -179.04999999999998,
           ...
                        179.05,  179.15000000000003,  179.25000000000006,
            179.35000000000002,  179.45000000000005,              179.55,
            179.65000000000003,  179.75000000000006,  179.85000000000002,
            179.95000000000005],
          dtype='float64', name='longitude', length=3600))

• Attributes: (1)

  nodata :

  -1.2676506002282294e+30

data.hvplot("longitude", "latitude")

Upload this notebook

You can upload the notebook to anyplace you like, but one of the easiest ones is notebook sharing space. Just change the following cell from “Raw” to “Code”, run it and copy the output link.

Before uploading make sure: 1) you have not hard-coded any secrets or access keys. 2) you have saved this notebook. Hint (ctrl+s) will do it

!nbss-upload new-item.ipynb