GHSL — built-up growth (multi-epoch + aggregate)¶

A GHSL strength is its multi-temporal record. Request GHS-BUILT-S at two epochs (2000 and 2020) and reduce across them with aggregate= to map the maximum built-up surface — a simple growth proxy.

Setup¶

The imports up front: tempfile for a scratch download directory, EarthLens for the unified entry point, and AggregationConfig to describe the across-epoch reduction. The plotting imports (matplotlib, numpy, pyramids' Dataset) come along here too so the later cells stay focused on one task each.

In [ ]:

import tempfile

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

from earthlens import EarthLens
from earthlens.aggregate import AggregationConfig

import tempfile import matplotlib.pyplot as plt import numpy as np from pyramids.dataset import Dataset from earthlens import EarthLens from earthlens.aggregate import AggregationConfig

Build the request¶

GHSL is a raster backend, so download() writes GeoTIFFs to disk. We point it at a temporary directory and describe the request: GHS-BUILT-S over a small AOI around Lisbon, at the 2000 and 2020 epochs.

In [ ]:

out = tempfile.mkdtemp()
growth = EarthLens(
    data_source='ghsl',
    variables=['GHS_BUILT_S'],
    start='2000-01-01',
    end='2020-12-31',
    epochs=[2000, 2020],
    aoi=[-9.20, 38.68, -9.10, 38.78],
    path=out,
)

out = tempfile.mkdtemp() growth = EarthLens( data_source='ghsl', variables=['GHS_BUILT_S'], start='2000-01-01', end='2020-12-31', epochs=[2000, 2020], aoi=[-9.20, 38.68, -9.10, 38.78], path=out, )

Download and aggregate across epochs¶

download() fetches both epochs and chains the aggregation in the same call: AggregationConfig with op='max' over a '100YS' window collapses the two epochs into a single across-epoch maximum raster. The call returns the list of written paths.

In [ ]:

rasters = growth.download(
    progress_bar=False,
    aggregate=AggregationConfig(freq='100YS', op='max'),
)
rasters

rasters = growth.download( progress_bar=False, aggregate=AggregationConfig(freq='100YS', op='max'), ) rasters

Plot the across-epoch maximum¶

Read the single written raster back with pyramids' Dataset, mask the nodata fill (negative values), and render the maximum built-up surface for the Lisbon AOI.

In [ ]:

arr = np.asarray(Dataset.read_file(str(rasters[0])).read_array(), dtype='float32')
arr = arr[0] if arr.ndim == 3 else arr
arr = np.where(arr < 0, np.nan, arr)

plt.figure(figsize=(6, 5))
plt.imshow(arr, cmap='inferno')
plt.colorbar(label='max built m2')
plt.title('GHS-BUILT-S 2000-2020 max — Lisbon')
plt.axis('off')
plt.show()

arr = np.asarray(Dataset.read_file(str(rasters[0])).read_array(), dtype='float32') arr = arr[0] if arr.ndim == 3 else arr arr = np.where(arr < 0, np.nan, arr) plt.figure(figsize=(6, 5)) plt.imshow(arr, cmap='inferno') plt.colorbar(label='max built m2') plt.title('GHS-BUILT-S 2000-2020 max — Lisbon') plt.axis('off') plt.show()