Flood Fill, in depth: connectivity, stack-safe fills, scanline spans, tolerance and region labeling

How a scanline fill moves

First principles: reachability under a predicate

Strip away the pixels and a flood fill is a reachability question. The cells are vertices. Two cells are joined by an edge if they are neighbours and both satisfy a predicate, usually "has the same colour as the seed". The fill visits the connected component that contains the seed and applies an action to every member. Depth-first, breadth-first and span-based fills differ in order, memory and speed, never in which cells end up filled. See BFS and DFS, in depth for the general traversal theory.

Two choices define the result. The first is the predicate: exact colour match, colour within a tolerance, or "not a boundary colour" (a boundary fill, which stops at a specific outline colour rather than at any change). The second is connectivity: 4-connected neighbours share an edge (up, down, left, right); 8-connected neighbours also include the four diagonals.

Connectivity matters more than it looks. Draw a closed outline with a one-pixel diagonal line, so that consecutive outline pixels touch only at corners. A 4-connected fill inside it stays inside, because it cannot step diagonally through the corner gap. An 8-connected fill leaks straight out. The rule from digital topology is to use opposite connectivities for a region and its boundary: if outlines are drawn 8-connected, as most line-drawing algorithms do, fill with 4-connectivity.

The recursive version and why it crashes

The version in most tutorials is four lines of recursion:

def fill(img, r, c, old, new):
    if r < 0 or r >= len(img) or c < 0 or c >= len(img[0]):
        return
    if img[r][c] != old:
        return
    img[r][c] = new
    fill(img, r + 1, c, old, new); fill(img, r - 1, c, old, new)
    fill(img, r, c + 1, old, new); fill(img, r, c - 1, old, new)

It is correct and it is a time bomb. Recursion depth equals the length of the current path, and a path through a region can visit every cell. On a 1000 by 1000 canvas that can approach a million frames. CPython's default recursion limit is 1000, so the fill raises RecursionError on modest regions, and raising the limit trades the exception for a segfault. The JVM throws StackOverflowError at a depth set by the thread's stack size. C and C++ simply crash.

There is also a quiet bug: if new == old, the colour check never stops the walk, and the fill revisits cells forever (or until the stack dies). Every production fill begins with if new == old: return, or uses a separate visited bitmap.

An explicit stack: the safe baseline

The fix is to keep the frontier on the heap. Push the seed, pop a cell, fill it and push its matching neighbours. Mark a cell when you push it, not when you pop it, or the same cell can enter the stack up to four times.

def flood_fill(img, sr, sc, new):
    """4-connected fill, in place. img is a list of lists (or a 2-D array)."""
    h, w = len(img), len(img[0])
    old = img[sr][sc]
    if old == new:
        return 0
    img[sr][sc] = new                 # mark on push
    stack, filled = [(sr, sc)], 1
    while stack:
        r, c = stack.pop()
        for nr, nc in ((r + 1, c), (r - 1, c), (r, c + 1), (r, c - 1)):
            if 0 <= nr < h and 0 <= nc < w and img[nr][nc] == old:
                img[nr][nc] = new
                stack.append((nr, nc))
                filled += 1
    return filled

A collections.deque with popleft() gives a breadth-first fill with the identical result. Writing the new colour doubles as the visited mark, so no extra memory is needed. When the predicate is not plain equality (tolerance fills below), you need a separate visited array, because a filled cell might still satisfy the predicate. For more on converting recursion to an explicit stack, see Iterative DFS, in depth.

Scanline fill: one push per run

The per-cell fill touches each pixel and checks its four neighbours, so it does about four reads per pixel and pushes once per pixel. Hopping between rows also hurts cache locality. The scanline fill uses the fact that regions are made of horizontal runs. Pop a seed, extend left and right to fill the whole run in that row, then scan the rows above and below along the run and push one seed for each separate stretch of matching cells. Only one push happens per run, and the inner loops walk contiguous memory.

def scanline_fill(img, sr, sc, new):
    h, w = len(img), len(img[0])
    old = img[sr][sc]
    if old == new:
        return
    stack = [(sr, sc)]
    while stack:
        r, c = stack.pop()
        if img[r][c] != old:          # filled since it was pushed
            continue
        left = c
        while left > 0 and img[r][left - 1] == old:
            left -= 1
        right = c
        while right < w - 1 and img[r][right + 1] == old:
            right += 1
        for x in range(left, right + 1):
            img[r][x] = new
        for nr in (r - 1, r + 1):     # seed each new run above and below
            if 0 <= nr < h:
                in_run = False
                for x in range(left, right + 1):
                    if img[nr][x] == old:
                        if not in_run:
                            stack.append((nr, x))
                            in_run = True
                    else:
                        in_run = False

The continue at the top matters: two parents can push seeds inside the same run, and the second one must be skipped rather than re-filled. For 8-connectivity, widen the scan of the neighbouring rows by one cell at each end, from left - 1 to right + 1, so diagonal contacts are seeded. Heckbert's version in Graphics Gems (1990) also skips re-scanning the parent row.

Worked example

Follow the code above on the grid in the diagram, seeded at row 2, column 1, with the stack popping its last entry. Pop (2,1): the run extends to columns 1 to 4, which become span 1. Scanning row 1 over those columns finds one run (columns 1 to 3; column 4 is a wall) and pushes (1,1). Scanning row 3 finds columns 3 and 4 open and pushes (3,3). Pop (3,3): span 2 covers columns 3 to 6. In row 2, columns 3 and 4 are already filled and column 5 is a wall, so only column 6 is pushed; in row 4, column 3 is pushed.

Pop (4,3): it extends left to column 1 (span 3) and finds nothing new. Pop (2,6): a one-cell span 4 that pushes (1,6). Pop (1,6): span 5 covers columns 6 to 8 and pushes (2,8). The right-hand corridor then unrolls downward as spans 6, 7 and 8: (2,8), (3,8), and row 4 columns 7 to 8. Finally the stack returns to the very first push, (1,1), which becomes span 9. Every one of the 22 open cells is filled with nine pushes, nine pops and no wasted pops. The per-cell fill would push 22 times and make 88 neighbour checks, and many of those checks jump between rows.

On a 3840 by 2160 frame (8.3 million pixels) the gap is larger: millions of pushes for the per-cell fill, a few thousand for the scanline fill on a convex region.

Tolerance fills that do not drift

Photographs have no flat regions, so paint tools use a tolerance: a cell matches if its colour is within distance t of a reference. The critical choice is the reference. Comparing each cell to its neighbour lets the fill creep along a smooth gradient, a tiny step at a time, until it has eaten a whole sky and the horizon. Compare every cell to the seed colour instead, and the region is bounded by a fixed ball in colour space.

def within(a, b, t):
    # a, b are (r, g, b) tuples in 0..255. Max-channel distance is cheap and predictable;
    # Euclidean RGB is common; distance in a perceptual space (CIELAB) matches the eye best.
    return max(abs(a[0] - b[0]), abs(a[1] - b[1]), abs(a[2] - b[2])) <= t

def tolerance_fill(img, sr, sc, new, t):
    h, w = len(img), len(img[0])
    seed = img[sr][sc]
    seen = bytearray(h * w)           # needed: a filled cell may still be "within" the seed
    seen[sr * w + sc] = 1
    stack = [(sr, sc)]
    while stack:
        r, c = stack.pop()
        img[r][c] = new
        for nr, nc in ((r + 1, c), (r - 1, c), (r, c + 1), (r, c - 1)):
            if 0 <= nr < h and 0 <= nc < w and not seen[nr * w + nc] \
                    and within(img[nr][nc], seed, t):
                seen[nr * w + nc] = 1
                stack.append((nr, nc))

Anti-aliased outlines leave a ring of blended pixels: a strict fill leaves a halo, a generous one eats the line. Editors build a soft mask whose coverage falls with colour distance and composite through it. Decide whether alpha takes part, since fully transparent pixels often carry arbitrary RGB values.

Labeling every region at once

Sometimes you need every region at once: count the cells in a microscope slide, label land masses on a map, or give each object in a binary mask an id. A flood fill from every unlabeled foreground pixel does this in linear time, since each pixel is labeled once.

The classical alternative is two-pass connected-component labeling. Pass one scans rows in order and gives each foreground pixel the label of its already-visited neighbours (left and up for 4-connectivity). When two different labels meet, it records that they are equivalent in a disjoint-set structure. Pass two replaces each label with its set representative. Access is sequential, there is no stack, and it extends to images streamed row by row. See Union-Find, in depth for the disjoint-set operations.

def label(mask):                      # mask: list of lists of 0/1
    h, w = len(mask), len(mask[0])
    lab = [[0] * w for _ in range(h)]
    parent = [0]
    def find(x):
        while parent[x] != x:
            parent[x] = parent[parent[x]]
            x = parent[x]
        return x
    for r in range(h):
        for c in range(w):
            if not mask[r][c]:
                continue
            up = lab[r - 1][c] if r else 0
            left = lab[r][c - 1] if c else 0
            if up and left:
                a, b = find(up), find(left)
                lab[r][c] = min(a, b)
                parent[max(a, b)] = min(a, b)
            elif up or left:
                lab[r][c] = up or left
            else:
                parent.append(len(parent))
                lab[r][c] = len(parent) - 1
    for r in range(h):
        for c in range(w):
            if lab[r][c]:
                lab[r][c] = find(lab[r][c])
    return lab

Labels after pass two are representatives, not consecutive numbers; relabel with a dictionary if you need ids from 1 to k. Libraries already do this well: scipy.ndimage.label and OpenCV's connectedComponents and connectedComponentsWithStats both accept a connectivity setting and run in optimized C or C++. GPU implementations use parallel union-find or label propagation instead.

Failure modes

Trade-offs

What to do next

1. Search your codebase for recursive fills and replace them with the explicit-stack version; add a test with a 2000 by 2000 open grid and a spiral so a regression shows up as a test failure, not a crash.
2. Add the new == old guard and a test that fills a region with its own colour.
3. Write down the connectivity your outlines use and pick the opposite for fills; add a test with a diagonal-only boundary.
4. If you fill photographs, compare against the seed colour, expose the tolerance, and decide how alpha takes part.
5. For single fills on large images, implement the scanline fill and measure pushes and time against the per-cell version on your real data.
6. For whole-image labeling, use scipy.ndimage.label or OpenCV before writing your own; write your own two-pass labeler only if you need streaming or custom equivalence rules.
7. Continue with BFS on grids for distances and multi-source floods.