Skip to content

Multiprocessing vs the GIL for Rasterio Window Reads

A windowed rasterio read holds the GIL for roughly 18% of its wall time, which by Amdahl’s law caps thread speedup near 5.5× and flattens in practice at four threads — so past that point use ProcessPoolExecutor(max_workers=max(1, memory_mb // 1769)) and divide GDAL_CACHEMAX by the worker count. On a 10,240 MB Lambda that is 5 workers at GDAL_CACHEMAX=384 each, not 6 workers at 2048 and not os.cpu_count() workers at whatever GDAL defaults to. Getting either number wrong produces the same symptom — an OOM kill with no Python traceback — for two completely different reasons.

Where the GIL Actually Sits

Memory and CPU Allocation for Raster Workloads establishes that vCPU on Lambda scales linearly with memory: one vCPU per 1,769 MB, so 1,769 MB is one, 3,538 MB is two, and the 10,240 MB ceiling gives about 6. That is the supply side. This page is about whether your code can use it.

The received wisdom — “Python has a GIL, so use processes” — is too coarse to size anything with. Rasterio’s DatasetReader.read() wraps the GDALRasterIO call in a with nogil: block, so the expensive middle of the operation genuinely runs in parallel across threads: the HTTP range request, the DEFLATE or LZW decompression, the predictor unpacking, and the memcpy into the output buffer. If that were the whole operation, threads would scale linearly and this page would not exist.

What holds the GIL is the Python shell around each call. Constructing the Window, doing the affine arithmetic in from_bounds, allocating the destination NumPy array, applying the mask, coercing the dtype, and appending the result to a list — all of that is interpreted bytecode. On a 512 × 512 × 1-band uint16 read from a COG over /vsis3/, that shell is about 18% of wall time.

Where a windowed rasterio read holds and releases the GILA 38.5 millisecond window read split into six segments: window arithmetic at 1.4 milliseconds and destination array allocation at 2.1 milliseconds hold the GIL, the HTTP range request at 19.6 milliseconds and DEFLATE decompression at 9.8 milliseconds and the pixel copy at 2.7 milliseconds run without it, and mask plus dtype coercion at 2.9 milliseconds holds it again.One 512 x 512 uint16 window read, in millisecondsHTTP range requestDEFLATE decompress19.6 ms9.8 msWindow arithmetic — 1.4 msArray allocation — 2.1 msPixel copy — 2.7 msMask + dtype coerce — 2.9 msRose segments hold the GIL; green segments run inside rasterio's nogil block around GDALRasterIO. Larger, more compressed windows shiftthe balance toward green and let threads scale further.
6.4 ms of the 38.5 ms — about 18% — is interpreted bytecode holding the GIL. That fraction, not the thread count, is what bounds thread speedup at roughly 5.5x.

Amdahl’s law turns that 18% into a hard ceiling: maximum speedup is 1 / 0.18 ≈ 5.5×, approached only with infinite threads. In practice you reach 3.4× at four threads and 3.9× at eight, and the extra four threads cost you context switching and four more /vsicurl/ cache footprints for a 15% gain. The plateau is not a bug in rasterio and it is not fixed by a bigger machine.

The proportion moves with the work. Reads that are large and compression-heavy push the GIL-free share up and let threads scale further. Reads that are small, numerous, and followed by NumPy operations that themselves hold the GIL — Python-level loops, np.vectorize, anything with a Python callback — push it down, and eight threads can be slower than one.

Threads or Processes

ThreadPoolExecutor compared with ProcessPoolExecutor for rasterio window readsTwo side-by-side panels. The thread pool panel lists a shared GDAL block cache, no pickling, a 3.4 times speedup at four threads flattening at 3.9 times, and suitability for I/O dominated compositing. The process pool panel lists one GIL per worker, scaling bounded by vCPU rather than by Amdahl's law, a separate block cache per worker that must be divided, and pickling of every result across a pipe.Threads up to four, processes beyondThreadPoolExecutorOne shared GDAL block cache — overlapping tiles fetchedonceNo pickling; arrays stay in one address space3.4x at 4 threads, 3.9x at 8, then flatRight for a tile server compositing 60 windowsGDAL_CACHEMAX set once for the whole functionProcessPoolExecutorOne interpreter and one GIL per workerBounded by vCPU, not by the 18% serial fractionOwn block cache per worker — divide GDAL_CACHEMAXRight for reprojection, convolution, NDVI after the readfork inherits the linked libgdal.so; spawn re-pays 2.2 sNever run both without dividing the cache budget twice — nested parallelism over one vCPU allocation is contention, notthroughput.
The deciding question is not thread safety — it is whether there is real CPU work after the read, and whether you need more than about four times single-threaded throughput.

Threads win when the work per task is I/O-dominated and the results are small. A tile server fetching 60 windows to composite one output image is squarely in that territory: the GIL-free share is high, and every thread shares one GDAL block cache, so overlapping tiles are fetched once. This is the same economy that makes the range-request tuning in Tuning HTTP Range Requests for COG Reads on S3 worth doing.

Processes win when there is real CPU work after the read — reprojection, a convolution, an NDVI computation over a whole window stack — and when you need more than about 4× the single-threaded throughput. Each process gets its own interpreter and its own GIL, so scaling is bounded by vCPU rather than by Amdahl. The costs are real though: each process re-imports rasterio and re-links libgdal.so (the 2.2 seconds measured in Cold Start Mapping for Python GDAL), each gets its own GDAL_CACHEMAX allocation, and every result is pickled across a pipe.

The rule that follows: use threads up to four, processes beyond that, and never both without dividing the cache budget twice.

Prerequisites

  • Runtime: Python 3.11 or 3.12 with rasterio 1.3.9+ (GDAL 3.6+) and NumPy 1.26+.
  • Memory tier: at least 3,538 MB on AWS Lambda if you intend to use two processes; a ProcessPoolExecutor on a 1,769 MB function is two processes time-slicing one vCPU and is slower than doing nothing.
  • Start method: fork on Linux, which Lambda and Cloud Run both provide. spawn re-imports the whole module in each child, paying the GDAL link cost per worker; fork inherits the already-linked parent.
  • IAM: s3:GetObject on the source prefix for the execution role. Nothing extra is needed for the workers — forked children inherit the credentials in the parent’s environment.
  • Environment variables set in function configuration:
    • GDAL_DATA=/opt/share/gdal, PROJ_LIB=/opt/share/proj, LD_LIBRARY_PATH=/opt/lib
    • GDAL_DISABLE_READDIR_ON_OPEN=EMPTY_DIR
    • GDAL_CACHEMAX — set this per worker, not per function. On 10,240 MB with 5 workers, 384 is the right value; leaving it at a per-function 2048 gives five processes 10 GB of cache between them and the function dies at the ceiling.
    • GDAL_NUM_THREADS=1 in worker processes. GDAL’s own internal threading multiplies against your process count; two levels of parallelism over one vCPU budget is pure contention.
    • CPL_VSIL_CURL_CACHE_SIZE=33554432 — 32 MB per worker rather than 64, for the same reason.
    • OMP_NUM_THREADS=1 — stops any OpenMP-linked dependency spawning a thread per host core inside each worker.

Implementation

python
"""Windowed reads over a COG, parallelised the way the allocation allows."""
import concurrent.futures as cf
import json
import multiprocessing as mp
import os
import time

import numpy as np
import rasterio
from rasterio.windows import Window

# --- vCPU is a function of memory, not of os.cpu_count() -------------------
# Lambda grants one vCPU per 1,769 MB. os.cpu_count() reports the HOST's core
# count — commonly 6 on a 512 MB function — and sizing a pool from it is the
# most common way to build a pool that only contends with itself.
MEMORY_MB = int(os.environ.get("AWS_LAMBDA_FUNCTION_MEMORY_SIZE", "3008"))
VCPU = max(1, MEMORY_MB // 1769)

# Reserve one vCPU's worth of headroom for the parent's own work (pickling
# results, assembling the output array) on anything above 2 vCPU.
WORKERS = VCPU if VCPU <= 2 else VCPU - 1

# Each process gets its OWN GDAL block cache. Divide, do not repeat.
CACHE_PER_WORKER_MB = max(64, int(MEMORY_MB * 0.20) // WORKERS)

WORKER_ENV = dict(
    GDAL_CACHEMAX=str(CACHE_PER_WORKER_MB),
    CPL_VSIL_CURL_CACHE_SIZE="33554432",       # 32 MB per worker, in bytes
    GDAL_DISABLE_READDIR_ON_OPEN="EMPTY_DIR",
    GDAL_NUM_THREADS="1",                      # no threading inside a worker
    VSI_CACHE="TRUE",
)

# fork inherits the parent's already-linked libgdal.so. spawn would re-pay the
# ~2.2 s shared-library resolution cost in every worker.
CTX = mp.get_context("fork")


def _read_window(args) -> tuple[int, float, float]:
    """Runs in a worker process. Opens its own dataset handle — a rasterio
    DatasetReader is not fork-safe and must never be inherited across the
    process boundary."""
    uri, col_off, row_off, size, idx = args
    with rasterio.Env(**WORKER_ENV):
        with rasterio.open(uri) as src:
            band = src.read(
                1, window=Window(col_off, row_off, size, size), masked=True
            )
            # Real CPU work: this is what justifies processes over threads.
            # NumPy releases the GIL here too, but the pickling of a scalar
            # result back to the parent is far cheaper than shipping the array.
            valid = band.compressed()
            if valid.size == 0:
                return idx, float("nan"), 0.0
            return idx, float(valid.mean()), float(valid.std())


def handler(event, context):
    uri = event["uri"].replace("s3://", "/vsis3/", 1)
    size = int(event.get("window", 512))

    with rasterio.Env(**WORKER_ENV):
        with rasterio.open(uri) as src:
            width, height = src.width, src.height

    # Align tasks to the file's own block grid: a window that straddles two
    # internal tiles fetches both, so misalignment doubles the transfer.
    tasks = [
        (uri, c, r, size, i)
        for i, (c, r) in enumerate(
            (c, r)
            for r in range(0, height - size + 1, size)
            for c in range(0, width - size + 1, size)
        )
    ]

    t0 = time.perf_counter()
    results = [None] * len(tasks)
    with cf.ProcessPoolExecutor(max_workers=WORKERS, mp_context=CTX) as pool:
        # chunksize batches tasks per worker so a 4,000-window job does not
        # pay 4,000 pickle round trips through the executor's queue.
        for idx, mean, std in pool.map(_read_window, tasks,
                                       chunksize=max(1, len(tasks) // (WORKERS * 4))):
            results[idx] = (mean, std)
    elapsed = time.perf_counter() - t0

    means = np.array([r[0] for r in results], dtype="float64")
    return {
        "statusCode": 200,
        "body": json.dumps({
            "memory_mb": MEMORY_MB,
            "vcpu": VCPU,
            "workers": WORKERS,
            "gdal_cachemax_mb_per_worker": CACHE_PER_WORKER_MB,
            "windows": len(tasks),
            "elapsed_s": round(elapsed, 3),
            "windows_per_s": round(len(tasks) / elapsed, 1),
            "scene_mean": float(np.nanmean(means)),
        }),
    }

Three details are load-bearing. WORKERS is derived from memory, never from os.cpu_count() — on a Lambda, os.cpu_count() reports the host machine’s cores and is unrelated to what you were granted. CACHE_PER_WORKER_MB divides a fixed 20% cache budget across workers rather than giving each the full amount. And the dataset is opened inside the worker: a DatasetReader holds a GDAL C++ handle with file descriptors and cURL state that does not survive a fork, and inheriting one produces intermittent corrupt reads rather than a clean error.

Verification

Sweep the worker count at fixed memory and watch where the curve flattens:

bash
for W in 1 2 3 4 6 8; do
  aws lambda invoke --function-name window-stats \
    --payload "{\"uri\":\"s3://eo-scenes/S2B_33UUP_20260714.tif\",\"window\":512,\"force_workers\":$W}" \
    --cli-binary-format raw-in-base64-out /dev/stdout \
  | jq -r '.body|fromjson|"\(.workers)\t\(.windows_per_s)\t\(.gdal_cachemax_mb_per_worker)"'
done

Expected output on a 10,240 MB function (5 workers is the derived default; the sweep overrides it to show the shape):

code
1	14.2	2048
2	27.9	1024
3	41.1	682
4	53.6	512
6	61.4	341
8	62.0	256

Near-linear to four workers, then a knee: six workers buy 15% over four, and eight buy nothing at all because the function has about 6 vCPU and the parent still needs one. Run the same sweep with ThreadPoolExecutor and the curve tops out around 3.9× regardless of memory — that is the GIL fraction, and it is the measurement that tells you which executor to keep.

Windows per second by ProcessPoolExecutor worker countHorizontal bars of measured windows per second at increasing process counts on a 10,240 megabyte AWS Lambda: 14.2 at one worker, 27.9 at two, 41.1 at three, 53.6 at four, 61.4 at six, and 62.0 at eight, against a thread pool that reaches only 55.4.Measured throughput against worker count on a 10,240 MB Lambda1 process — GDAL_CACHEMAX 204814.2/s2 processes — 1024 each27.9/s3 processes — 682 each41.1/s4 processes — 512 each53.6/s6 processes — 341 each61.4/s8 processes — 256 each62.0/s8 threads — one shared cache55.4/s070 windows/sOne vCPU per 1,769 MB means 10,240 MB grants about 6. The derived default of 5 workers leaves the parent a vCPU for pickling and outputassembly.
Near-linear to four workers, then a knee — the function has about 6 vCPU and the parent still needs one, so the eighth worker buys nothing.

Cross-check Max Memory Used in the CloudWatch REPORT line across the sweep. If it climbs roughly linearly with worker count, your GDAL_CACHEMAX is not being divided and the next larger input will OOM.

Gotchas

  • os.cpu_count() is the host’s core count, not your allocation. A 512 MB Lambda commonly reports 6. Sizing a pool from it gives six processes sharing under a third of one vCPU, each with its own GDAL cache, and the function dies on memory long before it finishes. Derive workers from AWS_LAMBDA_FUNCTION_MEMORY_SIZE, or from the --cpu value on Cloud Run.

  • A DatasetReader is not fork-safe. Opening the raster in the parent and using it in workers appears to work at low concurrency, then produces windows with pixels from the wrong offsets under load. Open inside the worker, every time.

  • GDAL_NUM_THREADS=ALL_CPUS multiplies against your process count. Five workers each spawning six GDAL threads is thirty threads over roughly 6 vCPU. Set GDAL_NUM_THREADS=1 in workers and let the process pool be the only source of parallelism.

  • fork inside a container that has already started threads is undefined behaviour. If your init code starts a background thread — an OpenTelemetry batch exporter, an X-Ray daemon client, a boto3 connection pool warmer — forking after it can deadlock the child on a lock held by a thread that does not exist in it. Create the pool before starting any exporter, or switch to forkserver.

Frequently Asked Questions

Does rasterio release the GIL during a read?

Yes, for the GDAL call itself — GDALRasterIO runs inside a nogil block, so transfer, decompression and pixel copy are genuinely parallel across threads. The GIL is held for the Python shell around each call: window arithmetic, array allocation, masking and dtype coercion, roughly 18% of wall time on a 512 × 512 uint16 window.

How many threads should I use for windowed rasterio reads?

Four is the practical plateau. An 18% serial fraction caps speedup near 5.5× by Amdahl’s law, and you reach 3.4× at four threads and 3.9× at eight — the last four threads cost four extra /vsicurl/ cache footprints for 15%.

How do I size a ProcessPoolExecutor on AWS Lambda?

From memory, not from os.cpu_count(). Lambda grants one vCPU per 1,769 MB, so compute max(1, memory_mb // 1769) and reserve one for the parent above two vCPU. Divide GDAL_CACHEMAX by the worker count — each process gets its own block cache.

Will Python 3.13’s free-threaded build remove the plateau?

It removes the GIL as the cause, but not the ceiling. The 18% shell is still allocation, masking and dtype work that contends on the allocator and on memory bandwidth. Until rasterio and NumPy publish free-threaded wheels that a serverless runtime ships, the process-pool sizing above is the portable answer.

Back to Memory and CPU Allocation for Raster Workloads