Order Service Example

Order Service Example#

This notebook will demonstrate how to submit a request to the Hydrosat Fusion Hub Order Service and load the data.

import json
import requests

import os
import base64
from pystac_client import Client
import time
import urllib
import pandas as pd

This next cell opens a file creds.json which you will need to create in the same directory as the notebook. The format of the file should be:

{
"username":"your_username",
"password":"your_password"
}

with open('creds.json') as f:
    creds = json.loads(f.read())

and you have updated with your username and password. The next cells extract your credentials to pass to the order service endpoint.

userpass = f"{creds['username']}:{creds['password']}"
b64 = base64.b64encode(userpass.encode()).decode()
headers = {'Authorization':'Basic ' + b64}

cat_url = 'https://fusion-stac.hydrosat.com/'
catalog = Client.open(cat_url, headers)

Let’s grab some items that will reflect parameters for a request to the order service so we can grab the appropriate coordinate reference system to pass to the order service payload.

bbox = [-120.0373077392578, 36.16836821871061, -119.99130249023436, 36.19109202182454]
start_date ="2021-08-07T00:00:00Z"
end_date = "2021-08-17T23:59:59Z"    
collections = ["starfm_predictions_modis_landsat", "pydms_sharpened_landsat_hires" ]

search = catalog.search(
    collections = collections,
    bbox = bbox,
    datetime = [start_date, end_date],
    max_items = 500
)

items = list(search.items())
print(f"number of items returned: {len(items)}")

number of items returned: 10

# check the coordinate system for one of the items:
epsg = items[0].to_dict()['assets']['lst']['proj:epsg']
print(epsg)

We will use the CRS printed out above, as well as the start_date, end_date, bbox, and the geometry information from the bbox, to create a payload which needs to be passed the order service endpoint. Change the start_date, end_date, projection, and assets to match what is needed for the order. Please visit our documentation page for additional information on the order service.

Note #1: for the order service, start date needs to be formatted as “YYYY-MM-DD”, so we will use string formatting to extract that information from the format we passed to the STAC API.

Note #2: To get the geometry coordinates from the bbox list of coordinates, we will use shapely to convert that list to a polygon and get the coordinates from the polygon boundary. We will use that polygon to see where on the map we are, too.

from shapely.geometry import shape, box
import geopandas as gpd
import folium 
from pyproj.crs import CRS

polygon = box(*bbox)
coords = list(polygon.boundary.coords)

# Use WGS 84 (epsg:4326) as the geographic coordinate system
df = gpd.GeoDataFrame({'geometry': [polygon]}, crs=CRS(4326))

p_geom = polygon.centroid
m = folium.Map(location=[p_geom.y, p_geom.x], zoom_start=13, tiles='CartoDB positron')

# add the polygon and centroid
for _, r in df.iterrows():
    # Without simplifying the representation of each polygon,
    # the map might not be displayed
    #sim_geo = gpd.GeoSeries(r['geometry']).simplify(tolerance=0.001)
    sim_geo = r[0].simplify(tolerance=0.001)
    geo_j = gpd.GeoSeries(sim_geo).to_json()
    geo_j = folium.GeoJson(data=geo_j,
                           style_function=lambda x: {'fillColor': 'orange'})
    
    geo_j.add_to(m)
    
    lat = sim_geo.centroid.y
    lon = sim_geo.centroid.x
    folium.Marker(location=[lat, lon]).add_to(m)
    
m

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

now we construct the payload as described

payload = {
    "start_date": start_date.split('T')[0],
    "end_date": end_date.split('T')[0],
    "projection": f"epsg:{epsg}",
    "resolution": 20.0,
    "output_format": "GeoTIFF",
    "collections": collections,
    "assets": [
      "lst",
      "combined_qa"
    ],
    "bbox": bbox,
    "geometry": {
      "type": "Polygon",
      "coordinates": [coords]
    }
}

order_url = 'https://fusion-stac.hydrosat.com/orders'

p = requests.post(order_url, headers=headers, json=payload)

p.json()

{'orderId': '3ab7b536-90ef-4369-a32b-af15b5aafe10',
 'url': 'https://fusion-stac.hydrosat.com/orders/3ab7b536-90ef-4369-a32b-af15b5aafe10'}

sometimes the order service takes a bit to start, so we will check for a status code of 500 and once the order service begins processing, the status will change to pending. Once it is finished, we will download the data to a local folder called ‘my_order’ (change it to whatever is appropriate).

In the order payload, we specified the combined_qa asset to be ordered in addition to lst. This example only downloads the lst asset for each item. To download the combined_qa asset, swap out ‘lst’ for ‘combined_qa’ in the assemble the assets and urls section of the code in this block.

# ping the order status until it completes
while (requests.get(p.json()['url'], headers=headers).status_code == 500) or (requests.get(p.json()['url'], headers=headers).json()['status'] == 'pending'):
    time.sleep(15)
    pass

# get the compeleted order service response, which will contain URLs for the items to download
g = requests.get(p.json()['url'], headers=headers)
res = g.json()

# specify a folder to hold the data
local_folder = 'my_order'

# create it if it does not exist
if not os.path.exists(local_folder):
    os.makedirs(local_folder)
    
# assemble the assets and urls
dl_assets = [r['assets']['lst']['href'] for r in res['output']['features']]
dl_urls =[d.split('?')[0] for d in dl_assets]

# illegal character set to check against
sets = [':','*','?','"','<','>','|']

# for each of the items, download it to the folder.
for dl_url, dl_asset in zip(dl_urls, dl_assets):
    if local_folder is not None:
        outfile = os.path.join(local_folder, os.path.basename(dl_url)) #+ '.tif')
        for char in outfile:
            if char in sets:
                outfile = outfile.replace(char,'_')
    else:
        outfile = os.path.basename(dl_url) #+ '.tif'
        for char in outfile:
            if char in sets:
                outfile = outfile.replace(char,'_')

    # download file if it doesn't exist. If you wish to overwrite, comment this if-else statement and uncomment the
    # line below it.
    if not os.path.exists(outfile):
        urllib.request.urlretrieve(dl_asset, outfile)
    else:
        pass
    
    #urllib.request.urlretrieve(dl_asset, outfile)

Once the data is downloaded, we can load it as a data cube using rioxarray. The dataset will not be time-aware immediately, so we will assign the time dimension of each file using information from the filename.

import rioxarray as rxr
import xarray as xr
from glob import glob

all_tifs = glob(f"{local_folder}/*.tif")
order_ds = xr.concat([rxr.open_rasterio(t) for t in all_tifs], dim='time')

# print out the shape of the data retrieved from the order
order_ds.shape

(10, 1, 134, 212)

# now lets plot them. If your order is large, you will want to subset the time dimension
# by calling order_ds.isel(time=[0,1,2,3])... see xarray documentation for more use of selecting data
order_ds.plot(x='x', y='y', col='time', col_wrap=5)

<xarray.plot.facetgrid.FacetGrid at 0x7f45be197250>

../_images/319c275c4749072122144913f43b3c2cd8f6d824d1738d1c028ecbbfef279887.png