mdsumner

vrt = gdalxarray.warp("/vsicurl/https://projects.pawsey.org.au/idea-gebco-tif/GEBCO_2024.tif", 
     cutlineDSName =  "/vsicurl/https://github.com/mdsumner/geoboundaries/releases/download/latest/geoBoundariesCGAZ_ADM0.parquet",      
     cutlineSQL = "SELECT shapeGroup FROM geoBoundariesCGAZ_ADM0 WHERE shapeGroup IN ('AUS')",  
     cropToCutline = True, resolution = 0.1)
ds = xarray.open_dataset(vrt, engine = "gdalxarray", multidim = False)
ds
<xarray.Dataset> Size: 858kB
Dimensions: (band: 1, y: 447, x: 947)

mdsumner

## list the structure elements (Icechunk) gdal vsi ls  "/vsis3/earthmover-icechunk-era5/icechunkV2" --config AWS_NO_SIGN_REQUEST YES --config AWS_REGION us-east-1
repo
chunks
manifests
overwritten
snapshots
transactions

mdsumner

gdal mdim info "/vsis3/earthmover-icechunk-era5/icechunkV2" --config AWS_NO_SIGN_REQUEST YES --config AWS_REGION us-east-1 | head -n 120
ERROR 6: Unsupported codec: numcodecs.pcodec
{
  "type": "group",
  "driver": "Icechunk",
  "name": "/",
  "groups": {
    "500hPa": {
 "groups": {

mdsumner

from osgeo import gdal
gdal.SetConfigOption("AWS_NO_SIGN_REQUEST", "YES"); gdal.SetConfigOption("AWS_REGION", "us-west-2")

import gdalxarray
from gdalxarray import GDALBackendEntrypoint
backend = GDALBackendEntrypoint()
xds = backend.open_dataset(
     "/vsis3/dynamical-ecmwf-aifs-single/ecmwf-aifs-single-forecast/v0.1.0.icechunk"
   )

mdsumner

  xyll <- expand.grid(x= seq(100, 150, length.out = 15), y = seq(-70, -40, length.out = 12))

library(tissot)
ti <- tissot(xyll, "EPSG:3031")
atan2(ti$dy_dlam, ti$dx_dlam) 
#>   [1] -1.745330 -1.807663 -1.869996 -1.932330 -1.994663 -2.056996 -2.119329
#>   [8] -2.181662 -2.243996 -2.306329 -2.368662 -2.430995 -2.493328 -2.555662
#>  [15] -2.617995 -1.745330 -1.807663 -1.869996 -1.932330 -1.994663 -2.056996
#>  [22] -2.119329 -2.181662 -2.243996 -2.306329 -2.368662 -2.430995 -2.493328

mdsumner

finally, do the basic thing - this is like one line of basic generic config and code in GDAL ...

import fsspec
import xarray as xr
import zarr

pawsey_s3_config = {
     "anon": True,
     "client_kwargs": {

mdsumner

f <- c("/rdsi/PUBLIC/raad/data/data.marine.copernicus.eu/SEALEVEL_GLO_PHY_L4_NRT_008_046/cmems_obs-sl_glo_phy-ssh_nrt_allsat-l4-duacs-0.25deg_P1D_202311/2022/01/nrt_global_allsat_phy_l4_20220111_20220111.nc", 
  "/rdsi/PUBLIC/raad/data/data.marine.copernicus.eu/SEALEVEL_GLO_PHY_L4_NRT_008_046/cmems_obs-sl_glo_phy-ssh_nrt_allsat-l4-duacs-0.25deg_P1D_202311/2022/01/nrt_global_allsat_phy_l4_20220112_20220118.nc",
"/rdsi/PUBLIC/raad/data/data.marine.copernicus.eu/SEALEVEL_GLO_PHY_L4_NRT_008_046/cmems_obs-sl_glo_phy-ssh_nrt_allsat-l4-duacs-0.25deg_P1D_202311/2022/01/nrt_global_allsat_phy_l4_20220113_20220116.nc")

fs::file_size(f)
#9.05M 7.49M 9.03M

mdsumner

This local example in the Usage for VirtualiZarr is (still) broken:

https://virtualizarr.readthedocs.io/en/latest/usage.html#__tabbed_1_7

I personally find it impossible to remember all this rigmarole (what's wrong with a file path?), so I need a reference example that actually works:

## local file 'data.nc' in your wd e.g. 
# wget -O data.nc https://www.ncei.noaa.gov/data/sea-surface-temperature-optimum-interpolation/v2.1/access/avhrr/198109/oisst-avhrr-v02r01.19810901.nc

mdsumner

+proj=stere +a=6378137.0 +b=6356752.3 +lat_0=-90 +lon_0=64 219383.7 309876.1 2473114 2944041

mdsumner

Bare minimum capability and library alignment for using GDAL R and Python together, using mdim or classic raster:

GDAL R and Python alignment

Single GDAL, single PROJ, single GEOS, R packages and Python packages linking the same libraries, reticulate bridging cleanly between them, modern numpy, modern Python, modern R, current GDAL master available or latest release.