Skip to content

Instantly share code, notes, and snippets.

View mdsumner's full-sized avatar

Michael Sumner mdsumner

564 followers · 1.5k following
  • Integrated Digital East Antarctica, Australian Antarctic Division
  • Hobart, Australia
View GitHub Profile
@mdsumner
mdsumner / reproj_array.ipynb
Last active July 4, 2025 12:32
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
@mdsumner
mdsumner / pick_ht.md
Last active July 4, 2025 06:13

https://bsky.app/profile/mdsumner.bsky.social/post/3lt4lhylxhs2v

info <- vapour::vapour_raster_info(dsn <- "/vsicurl/https://projects.pawsey.org.au/idea-gebco-tif/GEBCO_2024.tif")
## remotes::install_github("hypertidy/grout")
g <- grout::grout(info$dimension, info$extent, blocksize = info$block)
idx <- grout::tile_index(g)
options(parallelly.fork.enable = TRUE, future.rng.onMisuse = "ignore")
library(furrr); plan(multicore)
@mdsumner
mdsumner / erddap_vs_rangedownload.md
Last active July 3, 2025 20:31

ERDDAP doesn't support range downloading

export GDAL_DISABLE_READDIR_ON_OPEN=EMPTY_DIR

gdalinfo /vsicurl/https://coastwatch.pfeg.noaa.gov/erddap/files/jplMURSST41/20250701090000-JPL-L4_GHRSST-SSTfnd-MUR-GLOB-v02.0-fv04.1.nc
ERROR 1: 416: Range downloading not supported by this server: Error {
    code=416;
    message="Requested Range Not Satisfiable: REQUESTED_RANGE_NOT_SATISFIABLE: Don't try to connect to .nc or .hdf files on ERDDAP's /files/ system as if they were local files. It is horribly inefficient and often causes other problems. Instead: a) Use (OPeN)DAP client software to connect to ERDDAP's DAP services for this dataset (which have /griddap/ or /tabledap/ in the URL). That's what DAP is for. b) Or, use the dataset's Data Access Form to request a subset of data. c) Or, if you need the entire file or repeated access over a long period of time, use curl, wget, or your browser to download the entire file, then access the data from your local copy of the file.";
}
@mdsumner
mdsumner / vzarr_thredds.md
Last active July 1, 2025 05:22

This works to load but we can't sel() it sensibly, any ideas?

import virtualizarr
#virtualizarr.__version__
#'1.3.3.dev81+ga5d04d7'

from obstore.store import HTTPStore
from virtualizarr.parsers import HDFParser
@mdsumner
mdsumner / pytest_gdal_trouble.md
Created June 27, 2025 00:12 
pytest autotest/gcore/geoloc.py::test_geoloc_normalize_longitude

ImportError while loading conftest '/perm_storage/home/mdsumner/gdal/autotest/conftest.py'.
autotest/conftest.py:7: in <module>
    from osgeo import gdal, ogr, osr
/usr/lib/python3/dist-packages/osgeo/__init__.py:35: in <module>
 _gdal = swig_import_helper()
@mdsumner
mdsumner / odc_geobox.md
Last active June 26, 2025 04:03

Here's my take on defining a "geobox" at a longlat point, where I only need a longitude, latitude, distance, and n-pixels.

I use a local projection to ease the geographic extent definition, then transform that to (auto) UTM, then get the longlat bbox (for STAC query), and the UTM GeoBox (for odc to render to).

from odc.geo import BoundingBox
from odc.geo.geobox import GeoBox
import math
ltln = (-54.50,158.93)
@mdsumner
mdsumner / pq_mpic.md
Last active June 25, 2025 05:45 
pq <- function(x, ex = NULL) {
  env <- Sys.getenv("AWS_NO_SIGN_REQUEST")
  Sys.setenv("AWS_NO_SIGN_REQUEST" = "YES")
  on.exit(Sys.setenv("AWS_NO_SIGN_REQUEST" = env))
  if (!is.null(ex)) {
    x <- terra::crop(x, ex)
  } else {
@mdsumner
mdsumner / bike_cog_shedding.md
Created June 24, 2025 18:06 
## python takes the same amount of time about 2 seconds
#from osgeo import gdal
#gdal.UseExceptions()

#import time
#t0 = time.time()
#ds = gdal.Open("TCI.tif")
#d = ds.ReadRaster()
#t1 = time.time()
@mdsumner
mdsumner / Sentinel_recents.md
Created June 24, 2025 02:47

Doing a Sentinel query of the last 120 days, the availability of imagery is very location dependent.

Scullins Monument

image

Hurd Point

@mdsumner
mdsumner / GDETM30.md
Last active June 23, 2025 05:35