Choosing the Unit of Work in a Batch Job: File, Layer, Feature or Tile

Problem statement

A batch job buffers every country in the world, one file per country — Natural Earth's 258 of them (the setup is under The data used on this page):

import time

import geopandas as gpd

timings = {}
for path in sorted(Path("data/countries").glob("*.gpkg")):
    g = gpd.read_file(path)
    t0 = time.perf_counter()
    g.geometry.buffer(5_000)
    timings[path.stem] = time.perf_counter() - t0
cost = pd.Series(timings).sort_values()
print(cost.head(3).round(4).to_string())
print("...")
print(cost.tail(3).round(3).to_string())
print(f"total {cost.sum():.2f} s — the largest file is {cost.max() / cost.sum():.0%} of it")

BRT    0.0001
BRI    0.0001
CSI    0.0001
...
USA    0.119
RUS    0.133
CAN    0.301
total 1.56 s — the largest file is 19% of it

Nothing is broken, and at these sizes nothing is slow — but the shape is the one that ruins big jobs. The unit of work is "one file", and one file holds a fifth of the work — thousands of times the median. Parallelism cannot help much — however many workers there are, one of them is busy with Canada. Resumability cannot help — the job either finishes that item or has nothing to resume from. Scale the same shape up to hours instead of seconds, and progress reporting says "257 of 258" for most of the night.

The unit of work — what one iteration of the batch handles — is a design decision that gets made by default, usually as "whatever the files happen to be". It determines how well the job parallelises, how finely it resumes, how useful the progress bar is, and whether one bad input costs you one item or the whole run.

Quick answer

Four candidates, and the right one depends on how evenly the work divides:

# the diagnostic that makes the choice for you
from pathlib import Path

import pandas as pd

for folder in ("data/countries", "data/ga_counties"):
    sizes = pd.Series({p.name: p.stat().st_size for p in Path(folder).glob("*.gpkg")})
    print(folder)
    print(f"  files      {len(sizes)}")
    print(f"  total      {sizes.sum() / 1e6:.1f} MB")
    print(f"  median     {sizes.median() / 1e3:.0f} kB")
    print(f"  largest    {sizes.max() / 1e3:.0f} kB  ({sizes.idxmax()})")
    print(f"  skew       {sizes.max() / sizes.median():.0f}× median")
    print(f"  top 1% is  {sizes.nlargest(max(1, len(sizes) // 100)).sum() / sizes.sum():.0%} of the bytes")

data/countries
  files      258
  total      33.8 MB
  median     106 kB
  largest    1192 kB  (CAN.gpkg)
  skew       11× median
  top 1% is  6% of the bytes
data/ga_counties
  files      159
  total      347.8 MB
  median     1868 kB
  largest    9978 kB  (121.gpkg)
  skew       5× median
  top 1% is  3% of the bytes

A skew above about 20× means the file is probably the wrong unit — and by that measure both folders pass. For Georgia's counties that verdict is right: the timings agree that no county dominates. For the countries it is wrong. Every country file carries the same GeoPackage overhead, so the tiny ones are not tiny on disk and the byte skew is only 11×; the work skew, from the timings above, is in the thousands. File size is a proxy, and here a poor one. When an item's cost is not proportional to its bytes, measure the cost.

The four units

File, layer, tile or feature chunk

The four properties a unit determines

Whatever you pick, it fixes these four at once:

# 1. Parallelism ceiling — wall clock cannot go below the longest single item
longest_item_seconds     # Canada's buffer, in the opening example
# no number of workers makes the job faster than that

# 2. Resume granularity — a crash loses at most one unit of work
#    file unit:  the whole of the largest item
#    tile unit:  one cell of it

# 3. Progress fidelity — "257/258" is only meaningful if items are comparable
#    with hundreds-fold skew, percent-complete is a lie

# 4. Failure blast radius — one corrupt input costs exactly one unit
#    file unit:  lose a whole region
#    tile unit:  lose one cell, and know which

That is why this is worth deciding rather than inheriting. Every operational property people later try to bolt on — parallelism, resumability, progress, quarantine — is capped by a choice made in the first loop that was written.

File: the default, right more often than not

for path in sorted(src.rglob("*.gpkg")):
    process(path)

Right when the files are roughly comparable — county extracts of similar density, daily deliveries, tiles somebody else already made. It is the simplest to implement, the easiest to resume (an output file either exists or does not), and the most natural to report on.

Check before committing to it:

def file_unit_is_ok(paths, max_skew=20):
    sizes = pd.Series({p: p.stat().st_size for p in paths})
    skew = sizes.max() / sizes.median()
    return bool(skew < max_skew), round(float(skew), 1)


for folder in ("data/ga_counties", "data/countries"):
    print(folder, file_unit_is_ok(sorted(Path(folder).glob("*.gpkg"))))

data/ga_counties (True, 5.3)
data/countries (True, 11.2)

File size is a proxy for work, not a measurement of it. A raster of uniform sea is large and trivial; a dense urban vector file is small and expensive. Where the operation cost is not proportional to bytes, measure feature counts instead:

import fiona
counts = {p: len(fiona.open(p)) for p in paths}      # reads the header, not the data

Layer: when a container holds many independent things

import fiona

def discover_layers(paths):
    for path in paths:
        for layer in fiona.listlayers(path):
            yield (path, layer)

for path, layer in discover_layers(paths):
    gdf = gpd.read_file(path, layer=layer)
    ...

A GeoPackage holding 40 layers processed as one unit is one item that takes 40× as long as it needs to. Splitting to the layer level costs one line in discovery and immediately buys 40× the parallelism.

The catch: writing to one GeoPackage from several processes at once will corrupt it. SQLite allows one writer. Either write one output container per input, or collect results and write serially at the end.

Tile: when the data is spatially uneven

This is the answer to the opening problem. Instead of one item per file, one item per spatial cell — and the cells are the same size regardless of how much data falls in them.

import geopandas as gpd
from shapely.geometry import box

def make_tiles(bounds, size):
    minx, miny, maxx, maxy = bounds
    cells = []
    y = miny
    while y < maxy:
        x = minx
        while x < maxx:
            cells.append(box(x, y, min(x + size, maxx), min(y + size, maxy)))
            x += size
        y += size
    return gpd.GeoSeries(cells)

canada = gpd.read_file("data/countries/CAN.gpkg")
tiles = make_tiles(canada.total_bounds, 1_000_000)      # 1,000 km cells
print(f"{len(tiles)} tiles, {int(tiles.intersects(canada.geometry.iloc[0]).sum())} of them touching Canada")

28 tiles, 23 of them touching Canada

Tiles equalise the area per item, not the work per item — a coastal tile of Canada holds far more vertices than a prairie one. But the largest item shrinks to a fraction of what it was, and Example 2 measures what that does to the wall clock. See how to split a large layer into tiles.

Tiling has one hard requirement: the operation must be local. Reprojecting, clipping, buffering and attribute work are all local — a feature's result depends only on that feature. Dissolving, aggregating to a region, or finding nearest neighbours are not: a feature near a tile edge needs data from the neighbouring tile.

# local: safe to tile
gdf.to_crs(27700)
gdf.geometry.buffer(50)
gdf.assign(area=gdf.geometry.area)

# not local: tiling changes the answer
gdf.dissolve(by="ward")                     # a ward spanning two tiles dissolves twice
gpd.sjoin_nearest(points, hospitals)        # the nearest hospital may be in the next tile

The standard workaround is to read a buffered window and write only the core:

def process_tile(tile, source, buffer_m=500):
    window = tile.buffer(buffer_m)
    gdf = gpd.read_file(source, bbox=window.bounds)   # neighbours included
    out = expensive_operation(gdf)
    return out[out.geometry.representative_point().within(tile)]   # core only

Feature chunk: when rows are uniform and independent

CHUNK = 50_000
for start in range(0, len(gdf), CHUNK):
    process(gdf.iloc[start:start + CHUNK])

The most even division available, and the least spatially aware. Right for attribute work, geocoding, per-row API calls and writing to a database. Wrong for anything where neighbouring features interact, because an arbitrary chunk boundary cuts through the middle of a neighbourhood.

For data too large to load at all, chunk at read time rather than after:

import pyogrio
for chunk in pyogrio.read_dataframe(path, return_fids=True, batch_size=50_000):
    process(chunk)

Mixing units is normal

Real jobs often use two:

def discover(src):
    """File-level items normally; tile-level for the few that are too big."""
    for path in sorted(src.rglob("*.gpkg")):
        size = path.stat().st_size
        if size < 500_000_000:
            yield {"kind": "file", "path": path}
        else:
            for tile in tiles_for(path, size=10_000):
                yield {"kind": "tile", "path": path, "bbox": tile.bounds}

Because discovery produces a work list rather than a loop, this costs nothing structurally — apply_one dispatches on kind, and everything downstream is unchanged. That is the payoff of the discover/apply/report shape.

Georgia counties against world countries: measuring skew

The data used on this page

Two real folders with opposite shapes. data/ga_counties/ is Georgia's 232,717 2020 census blocks from the US Census Bureau, one GeoPackage per county — 159 files of broadly similar size. data/countries/ is Natural Earth's 1:10m countries, one file per country — 258 files, one of which is Canada. Both are reprojected to equal-area CRS so buffers and areas are in metres.

from pathlib import Path

import geopandas as gpd
import pandas as pd

blocks = gpd.read_file("https://www2.census.gov/geo/tiger/TIGER2020/TABBLOCK20/tl_2020_13_tabblock20.zip",
                       columns=["COUNTYFP20", "TRACTCE20", "GEOID20", "POP20"]).to_crs(5070)
Path("data/ga_counties").mkdir(parents=True, exist_ok=True)
for county, rows in blocks.groupby("COUNTYFP20"):
    rows.to_file(f"data/ga_counties/{county}.gpkg")

world = gpd.read_file("https://naciscdn.org/naturalearth/10m/cultural/ne_10m_admin_0_countries.zip",
                      columns=["ADM0_A3", "NAME"]).to_crs(8857)                # Equal Earth
Path("data/countries").mkdir(parents=True, exist_ok=True)
for row in world.itertuples():
    world[world["ADM0_A3"] == row.ADM0_A3].to_file(f"data/countries/{row.ADM0_A3}.gpkg")
print(f"{len(list(Path('data/ga_counties').glob('*.gpkg')))} county files, "
      f"{len(list(Path('data/countries').glob('*.gpkg')))} country files")

159 county files, 258 country files

Example 1: measuring the skew properly

import fiona, pandas as pd
from pathlib import Path

def work_profile(paths, by="features"):
    """Estimate cost per item, so the unit can be chosen on evidence."""
    rows = []
    for p in paths:
        try:
            with fiona.open(p) as src:
                rows.append({"path": str(p), "bytes": p.stat().st_size,
                             "features": len(src), "layers": 1})
        except Exception as exc:
            rows.append({"path": str(p), "bytes": p.stat().st_size,
                         "features": None, "error": str(exc)})
    df = pd.DataFrame(rows)
    cost = df[by].dropna()
    return {
        "items": len(df),
        "median": float(cost.median()),
        "p95": float(cost.quantile(0.95)),
        "max": float(cost.max()),
        "skew_max_over_median": round(float(cost.max() / max(cost.median(), 1)), 1),
        "top_1pct_share": round(cost.nlargest(max(1, len(cost)//100)).sum() / cost.sum(), 3),
        "recommendation": (
            "file is fine" if cost.max() / max(cost.median(), 1) < 20
            else "split the largest items — tile or layer"
        ),
    }


for folder in ("data/ga_counties", "data/countries"):
    profile = work_profile(sorted(Path(folder).glob("*.gpkg")), by="features")
    print(folder, {k: profile[k] for k in ("items", "median", "max", "skew_max_over_median", "recommendation")})

data/ga_counties {'items': 159, 'median': 1092.0, 'max': 11748.0, 'skew_max_over_median': 10.8, 'recommendation': 'file is fine'}
data/countries {'items': 258, 'median': 1.0, 'max': 1.0, 'skew_max_over_median': 1.0, 'recommendation': 'file is fine'}

len(src) on a Fiona handle reads the feature count from the header without reading geometry, so this profiles hundreds of files in seconds. It also shows the limit of counting features: every country file holds exactly one, so by features the world looks perfectly even and the profile recommends the file unit — while the timing at the top of the page says one of those features is a quarter of the work. Profile by what drives the cost: bytes or vertices for geometry operations, rows for attribute work.

Example 2: what the skew costs in wall clock

def wall_clock(costs, workers):
    """Longest-processing-time-first scheduling — a good estimate of a real pool."""
    lanes = [0.0] * workers
    for c in sorted(costs, reverse=True):
        i = lanes.index(min(lanes))
        lanes[i] += c
    return max(lanes)

import shapely

costs = list(cost)                                     # the opening example's measured seconds

# the tile unit: split the three most expensive countries into 1,000 km cells and time each piece
tiled = dict(timings)
for code in cost.index[-3:]:
    geom = gpd.read_file(f"data/countries/{code}.gpkg").geometry.iloc[0]
    del tiled[code]
    for k, cell in enumerate(make_tiles(geom.bounds, 1_000_000)):
        piece = shapely.intersection(geom, cell)
        if not piece.is_empty:
            t0 = time.perf_counter()
            piece.buffer(5_000)
            tiled[f"{code}_{k}"] = time.perf_counter() - t0

print(f"{'workers':>8} {'file unit':>10} {'tile unit':>10}")
for w in (1, 4, 8, 16, 32):
    print(f"{w:>8} {wall_clock(costs, w):>9.2f}s {wall_clock(list(tiled.values()), w):>9.2f}s")
print(f"longest item: {max(costs):.2f} s as a file, {max(tiled.values()):.2f} s as a tile")

 workers  file unit  tile unit
       1      1.56s      1.27s
       4      0.39s      0.32s
       8      0.30s      0.16s
      16      0.30s      0.08s
      32      0.30s      0.06s
longest item: 0.30 s as a file, 0.06 s as a tile

With files as the unit, the wall clock stops improving at a handful of workers: after that the floor is the single longest item, and no scheduler can go below it. Tiling the three biggest countries lowers the floor, and more workers help again. These are seconds because the example is small; the shape is identical when each second is an hour. This table is the argument for changing the unit, in a form that survives a conversation with someone who wants to add more workers instead.

Example 3: the unit as an explicit parameter

from enum import Enum

class Unit(str, Enum):
    FILE = "file"
    LAYER = "layer"
    TILE = "tile"

def discover(src, unit: Unit = Unit.FILE, tile_size=10_000):
    paths = sorted(Path(src).rglob("*.gpkg"))
    if unit is Unit.FILE:
        return [{"kind": "file", "path": p} for p in paths]
    if unit is Unit.LAYER:
        return [{"kind": "layer", "path": p, "layer": l}
                for p in paths for l in fiona.listlayers(p)]
    return [{"kind": "tile", "path": p, "bbox": t.bounds}
            for p in paths for t in tiles_for(p, tile_size)]

python -m pipeline run --unit file      # normal nightly
python -m pipeline run --unit tile      # the annual full rebuild

Making it a flag means the choice can be revisited when the data changes shape, without a rewrite.

The longest unit sets the floor

The governing fact is that the longest single unit of work sets a floor on wall-clock time, and no amount of parallelism goes below it. That floor is chosen when you choose the unit, before any code is written, and it is invisible until the day a big input arrives.

This is why "add more workers" so often disappoints. Parallelism divides the total work among workers; it cannot divide one item. A job whose largest item is 62% of the total is 62% serial no matter what, which is Amdahl's law wearing GIS clothes.

The same choice quietly sets three other properties:

Resume granularity. A crash loses the item in flight. With file units and one item that is a fifth of the job, a crash near its end loses a fifth of the night. With tiles, one cell's worth. Resumability is often described as a feature you add; it is mostly a consequence of how finely the work is divided. See how to build a resumable batch GIS job.

Progress honesty. "129 of 258" implies 50% done. When one item is thousands of times the median, it might mean 5% or 95%. Progress reporting is only meaningful when items are comparable, which is a property of the unit rather than of the progress bar. See how to add a progress bar and ETA.

Blast radius. One corrupt input costs exactly one unit. Quarantining a file means losing a region; quarantining a tile means losing a 10 km cell, and knowing precisely which. See failure policy in batch processing.

The constraint that limits the choice is locality. Tiles and chunks are only correct when a feature's result depends on that feature alone. The moment the operation needs neighbours — dissolve, nearest, aggregate-to-region — an arbitrary boundary changes the answer, and the buffered-window trick is the price of admission: read wider than you write, and keep only the core.

Size as a proxy, shared GeoPackages and straddling features

• File size is a proxy, not a measurement. A sparse raster is large and cheap; a dense vector file is small and expensive. Profile feature counts where they diverge.
• Concurrent writes to one GeoPackage corrupt it. SQLite allows a single writer; one output per worker, merged afterwards, is the safe pattern.
• Tiles double-count features that straddle a boundary. Assign each to exactly one tile by its representative point.
• Buffered reads cost time proportional to the buffer. A 500 m buffer on 1 km tiles reads far more than it writes.
• fiona.open() gives feature counts cheaply; gpd.read_file() does not — it reads everything.
• Chunking by row breaks spatial locality, so any operation using a spatial index inside a chunk gets slower, not faster.
• A mixed-unit work list needs a kind field so results can be reassembled correctly.
• Very small units have overhead of their own — process startup, file open, index build. Below about a second per item, the overhead dominates.

Internal links

• The anatomy of a batch job: discover, apply, report — where the unit is decided
• How to split a large layer into tiles for batch processing — the tile unit in practice
• How to speed up batch GIS jobs with parallel processing — why workers alone do not fix skew
• How to build a resumable batch GIS job in Python — resume granularity
• How to add a progress bar and ETA to a long GIS batch job — why comparable items matter
• Failure policy in batch processing — blast radius
• How to process a very large GeoPackage in chunks — the feature-chunk unit
• How to build an inventory of a GIS data folder in Python — the profiling this page starts from

FAQ

How do I know the file is the wrong unit?

Divide the largest item's cost by the median. Above about 20× the file is the wrong unit; above 100× no amount of parallelism will help.

Will more workers fix one enormous file?

No. The longest single item is a floor on wall-clock time. In Example 2, eight, sixteen and thirty-two workers all took the same 0.30 seconds with country files as the unit, because Canada alone took that long; tiling the three largest countries let the same thirty-two workers finish in a fifth of the time.

When is tiling wrong?

When the operation needs neighbouring data — dissolve, aggregate-to-region, nearest neighbour. Read a buffered window and write only the core, or use a different unit.

Can I mix units in one job?

Yes, and it is common: files normally, tiles for the few that are too big. Because discovery produces a list, the rest of the job does not change.

Should the unit match the output files?

Not necessarily. Tiled processing writing one file per region is normal — collect the tile results and write once at the end.

How small should a unit be?

Large enough that per-item overhead — process start, file open, index build — is a small fraction of the work. Below roughly a second per item, the overhead starts to dominate.

Does the unit affect memory?

Directly. Peak memory is roughly the largest single unit, so splitting a huge file into tiles is often the fix for a job that dies with a MemoryError.