🧭 Concept Performance & scale Batch Processing · step 6 of 26 Updated

Threads, Processes and the GIL: Parallelism in Python GIS Explained

Problem statement

Eight cores, one busy:

from concurrent.futures import ThreadPoolExecutor

with ThreadPoolExecutor(max_workers=8) as pool:
    results = list(pool.map(buffer_and_dissolve, files))
# wall clock: 402 s — the same as the serial loop

Switch to processes and it gets worse:

from concurrent.futures import ProcessPoolExecutor

with ProcessPoolExecutor(max_workers=8) as pool:
    results = list(pool.map(buffer_and_dissolve, gdf_chunks))
# wall clock: 511 s — now with 8× the memory

Both outcomes are predictable once you know two things: what Python's Global Interpreter Lock actually locks, and where your GIS work is really spending its time. Geospatial Python is unusual here — much of the heavy lifting happens in C libraries that release the lock, so the usual "threads are useless in Python" advice is wrong about half the time.

Quick answer

Choose by what the work is waiting for:

  1. I/O-bound (reading files, network, database) → threads; GDAL releases the GIL during I/O
  2. CPU-bound in C (GEOS geometry ops, NumPy, rasterio) → threads often work, because those libraries release the GIL
  3. CPU-bound in Python (per-feature loops, apply, dict building) → processes, or rewrite it vectorised
  4. Whole-file batches → processes, one file per worker, no shared state
  5. Always measure; parallel overhead can exceed the saving
from concurrent.futures import ThreadPoolExecutor, ProcessPoolExecutor
import geopandas as gpd

# I/O-bound: many small reads → threads
def read_info(path):
    import pyogrio
    return pyogrio.read_info(path)

with ThreadPoolExecutor(max_workers=16) as pool:
    infos = list(pool.map(read_info, paths))

# whole-file CPU work → processes (module-level function, paths not objects)
def process_one(path_str: str) -> str:
    gdf = gpd.read_file(path_str)
    out = gdf.to_crs(gdf.estimate_utm_crs()).buffer(25)
    dest = path_str.replace("raw", "out")
    gpd.GeoDataFrame(geometry=out).to_file(dest, driver="GPKG")
    return dest

if __name__ == "__main__":
    with ProcessPoolExecutor(max_workers=4) as pool:
        list(pool.map(process_one, [str(p) for p in paths]))

The if __name__ == "__main__" guard is not optional for processes, and passing paths rather than GeoDataFrames is what keeps the transfer cost near zero.

What the GIL actually blocks

Panels showing Python bytecode holding the GIL versus C extensions releasing it during work.
The lock protects the interpreter, not your computation — and GEOS, GDAL and NumPy let go of it.

Step-by-step solution

Grid mapping workload type to threads, processes, vectorisation or a database.
Four kinds of work, four different answers — and one of them is “do not parallelise”.

What the GIL is

The Global Interpreter Lock is a mutex inside CPython that ensures only one thread executes Python bytecode at a time. It exists because CPython's memory management — reference counting on every object — is not thread-safe, and a lock around the interpreter was far simpler and faster for single-threaded code than fine-grained locking everywhere.

The crucial detail is what it does not cover. A C extension can release the GIL while it works and reacquire it before touching Python objects. During that window other threads run freely.

# rough shape of a GIL-releasing C function
# Py_BEGIN_ALLOW_THREADS
#     GEOSIntersection_r(handle, a, b);      /* no Python objects touched */
# Py_END_ALLOW_THREADS

Shapely 2 releases the GIL around GEOS calls. GDAL releases it during I/O. NumPy releases it for most array operations. That is why threads are genuinely useful in a geospatial workload in a way they are not in, say, a pure-Python parser.

Find out what your work is actually doing 

import time
import geopandas as gpd

def profile_step(fn, *args, label=""):
    t0 = time.perf_counter()
    cpu0 = time.process_time()
    result = fn(*args)
    wall = time.perf_counter() - t0
    cpu = time.process_time() - cpu0
    kind = "CPU-bound" if cpu / wall > 0.8 else "I/O-bound or waiting"
    print(f"{label:28} wall {wall:6.2f}s  cpu {cpu:6.2f}s  → {kind}")
    return result

gdf = profile_step(gpd.read_file, "data/raw/parcels.gpkg", label="read_file")
profile_step(lambda g: g.to_crs(3857), gdf, label="to_crs")
profile_step(lambda g: g.buffer(25), gdf.geometry, label="buffer (GEOS)")
profile_step(lambda g: g.apply(lambda x: x.area), gdf.geometry, label="apply (Python loop)")

If CPU time is close to wall time, the work is computing; if it is much lower, it is waiting. That single ratio tells you which parallelism to reach for.

Threads for I/O 

from concurrent.futures import ThreadPoolExecutor, as_completed
import pyogrio

def inventory(paths, workers=16):
    rows = []
    with ThreadPoolExecutor(max_workers=workers) as pool:
        futures = {pool.submit(pyogrio.read_info, p): p for p in paths}
        for future in as_completed(futures):
            path = futures[future]
            try:
                info = future.result()
                rows.append({"path": str(path), "features": info["features"],
                             "crs": info["crs"]})
            except Exception as exc:
                rows.append({"path": str(path), "error": str(exc)})
    return rows

Reading metadata for 400 files over a network share is almost entirely waiting. Sixteen threads turn four minutes into fifteen seconds, and there is no pickling, no memory multiplication and no __main__ guard needed.

Threads for GEOS-bound geometry work 

from concurrent.futures import ThreadPoolExecutor
import numpy as np
import geopandas as gpd

def buffer_chunk(geoms):
    return geoms.buffer(25)          # GEOS releases the GIL

gdf = gpd.read_file("data/raw/parcels.gpkg").to_crs(27700)
chunks = np.array_split(gdf.geometry, 8)

with ThreadPoolExecutor(max_workers=8) as pool:
    parts = list(pool.map(buffer_chunk, chunks))

buffered = gpd.GeoSeries(np.concatenate([p.values for p in parts]), crs=gdf.crs)

This genuinely uses multiple cores, because each thread spends its time inside GEOS with the lock released. It also shares memory, so there is no per-worker copy of the data — the big advantage over processes.

Processes for Python-bound work 

from concurrent.futures import ProcessPoolExecutor
from pathlib import Path
import geopandas as gpd

def process_file(path_str: str) -> dict:
    path = Path(path_str)
    gdf = gpd.read_file(path)
    metric = gdf.to_crs(gdf.estimate_utm_crs())
    metric["geometry"] = metric.geometry.buffer(25)
    dissolved = metric.dissolve(by="class")
    dest = Path("data/out") / f"{path.stem}.gpkg"
    dissolved.to_file(dest, driver="GPKG")
    return {"file": path.name, "features": len(dissolved)}

if __name__ == "__main__":
    paths = [str(p) for p in sorted(Path("data/raw").glob("*.gpkg"))]
    with ProcessPoolExecutor(max_workers=4) as pool:
        for result in pool.map(process_file, paths):
            print(result)

Whole-file batches are the ideal process workload: each worker is independent, arguments are small strings, and results are small dicts. The costs — start-up, pickling, memory per worker — are amortised over minutes of work per file.

The trap: oversubscription 

import os

# GDAL, NumPy's BLAS and PROJ each start their own thread pools
for var in ("OMP_NUM_THREADS", "OPENBLAS_NUM_THREADS", "MKL_NUM_THREADS",
            "GDAL_NUM_THREADS", "NUMEXPR_NUM_THREADS"):
    os.environ.setdefault(var, "1")

import geopandas as gpd     # import after setting the variables

Four worker processes, each spawning eight BLAS threads, is thirty-two threads competing for eight cores. They spend their time context-switching, and the parallel version is slower than the serial one. Setting the native thread counts to 1 in worker processes is the standard fix — and the most common explanation for "parallel made it worse".

The other trap: transferring data 

import geopandas as gpd
from concurrent.futures import ProcessPoolExecutor

big = gpd.read_file("data/raw/buildings.gpkg")      # 2 GB in memory

# bad: pickles the whole frame into every worker
with ProcessPoolExecutor(4) as pool:
    pool.map(some_function, [big] * 4)

# good: workers read their own slice
def work_on_window(args):
    path, start, count = args
    part = gpd.read_file(path, engine="pyogrio",
                         skip_features=start, max_features=count)
    return len(part)

if __name__ == "__main__":
    windows = [("data/raw/buildings.gpkg", s, 100_000) for s in range(0, 800_000, 100_000)]
    with ProcessPoolExecutor(4) as pool:
        print(sum(pool.map(work_on_window, windows)))

Pickling a GeoDataFrame serialises every geometry to WKB, copies it through a pipe, and rebuilds it on the other side — frequently more expensive than the work itself.

Sometimes the answer is not parallelism 

import time
import geopandas as gpd
import numpy as np

gdf = gpd.read_file("data/raw/parcels.gpkg").to_crs(27700)

t0 = time.perf_counter()
areas_loop = gdf.geometry.apply(lambda g: g.area)          # Python per feature
loop_s = time.perf_counter() - t0

t0 = time.perf_counter()
areas_vec = gdf.geometry.area                              # vectorised in GEOS
vec_s = time.perf_counter() - t0

print(f"apply      : {loop_s:6.2f} s")
print(f"vectorised : {vec_s:6.2f} s  ({loop_s/vec_s:.0f}× faster, no parallelism at all)")

Vectorising usually beats parallelising a slow implementation, and it is simpler. Reach for parallelism after the obvious single-threaded wins, not before.

And sometimes the answer is a database 

import geopandas as gpd
from sqlalchemy import create_engine

engine = create_engine("postgresql+psycopg://[email protected]/gis")

# PostGIS runs this with its own parallel workers, against indexed data
counts = gpd.pd.read_sql("""
    SELECT w.ward_name, count(b.id) AS buildings
    FROM wards w LEFT JOIN buildings b ON ST_Within(b.geom, w.geom)
    GROUP BY w.ward_name
""", engine)

A spatial join across millions of rows is a solved problem in PostGIS, with indexes and a parallel query planner. Moving the data to Python to parallelise it by hand is usually the slower path.

Code examples

Example 1: pick the executor from the workload 

from concurrent.futures import ThreadPoolExecutor, ProcessPoolExecutor
import os

def run_parallel(fn, items, kind="auto", workers=None, sample_ratio=None):
    """kind: 'thread' | 'process' | 'auto' (decide from a measured CPU ratio)."""
    workers = workers or (os.cpu_count() or 4)

    if kind == "auto":
        import time
        t0, c0 = time.perf_counter(), time.process_time()
        fn(items[0])
        wall, cpu = time.perf_counter() - t0, time.process_time() - c0
        ratio = cpu / wall if wall else 1
        kind = "process" if ratio > 0.8 else "thread"
        print(f"cpu/wall = {ratio:.2f} → using {kind}s")

    Executor = ProcessPoolExecutor if kind == "process" else ThreadPoolExecutor
    if kind == "thread":
        workers = min(workers * 2, 32)          # I/O tolerates oversubscription

    with Executor(max_workers=workers) as pool:
        return list(pool.map(fn, items))

The measured cpu/wall ratio on one item is a surprisingly good decision rule, and it costs one extra call.

Example 2: worker initialiser that tames native threads 

from concurrent.futures import ProcessPoolExecutor
import os

def init_worker():
    for var in ("OMP_NUM_THREADS", "OPENBLAS_NUM_THREADS", "MKL_NUM_THREADS",
                "GDAL_NUM_THREADS"):
        os.environ[var] = "1"

if __name__ == "__main__":
    with ProcessPoolExecutor(max_workers=4, initializer=init_worker,
                             max_tasks_per_child=25) as pool:
        list(pool.map(process_file, paths))

max_tasks_per_child recycles workers, which caps any per-worker memory growth — the two settings that make a long parallel batch survivable.

Example 3: a scaling test before committing 

import time
from concurrent.futures import ThreadPoolExecutor

def scaling_test(fn, items, max_workers=8):
    baseline = None
    for n in [1, 2, 4, max_workers]:
        t0 = time.perf_counter()
        with ThreadPoolExecutor(max_workers=n) as pool:
            list(pool.map(fn, items))
        elapsed = time.perf_counter() - t0
        baseline = baseline or elapsed
        print(f"{n:>2} workers: {elapsed:6.2f} s  speed-up {baseline/elapsed:4.1f}×  "
              f"efficiency {baseline/elapsed/n:4.0%}")

scaling_test(buffer_chunk, chunks)

Efficiency below about 50% means the parallelism is not paying for itself — usually contention, oversubscription or transfer cost.

Example 4: Dask, when the dataset does not fit 

import dask_geopandas as dgpd
import geopandas as gpd

ddf = dgpd.read_parquet("data/out/buildings/", npartitions=16)
ddf = ddf.to_crs(27700)
ddf["area_m2"] = ddf.geometry.area

result = ddf.groupby("class")["area_m2"].sum().compute()
print(result)

dask-geopandas partitions the data, builds a task graph and runs it across threads or processes with spatial partitioning built in. It is the right tool when the dataset exceeds memory and the operations are the standard ones; hand-rolled parallelism remains simpler for whole-file batches.

Explanation

Two ideas explain almost every parallel-GIS outcome: what the GIL protects, and what parallelism costs.

Bar chart of the cost components of process versus thread parallelism.
Processes buy independence and pay for it in start-up, memory and data transfer.

The GIL serialises the interpreter, not the machine. Any C extension that is not touching Python objects can release it, and the geospatial stack does so extensively: Shapely 2 around GEOS calls, GDAL during I/O and format decoding, NumPy for array arithmetic. Consequently a workload dominated by buffer, intersection, to_crs or file reading really does scale with threads, which is the opposite of the folk wisdom about Python. A workload dominated by per-feature Python — apply, dict building, string handling — does not, because that code holds the lock the entire time.

Processes sidestep the lock by having one interpreter each, and they pay for it three times: start-up (tens to hundreds of milliseconds each, more with spawn), memory (each worker holds its own copy of everything it touches), and transfer (arguments and results are pickled). For file-sized units of work those costs vanish into the noise. For fine-grained work they dominate, which is exactly why the naïve ProcessPoolExecutor version of a chunked geometry operation can be slower than the serial loop.

Oversubscription is the third factor and the least obvious. GDAL, BLAS and PROJ each maintain their own thread pools sized to the machine, so four processes can quietly become thirty-two threads. Cores then spend their time switching contexts and invalidating caches. Setting OMP_NUM_THREADS=1 and its siblings inside workers is not a micro-optimisation; it is often the difference between a speed-up and a slow-down.

Which leads to the practical ordering. First make the single-threaded version fast: vectorise instead of looping, narrow the read with columns= and where=, use the spatial index. Then, if it is still too slow, measure the CPU-to-wall ratio to see what you are waiting for. Use threads for I/O and GEOS-heavy work, processes for whole-file batches, a database for large joins and aggregations, and Dask when the data does not fit at all. And measure the result — parallel code that is not faster is just more complicated.

Edge cases or notes

  • fork is unsafe with GDAL and GEOS: Native libraries can hold locks owned by threads that do not survive the fork. Use spawn.
  • The __main__ guard is mandatory: With spawn, children re-import your module; without the guard they re-run it.
  • Free-threaded Python (3.13t+): Experimental builds remove the GIL, but the geospatial wheels are not there yet. Do not plan around it.
  • asyncio does not help CPU work: It is a concurrency model for I/O; geometry operations still block the loop.
  • Threads and PyQGIS do not mix: Qt objects belong to their creating thread. Use one process per task instead.
  • Shared memory is possible: multiprocessing.shared_memory avoids copying large NumPy arrays, but geometry arrays need re-wrapping — usually not worth it.
  • Reading is parallel-safe, writing usually is not: Several processes writing one GeoPackage will corrupt it.

FAQ

Does the GIL make threads useless in Python GIS?

No. GEOS, GDAL and NumPy release the GIL while they work, so geometry operations and file I/O do scale with threads. Only per-feature Python code is serialised.

When should I use processes instead of threads?

When the work is dominated by Python-level code, or when each unit is a whole file. Processes cost start-up, memory and pickling, so the unit of work needs to be large enough to absorb that.

Why did parallel processing make my job slower?

Most often oversubscription — each worker starting its own native thread pool — or the cost of pickling GeoDataFrames between processes. Set OMP_NUM_THREADS=1 in workers and pass paths, not data.

How many workers should I use?

For CPU work, about the number of physical cores, and check that memory allows it. For I/O work, more than the core count is fine — 16 to 32 threads is common for reading many small files.

Is asyncio useful for GIS work?

Only for network I/O — fetching many tiles or API responses. It does nothing for geometry computation, which blocks the event loop.

What about Dask?

dask-geopandas is the right tool when the dataset does not fit in memory and you want the standard operations parallelised for you. For whole-file batches, a process pool is simpler.

How do I know whether my code is CPU-bound?

Compare time.process_time() with time.perf_counter() over the same call. A ratio near 1 means CPU-bound; much less means it is waiting on I/O.