ODC-geo with H5NetCDF and earthaccess

import earthaccess
import xarray as xr
from odc.geo.geobox import GeoBox
from odc.geo.xr import xr_reproject

def warp_resample():
    bucket = "podaac-ops-cumulus-protected"
    input_uri = "MUR-JPL-L4-GLOB-v4.1/20020601090000-JPL-L4_GHRSST-SSTfnd-MUR-GLOB-v02.0-fv04.1.nc"
    variable = "analysed_sst"
    src = f"s3://{bucket}/{input_uri}"
    dstSRS = "EPSG:3857"
    srcSRS = "EPSG:4326"
    width = height = 256
    te = [
        -20037508.342789244,
        -20037508.342789244,
        20037508.342789244,
        20037508.342789244,
    ]
    earthaccess.login()
    fs = earthaccess.get_s3fs_session(daac="PODAAC")
    fsspec_caching = {
        "cache_type": "none",
    }
    gbox = GeoBox.from_bbox(te, dstSRS, shape=(height, width))
    f = fs.open(src, **fsspec_caching)
    da = xr.open_dataset(f, engine="h5netcdf")[variable]
    da = da.odc.assign_crs(srcSRS)
    return xr_reproject(da, gbox)

if __name__ == "__main__":
    da = warp_resample()