Climate DT Parameter Series Plot- Data Access using DEDL HDA

import destinelab as deauth
import json
import datetime
import importlib.metadata

import requests
from requests.adapters import HTTPAdapter
from urllib3.util.retry import Retry
import os
from getpass import getpass
from tqdm import tqdm
import time
from datetime import datetime
from urllib.parse import unquote
from time import sleep
from IPython.display import JSON
import ipywidgets as w

DESP_USERNAME = input("Please input your DESP username: ")
DESP_PASSWORD = getpass("Please input your DESP password: ")

auth = deauth.AuthHandler(DESP_USERNAME, DESP_PASSWORD)
access_token = auth.get_token()
if access_token is not None:
    print("DEDL/DESP Access Token Obtained Successfully")
else:
    print("Failed to Obtain DEDL/DESP Access Token")

auth_headers = {"Authorization": f"Bearer {access_token}"}

import importlib
installed_version = importlib.metadata.version("destinelab")
version_number = installed_version.split('.')[1]
if((int(version_number) >= 8 and float(installed_version) < 1) or float(installed_version) >= 1):
    auth.is_DTaccess_allowed(access_token)

HDA_STAC_ENDPOINT="https://hda.data.destination-earth.eu/stac/v2"
print("STAC endpoint: ", HDA_STAC_ENDPOINT)

HDA_DISCOVERY_ENDPOINT = HDA_STAC_ENDPOINT+'/collections'
print("HDA discovery endpoint: ", HDA_DISCOVERY_ENDPOINT)

discovery_json=(requests.get(HDA_DISCOVERY_ENDPOINT,params = {"q": '"Climate Change Adaptation Digital Twin","Future Projection","IFS-NEMO"'}).json())

#print("Result from the free text search on the HDA Discovery API : ")
#JSON(discovery_json)

print("The discovery result give us:\nthe collection ID : ", discovery_json["collections"][0].get("id"))
print("\nIts time extension : ", discovery_json["collections"][0].get("extent").get("temporal").get("interval"))
#print("\nThe available parameters: ", discovery_json["collections"][0].get("cube:variables").keys())
print("\nThe available parameters: ")
print("(Type in the text box to narrow the list. When only one item remains, its details will appear)")


keys = sorted(discovery_json["collections"][0]["cube:variables"])
txt = w.Text(description="")
out = w.Output()
def run(_):
    q = txt.value.lower()
    with out:
        out.clear_output(); 
        matches = [k for k in keys if q in k.lower()]
        for k in matches:
            print(k)
            if len(matches)==1:            
                var = discovery_json["collections"][0]["cube:variables"][k]
                print(var["description"], "\n")
                print(json.dumps(var["attrs"], indent=2))
                #print(json.dumps(discovery_json["collections"][0]["cube:variables"][k], indent=2))

                
txt.observe(run, names="value"); display(txt, out); run(None)

#print("\nThe information related to each parameter, e.g. the 2 metre temperature: ")
#JSON(discovery_json["collections"][0].get("cube:variables").get("2_metre_temperature(sfc)"))

COLLECTION_ID="EO.ECMWF.DAT.D1.DT_CLIMATE.G1.SCENARIOMIP_SSP3-7.0_IFS-NEMO.R1"

datechoice = "20200731"

filter_params = {
        "resolution": "high",      # standard/ high 
        "time": "0000",            # choose the hourly slot(s)
        "type": "fc",              # fixed forecasted fields
        "levtype": "sfc",          # Surface fields (levtype=sfc),
        "param": "167",            # 2m Temperature parameter
        "date":datechoice
    }

# Initialize a list to store filenames
filenames = []

# Define start and end years
start_year = 2020
end_year = 2039

#timeout and step for polling (sec)
TIMEOUT = 300
STEP = 1
ONLINE_STATUS = "online"

# Loop 
for year in range(start_year, end_year+1,4):
    # Create a datetime object 
    obsdate = datetime(year, 7, 31)
    datechoice = obsdate.strftime("%Y%m%d")
    filter_params["date"]=datechoice
    response = requests.post(f"{HDA_STAC_ENDPOINT}/collections/{COLLECTION_ID}/order", json=filter_params, headers=auth_headers)

    if response.status_code != 200:
        print(response.content)
    response.raise_for_status()

    ordered_item = response.json()

    product_id = ordered_item["id"]
    storage_tier = ordered_item["properties"].get("storage:tier", "online")
    order_status = ordered_item["properties"].get("order:status", "unknown")
    federation_backend = ordered_item["properties"].get("federation:backends", [None])[0]

    print(f"Product ordered: {product_id}")
    print(f"Provider: {federation_backend}")
    print(f"Storage tier: {storage_tier} (product must have storage tier \"online\" to be downloadable)")
    print(f"Order status: {order_status}")      

    self_url = f"{HDA_STAC_ENDPOINT}/collections/{COLLECTION_ID}/items/{product_id}"
    item = {}

    for i in range(0, TIMEOUT, STEP):
        print(f"Polling {i + 1}/{TIMEOUT // STEP}")

        response = requests.get(self_url, headers=auth_headers)
        if response.status_code != 200:
            print(response.content)
        response.raise_for_status()
        item = response.json()

        storage_tier = item["properties"].get("storage:tier", ONLINE_STATUS)

        if storage_tier == ONLINE_STATUS:
            download_url = item["assets"]["downloadLink"]["href"]
            print("Product is ready to be downloaded.")
            print(f"Asset URL: {download_url}")
            break
        sleep(STEP)
    else:
        order_status = item["properties"].get("order:status", "unknown")
        print(f"We could not download the product after {TIMEOUT // STEP} tries. Current order status is {order_status}")
    
    response = requests.get(download_url, stream=True, headers=auth_headers)
    response.raise_for_status()

    content_disposition = response.headers.get('Content-Disposition')
    total_size = int(response.headers.get("content-length", 0))
    if content_disposition:
        filename = content_disposition.split('filename=')[1].strip('"')
        filename = unquote(filename)
    else:
        filename = os.path.basename(url)

    # Open a local file in binary write mode and write the content
    print(f"downloading {filename}")

    with tqdm(total=total_size, unit="B", unit_scale=True) as progress_bar:
        with open(filename, 'wb') as f:
            for data in response.iter_content(1024):
                progress_bar.update(len(data))
                f.write(data)
               
    # Add the filename to the list
    filenames.append(filename)

import earthkit.data
import earthkit.plots

# Iterate over filenames
data=[]
for filename in filenames:
    data.append( earthkit.data.from_source("file", filename))

STYLE = earthkit.plots.styles.Style(
    colors="Spectral_r",
    levels=range(-20, 45),
    units="celsius",
    # Extend the colorbar at both ends
    extend="both",
)

figure = earthkit.plots.Figure(rows=3, columns=2, size=(10, 10))

for i, year in enumerate(range(start_year, end_year+1,4)):
    subplot = figure.add_map()
    subplot.contourf(data[i], style=STYLE)
    subplot.title(f"{year} - Climate DT - IFS NEMO {{short_name!u}} {{time:%d %B %Y}}")


figure.title("{variable_name} in {time:%B %d} of {time:%Y}", fontsize=14)

figure.legend(label="{variable_name!l} ({units})")

figure.show()

for filename in filenames:
    os.remove(filename)