diff --git a/.github/workflows/javascript-test.yml b/.github/workflows/javascript-test.yml
index a6f3dd0..e51fba5 100644
--- a/.github/workflows/javascript-test.yml
+++ b/.github/workflows/javascript-test.yml
@@ -50,7 +50,7 @@ jobs:
 
       - name: Build WASM
         working-directory: .
-        run: wasm-pack build packages/rust/wasm --target bundler --out-dir ../../javascript/src/wasm-core
+        run: wasm-pack build packages/rust/wasm --target web --out-dir ../../javascript/src/wasm-core
 
       - uses: oven-sh/setup-bun@v2
 
diff --git a/.github/workflows/release.yml b/.github/workflows/release.yml
index 304439e..a3dbcb5 100644
--- a/.github/workflows/release.yml
+++ b/.github/workflows/release.yml
@@ -149,7 +149,7 @@ jobs:
 
       - name: Build WASM
         working-directory: ${{ github.workspace }}
-        run: wasm-pack build packages/rust/wasm --target bundler --out-dir ../../javascript/src/wasm-core
+        run: wasm-pack build packages/rust/wasm --target web --out-dir ../../javascript/src/wasm-core
 
       - uses: actions/setup-node@v4
         with:
diff --git a/AGENTS.md b/AGENTS.md
index 3e280f7..0ff45c2 100644
--- a/AGENTS.md
+++ b/AGENTS.md
@@ -34,11 +34,11 @@ uv run prek run --all-files --config .pre-commit-config.yaml   # lint (ruff, myp
 
 ### JavaScript (`packages/javascript/`)
 
-The WASM bundle must be built into `src/wasm-core/` before anything else works:
+The WASM bundle must be built into `src/wasm-core/` before anything else works. Use `--target web`: `core.ts` inlines the `.wasm` bytes (tsup binary loader) and initialises with the generated `initSync`, so it never touches wasm-bindgen's internal glue.
 
 ```bash
 # from repo root; requires wasm-pack and the wasm32-unknown-unknown target
-wasm-pack build packages/rust/wasm --target bundler --out-dir ../../javascript/src/wasm-core
+wasm-pack build packages/rust/wasm --target web --out-dir ../../javascript/src/wasm-core
 
 cd packages/javascript
 bun install
diff --git a/packages/javascript/README.md b/packages/javascript/README.md
index 0d4587d..48565c8 100644
--- a/packages/javascript/README.md
+++ b/packages/javascript/README.md
@@ -13,6 +13,7 @@ High-quality dithering algorithms for e-paper/e-ink displays, powered by a Rust/
 - **OKLab Color Matching**: Weighted Cartesian OKLab — preserves hue without the achromatic-attractor bug that plagues LCH-weighted approaches
 - **Pre-dither Adjustments**: Per-image exposure, saturation, shadows, highlights, dynamic-range compression, and gamut compression — all orthogonal knobs
 - **Serpentine Scanning**: Alternates row direction to eliminate directional artifacts
+- **DBS Refinement**: Optional Direct Binary Search pass that optimizes the dithered pixels against a model of the eye
 - **Universal**: Works in browser (Canvas API) and Node.js (≥18)
 - **Zero Dependencies**: WASM binary bundled inline, no image library required
 
@@ -124,8 +125,29 @@ ditherImage(image: ImageBuffer, palette: ColorScheme | ColorPalette, options?: D
 | `highlights` | `number` | `0.0` | Highlight compression strength (S-curve upper half). `0.0` = off, `1.0` = strong |
 | `tone` | `number \| 'auto' \| 'off'` | `0.0` | Dynamic range compression. `0.0`/`'off'` = disabled; `'auto'` = histogram-based; numeric = fixed strength. Ignored for `ColorScheme` |
 | `gamut` | `number \| 'auto' \| 'off'` | `0.0` | Pre-dither gamut compression. `0.0`/`'off'` = disabled; `'auto'` = activate when image exceeds palette gamut; numeric = fixed. Ignored for `ColorScheme` |
+| `dbs` | `boolean \| DbsParams` | `false` | Direct Binary Search refinement after dithering. `true` = default parameters. See below |
 
-Pre-processing pipeline: `exposure → saturation → shadows/highlights → tone → gamut → dither`. Each step is a no-op at its identity value.
+Pipeline: `exposure → saturation → shadows/highlights → tone → gamut → dither → dbs`. Each step is a no-op at its identity value.
+
+#### DBS refinement
+
+Direct Binary Search (DBS) improves the dithered result directly: it repeatedly changes pixels to other inks, or swaps them with neighbours, whenever that lowers the error *as perceived* through a model of the eye — sharp for brightness, blurrier for colour, scaled to viewing distance and panel resolution. It is much slower than plain dithering (seconds on a full panel in a single-threaded browser) and blocks the calling thread, so run it in a Web Worker.
+
+```typescript
+// Default geometry: 40 cm from a ~127 ppi panel (7.3" 800×480)
+ditherImage(imageBuffer, SPECTRA_7_3_6COLOR, { dbs: true });
+
+// Match your setup
+ditherImage(imageBuffer, SPECTRA_7_3_6COLOR, { dbs: { viewingDistanceCm: 150, ppi: 150, maxPasses: 10 } });
+```
+
+| `DbsParams` field | Default | Description |
+|---|---|---|
+| `viewingDistanceCm` | `40` | Distance between viewer and panel. Larger values blur more |
+| `ppi` | `127` | Panel pixel density |
+| `maxPasses` | `10` | Upper bound on refinement passes; stops early once nothing improves |
+
+Pixels pinned to exact display colors (see below) are never changed by DBS.
 
 `DitherMode.NONE` performs direct nearest-color mapping without error diffusion or ordered dithering. It is intended for already-quantized graphics, not continuous-tone photos: because there is no error diffusion, on limited palettes (especially BWR) a continuous-tone image or a large flat mid-tone area can map to an unexpected ink — for example, a solid mid-gray region can render as solid red. Use an error-diffusion mode (e.g. `FLOYD_STEINBERG`, `BURKES`) for photographic input. Built-in measured palettes carry their canonical firmware `scheme`, so pure display colors map to the corresponding firmware palette index even when measured RGB values are used for matching.
 
@@ -210,7 +232,7 @@ Features: drag & drop or paste from clipboard, live re-render on every setting c
 bun install
 
 # When Rust source changes, rebuild the WASM (from repo root):
-wasm-pack build packages/rust/wasm --target bundler --out-dir ../../javascript/src/wasm-core
+wasm-pack build packages/rust/wasm --target web --out-dir ../../javascript/src/wasm-core
 
 bun run test        # vitest
 bun run build       # tsup → dist/
diff --git a/packages/javascript/src/core.ts b/packages/javascript/src/core.ts
index ff3c83b..294bddc 100644
--- a/packages/javascript/src/core.ts
+++ b/packages/javascript/src/core.ts
@@ -2,30 +2,35 @@ import type { ImageBuffer, PaletteImageBuffer, ColorPalette } from './types';
 import { DitherMode } from './enums';
 import { ColorScheme, getPalette } from './palettes';
 import {
+  initSync,
   dither_image as wasmDitherImage,
   composite_rgba as wasmCompositeRgba,
-  __wbg_set_wasm,
-  __wbindgen_init_externref_table,
-  __wbindgen_cast_0000000000000001,
-} from './wasm-core/epaper_dithering_wasm_bg.js';
+} from './wasm-core/epaper_dithering_wasm.js';
 import wasmBytes from './wasm-core/epaper_dithering_wasm_bg.wasm';
 
-// Synchronous WASM initialization — runs once at module load time.
-const wasmModule = new WebAssembly.Module(wasmBytes as unknown as ArrayBuffer);
-const wasmInstance = new WebAssembly.Instance(wasmModule, {
-  './epaper_dithering_wasm_bg.js': {
-    __wbindgen_init_externref_table,
-    __wbindgen_cast_0000000000000001,
-  },
-});
-__wbg_set_wasm(wasmInstance.exports);
-(wasmInstance.exports as Record<string, () => void>).__wbindgen_start?.();
+// Synchronous WASM initialization from the inlined binary, once at module load time, via
+// wasm-bindgen's own `initSync` (built with `wasm-pack --target web`).
+initSync({ module: wasmBytes });
+
+/**
+ * Direct Binary Search refinement settings. DBS improves a dithered image by repeatedly
+ * changing pixels to other inks or swapping them with neighbours whenever that lowers the
+ * error as perceived through a model of the eye at the given viewing geometry.
+ */
+export interface DbsParams {
+  /** Distance between viewer and panel in cm. Larger values blur more. Default: `40`. */
+  viewingDistanceCm?: number;
+  /** Panel pixel density. A 7.3" 800×480 panel is about 127 ppi. Default: `127`. */
+  ppi?: number;
+  /** Upper bound on full-image passes; stops early once nothing improves. Default: `10`. */
+  maxPasses?: number;
+}
 
 /**
  * Options for `ditherImage`. All fields are optional — defaults are sensible.
  *
- * Pre-processing pipeline (each step is a no-op at its identity value):
- * `exposure → saturation → shadows/highlights → tone → gamut → dither`.
+ * Pipeline (each step is a no-op at its identity value):
+ * `exposure → saturation → shadows/highlights → tone → gamut → dither → dbs`.
  */
 export interface DitherOptions {
   /** Dithering algorithm. Default: `DitherMode.BURKES`. */
@@ -44,6 +49,12 @@ export interface DitherOptions {
   tone?: number | 'auto' | 'off';
   /** Gamut compression: `0.0`/`'off'` disables, `'auto'` opts in. Default: `0.0`. */
   gamut?: number | 'auto' | 'off';
+  /**
+   * Direct Binary Search refinement after dithering: `false` disables, `true` uses default
+   * {@link DbsParams}. Slow — seconds on a full panel, and it blocks the calling thread, so
+   * run it in a Web Worker in the browser. Default: `false`.
+   */
+  dbs?: boolean | DbsParams;
 }
 
 /**
@@ -68,7 +79,9 @@ export function ditherImage(
     highlights = 0.0,
     tone = 0.0,
     gamut = 0.0,
+    dbs = false,
   } = options;
+  const dbsArgs = parseDbs(dbs);
 
   const expectedLength = image.width * image.height * 4;
   if (image.data.length !== expectedLength) {
@@ -114,11 +127,31 @@ export function ditherImage(
     mode as number, serpentine,
     exposure, saturation, shadows, highlights,
     toneArg, gamutArg,
+    ...dbsArgs,
   );
 
   return { width: image.width, height: image.height, indices, palette: outputColors };
 }
 
+/**
+ * Map the `dbs` option onto the binding's `(dbs_passes, viewing_distance_cm, ppi)` arguments;
+ * `dbs_passes = undefined` disables DBS. Validated here because wasm-bindgen silently coerces
+ * numbers into `u32` (e.g. `-1` → 4294967295).
+ */
+function parseDbs(dbs: boolean | DbsParams): [number | undefined, number, number] {
+  if (dbs === false) return [undefined, 40, 127];
+  const { viewingDistanceCm = 40, ppi = 127, maxPasses = 10 } = dbs === true ? {} : dbs;
+  for (const [name, value] of [['viewingDistanceCm', viewingDistanceCm], ['ppi', ppi]] as const) {
+    if (!(Number.isFinite(value) && value > 0)) {
+      throw new Error(`dbs.${name} must be finite and greater than 0, got ${value}`);
+    }
+  }
+  if (!(Number.isInteger(maxPasses) && maxPasses >= 0 && maxPasses <= 0xffffffff)) {
+    throw new Error(`dbs.maxPasses must be a non-negative integer, got ${maxPasses}`);
+  }
+  return [maxPasses, viewingDistanceCm, ppi];
+}
+
 /** Map `'auto'` → `undefined` (Rust None = auto), `'off'` → 0.0, number → pass through. */
 function parseCompression(v: number | 'auto' | 'off'): number | undefined {
   if (v === 'auto') return undefined;
diff --git a/packages/javascript/src/index.ts b/packages/javascript/src/index.ts
index 16166c0..33a452d 100644
--- a/packages/javascript/src/index.ts
+++ b/packages/javascript/src/index.ts
@@ -1,5 +1,5 @@
 export { ditherImage } from './core';
-export type { DitherOptions } from './core';
+export type { DbsParams, DitherOptions } from './core';
 export { DitherMode } from './enums';
 export {
   ColorScheme,
diff --git a/packages/javascript/src/wasm.d.ts b/packages/javascript/src/wasm.d.ts
index ad2e62b..a6219fd 100644
--- a/packages/javascript/src/wasm.d.ts
+++ b/packages/javascript/src/wasm.d.ts
@@ -3,24 +3,3 @@ declare module '*.wasm' {
   const bytes: Uint8Array;
   export default bytes;
 }
-
-/**
- * Hand-written types for the wasm-bindgen-generated `_bg.js` shim. wasm-pack only
- * emits `.d.ts` for the main module; the `_bg.js` exports we need for synchronous
- * initialization are declared here.
- */
-declare module '*epaper_dithering_wasm_bg.js' {
-  /** Throws if `rgba.length` is not a multiple of 4. */
-  export function composite_rgba(rgba: Uint8Array): Uint8Array;
-  export function dither_image(
-    pixels: Uint8Array, width: number,
-    scheme_id: number | null | undefined, palette_bytes: Uint8Array, accent_idx: number,
-    mode_id: number, serpentine: boolean,
-    exposure: number, saturation: number, shadows: number, highlights: number,
-    tone?: number | null, gamut?: number | null,
-  ): Uint8Array;
-  export function measured_palettes(): string;
-  export function __wbg_set_wasm(wasm: unknown): void;
-  export function __wbindgen_init_externref_table(): void;
-  export function __wbindgen_cast_0000000000000001(arg0: number, arg1: number): unknown;
-}
diff --git a/packages/javascript/tests/dithering.test.ts b/packages/javascript/tests/dithering.test.ts
index 6ff5df1..7cb00a8 100644
--- a/packages/javascript/tests/dithering.test.ts
+++ b/packages/javascript/tests/dithering.test.ts
@@ -8,9 +8,10 @@ import {
   fromValue,
   SPECTRA_7_3_6COLOR,
 } from '../src';
+import type { DbsParams, DitherOptions } from '../src';
 import { createTestImage, createGradient, createTransparentTestImage } from './fixtures';
-// Imported after '../src' so that core.ts has already run __wbg_set_wasm.
-import { composite_rgba as compositeRgba } from '../src/wasm-core/epaper_dithering_wasm_bg.js';
+// Imported after '../src' so that core.ts has already run initSync.
+import { composite_rgba as compositeRgba } from '../src/wasm-core/epaper_dithering_wasm.js';
 
 describe('Dithering Algorithms', () => {
   it.each(Object.values(DitherMode).filter((v) => typeof v === 'number'))(
@@ -544,3 +545,46 @@ describe('Dizzy dithering', () => {
     expect(Array.from(a.indices)).toEqual(Array.from(b.indices));
   });
 });
+
+describe('DBS refinement', () => {
+  // Grayscale gradient on BWR: enough structure that DBS changes the halftone.
+  const image = createGradient(40, 30);
+  const indices = (options: DitherOptions = {}) => Array.from(ditherImage(image, ColorScheme.BWR, options).indices);
+
+  it('changes the output for a color scheme', () => {
+    expect(indices({ dbs: true })).not.toEqual(indices());
+  });
+
+  it('changes the output for a measured palette', () => {
+    const plain = ditherImage(image, SPECTRA_7_3_6COLOR).indices;
+    const refined = ditherImage(image, SPECTRA_7_3_6COLOR, { dbs: true }).indices;
+    expect(Array.from(refined)).not.toEqual(Array.from(plain));
+  });
+
+  it('true means default parameters', () => {
+    expect(indices({ dbs: true })).toEqual(indices({ dbs: {} }));
+    expect(indices({ dbs: true })).toEqual(indices({ dbs: { viewingDistanceCm: 40, ppi: 127, maxPasses: 10 } }));
+  });
+
+  it('false and omitted match plain dithering', () => {
+    expect(indices({ dbs: false })).toEqual(indices());
+  });
+
+  it('passes viewing geometry through', () => {
+    expect(indices({ dbs: { viewingDistanceCm: 15 } })).not.toEqual(indices({ dbs: { viewingDistanceCm: 150 } }));
+  });
+
+  it('zero passes matches plain dithering', () => {
+    expect(indices({ dbs: { maxPasses: 0 } })).toEqual(indices());
+  });
+
+  it.each([
+    [{ ppi: 0 }, /ppi/],
+    [{ ppi: Number.NaN }, /ppi/],
+    [{ viewingDistanceCm: -1 }, /viewingDistanceCm/],
+    [{ maxPasses: -1 }, /maxPasses/],
+    [{ maxPasses: 1.5 }, /maxPasses/],
+  ] as [DbsParams, RegExp][])('rejects invalid parameters %o', (dbs, message) => {
+    expect(() => ditherImage(image, ColorScheme.BWR, { dbs })).toThrow(message);
+  });
+});
diff --git a/packages/python/README.md b/packages/python/README.md
index 401c09e..f04aaaa 100644
--- a/packages/python/README.md
+++ b/packages/python/README.md
@@ -28,6 +28,7 @@ pip install epaper-dithering
 - **10 Color Schemes**: Support for mono, 3-color, 4-color, 6-color, 7-color, and grayscale e-paper displays
 - **Pre-dither Adjustments**: Per-image exposure, saturation, shadows, highlights, dynamic-range compression, and gamut compression — all orthogonal knobs you can mix freely
 - **Serpentine Scanning**: Reduces directional artifacts in error diffusion (enabled by default)
+- **DBS Refinement**: Optional Direct Binary Search pass that optimizes the dithered pixels against a model of the eye (slower: ~0.3 s per 800×480 panel on a multi-core machine)
 - **RGBA Support**: Automatic compositing on white background for transparent images
 
 ## Quick Start
@@ -60,6 +61,7 @@ dither_image(
     highlights=0.0,            # highlight compression, S-curve upper half
     tone="auto",               # dynamic-range compression: "auto" | 0.0–1.0
     gamut="auto",              # gamut compression: "auto" | 0.0–1.0
+    dbs=None,                  # DBS refinement: None | True | DbsParams(...)
 )
 ```
 
@@ -185,6 +187,22 @@ result = dither_image(img, SPECTRA_7_3_6COLOR, mode=DitherMode.BURKES, gamut="of
 
 Note: `gamut` also has no effect for theoretical `ColorScheme` palettes.
 
+#### DBS Refinement
+
+Direct Binary Search (DBS) runs after dithering and improves the result directly: it repeatedly changes pixels to other inks, or swaps them with neighbours, whenever that lowers the error *as perceived* through a model of the eye — sharp for brightness, blurrier for colour, scaled to your viewing distance and panel resolution. It typically cuts visible tone and colour error noticeably, at the cost of runtime (about 0.3 s per 800×480 panel on a multi-core machine, versus ~30 ms for plain dithering). The dither call releases the GIL, so other threads keep running meanwhile.
+
+```python
+from epaper_dithering import DbsParams
+
+# Default viewing geometry: 40 cm from a ~127 ppi panel (7.3" 800×480)
+result = dither_image(img, SPECTRA_7_3_6COLOR, dbs=True)
+
+# Match your setup: a wall-mounted 13.3" panel viewed from 1.5 m
+result = dither_image(img, SPECTRA_7_3_6COLOR, dbs=DbsParams(viewing_distance_cm=150, ppi=150))
+```
+
+`max_passes` (default 10) caps the number of refinement passes; most of the improvement happens in the first few. Pixels pinned to exact display colors (see below) are never changed by DBS.
+
 `DitherMode.NONE` performs direct nearest-color mapping without error diffusion or ordered dithering. It is intended for already-quantized graphics, not continuous-tone photos: because there is no error diffusion, on limited palettes (especially BWR) a continuous-tone image or a large flat mid-tone area can map to an unexpected ink — for example, a solid mid-gray region can render as solid red. Use an error-diffusion mode (e.g. `FLOYD_STEINBERG`, `BURKES`) for photographic input. For built-in measured palettes, pure canonical display colors such as `(255, 0, 0)` map directly to the corresponding firmware palette index even though matching uses measured display RGB values.
 
 For built-in measured palettes, exact canonical display colors are also protected in ordered and error-diffusion modes when pre-processing is off: an image made entirely of display colors is returned as a direct palette-index map, and exact display-color pixels inside a mixed image keep their canonical index instead of being rematched to the measured RGB palette. Pre-processing runs before that exact-pixel check, so explicit `tone="auto"`/`gamut="auto"` or other adjustments may intentionally alter those pixels first.
@@ -204,8 +222,7 @@ result = dither_image(img, SPECTRA_7_3_6COLOR, shadows=0.5)
 result = dither_image(img, SPECTRA_7_3_6COLOR, highlights=0.7)
 
 # Combine for a "vivid photo" look
-result = dither_image(img, SPECTRA_7_3_6COLOR,
-                      exposure=1.1, saturation=1.3, shadows=0.3, highlights=0.5)
+result = dither_image(img, SPECTRA_7_3_6COLOR, exposure=1.1, saturation=1.3, shadows=0.3, highlights=0.5)
 ```
 
 Pipeline order: `exposure → saturation → shadows/highlights → tone → gamut → dither`.
@@ -266,11 +283,11 @@ from epaper_dithering import dither_image, ColorPalette, DitherMode
 # Your measured RGB values
 my_display = ColorPalette(
     colors={
-        'black': (5, 5, 5),           # Measured from your display
-        'white': (185, 190, 180),     # Much darker than (255,255,255)
-        'red': (120, 15, 5),          # Much darker than (255,0,0)
+        "black": (5, 5, 5),  # Measured from your display
+        "white": (185, 190, 180),  # Much darker than (255,255,255)
+        "red": (120, 15, 5),  # Much darker than (255,0,0)
     },
-    accent='red'
+    accent="red",
 )
 
 # Use it directly
diff --git a/packages/python/docs/CALIBRATION.md b/packages/python/docs/CALIBRATION.md
index 7861ce6..105eddf 100644
--- a/packages/python/docs/CALIBRATION.md
+++ b/packages/python/docs/CALIBRATION.md
@@ -141,14 +141,14 @@ from epaper_dithering import ColorPalette
 
 my_display = ColorPalette(
     colors={
-        'black':  (26, 13, 35),
-        'white':  (185, 202, 205),
-        'yellow': (202, 184, 0),
-        'red':    (121, 9, 0),
-        'blue':   (0, 69, 139),
-        'green':  (40, 82, 57),
+        "black": (26, 13, 35),
+        "white": (185, 202, 205),
+        "yellow": (202, 184, 0),
+        "red": (121, 9, 0),
+        "blue": (0, 69, 139),
+        "green": (40, 82, 57),
     },
-    accent='red'
+    accent="red",
 )
 ```
 
@@ -186,11 +186,11 @@ To contribute your measurements:
    ```python
    MY_DISPLAY_BWR = ColorPalette(
        colors={
-           'black': (5, 5, 5),
-           'white': (185, 190, 180),
-           'red': (120, 15, 5),
+           "black": (5, 5, 5),
+           "white": (185, 190, 180),
+           "red": (120, 15, 5),
        },
-       accent='red'
+       accent="red",
    )
    ```
 
diff --git a/packages/python/src/epaper_dithering/__init__.py b/packages/python/src/epaper_dithering/__init__.py
index 76aa422..e941d0e 100644
--- a/packages/python/src/epaper_dithering/__init__.py
+++ b/packages/python/src/epaper_dithering/__init__.py
@@ -4,7 +4,7 @@
 for limited-color e-paper/e-ink displays.
 """
 
-from .core import dither_image
+from .core import DbsParams, dither_image
 from .enums import DitherMode
 from .palettes import (
     BWRY_3_97,
@@ -23,6 +23,7 @@
 
 __all__ = [
     "dither_image",
+    "DbsParams",
     "DitherMode",
     "ColorPalette",
     "ColorScheme",
diff --git a/packages/python/src/epaper_dithering/_rs.pyi b/packages/python/src/epaper_dithering/_rs.pyi
index 9f9b0b4..b7e0cf6 100644
--- a/packages/python/src/epaper_dithering/_rs.pyi
+++ b/packages/python/src/epaper_dithering/_rs.pyi
@@ -14,6 +14,9 @@ def dither_image(
     highlights: float = ...,
     tone: float | None = ...,
     gamut: float | None = ...,
+    dbs_passes: int | None = ...,
+    dbs_viewing_distance_cm: float = ...,
+    dbs_ppi: float = ...,
 ) -> bytes: ...
 def composite_rgba(rgba: bytes) -> bytes: ...
 def measured_palettes() -> list[tuple[str, list[int], list[str], int, int]]: ...
diff --git a/packages/python/src/epaper_dithering/core.py b/packages/python/src/epaper_dithering/core.py
index 14f6fb5..3006c63 100644
--- a/packages/python/src/epaper_dithering/core.py
+++ b/packages/python/src/epaper_dithering/core.py
@@ -3,6 +3,8 @@
 from __future__ import annotations
 
 import logging
+import math
+from dataclasses import dataclass
 from typing import cast
 
 from PIL import Image
@@ -15,6 +17,49 @@
 _LOGGER = logging.getLogger(__name__)
 
 
+@dataclass(frozen=True)
+class DbsParams:
+    """Direct Binary Search refinement settings.
+
+    DBS improves a dithered image by repeatedly changing pixels to other inks or swapping
+    them with neighbours whenever that lowers the error as perceived through a model of the
+    eye at the given viewing geometry. It is slow (seconds on a full panel), so it is off
+    unless requested.
+
+    Attributes:
+        viewing_distance_cm: Distance between viewer and panel. Larger values blur more.
+        ppi: Panel pixel density. A 7.3" 800x480 panel is about 127 ppi.
+        max_passes: Upper bound on full-image passes; stops early once nothing improves.
+    """
+
+    viewing_distance_cm: float = 40.0
+    ppi: float = 127.0
+    max_passes: int = 10
+
+    def __post_init__(self) -> None:
+        for name in ("viewing_distance_cm", "ppi"):
+            value = getattr(self, name)
+            if not (math.isfinite(value) and value > 0):
+                raise ValueError(f"{name} must be finite and greater than 0, got {value}")
+        if self.max_passes < 0:
+            raise ValueError(f"max_passes must be >= 0, got {self.max_passes}")
+
+
+def _dbs_kwargs(dbs: DbsParams | bool | None) -> dict[str, object]:
+    """Map the `dbs` argument onto the binding's flat `dbs_*` keyword arguments."""
+    if dbs is None or dbs is False:
+        return {}
+    if dbs is True:
+        dbs = DbsParams()
+    if not isinstance(dbs, DbsParams):
+        raise TypeError(f"dbs must be DbsParams, bool, or None, got {type(dbs).__name__}")
+    return {
+        "dbs_passes": dbs.max_passes,
+        "dbs_viewing_distance_cm": dbs.viewing_distance_cm,
+        "dbs_ppi": dbs.ppi,
+    }
+
+
 def _to_rgb_bytes(image: Image.Image) -> tuple[bytes, int, int]:
     """Convert PIL image to flat RGB bytes. Composites RGBA on white."""
     if image.mode == "RGBA":
@@ -48,6 +93,7 @@ def dither_image(  # pylint: disable=too-many-arguments
     highlights: float = 0.0,
     tone: float | str = 0.0,
     gamut: float | str = 0.0,
+    dbs: DbsParams | bool | None = None,
 ) -> Image.Image:
     """Apply dithering to an image for e-paper display.
 
@@ -69,6 +115,9 @@ def dither_image(  # pylint: disable=too-many-arguments
             to display range, 0.0–1.0 = fixed strength. Only meaningful for measured palettes.
         gamut: Gamut compression for out-of-gamut pixels. 0.0 = off, "auto" = full
             strength on out-of-gamut pixels (smoothstep), 0.0–1.0 = fixed strength.
+        dbs: Direct Binary Search refinement after dithering. `None`/`False` = off,
+            `True` = default `DbsParams()`, or a `DbsParams` for custom viewing geometry.
+            Slow: seconds on a full panel.
 
     Returns:
         Dithered palette-mode (`"P"`) PIL Image matching the color scheme.
@@ -92,6 +141,7 @@ def dither_image(  # pylint: disable=too-many-arguments
         "highlights": highlights,
         "tone": _compression(tone),
         "gamut": _compression(gamut),
+        **_dbs_kwargs(dbs),
     }
 
     if isinstance(color_scheme, ColorScheme):
diff --git a/packages/python/src/lib.rs b/packages/python/src/lib.rs
index 98d396a..32b97a9 100644
--- a/packages/python/src/lib.rs
+++ b/packages/python/src/lib.rs
@@ -1,6 +1,7 @@
 use epaper_dithering_core::{
     color_space_lab::rgb_to_oklab,
     composite::composite_rgba_on_white,
+    dbs::DbsParams,
     dither, dither_with_canonical, DitherConfig,
     enums::{DitherMode, GamutCompression, ToneCompression},
     measured_palettes::CATALOG,
@@ -48,10 +49,21 @@ fn parse_gamut(v: Option<f64>) -> GamutCompression {
     }
 }
 
+/// `dbs_passes = None` disables DBS refinement; otherwise build and validate its parameters.
+fn parse_dbs(passes: Option<u32>, viewing_distance_cm: f64, ppi: f64) -> PyResult<Option<DbsParams>> {
+    passes
+        .map(|max_passes| {
+            let params = DbsParams { viewing_distance_cm, ppi, max_passes };
+            params.validate().map(|()| params).map_err(|e| PyValueError::new_err(e.to_string()))
+        })
+        .transpose()
+}
+
 /// Dither a flat RGB image.
 ///
 /// Pass either `scheme_id` (idealized color scheme) or `palette_bytes` + `accent_idx`
-/// (measured palette). If both are given, `palette_bytes` wins.
+/// (measured palette). If both are given, `palette_bytes` wins. `dbs_passes` enables
+/// Direct Binary Search refinement with up to that many passes.
 #[pyfunction]
 #[pyo3(signature = (
     pixels, width, height, *,
@@ -59,6 +71,7 @@ fn parse_gamut(v: Option<f64>) -> GamutCompression {
     mode_id=1, serpentine=true,
     exposure=1.0, saturation=1.0, shadows=0.0, highlights=0.0,
     tone=0.0, gamut=0.0,
+    dbs_passes=None, dbs_viewing_distance_cm=40.0, dbs_ppi=127.0,
 ))]
 #[allow(clippy::too_many_arguments)]
 fn dither_image(
@@ -77,6 +90,9 @@ fn dither_image(
     highlights: f64,
     tone: Option<f64>,
     gamut: Option<f64>,
+    dbs_passes: Option<u32>,
+    dbs_viewing_distance_cm: f64,
+    dbs_ppi: f64,
 ) -> PyResult<Vec<u8>> {
     validate_image(&pixels, width)?;
     // Additionally validate against the caller's height instead of silently
@@ -105,6 +121,7 @@ fn dither_image(
         highlights,
         tone:  parse_tone(tone),
         gamut: parse_gamut(gamut),
+        dbs:   parse_dbs(dbs_passes, dbs_viewing_distance_cm, dbs_ppi)?,
     };
 
     // Validation and `Palette` construction stay under the GIL; only the
diff --git a/packages/python/tests/test_dithering.py b/packages/python/tests/test_dithering.py
index 04b8ddb..4028ca9 100644
--- a/packages/python/tests/test_dithering.py
+++ b/packages/python/tests/test_dithering.py
@@ -557,3 +557,89 @@ def test_is_deterministic(self):
         a = dither_image(img, ColorScheme.BWR, mode=DitherMode.DIZZY)
         b = dither_image(img, ColorScheme.BWR, mode=DitherMode.DIZZY)
         assert list(a.getdata()) == list(b.getdata())
+
+
+def _color_gradient(width: int = 40, height: int = 30) -> Image.Image:
+    """Smooth colour gradient: enough structure that DBS changes the halftone."""
+    img = Image.new("RGB", (width, height))
+    img.putdata(
+        [
+            (x * 255 // (width - 1), y * 255 // (height - 1), 255 - x * 255 // (width - 1))
+            for y in range(height)
+            for x in range(width)
+        ]
+    )
+    return img
+
+
+class TestDbs:
+    """Direct Binary Search refinement (`dbs=`)."""
+
+    def test_changes_output_for_measured_palette(self):
+        from epaper_dithering import SPECTRA_7_3_6COLOR
+
+        img = _color_gradient()
+        plain = dither_image(img, SPECTRA_7_3_6COLOR)
+        refined = dither_image(img, SPECTRA_7_3_6COLOR, dbs=True)
+
+        assert refined.mode == "P"
+        assert refined.size == plain.size
+        assert refined.tobytes() != plain.tobytes()
+
+    def test_changes_output_for_color_scheme(self):
+        img = _color_gradient()
+        plain = dither_image(img, ColorScheme.BWR)
+        refined = dither_image(img, ColorScheme.BWR, dbs=True)
+        assert refined.tobytes() != plain.tobytes()
+
+    def test_true_means_default_params(self):
+        from epaper_dithering import DbsParams
+
+        img = _color_gradient()
+        assert (
+            dither_image(img, ColorScheme.BWR, dbs=True).tobytes()
+            == dither_image(img, ColorScheme.BWR, dbs=DbsParams()).tobytes()
+        )
+
+    @pytest.mark.parametrize("off", [None, False])
+    def test_off_values_match_plain_dither(self, off):
+        img = _color_gradient()
+        assert dither_image(img, ColorScheme.BWR, dbs=off).tobytes() == dither_image(img, ColorScheme.BWR).tobytes()
+
+    def test_params_are_passed_through(self):
+        """Viewing geometry sets the eye-model blur width, so it must change the result."""
+        from epaper_dithering import DbsParams
+
+        img = _color_gradient()
+        near = dither_image(img, ColorScheme.BWR, dbs=DbsParams(viewing_distance_cm=15.0))
+        far = dither_image(img, ColorScheme.BWR, dbs=DbsParams(viewing_distance_cm=150.0))
+        assert near.tobytes() != far.tobytes()
+
+    def test_zero_passes_matches_plain_dither(self):
+        from epaper_dithering import DbsParams
+
+        img = _color_gradient()
+        assert (
+            dither_image(img, ColorScheme.BWR, dbs=DbsParams(max_passes=0)).tobytes()
+            == dither_image(img, ColorScheme.BWR).tobytes()
+        )
+
+    @pytest.mark.parametrize(
+        "kwargs",
+        [{"ppi": 0.0}, {"ppi": float("nan")}, {"viewing_distance_cm": -1.0}, {"max_passes": -1}],
+    )
+    def test_invalid_params_raise_value_error(self, kwargs):
+        from epaper_dithering import DbsParams
+
+        with pytest.raises(ValueError):
+            DbsParams(**kwargs)
+
+    def test_rejects_non_params_values(self):
+        with pytest.raises(TypeError):
+            dither_image(_color_gradient(), ColorScheme.BWR, dbs=10)  # type: ignore[arg-type]
+
+    def test_binding_validates_geometry(self):
+        """The Rust binding reports bad geometry as ValueError instead of panicking."""
+        pixels, width, height = _to_rgb_bytes(_color_gradient(8, 8))
+        with pytest.raises(ValueError, match="ppi"):
+            _rs.dither_image(pixels, width, height, scheme_id=1, dbs_passes=10, dbs_ppi=0.0)
diff --git a/packages/rust/core/examples/dbs_compare.rs b/packages/rust/core/examples/dbs_compare.rs
new file mode 100644
index 0000000..952d5cf
--- /dev/null
+++ b/packages/rust/core/examples/dbs_compare.rs
@@ -0,0 +1,277 @@
+//! Visual comparison of plain error diffusion against the same output refined by DBS.
+//!
+//! For each fixture image (resized to fill the panel), renders every variant in measured ink
+//! colours and writes PNGs plus an `index.html` for side-by-side viewing:
+//!
+//!   <out>/<image>/source.png, <variant>.png  full-size renders
+//!   <out>/<image>/strip.png                  source | variants, side by side
+//!   <out>/<image>/zoom.png                   centre crop at ZOOM× nearest-neighbour, same order
+//!
+//! Metrics: mean OKLab dE over 4×4 block averages (neutral — independent of the DBS eye model)
+//! and DBS perceived error (the quantity DBS minimises, so DBS wins it by construction).
+//!
+//!   cargo run --release --example dbs_compare -- [options] [image files…]
+//!
+//! Options (defaults in brackets):
+//!   --palette spectra6|mono   [spectra6]   --passes N        [10]
+//!   --size WxH                [800x480]    --distance CM     [40]
+//!   --ppi P                   [127]        --out DIR         [target/dbs_compare]
+//!   --variants a,b,…          [all]        --baseline DIR    compare renders against an
+//!                                                             earlier --out (% pixels differing)
+//!
+//! Without image arguments, every image in tests/fixtures/images is used.
+
+use std::path::{Path, PathBuf};
+use std::time::Instant;
+
+use epaper_dithering_core::{
+    color_space::srgb_channel_to_linear,
+    color_space_lab::rgb_to_oklab,
+    dbs::{perceived_error, DbsParams},
+    dither, dither_with_canonical,
+    enums::{DitherMode, GamutCompression, ToneCompression},
+    measured_palettes::SPECTRA_7_3_6COLOR,
+    palettes::{ColorScheme, Palette},
+    types::ImageBuffer,
+    DitherConfig,
+};
+use image::{imageops, RgbImage};
+
+const BLOCK: usize = 4;
+const ZOOM: u32 = 4;
+const CROP: (u32, u32) = (160, 120);
+
+struct Args {
+    measured: bool,
+    width: u32,
+    height: u32,
+    params: DbsParams,
+    out: PathBuf,
+    variants: Vec<String>,
+    baseline: Option<PathBuf>,
+    images: Vec<PathBuf>,
+}
+
+fn parse_args() -> Args {
+    let mut a = Args {
+        measured: true,
+        width: 800,
+        height: 480,
+        params: DbsParams::default(),
+        out: Path::new(env!("CARGO_MANIFEST_DIR")).join("../target/dbs_compare"),
+        variants: Vec::new(),
+        baseline: None,
+        images: Vec::new(),
+    };
+    let mut it = std::env::args().skip(1);
+    while let Some(arg) = it.next() {
+        let mut val = || it.next().unwrap_or_else(|| panic!("{arg} needs a value"));
+        match arg.as_str() {
+            "--palette" => a.measured = match val().as_str() {
+                "spectra6" => true,
+                "mono" => false,
+                p => panic!("unknown palette {p} (spectra6|mono)"),
+            },
+            "--size" => {
+                let v = val();
+                let (w, h) = v.split_once('x').expect("--size WxH");
+                (a.width, a.height) = (w.parse().expect("width"), h.parse().expect("height"));
+            }
+            "--passes" => a.params.max_passes = val().parse().expect("--passes N"),
+            "--distance" => a.params.viewing_distance_cm = val().parse().expect("--distance CM"),
+            "--ppi" => a.params.ppi = val().parse().expect("--ppi P"),
+            "--out" => a.out = val().into(),
+            "--variants" => a.variants = val().split(',').map(str::to_owned).collect(),
+            "--baseline" => a.baseline = Some(val().into()),
+            _ if arg.starts_with("--") => panic!("unknown option {arg}"),
+            _ => a.images.push(arg.into()),
+        }
+    }
+    if a.images.is_empty() {
+        let dir = Path::new(env!("CARGO_MANIFEST_DIR")).join("tests/fixtures/images");
+        a.images = std::fs::read_dir(&dir)
+            .expect("fixtures/images")
+            .filter_map(|e| Some(e.ok()?.path()))
+            .filter(|p| p.is_file() && matches!(p.extension().and_then(|e| e.to_str()), Some("png" | "jpg" | "jpeg")))
+            .collect();
+        a.images.sort();
+    }
+    a
+}
+
+/// Mean OKLab dE between two flat RGB buffers, averaged over BLOCK × BLOCK tiles.
+fn block_delta_e(a: &[u8], b: &[u8], width: usize, height: usize) -> f64 {
+    let mut total = 0.0;
+    let mut blocks = 0usize;
+    for by in (0..height).step_by(BLOCK) {
+        for bx in (0..width).step_by(BLOCK) {
+            let (mut sa, mut sb, mut n) = ([0.0; 3], [0.0; 3], 0.0);
+            for y in by..(by + BLOCK).min(height) {
+                for x in bx..(bx + BLOCK).min(width) {
+                    let i = (y * width + x) * 3;
+                    for c in 0..3 {
+                        sa[c] += srgb_channel_to_linear(a[i + c]);
+                        sb[c] += srgb_channel_to_linear(b[i + c]);
+                    }
+                    n += 1.0;
+                }
+            }
+            let la = rgb_to_oklab(sa[0] / n, sa[1] / n, sa[2] / n);
+            let lb = rgb_to_oklab(sb[0] / n, sb[1] / n, sb[2] / n);
+            total += ((la.l - lb.l).powi(2) + (la.a - lb.a).powi(2) + (la.b - lb.b).powi(2)).sqrt();
+            blocks += 1;
+        }
+    }
+    total / blocks as f64
+}
+
+fn render(indices: &[u8], palette: &Palette, w: u32, h: u32) -> RgbImage {
+    let rgb = indices.iter().flat_map(|&i| palette.colors[i as usize]).collect();
+    RgbImage::from_raw(w, h, rgb).expect("render buffer size")
+}
+
+/// Images side by side with a 4 px white gutter.
+fn hstack(images: &[RgbImage]) -> RgbImage {
+    const GAP: u32 = 4;
+    let h = images.iter().map(RgbImage::height).max().unwrap_or(0);
+    let w = images.iter().map(|i| i.width() + GAP).sum::<u32>() - GAP;
+    let mut out = RgbImage::from_pixel(w, h, image::Rgb([255, 255, 255]));
+    let mut x = 0;
+    for img in images {
+        imageops::replace(&mut out, img, x.into(), 0);
+        x += img.width() + GAP;
+    }
+    out
+}
+
+fn zoom_crop(img: &RgbImage) -> RgbImage {
+    let (cw, ch) = (CROP.0.min(img.width()), CROP.1.min(img.height()));
+    let crop = imageops::crop_imm(img, (img.width() - cw) / 2, (img.height() - ch) / 2, cw, ch).to_image();
+    imageops::resize(&crop, cw * ZOOM, ch * ZOOM, imageops::FilterType::Nearest)
+}
+
+struct Variant {
+    name: &'static str,
+    mode: DitherMode,
+    dbs: bool,
+}
+
+const VARIANTS: [Variant; 4] = [
+    Variant { name: "burkes", mode: DitherMode::Burkes, dbs: false },
+    Variant { name: "burkes_dbs", mode: DitherMode::Burkes, dbs: true },
+    Variant { name: "dizzy", mode: DitherMode::Dizzy, dbs: false },
+    Variant { name: "dizzy_dbs", mode: DitherMode::Dizzy, dbs: true },
+];
+
+fn main() {
+    let args = parse_args();
+    let palette: &Palette = if args.measured { &SPECTRA_7_3_6COLOR } else { ColorScheme::Mono.palette() };
+    let (w, h) = (args.width as usize, args.height as usize);
+    std::fs::create_dir_all(&args.out).expect("create output dir");
+
+    println!(
+        "palette {}, {}x{}, {} cm, {} ppi, ≤{} passes → {}\n",
+        if args.measured { "spectra6 (measured)" } else { "mono" },
+        w, h, args.params.viewing_distance_cm, args.params.ppi, args.params.max_passes,
+        args.out.display()
+    );
+    let variants: Vec<&Variant> = VARIANTS
+        .iter()
+        .filter(|v| args.variants.is_empty() || args.variants.iter().any(|n| n == v.name))
+        .collect();
+    assert!(!variants.is_empty(), "no variant matches {:?}", args.variants);
+    // Per variant: (sum block dE, sum perceived E, total seconds, differing px, compared px).
+    let mut totals = vec![(0.0, 0.0, 0.0, 0usize, 0usize); variants.len()];
+
+    println!(
+        "{:<22} {:<12} {:>10} {:>14} {:>10} {:>9}",
+        "image", "variant", "block dE", "perceived E", "time", "Δpx %"
+    );
+
+    let mut html = String::from(
+        "<!doctype html><meta charset=utf-8><title>DBS comparison</title>\
+         <style>body{font:14px system-ui;margin:16px;background:#888}img{max-width:100%;image-rendering:pixelated}\
+         h2{margin:24px 0 4px}p{margin:0 0 8px}</style>",
+    );
+    html += &format!(
+        "<h1>DBS comparison</h1><p>{}, {w}×{h}, {} cm, {} ppi, ≤{} passes. Order: source, {}.</p>",
+        if args.measured { "Spectra 6 measured" } else { "Mono" },
+        args.params.viewing_distance_cm, args.params.ppi, args.params.max_passes,
+        variants.iter().map(|v| v.name).collect::<Vec<_>>().join(", ")
+    );
+
+    for path in &args.images {
+        let stem = path.file_stem().and_then(|s| s.to_str()).expect("image name");
+        let src = image::open(path)
+            .unwrap_or_else(|e| panic!("failed to load {}: {e}", path.display()))
+            .resize_to_fill(args.width, args.height, imageops::FilterType::Lanczos3)
+            .to_rgb8();
+        let dir = args.out.join(stem);
+        std::fs::create_dir_all(&dir).expect("create image dir");
+        src.save(dir.join("source.png")).expect("save source");
+
+        let buf = ImageBuffer::new(src.as_raw(), w);
+        let mut renders = vec![src.clone()];
+        let mut rows = String::new();
+
+        for (vi, v) in variants.iter().enumerate() {
+            let cfg = DitherConfig {
+                mode: v.mode,
+                dbs: v.dbs.then_some(args.params),
+                ..Default::default()
+            };
+            let t = Instant::now();
+            let indices = if args.measured {
+                let cfg = DitherConfig { tone: ToneCompression::Auto, gamut: GamutCompression::Auto, ..cfg };
+                dither_with_canonical(&buf, palette, ColorScheme::Bwgbry.palette(), cfg)
+            } else {
+                dither(&buf, palette, cfg)
+            };
+            let elapsed = t.elapsed();
+
+            let img = render(&indices, palette, args.width, args.height);
+            let de = block_delta_e(src.as_raw(), img.as_raw(), w, h);
+            let pe = perceived_error(src.as_raw(), w, h, palette, &indices, &args.params);
+            let diff = args.baseline.as_ref().and_then(|b| {
+                let base = image::open(b.join(stem).join(format!("{}.png", v.name))).ok()?.to_rgb8();
+                (base.dimensions() == img.dimensions())
+                    .then(|| base.pixels().zip(img.pixels()).filter(|(a, b)| a != b).count())
+            });
+            let t = &mut totals[vi];
+            (t.0, t.1, t.2) = (t.0 + de, t.1 + pe, t.2 + elapsed.as_secs_f64());
+            if let Some(d) = diff {
+                (t.3, t.4) = (t.3 + d, t.4 + w * h);
+            }
+            let diff_col = diff.map_or("-".into(), |d| format!("{:.3}", 100.0 * d as f64 / (w * h) as f64));
+            println!(
+                "{stem:<22} {:<12} {de:>10.4} {pe:>14.4e} {:>9.2}s {diff_col:>9}",
+                v.name,
+                elapsed.as_secs_f64()
+            );
+            rows += &format!("<tr><td>{}<td>{de:.4}<td>{pe:.4e}<td>{:.2}s", v.name, elapsed.as_secs_f64());
+
+            img.save(dir.join(format!("{}.png", v.name))).expect("save variant");
+            renders.push(img);
+        }
+
+        hstack(&renders).save(dir.join("strip.png")).expect("save strip");
+        hstack(&renders.iter().map(zoom_crop).collect::<Vec<_>>())
+            .save(dir.join("zoom.png"))
+            .expect("save zoom");
+
+        html += &format!(
+            "<h2>{stem}</h2><table><tr><th>variant<th>block dE<th>perceived E<th>time{rows}</table>\
+             <p><img src=\"{stem}/strip.png\"></p><p><img src=\"{stem}/zoom.png\"></p>"
+        );
+    }
+
+    println!("\n{:<22} {:<12} {:>10} {:>14} {:>10} {:>9}", "TOTAL", "variant", "mean dE", "mean E", "time", "Δpx %");
+    let n = args.images.len() as f64;
+    for (v, t) in variants.iter().zip(&totals) {
+        let diff_col = if t.4 > 0 { format!("{:.3}", 100.0 * t.3 as f64 / t.4 as f64) } else { "-".into() };
+        println!("{:<22} {:<12} {:>10.4} {:>14.4e} {:>9.2}s {diff_col:>9}", "", v.name, t.0 / n, t.1 / n, t.2);
+    }
+
+    std::fs::write(args.out.join("index.html"), html).expect("write index.html");
+    println!("\nopen {}", args.out.join("index.html").display());
+}
diff --git a/packages/rust/core/src/composite.rs b/packages/rust/core/src/composite.rs
index 178f29f..3ecfc2d 100644
--- a/packages/rust/core/src/composite.rs
+++ b/packages/rust/core/src/composite.rs
@@ -39,7 +39,7 @@ pub fn composite_rgba_on_white(rgba: &[u8]) -> Result<Vec<u8>, InvalidRgbaLength
         return Err(InvalidRgbaLength { len: rgba.len() });
     }
     let mut rgb = Vec::with_capacity(rgba.len() / 4 * 3);
-    for px in rgba.chunks_exact(4) {
+    for px in rgba.as_chunks::<4>().0 {
         let a = px[3] as f64 / 255.0;
         let inv = 1.0 - a;
         rgb.push((px[0] as f64 * a + 255.0 * inv).round() as u8);
diff --git a/packages/rust/core/src/dbs.rs b/packages/rust/core/src/dbs.rs
new file mode 100644
index 0000000..7a34b5d
--- /dev/null
+++ b/packages/rust/core/src/dbs.rs
@@ -0,0 +1,751 @@
+//! Colour Direct Binary Search (DBS): iterative refinement of an already-dithered image.
+//!
+//! Minimises perceived error `E = Σ_ch Σ_m ((p_ch ∗ e_ch)(m))²`, where `e = halftone − target`
+//! in YyCxCz and `p_ch` is a per-channel eye-model blur: narrow for luminance, wider for
+//! chrominance. YyCxCz is a linear transform of XYZ (linearised CIELAB), so blurring the error
+//! averages light the way the eye does — unlike OKLab, which `color_space_lab` uses for
+//! per-pixel matching.
+//!
+//! Each pass visits every pixel (raster order within tiles; tiles run in parallel, see
+//! `dbs_refine`) and tries changing it to every other ink and swapping it with each differing
+//! 8-neighbour, applying the move with the largest decrease in `E`. Trial cost is O(1) via the cached correlation `c_pe = c_pp ⋆ e`, where
+//! `c_pp = p ⋆ p` is the blur's autocorrelation; an accepted move updates `c_pe` over the
+//! footprint of `c_pp` only.
+//!
+//! Eye-model blurs are stored as full 2-D tables and `c_pp` is derived from them generically,
+//! so a non-separable model (e.g. the exponential CSF's spatial kernel) can be swapped in.
+//!
+//! References: Analoui & Allebach (1992); Lieberman & Allebach (1997); Agar & Allebach (2005);
+//! Flohr, Kolpatzik et al. (1993, YyCxCz); Kolpatzik & Bouman (1992, CSF constants);
+//! Mullen (1985). Ported from OpenDithering's `dbs.ts` (MIT, © Guy Sie).
+
+use rayon::prelude::*;
+
+use crate::color_space::srgb_channel_to_linear;
+use crate::palettes::Palette;
+
+// Kolpatzik & Bouman (1992) exponential CSF decay constants, W(f) ∝ exp(−α·f), f in cycles/degree.
+// Luminance: Näsänen's α(L) = 1 / (0.525·ln L + 3.91) at L = 11 cd/m² ≈ 0.193.
+// Chrominance: fitted to Mullen (1985), α = 0.419 for both opponent channels.
+const ALPHA_LUM: f64 = 0.193;
+const ALPHA_CHROMA: f64 = 0.419;
+
+/// Linear sRGB → CIE XYZ (D65), IEC 61966-2-1.
+const SRGB_TO_XYZ: [[f64; 3]; 3] = [
+    [0.412_456_4, 0.357_576_1, 0.180_437_5],
+    [0.212_672_9, 0.715_152_2, 0.072_175_0],
+    [0.019_333_9, 0.119_192_0, 0.950_304_1],
+];
+
+/// sRGB bytes → YyCxCz: `Yy = 116·Y/Yn`, `Cx = 500·(X/Xn − Y/Yn)`, `Cz = 200·(Y/Yn − Z/Zn)`.
+/// The white point is sRGB white, so (255, 255, 255) maps to zero chroma. Measured palettes
+/// are normalised to paper white, which therefore lands at the same point.
+fn srgb_to_yycxcz(rgb: [u8; 3]) -> [f64; 3] {
+    let lin = rgb.map(srgb_channel_to_linear);
+    let [x, y, z] = SRGB_TO_XYZ.map(|row| row[0] * lin[0] + row[1] * lin[1] + row[2] * lin[2]);
+    let [xn, yn, zn] = SRGB_TO_XYZ.map(|row| row[0] + row[1] + row[2]);
+    [116.0 * y / yn, 500.0 * (x / xn - y / yn), 200.0 * (y / yn - z / zn)]
+}
+
+/// Square, centred 2-D table: `data[(dy + r)·size + (dx + r)]` for `|dx|, |dy| ≤ r`.
+#[derive(Clone)]
+struct Kernel2D {
+    radius: usize,
+    size: usize,
+    data: Vec<f64>,
+}
+
+impl Kernel2D {
+    #[inline]
+    fn at(&self, dx: isize, dy: isize) -> f64 {
+        let r = self.radius as isize;
+        self.data[((dy + r) * self.size as isize + dx + r) as usize]
+    }
+}
+
+/// Blur truncation in standard deviations. Update cost grows with the square of the
+/// autocorrelation radius (2× the blur radius), so this is the main speed/accuracy knob.
+const BLUR_TRUNCATION_SIGMAS: f64 = 2.5;
+
+/// Normalised (unit-sum) sampled Gaussian, truncated at `ceil(BLUR_TRUNCATION_SIGMAS·σ)`.
+fn gaussian_blur(sigma: f64) -> Kernel2D {
+    let r = (BLUR_TRUNCATION_SIGMAS * sigma).ceil().max(1.0) as usize;
+    let size = 2 * r + 1;
+    let g: Vec<f64> = (0..size)
+        .map(|i| {
+            let d = i as f64 - r as f64;
+            (-(d * d) / (2.0 * sigma * sigma)).exp()
+        })
+        .collect();
+    let sum: f64 = g.iter().sum();
+    let data = (0..size * size).map(|i| g[i / size] * g[i % size] / (sum * sum)).collect();
+    Kernel2D { radius: r, size, data }
+}
+
+/// `c(d) = Σ_u p(u)·p(u + d)`, with twice the radius of `p`.
+fn autocorrelate(p: &Kernel2D) -> Kernel2D {
+    let pr = p.radius as isize;
+    let radius = 2 * p.radius;
+    let size = 2 * radius + 1;
+    let mut data = vec![0.0; size * size];
+    for dy in -(radius as isize)..=radius as isize {
+        for dx in -(radius as isize)..=radius as isize {
+            let mut s = 0.0;
+            for uy in (-pr).max(-pr - dy)..=pr.min(pr - dy) {
+                for ux in (-pr).max(-pr - dx)..=pr.min(pr - dx) {
+                    s += p.at(ux, uy) * p.at(ux + dx, uy + dy);
+                }
+            }
+            data[(dy + radius as isize) as usize * size + (dx + radius as isize) as usize] = s;
+        }
+    }
+    Kernel2D { radius, size, data }
+}
+
+/// Eye-model blur width in pixels for a CSF decay constant at the given viewing geometry.
+fn eye_sigma_px(alpha: f64, params: &DbsParams) -> f64 {
+    let pixel_deg = (2.54 / params.ppi / params.viewing_distance_cm).to_degrees();
+    alpha / (2.0 * std::f64::consts::PI) / pixel_deg
+}
+
+/// Which eye-model blur each YyCxCz channel uses: Yy → luminance (0), Cx and Cz → chrominance (1).
+const CHANNEL_BLUR: [usize; 3] = [0, 1, 1];
+
+/// Eye-model blurs `[luminance, chrominance]`; index with `CHANNEL_BLUR`.
+fn eye_blurs(params: &DbsParams) -> [Kernel2D; 2] {
+    if let Err(e) = params.validate() {
+        panic!("{e}");
+    }
+    let lum = gaussian_blur(eye_sigma_px(ALPHA_LUM, params));
+    let chroma = gaussian_blur(eye_sigma_px(ALPHA_CHROMA, params));
+    [lum, chroma]
+}
+
+/// Palette inks and target pixels in YyCxCz.
+fn to_yycxcz(pixels: &[u8], width: usize, height: usize, palette: &Palette) -> (Vec<[f64; 3]>, Vec<[f64; 3]>) {
+    assert_eq!(pixels.len(), width * height * 3, "pixel buffer size mismatch");
+    let inks = palette.colors.iter().map(|&c| srgb_to_yycxcz(c)).collect();
+    let target = pixels.par_chunks_exact(3).map(|c| srgb_to_yycxcz([c[0], c[1], c[2]])).collect();
+    (inks, target)
+}
+
+/// `out(m) = Σ_d c(d)·e(m + d)` over the image (e is zero outside it).
+fn correlate(e: &[f64], width: usize, height: usize, c: &Kernel2D) -> Vec<f64> {
+    let r = c.radius as isize;
+    let (w, h) = (width as isize, height as isize);
+    let mut out = vec![0.0; e.len()];
+    out.par_chunks_mut(width).enumerate().for_each(|(y, row)| {
+        let y = y as isize;
+        for (x, o) in row.iter_mut().enumerate() {
+            let x = x as isize;
+            let mut s = 0.0;
+            for dy in (-r).max(-y)..=r.min(h - 1 - y) {
+                let base = ((y + dy) * w + x) as usize;
+                for dx in (-r).max(-x)..=r.min(w - 1 - x) {
+                    s += c.at(dx, dy) * e[(base as isize + dx) as usize];
+                }
+            }
+            *o = s;
+        }
+    });
+    out
+}
+
+/// 8-neighbourhood offsets `(dx, dy)` tried for swaps.
+const NEIGHBOURS: [(isize, isize); 8] =
+    [(-1, -1), (0, -1), (1, -1), (-1, 0), (1, 0), (-1, 1), (0, 1), (1, 1)];
+
+/// `(x + dx, y + dy)` if it lies inside a `width × height` image.
+fn neighbour(x: usize, y: usize, dx: isize, dy: isize, width: usize, height: usize) -> Option<(usize, usize)> {
+    let nx = x.checked_add_signed(dx).filter(|&v| v < width)?;
+    let ny = y.checked_add_signed(dy).filter(|&v| v < height)?;
+    Some((nx, ny))
+}
+
+/// Parameters for DBS refinement.
+#[derive(Debug, Clone, Copy, PartialEq)]
+pub struct DbsParams {
+    /// Distance between viewer and panel, in centimetres.
+    pub viewing_distance_cm: f64,
+    /// Panel pixel density, in pixels per inch.
+    pub ppi: f64,
+    /// Upper bound on full-image passes; refinement stops earlier once a pass changes nothing.
+    pub max_passes: u32,
+}
+
+impl Default for DbsParams {
+    /// 40 cm from a ~127 ppi panel (7.3" 800×480), at most 10 passes.
+    fn default() -> Self {
+        Self { viewing_distance_cm: 40.0, ppi: 127.0, max_passes: 10 }
+    }
+}
+
+/// A `DbsParams` field outside its valid range (finite and > 0).
+#[derive(Debug, Clone, Copy, PartialEq)]
+pub struct InvalidDbsParams {
+    pub field: &'static str,
+    pub value: f64,
+}
+
+impl std::fmt::Display for InvalidDbsParams {
+    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
+        write!(f, "DBS {} must be finite and greater than 0, got {}", self.field, self.value)
+    }
+}
+
+impl std::error::Error for InvalidDbsParams {}
+
+impl DbsParams {
+    /// Check that the viewing geometry is usable. FFI boundaries call this to report a clean
+    /// error; `dbs_refine` panics on invalid parameters.
+    pub fn validate(&self) -> Result<(), InvalidDbsParams> {
+        for (field, value) in [("ppi", self.ppi), ("viewing_distance_cm", self.viewing_distance_cm)] {
+            if !(value.is_finite() && value > 0.0) {
+                return Err(InvalidDbsParams { field, value });
+            }
+        }
+        Ok(())
+    }
+}
+
+/// Diagnostics from a `dbs_refine` run.
+#[derive(Debug, Clone, PartialEq)]
+pub struct DbsStats {
+    /// Passes actually run (≤ `max_passes`).
+    pub passes: u32,
+    /// Number of accepted changes in each pass.
+    pub accepted_per_pass: Vec<usize>,
+    /// Perceived error of the starting halftone.
+    pub initial_error: f64,
+    /// Perceived error of the refined halftone.
+    pub final_error: f64,
+}
+
+/// Calls `f(c_pe(x, y), k(x − mx, y − my))` over the footprint of `k` centred on `(mx, my)`,
+/// clipped to the image.
+fn sweep(
+    cpe: &mut [[f64; 3]],
+    width: usize,
+    height: usize,
+    k: &Kernel2D,
+    (mx, my): (usize, usize),
+    mut f: impl FnMut(&mut [f64; 3], f64),
+) {
+    let r = k.radius;
+    let (x0, x1) = (mx.saturating_sub(r), (mx + r).min(width - 1));
+    // Kernel column of x0 is x0 − mx + r.
+    let kx0 = x0 + r - mx;
+    for y in my.saturating_sub(r)..=(my + r).min(height - 1) {
+        let krow = (y + r - my) * k.size;
+        let kernel = &k.data[krow + kx0..=krow + kx0 + (x1 - x0)];
+        let row = &mut cpe[y * width + x0..=y * width + x1];
+        for (c, &kv) in row.iter_mut().zip(kernel) {
+            f(c, kv);
+        }
+    }
+}
+
+/// Read-only state shared by every tile.
+struct Shared<'a> {
+    width: usize,
+    height: usize,
+    inks: &'a [[f64; 3]],
+    cpp_lum: &'a Kernel2D,
+    cpp_chroma: &'a Kernel2D,
+    pinned: Option<&'a [bool]>,
+}
+
+impl Shared<'_> {
+    fn is_pinned(&self, i: usize) -> bool {
+        self.pinned.is_some_and(|p| p[i])
+    }
+}
+
+/// Inclusive-exclusive rectangle `[x0, x1) × [y0, y1)` in image coordinates.
+#[derive(Clone, Copy)]
+struct Rect {
+    x0: usize,
+    y0: usize,
+    x1: usize,
+    y1: usize,
+}
+
+impl Rect {
+    fn width(&self) -> usize {
+        self.x1 - self.x0
+    }
+}
+
+/// A tile's private working copy of `c_pe` and `indices` over `rect` (the tile plus a margin
+/// wide enough to hold every read and write its moves make).
+struct Region {
+    rect: Rect,
+    cpe: Vec<[f64; 3]>,
+    indices: Vec<u8>,
+}
+
+impl Region {
+    fn copy_from(rect: Rect, cpe: &[[f64; 3]], indices: &[u8], width: usize) -> Self {
+        let rows = rect.y0..rect.y1;
+        let span = |y: usize| y * width + rect.x0..y * width + rect.x1;
+        Self {
+            rect,
+            cpe: rows.clone().flat_map(|y| cpe[span(y)].iter().copied()).collect(),
+            indices: rows.flat_map(|y| indices[span(y)].iter().copied()).collect(),
+        }
+    }
+
+    fn copy_into(&self, cpe: &mut [[f64; 3]], indices: &mut [u8], width: usize) {
+        let w = self.rect.width();
+        for (ry, y) in (self.rect.y0..self.rect.y1).enumerate() {
+            let span = y * width + self.rect.x0..y * width + self.rect.x1;
+            cpe[span.clone()].copy_from_slice(&self.cpe[ry * w..(ry + 1) * w]);
+            indices[span].copy_from_slice(&self.indices[ry * w..(ry + 1) * w]);
+        }
+    }
+}
+
+/// One raster pass of DBS over the pixels of `tile`, working in `region`. Returns the number
+/// of accepted moves and the summed change in perceived error.
+fn refine_tile(sh: &Shared, region: &mut Region, tile: Rect) -> (usize, f64) {
+    const EPS: f64 = 1e-9;
+    let Region { rect, cpe, indices } = region;
+    let (rw, rh) = (rect.width(), rect.y1 - rect.y0);
+    let local = |x: usize, y: usize| (y - rect.y0) * rw + (x - rect.x0);
+    let cpp = [sh.cpp_lum, sh.cpp_chroma, sh.cpp_chroma];
+    let cpp0 = cpp.map(|k| k.at(0, 0));
+    let inks = sh.inks;
+
+    // c_pe += a_ch·c_pp_ch(· − m) for every channel. Cx and Cz share the chrominance kernel,
+    // so they are updated in a single sweep. The region margin covers the whole footprint, so
+    // clipping to the region only ever clips at true image edges.
+    let add_kernels = |cpe: &mut [[f64; 3]], (x, y): (usize, usize), a: [f64; 3]| {
+        let m = (x - rect.x0, y - rect.y0);
+        sweep(cpe, rw, rh, sh.cpp_lum, m, |c, k| c[0] += a[0] * k);
+        sweep(cpe, rw, rh, sh.cpp_chroma, m, |c, k| {
+            c[1] += a[1] * k;
+            c[2] += a[2] * k;
+        });
+    };
+
+    let (mut accepted, mut delta) = (0, 0.0);
+    for y in tile.y0..tile.y1 {
+        for x in tile.x0..tile.x1 {
+            if sh.is_pinned(y * sh.width + x) {
+                continue;
+            }
+            let m = local(x, y);
+            let k = indices[m] as usize;
+            let c = cpe[m];
+
+            enum Move { Toggle(usize), Swap { n: usize, nx: usize, ny: usize } }
+            let mut best_de = -EPS;
+            let mut best = None;
+
+            for j in (0..inks.len()).filter(|&j| j != k) {
+                let de: f64 = (0..3)
+                    .map(|ch| {
+                        let a = inks[j][ch] - inks[k][ch];
+                        cpp0[ch] * a * a + 2.0 * a * c[ch]
+                    })
+                    .sum();
+                if de < best_de {
+                    best_de = de;
+                    best = Some(Move::Toggle(j));
+                }
+            }
+
+            for (dx, dy) in NEIGHBOURS {
+                let Some((nx, ny)) = neighbour(x, y, dx, dy, sh.width, sh.height) else { continue };
+                let n = local(nx, ny);
+                let kn = indices[n] as usize;
+                if kn == k || sh.is_pinned(ny * sh.width + nx) {
+                    continue;
+                }
+                // m takes kn (Δ = a), n takes k (Δ = −a).
+                let de: f64 = (0..3)
+                    .map(|ch| {
+                        let a = inks[kn][ch] - inks[k][ch];
+                        2.0 * a * a * (cpp0[ch] - cpp[ch].at(dx, dy)) + 2.0 * a * (c[ch] - cpe[n][ch])
+                    })
+                    .sum();
+                if de < best_de {
+                    best_de = de;
+                    best = Some(Move::Swap { n, nx, ny });
+                }
+            }
+
+            let Some(mv) = best else { continue };
+            accepted += 1;
+            delta += best_de;
+            match mv {
+                Move::Toggle(j) => {
+                    add_kernels(cpe, (x, y), std::array::from_fn(|ch| inks[j][ch] - inks[k][ch]));
+                    indices[m] = j as u8;
+                }
+                Move::Swap { n, nx, ny } => {
+                    let kn = indices[n] as usize;
+                    let a: [f64; 3] = std::array::from_fn(|ch| inks[kn][ch] - inks[k][ch]);
+                    add_kernels(cpe, (x, y), a);
+                    add_kernels(cpe, (nx, ny), a.map(|v| -v));
+                    indices.swap(m, n);
+                }
+            }
+        }
+    }
+    (accepted, delta)
+}
+
+/// Smallest tile side. Larger tiles amortise the region copies; smaller ones expose more
+/// parallelism. The correctness bound (`2·margin + 1`) is enforced separately.
+const MIN_TILE: usize = 64;
+
+/// Refine palette `indices` in place to minimise perceived error against `pixels` (sRGB,
+/// `width × height × 3`), as seen through `palette` (pass the measured palette).
+///
+/// `pinned` pixels are never changed or swapped, and their target is taken to be their own
+/// ink, so they neither carry error nor pull neighbours to compensate for it.
+///
+/// Each pass runs in four phases over a 2×2 checkerboard of tiles. A move inside a tile reads
+/// within 1 px of it and writes within `margin = R + 1` of it (R = autocorrelation radius), so
+/// tiles of one phase — at least one tile side (≥ 2·margin + 1) apart — never touch each
+/// other's data and run in parallel on private copies. The result depends only on the fixed
+/// tile order, not on thread count or scheduling.
+///
+/// # Panics
+/// Panics if buffer lengths disagree with `width × height`, or on a non-positive `ppi` or
+/// viewing distance.
+pub fn dbs_refine(
+    pixels: &[u8],
+    width: usize,
+    height: usize,
+    palette: &Palette,
+    indices: &mut [u8],
+    pinned: Option<&[bool]>,
+    params: &DbsParams,
+) -> DbsStats {
+    let n_px = width * height;
+    assert_eq!(indices.len(), n_px, "index buffer size mismatch");
+    if let Some(p) = pinned {
+        assert_eq!(p.len(), n_px, "pinned mask size mismatch");
+    }
+
+    let (inks, mut target) = to_yycxcz(pixels, width, height, palette);
+    for (i, t) in target.iter_mut().enumerate() {
+        if pinned.is_some_and(|p| p[i]) {
+            *t = inks[indices[i] as usize];
+        }
+    }
+
+    let [cpp_lum, cpp_chroma] = eye_blurs(params).map(|p| autocorrelate(&p));
+    let cpp = [&cpp_lum, &cpp_chroma, &cpp_chroma];
+
+    // c_pe per pixel, channels interleaved so one footprint sweep touches each cache line once.
+    let per_channel: [Vec<f64>; 3] = std::array::from_fn(|ch| {
+        let e: Vec<f64> = (0..n_px).map(|i| inks[indices[i] as usize][ch] - target[i][ch]).collect();
+        correlate(&e, width, height, cpp[ch])
+    });
+    let mut cpe: Vec<[f64; 3]> = (0..n_px).map(|i| per_channel.each_ref().map(|c| c[i])).collect();
+    let mut err: f64 = (0..3)
+        .map(|ch| {
+            (0..n_px)
+                .map(|i| (inks[indices[i] as usize][ch] - target[i][ch]) * cpe[i][ch])
+                .sum::<f64>()
+        })
+        .sum();
+    let initial_error = err;
+
+    let sh = Shared { width, height, inks: &inks, cpp_lum: &cpp_lum, cpp_chroma: &cpp_chroma, pinned };
+    let margin = cpp_lum.radius.max(cpp_chroma.radius) + 1;
+    let side = MIN_TILE.max(2 * margin + 1);
+    let tiles = |phase_x: usize, phase_y: usize| {
+        let (cols, rows) = (width.div_ceil(side), height.div_ceil(side));
+        (0..rows)
+            .filter(move |ty| ty % 2 == phase_y)
+            .flat_map(move |ty| (0..cols).filter(move |tx| tx % 2 == phase_x).map(move |tx| (tx, ty)))
+            .map(|(tx, ty)| Rect {
+                x0: tx * side,
+                y0: ty * side,
+                x1: ((tx + 1) * side).min(width),
+                y1: ((ty + 1) * side).min(height),
+            })
+            .collect::<Vec<_>>()
+    };
+    let with_margin = |t: Rect| Rect {
+        x0: t.x0.saturating_sub(margin),
+        y0: t.y0.saturating_sub(margin),
+        x1: (t.x1 + margin).min(width),
+        y1: (t.y1 + margin).min(height),
+    };
+
+    let mut accepted_per_pass = Vec::new();
+    for _ in 0..params.max_passes {
+        let mut accepted = 0;
+        for (phase_x, phase_y) in [(0, 0), (1, 0), (0, 1), (1, 1)] {
+            let results: Vec<_> = tiles(phase_x, phase_y)
+                .into_par_iter()
+                .map(|tile| {
+                    let mut region = Region::copy_from(with_margin(tile), &cpe, indices, width);
+                    let (n, delta) = refine_tile(&sh, &mut region, tile);
+                    (region, n, delta)
+                })
+                .collect();
+            // Regions of one phase are disjoint; fold in tile order for a deterministic sum.
+            for (region, n, delta) in results {
+                region.copy_into(&mut cpe, indices, width);
+                accepted += n;
+                err += delta;
+            }
+        }
+        accepted_per_pass.push(accepted);
+        if accepted == 0 {
+            break;
+        }
+    }
+
+    DbsStats {
+        passes: accepted_per_pass.len() as u32,
+        accepted_per_pass,
+        initial_error,
+        final_error: err,
+    }
+}
+
+/// Perceived error `Σ_ch Σ_m ((p_ch ∗ e_ch)(m))²` of `indices` against `pixels`, computed
+/// directly by blurring the error over the whole plane (no correlation shortcut). This is the
+/// quantity `dbs_refine` minimises; useful for comparing halftones under the eye model.
+pub fn perceived_error(
+    pixels: &[u8],
+    width: usize,
+    height: usize,
+    palette: &Palette,
+    indices: &[u8],
+    params: &DbsParams,
+) -> f64 {
+    assert_eq!(indices.len(), width * height, "index buffer size mismatch");
+    let (inks, target) = to_yycxcz(pixels, width, height, palette);
+    let blurs = eye_blurs(params);
+    let (w, h) = (width as isize, height as isize);
+
+    (0..3)
+        .map(|ch| {
+            let p = &blurs[CHANNEL_BLUR[ch]];
+            let r = p.radius as isize;
+            let e: Vec<f64> = (0..indices.len()).map(|i| inks[indices[i] as usize][ch] - target[i][ch]).collect();
+            // The blurred error is non-zero up to r pixels outside the image.
+            (-r..h + r)
+                .into_par_iter()
+                .map(|oy| {
+                    let mut s = 0.0;
+                    for ox in -r..w + r {
+                        let mut v = 0.0;
+                        for y in (oy - r).max(0)..=(oy + r).min(h - 1) {
+                            for x in (ox - r).max(0)..=(ox + r).min(w - 1) {
+                                v += p.at(ox - x, oy - y) * e[(y * w + x) as usize];
+                            }
+                        }
+                        s += v * v;
+                    }
+                    s
+                })
+                .sum::<f64>()
+        })
+        .sum()
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use crate::algorithms::{error_diffusion_dither, BURKES};
+    use crate::measured_palettes::SPECTRA_7_3_6COLOR;
+    use crate::palettes::ColorScheme;
+    use approx::assert_relative_eq;
+
+    /// Smooth colour gradient: enough structure that Burkes leaves room for improvement.
+    fn gradient(w: usize, h: usize) -> Vec<u8> {
+        let mut px = Vec::with_capacity(w * h * 3);
+        for y in 0..h {
+            for x in 0..w {
+                px.push((x * 255 / (w - 1)) as u8);
+                px.push((y * 255 / (h - 1)) as u8);
+                px.push((255 - x * 255 / (w - 1)) as u8);
+            }
+        }
+        px
+    }
+
+    fn burkes(px: &[u8], w: usize, h: usize, palette: &Palette) -> Vec<u8> {
+        error_diffusion_dither(px, w, h, palette, &BURKES, true)
+    }
+
+    #[test]
+    fn validate_accepts_defaults() {
+        assert_eq!(DbsParams::default().validate(), Ok(()));
+    }
+
+    #[test]
+    fn validate_rejects_non_positive_or_non_finite_geometry() {
+        for (params, field) in [
+            (DbsParams { ppi: 0.0, ..Default::default() }, "ppi"),
+            (DbsParams { ppi: f64::NAN, ..Default::default() }, "ppi"),
+            (DbsParams { viewing_distance_cm: -1.0, ..Default::default() }, "viewing_distance_cm"),
+            (DbsParams { viewing_distance_cm: f64::INFINITY, ..Default::default() }, "viewing_distance_cm"),
+        ] {
+            let err = params.validate().unwrap_err();
+            assert_eq!(err.field, field);
+            assert!(err.to_string().contains(field), "message should name the field: {err}");
+        }
+    }
+
+    #[test]
+    fn perceived_error_is_zero_when_halftone_equals_target() {
+        let palette = &SPECTRA_7_3_6COLOR;
+        let indices: Vec<u8> = (0..24).map(|i| (i % palette.colors.len()) as u8).collect();
+        let px: Vec<u8> = indices.iter().flat_map(|&i| palette.colors[i as usize]).collect();
+        let e = perceived_error(&px, 6, 4, palette, &indices, &DbsParams::default());
+        assert!(e.abs() < 1e-9, "expected zero error, got {e}");
+    }
+
+    #[test]
+    fn perceived_error_is_positive_for_a_mismatched_halftone() {
+        let palette = ColorScheme::Mono.palette();
+        let px = vec![128u8; 8 * 8 * 3];
+        let indices = vec![0u8; 64]; // all black for a mid-grey target
+        assert!(perceived_error(&px, 8, 8, palette, &indices, &DbsParams::default()) > 0.0);
+    }
+
+    #[test]
+    fn refine_reduces_perceived_error() {
+        let (w, h) = (32, 24);
+        let px = gradient(w, h);
+        let palette = &SPECTRA_7_3_6COLOR;
+        let params = DbsParams::default();
+        let mut idx = burkes(&px, w, h, palette);
+        let before = perceived_error(&px, w, h, palette, &idx, &params);
+
+        dbs_refine(&px, w, h, palette, &mut idx, None, &params);
+
+        let after = perceived_error(&px, w, h, palette, &idx, &params);
+        assert!(after < before * 0.95, "DBS should cut error noticeably: {before} -> {after}");
+    }
+
+    /// The incremental bookkeeping (ΔE formulas, correlation updates) must agree with a
+    /// from-scratch blur-and-square of the final halftone. Catches wrong toggle/swap deltas.
+    ///
+    /// The 200×150 case spans several tiles, so it also checks that tile margins are wide
+    /// enough: an update clipped at a margin would silently leave `c_pe` wrong.
+    #[test]
+    fn reported_errors_match_direct_computation() {
+        for (w, h) in [(32, 24), (200, 150)] {
+            let px = gradient(w, h);
+            let palette = &SPECTRA_7_3_6COLOR;
+            let params = DbsParams::default();
+            let mut idx = burkes(&px, w, h, palette);
+            let before = perceived_error(&px, w, h, palette, &idx, &params);
+
+            let stats = dbs_refine(&px, w, h, palette, &mut idx, None, &params);
+
+            let after = perceived_error(&px, w, h, palette, &idx, &params);
+            assert_relative_eq!(stats.initial_error, before, max_relative = 1e-9);
+            assert_relative_eq!(stats.final_error, after, max_relative = 1e-9);
+        }
+    }
+
+    /// After convergence no single toggle or 8-neighbour swap may lower the error, judged by
+    /// brute-force recomputation — the defining property of a DBS local optimum.
+    #[test]
+    fn converged_result_is_a_local_minimum() {
+        let (w, h) = (10, 8);
+        let px = gradient(w, h);
+        let palette = &SPECTRA_7_3_6COLOR;
+        let params = DbsParams { max_passes: 200, ..Default::default() };
+        let mut idx = burkes(&px, w, h, palette);
+
+        let stats = dbs_refine(&px, w, h, palette, &mut idx, None, &params);
+        assert_eq!(stats.accepted_per_pass.last(), Some(&0), "should converge: {stats:?}");
+
+        let base = perceived_error(&px, w, h, palette, &idx, &params);
+        let tol = base * 1e-9 + 1e-9;
+        for m in 0..w * h {
+            for k in 0..palette.colors.len() as u8 {
+                let mut t = idx.clone();
+                t[m] = k;
+                let e = perceived_error(&px, w, h, palette, &t, &params);
+                assert!(e >= base - tol, "toggle at {m} to {k} lowers error {base} -> {e}");
+            }
+            for (dx, dy) in NEIGHBOURS {
+                let Some((nx, ny)) = neighbour(m % w, m / w, dx, dy, w, h) else { continue };
+                let mut t = idx.clone();
+                t.swap(m, ny * w + nx);
+                let e = perceived_error(&px, w, h, palette, &t, &params);
+                assert!(e >= base - tol, "swap {m} with ({nx},{ny}) lowers error {base} -> {e}");
+            }
+        }
+    }
+
+    #[test]
+    fn pinned_pixels_never_change() {
+        let (w, h) = (32, 24);
+        let px = gradient(w, h);
+        let palette = &SPECTRA_7_3_6COLOR;
+        let mut idx = burkes(&px, w, h, palette);
+        let pinned: Vec<bool> = (0..w * h).map(|i| i % 3 == 0).collect();
+        let original = idx.clone();
+
+        let stats = dbs_refine(&px, w, h, palette, &mut idx, Some(&pinned), &DbsParams::default());
+
+        assert!(stats.accepted_per_pass.iter().sum::<usize>() > 0, "unpinned pixels should still move");
+        for i in (0..w * h).filter(|&i| pinned[i]) {
+            assert_eq!(idx[i], original[i], "pinned pixel {i} changed");
+        }
+    }
+
+    /// Pins the colour-science decision that DBS measures error in linear light: a flat
+    /// sRGB 188 (linear ≈ 0.50) on black/white must come out ≈ 50% white ink, not the
+    /// ≈ 74% an sRGB-space average would give.
+    #[test]
+    fn flat_grey_reaches_linear_light_coverage() {
+        let (w, h) = (48, 48);
+        let px = vec![188u8; w * h * 3];
+        let palette = ColorScheme::Mono.palette();
+        let white = palette.colors.iter().position(|&c| c == [255, 255, 255]).unwrap() as u8;
+        let mut idx = burkes(&px, w, h, palette);
+
+        dbs_refine(&px, w, h, palette, &mut idx, None, &DbsParams { max_passes: 50, ..Default::default() });
+
+        let frac = idx.iter().filter(|&&i| i == white).count() as f64 / (w * h) as f64;
+        assert!((0.46..=0.54).contains(&frac), "white coverage {frac:.3}, expected ≈ 0.50");
+    }
+
+    /// Tiles run in parallel; the result must depend only on the fixed tile order, never on
+    /// thread count or scheduling. The image spans several tiles in each direction.
+    #[test]
+    fn result_is_independent_of_thread_count() {
+        let (w, h) = (200, 150);
+        let px = gradient(w, h);
+        let palette = &SPECTRA_7_3_6COLOR;
+        let init = burkes(&px, w, h, palette);
+        let run = |threads: usize| {
+            let pool = rayon::ThreadPoolBuilder::new().num_threads(threads).build().unwrap();
+            let mut idx = init.clone();
+            let stats = pool.install(|| dbs_refine(&px, w, h, palette, &mut idx, None, &DbsParams::default()));
+            (idx, stats)
+        };
+        let (one, s1) = run(1);
+        let (four, s4) = run(4);
+        assert_eq!(one, four);
+        assert_eq!(s1, s4);
+    }
+
+    #[test]
+    fn zero_passes_leaves_indices_untouched() {
+        let (w, h) = (16, 12);
+        let px = gradient(w, h);
+        let palette = &SPECTRA_7_3_6COLOR;
+        let mut idx = burkes(&px, w, h, palette);
+        let original = idx.clone();
+        let stats = dbs_refine(&px, w, h, palette, &mut idx, None, &DbsParams { max_passes: 0, ..Default::default() });
+        assert_eq!(idx, original);
+        assert_eq!(stats.passes, 0);
+    }
+}
diff --git a/packages/rust/core/src/lib.rs b/packages/rust/core/src/lib.rs
index 06d7918..330763a 100644
--- a/packages/rust/core/src/lib.rs
+++ b/packages/rust/core/src/lib.rs
@@ -2,6 +2,7 @@ pub mod algorithms;
 pub mod color_space;
 pub mod color_space_lab;
 pub mod composite;
+pub mod dbs;
 pub mod dizzy;
 pub mod enums;
 pub mod error;
@@ -12,6 +13,7 @@ pub mod tone_map;
 pub mod types;
 
 use crate::color_space::{linear_channel_to_srgb, srgb_channel_to_linear};
+use crate::dbs::DbsParams;
 use crate::enums::{DitherMode, GamutCompression, ToneCompression};
 use crate::palettes::Palette;
 use crate::types::ImageBuffer;
@@ -24,7 +26,7 @@ use crate::types::ImageBuffer;
 /// ```
 ///
 /// Pre-processing pipeline (applied in order, each step is a no-op at its identity value):
-/// `exposure → saturation → shadows/highlights → tone → gamut → dither`.
+/// `exposure → saturation → shadows/highlights → tone → gamut → dither → dbs`.
 #[derive(Debug, Clone, PartialEq)]
 pub struct DitherConfig {
     /// Dithering algorithm.
@@ -45,6 +47,9 @@ pub struct DitherConfig {
     pub tone: ToneCompression,
     /// Gamut compression. Auto = full strength on out-of-gamut; Fixed(s) = manual.
     pub gamut: GamutCompression,
+    /// Optional Direct Binary Search refinement of the dithered output (see `dbs`). Slow:
+    /// seconds on a full panel, so `None` by default.
+    pub dbs: Option<DbsParams>,
 }
 
 impl Default for DitherConfig {
@@ -58,6 +63,7 @@ impl Default for DitherConfig {
             highlights: 0.0,
             tone:       ToneCompression::Fixed(0.0),
             gamut:      GamutCompression::None,
+            dbs:        None,
         }
     }
 }
@@ -162,13 +168,16 @@ fn dither_impl(
         {
             return indices;
         }
-        return dispatch(img, p, canonical, config.mode, config.serpentine, pin_exact_pixels);
+        let indices = dispatch(img, p, canonical, config.mode, config.serpentine, pin_exact_pixels);
+        return refine(img, p, canonical, indices, config.dbs, pin_exact_pixels);
     }
 
     // Convert sRGB bytes → linear, apply pre-processing pipeline, convert back.
     let mut linear: Vec<[f64; 3]> = img
         .data
-        .chunks_exact(3)
+        .as_chunks::<3>()
+        .0
+        .iter()
         .map(|c| [
             srgb_channel_to_linear(c[0]),
             srgb_channel_to_linear(c[1]),
@@ -192,7 +201,31 @@ fn dither_impl(
         .collect();
 
     let processed_img = ImageBuffer::new(&processed, img.width);
-    dispatch(&processed_img, p, canonical, config.mode, config.serpentine, pin_exact_pixels)
+    let indices = dispatch(&processed_img, p, canonical, config.mode, config.serpentine, pin_exact_pixels);
+    refine(&processed_img, p, canonical, indices, config.dbs, pin_exact_pixels)
+}
+
+/// Final pipeline step: optional DBS refinement against the (pre-processed) target image.
+/// Pixels the dispatch pinned to an exact canonical colour stay pinned.
+fn refine(
+    img: &ImageBuffer,
+    p: &Palette,
+    canonical: &Palette,
+    mut indices: Vec<u8>,
+    params: Option<DbsParams>,
+    pin_exact_pixels: bool,
+) -> Vec<u8> {
+    let Some(params) = params else { return indices };
+    let pinned: Option<Vec<bool>> = pin_exact_pixels.then(|| {
+        img.data
+            .as_chunks::<3>()
+            .0
+            .iter()
+            .map(|rgb| algorithms::exact_palette_index(rgb, canonical).is_some())
+            .collect()
+    });
+    dbs::dbs_refine(img.data, img.width, img.height, p, &mut indices, pinned.as_deref(), &params);
+    indices
 }
 
 #[cfg(test)]
@@ -335,4 +368,36 @@ mod tests {
         );
         assert_eq!(output.len(), 4);
     }
+
+    #[test]
+    fn dbs_refinement_lowers_perceived_error() {
+        let pixels: Vec<u8> = (0..24 * 16)
+            .flat_map(|i| [(i % 24 * 10) as u8, (i / 24 * 16) as u8, 200 - (i % 24 * 8) as u8])
+            .collect();
+        let img = ImageBuffer::new(&pixels, 24);
+        let params = DbsParams::default();
+
+        let plain = dither(&img, &SPECTRA_7_3_6COLOR, DitherConfig::default());
+        let refined = dither(&img, &SPECTRA_7_3_6COLOR, DitherConfig { dbs: Some(params), ..Default::default() });
+
+        let err = |idx: &[u8]| dbs::perceived_error(&pixels, 24, 16, &SPECTRA_7_3_6COLOR, idx, &params);
+        assert!(err(&refined) < err(&plain), "{} !< {}", err(&refined), err(&plain));
+    }
+
+    #[test]
+    fn exact_canonical_pixels_stay_pinned_under_dbs() {
+        let mut image = pixels([128, 128, 128], 8);
+        image[0..3].copy_from_slice(&[0, 255, 0]);
+        image[9..12].copy_from_slice(&[0, 255, 0]);
+        let img = ImageBuffer::new(&image, 4);
+
+        let output = dither_with_canonical(
+            &img,
+            &SPECTRA_7_3_6COLOR,
+            ColorScheme::Bwgbry.palette(),
+            DitherConfig { dbs: Some(DbsParams::default()), ..Default::default() },
+        );
+        assert_eq!(output[0], 5);
+        assert_eq!(output[3], 5);
+    }
 }
diff --git a/packages/rust/ios/src/lib.rs b/packages/rust/ios/src/lib.rs
index ff0dd1a..5102a9d 100644
--- a/packages/rust/ios/src/lib.rs
+++ b/packages/rust/ios/src/lib.rs
@@ -67,7 +67,7 @@ fn build_palette(ptr: *const u8, len: usize, accent_idx: usize) -> Result<Option
     }
     // SAFETY: caller guarantees `ptr` is valid for `len` bytes; checked non-null above.
     let bytes = unsafe { std::slice::from_raw_parts(ptr, len) };
-    let colors: Vec<[u8; 3]> = bytes.chunks_exact(3).map(|c| [c[0], c[1], c[2]]).collect();
+    let colors: Vec<[u8; 3]> = bytes.as_chunks::<3>().0.to_vec();
     if colors.len() < 2 || accent_idx >= colors.len() {
         return Err(ED_ERR_BAD_PALETTE);
     }
diff --git a/packages/rust/wasm/src/lib.rs b/packages/rust/wasm/src/lib.rs
index 5bace9f..bd72bb1 100644
--- a/packages/rust/wasm/src/lib.rs
+++ b/packages/rust/wasm/src/lib.rs
@@ -1,5 +1,6 @@
 use epaper_dithering_core::{
     composite::composite_rgba_on_white,
+    dbs::DbsParams,
     dither, dither_with_canonical, DitherConfig,
     enums::{DitherMode, GamutCompression, ToneCompression},
     measured_palettes::CATALOG,
@@ -28,8 +29,19 @@ fn parse_gamut(v: Option<f64>) -> GamutCompression {
     }
 }
 
+/// `dbs_passes = None` disables DBS refinement; otherwise build and validate its parameters.
+fn parse_dbs(passes: Option<u32>, viewing_distance_cm: f64, ppi: f64) -> Result<Option<DbsParams>, JsValue> {
+    passes
+        .map(|max_passes| {
+            let params = DbsParams { viewing_distance_cm, ppi, max_passes };
+            params.validate().map(|()| params).map_err(|e| JsValue::from_str(&e.to_string()))
+        })
+        .transpose()
+}
+
 /// Dither a flat RGB image. Accepts either an idealized `scheme_id` or a measured
 /// `palette_bytes`/`accent_idx` pair; `palette_bytes` empty ⇒ use `scheme_id`.
+/// `dbs_passes` enables Direct Binary Search refinement with up to that many passes.
 ///
 /// Returns a `Uint8Array` of palette indices (one per pixel, length = width × height).
 #[wasm_bindgen]
@@ -48,6 +60,9 @@ pub fn dither_image(
     highlights: f64,
     tone: Option<f64>,
     gamut: Option<f64>,
+    dbs_passes: Option<u32>,
+    dbs_viewing_distance_cm: f64,
+    dbs_ppi: f64,
 ) -> Result<Vec<u8>, JsValue> {
     if width == 0 {
         return Err(JsValue::from_str("width must be greater than 0"));
@@ -70,6 +85,7 @@ pub fn dither_image(
         highlights,
         tone:  parse_tone(tone),
         gamut: parse_gamut(gamut),
+        dbs:   parse_dbs(dbs_passes, dbs_viewing_distance_cm, dbs_ppi)?,
     };
 
     if palette_bytes.is_empty() {