03. Build a Time Series of SWOT Simulated Data

Import libraries

In additional to libraries from the first tutorial, also import libraries needed to submit CMR requests.

import xarray as xr
import numpy as np
from IPython.display import display, JSON
from datetime import datetime, timedelta, time
import os

# highlight the harmony-py library
from harmony import BBox, Client, Collection, Request, Environment, LinkType 

import matplotlib.pyplot as plt
import cartopy.crs as ccrs

%matplotlib inline

# Additional libraries compared to the first tutorial.
from urllib import request
import json
import requests
import sys
import shutil
from urllib.parse import urlencode

Let’s start up the client from the harmony-py library and define the CMR url.

harmony_client = Client(env=Environment.PROD)
cmr_root = 'cmr.earthdata.nasa.gov'

KaRIn CalVal

Search by cycle and pass using CMR

CMR Search: Number of item returned is limited to 10,000 (or 1 million if targeting collections) https://cmr.earthdata.nasa.gov/search/site/docs/search/api.html#paging-details

provider = 'POCLOUD'
short_name = 'SWOT_SIMULATED_L2_KARIN_SSH_GLORYS_CALVAL_V1'

target_cycle = list(range(1,101))
target_pass = [2, 17]
page_size = 2000

granuleIDs = []
params_without_cycle = [
    ('page_size', page_size),
    ('sort_key', "-start_date"),
    ('provider', provider),
    ('ShortName', short_name),
    #('collection_concept_id', ccid), 
    #('token', token),
    #('cycle', targetcycle), # can only specify one cycle at a time when specifying passes 
    ('passes[0][pass]', target_pass[0]), 
    ('passes[1][pass]', target_pass[1]),
]

for v in target_cycle:
    params = [("cycle[]", v)]
    params.extend(params_without_cycle)
    # print(params)
    query = urlencode(params)
    cmr_url = f"https://{cmr_root}/search/granules.umm_json?{query}"
    # print(cmr_url)
    response = requests.get(cmr_url)
    response_body = response.json()
    
    for itm in response_body['items']:
        granuleIDs.append(itm['meta']['concept-id'])

len(granuleIDs) # Note the 200-granule limit

Then perform a spatial subset and download using Harmony

On the back end, the subsetting part of Harmony is powered by L2SS-py

# collection = Collection(id=ccid)
collection = Collection(id=short_name)

# start_day = datetime(2015,4,15,0,0,0)
# end_day = datetime(2015,4,17,0,0,0)

request = Request(
    collection=collection,
    spatial=BBox(-140, 20, -100, 50), # [20-50N], [-140W, -100W] CA Current
    #variables=[],
    # temporal={
    #     'start': start_day,
    #     'stop': end_day # goal: try up to 21 days at least,
    #},
    granule_id=granuleIDs,
)

request.is_valid()
print(harmony_client.request_as_curl(request))
job_id = harmony_client.submit(request)
print(f'Job ID: {job_id}')
harmony_client.status(job_id)
# results = harmony_client.result_urls(job_id, link_type=LinkType.s3)
# urls = list(results)
# print(urls)
# create a new folder to put the subsetted data in
os.makedirs("swot_ocean",exist_ok = True)
futures = harmony_client.download_all(job_id, directory='./swot_ocean', overwrite=True)
file_names = [f.result() for f in futures]
sorted(file_names)
ds = xr.open_mfdataset(sorted(file_names),combine='nested',concat_dim='num_lines')
ds
# Plot only a pair of passes at a time
i_time = np.arange(0,1725*2)+1725*90

fig = plt.figure(figsize=[11,7]) 
ax = plt.axes(projection=ccrs.PlateCarree())
ax.coastlines()
ax.set_extent([-150, -90, 10, 60])
plt.scatter(ds.longitude[i_time,:], ds.latitude[i_time,:], lw=1, c=ds.ssha_karin[i_time,:])
plt.colorbar(label='SSHA (m)')
plt.clim(-1,1)
plt.show()