Bathymetry — GEBCO vs ETOPO1 resolution¶

Pull the same area from both shipped global DEMs — GEBCO 2020 at 15 arc-seconds and ETOPO1 at 1 arc-minute — and compare them side by side. The Big Island of Hawaii gives a wide relief range (volcano summit above sea level down to deep flanks), so the roughly 4x resolution difference is easy to see.

In [ ]:

import tempfile

import matplotlib.pyplot as plt
import numpy as np
from pyramids.dataset import Dataset, GeoReference

from earthlens.core import EarthLens

AOI = [-156.4, 18.6, -154.6, 20.4]  # [west, south, east, north]

import tempfile import matplotlib.pyplot as plt import numpy as np from pyramids.dataset import Dataset, GeoReference from earthlens.core import EarthLens AOI = [-156.4, 18.6, -154.6, 20.4] # [west, south, east, north]

In [ ]:

def fetch(dataset):
    out = tempfile.mkdtemp()
    paths = EarthLens(
        data_source='bathymetry', dataset=dataset, aoi=AOI, path=out
    ).download(progress_bar=False)
    source = Dataset.read_file(paths[0])
    grid = np.asarray(source.read_array(), dtype='float32')
    grid = grid[0] if grid.ndim == 3 else grid
    return Dataset.from_array(
        np.where(grid == 32767, np.nan, grid),
        no_data_value=np.nan,
        geo_ref=GeoReference(geo=source.geotransform, epsg=source.epsg),
    )


gebco = fetch('gebco_2020')
etopo = fetch('etopo1_ice')
print('GEBCO 15-arcsec grid:', gebco.shape)
print('ETOPO1 1-arcmin grid:', etopo.shape)

def fetch(dataset): out = tempfile.mkdtemp() paths = EarthLens( data_source='bathymetry', dataset=dataset, aoi=AOI, path=out ).download(progress_bar=False) source = Dataset.read_file(paths[0]) grid = np.asarray(source.read_array(), dtype='float32') grid = grid[0] if grid.ndim == 3 else grid return Dataset.from_array( np.where(grid == 32767, np.nan, grid), no_data_value=np.nan, geo_ref=GeoReference(geo=source.geotransform, epsg=source.epsg), ) gebco = fetch('gebco_2020') etopo = fetch('etopo1_ice') print('GEBCO 15-arcsec grid:', gebco.shape) print('ETOPO1 1-arcmin grid:', etopo.shape)

Side by side¶

Same colour scale for a fair comparison; the finer GEBCO grid resolves rift zones and submarine canyons that the coarser ETOPO1 grid blurs.

In [ ]:

vmin, vmax = -5000, 4000
fig, axes = plt.subplots(1, 2, figsize=(11, 5))
panels = (
    (axes[0], gebco, f'GEBCO 2020 (15-arcsec) {gebco.shape}'),
    (axes[1], etopo, f'ETOPO1 (1-arcmin) {etopo.shape}'),
)
# Both panels share one scale, so they share one colour bar: two identical bars
# would imply the two maps are on different scales when the whole point is that
# they are not.
for ax, grid, title in panels:
    glyph = grid.plot(
        fig=fig,
        ax=ax,
        cmap='terrain',
        vmin=vmin,
        vmax=vmax,
        title=title,
        colorbar=False,
    )
fig.colorbar(glyph.im, ax=list(axes), label='elevation (m)', shrink=0.7)

vmin, vmax = -5000, 4000 fig, axes = plt.subplots(1, 2, figsize=(11, 5)) panels = ( (axes[0], gebco, f'GEBCO 2020 (15-arcsec) {gebco.shape}'), (axes[1], etopo, f'ETOPO1 (1-arcmin) {etopo.shape}'), ) # Both panels share one scale, so they share one colour bar: two identical bars # would imply the two maps are on different scales when the whole point is that # they are not. for ax, grid, title in panels: glyph = grid.plot( fig=fig, ax=ax, cmap='terrain', vmin=vmin, vmax=vmax, title=title, colorbar=False, ) fig.colorbar(glyph.im, ax=list(axes), label='elevation (m)', shrink=0.7)