README: correct the Recreate Plugin Environment location

It is a right-click context-menu item on the plugin row in
Preferences -> PCB Editor -> Action Plugins, not a control on the
Plugins page. Also document the manual venv-delete fallback.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>

.gitea/workflows/build-pcm.yml

@@ -35,11 +35,12 @@ jobs:
             dist/*.zip
             dist/metadata-registry.json
 
-      - name: Create release with the zip attached
+      - name: Create release with the zip, registry metadata and figures
         if: startsWith(github.ref, 'refs/tags/v')
         uses: akkuman/gitea-release-action@b8d9144f302c68610911db1aaf722708d5c02d94 # v1
         with:
           files: |
             dist/*.zip
             dist/metadata-registry.json
+            docs/img/*.png
           token: ${{ secrets.GITEA_TOKEN }}

.gitignore

@@ -3,3 +3,11 @@ __pycache__/
 *.pyc
 .pytest_cache/
 dist/
+
+# local AI-tooling artifacts, never publish
+.claude/
+.claude-flow/
+.swarm/
+.mcp.json
+CLAUDE.md
+ruvector.db

README.md

@@ -6,11 +6,23 @@ between two contacts, **single- or multi-layer**: the chosen net's fills
 solved as coupled finite-difference sheets linked by the net's **via
 and through-hole-pad barrels** (18 µm plating, configurable). At a user-set **frequency** the exact 1D foil/barrel
 skin-effect correction is applied (AC results are a rigorous lower
-bound — see *Model & limits*). Shows per-layer rasterized maps,
-potential, current density, and **power density**, reports **per-via
+bound; see *Model & limits*). Shows per-layer rasterized maps,
+potential, current density, and **power density**, and reports **per-via
 currents** (via ampacity!) and total dissipation at a **selectable test
 current**. PNGs + a text summary are saved per run.
 
+![Current density on a two-layer demo net](docs/img/demo-current.png)
+*Real output on a synthetic two-layer net: current from a soldered
+THT-pad contact (V+, injected at the drill-wall ring) squeezes past a
+notch in the F.Cu pour, transfers through the stitching-via field into
+the B.Cu pour and leaves at the V− lug. Per-via currents and the
+hottest via are reported.*
+
+![Potential on the two-layer demo net](docs/img/demo-potential.png)
+*The matching potential map with equipotential contour lines: they
+bunch where the field is strongest — nearly the whole 8.7 mV drop
+happens around the notch on F.Cu.*
+
 Uses the KiCad **IPC API** (`kicad-python` / `kipy`), not the deprecated
 SWIG API. Requires KiCad **10.0.1+**.
 
@@ -23,7 +35,7 @@ SWIG API. Requires KiCad **10.0.1+**.
    on Windows or `/usr/bin/python3` on Linux (after a 9→10 upgrade it
    can point at KiCad 9).
 3. **Deploy** (dev checkout; end users install the PCM zip instead, see
-   *Packaging*):
+   *Packaging / publishing*):
    ```powershell
    powershell -ExecutionPolicy Bypass -File deploy.ps1        # junction (dev)
    powershell -ExecutionPolicy Bypass -File deploy.ps1 -Mode Copy
@@ -35,18 +47,30 @@ SWIG API. Requires KiCad **10.0.1+**.
    ```
 4. **Restart KiCad**; first load builds the plugin venv (numpy, scipy,
    matplotlib, PySide6 — takes minutes; the Ω button appears when done).
-   If stuck: Preferences → Plugins → *Recreate Plugin Environment*.
+   If stuck: in the PCB editor, Preferences → *PCB Editor → Action
+   Plugins*, **right-click** the plugin's row → *Recreate Plugin
+   Environment* (context menu only — there is no button). Manual
+   equivalent: delete
+   `%LOCALAPPDATA%\kicad\10.0\python-environments\th.co.b4l.fill-resistance`
+   and restart KiCad.
 
 ## Usage
 
 1. Mark the current-injection terminals. Each terminal may have
-   **multiple parts** (all merged into one externally-bonded contact):
+   **multiple parts** (all merged into one externally bonded contact):
    - **V+ rectangles on `User.1`**, **V− rectangles on `User.2`**
      (marker layers, configurable via `ELECTRODE_POS_LAYER` /
      `ELECTRODE_NEG_LAYER`), any number per side, axis-aligned;
-   - **pads** (real copper shape; through-hole pad contacts all layers,
-     SMD pad its own layer) — selected pads fill a side that has no
-     rectangles;
+   - **pads and vias** (SMD pad: real copper shape on its own layer;
+     through-hole pads and vias become **barrel contacts**: the current
+     enters at the drill wall on every spanned layer, see below).
+     Selected pads/vias fill the side that has **no rectangles**, so
+     mixing both kinds is the everyday workflow: e.g. select **one
+     rectangle on `User.1`** (V+) **plus any number of pads / THT
+     holes** (Ctrl-click) — the pads together form the V− terminal
+     (a connector's pin group, a via cluster, …). All selected
+     pads/vias go to that one side; if both marker layers already
+     provide rectangles, selecting pads on top is an error;
    - legacy: exactly 2 selected contacts with no marker rectangles still
      works; empty selection scans the whole board's marker layers.
 2. **Select the contacts**, click the **Fill Resistance** Ω button.
@@ -55,13 +79,28 @@ SWIG API. Requires KiCad **10.0.1+**.
    ("All selected layers" = bolted-lug/through contact), the **test
    current**, and optionally a grid cell size. Multiple layers are coupled
    through the net's via/pad barrels automatically.
-4. Read R / voltage drop / total power in the figure titles and status
+4. Wait for the solve. Depending on board size, included layers, cell
+   size and your hardware it can take **considerable time** — large
+   multi-layer pours at fine cell sizes may run for minutes (on our
+   test setup a typical real-board run finishes in ≈ 8 s). Then read
+   R / voltage drop / total power in the figure titles and status
    bar. Outputs land in `<board dir>\fill_res_results\<timestamp>\`:
    per-layer `1_raster_map` / `2_potential` / `3_current_density` /
    `4_power_density` PNGs, `summary.txt` (incl. the busiest vias with
    per-via current and dissipation, and the **current through each
    injection area** — computed flux with the equipotential model,
    prescribed area share with the uniform model), `geometry_dump.json`.
+5. **Experimental — overlays inside KiCad** (dialog checkbox, default
+   off; KiCad ≥ 10.0.1): after the solve, the per-layer **|J| heatmaps
+   are pushed into the open board** as unlocked reference images on
+   `User.9`…`User.12` (`OVERLAY_LAYERS`; enable them in Board Setup),
+   copper layers mapped in stackup order, top first. Toggle them in the
+   Appearance panel like any layer; opaque over copper, transparent
+   elsewhere, cold end lifted so it stays visible on the dark canvas.
+   Reference images never plot to gerbers. Every push **replaces all
+   reference images on those layers**, so don't store unrelated images
+   there. Also available headless:
+   `python tools/kicad_heatmap_overlay.py --net X --amps 10`.
 
 ## Model & limits
 
@@ -76,20 +115,56 @@ SWIG API. Requires KiCad **10.0.1+**.
   capped" checkbox** (default on, `VIAS_CAPPED`) the mouth carries a
   thin copper cap (`CAP_PLATING_UM = 15`, fab spec) on the **outer**
   layers and is an open hole on inner layers; unchecked, mouths are open
-  holes everywhere. Layer-to-layer the cap never matters at DC (it is in
-  parallel with the annular-ring contact, not in series) — the checkbox
+  holes everywhere. The fab caps only small vias: drills above the
+  dialog's **"capped up to drill"** threshold (default
+  `CAP_MAX_DRILL_MM = 0.5`) keep open mouths even with capping
+  selected. Layer-to-layer the cap never matters at DC (it is in
+  parallel with the annular-ring contact, not in series), so the checkbox
   only affects in-plane conduction across outer-layer mouths. Sub-cell
   mouths scale their cells' sheet conductance by the true covered
   fraction (4×4 supersampling), so coarse grids see the correct small
-  perturbation instead of a whole-cell hole. THT-pad copper and drills
-  remain outside the model; at f > 0 the thickness scaling is applied
-  multiplicatively to the skin-corrected sheet conductance
+  perturbation instead of a whole-cell hole. Barrels are gathered in
+  **single-layer runs too** (drill mouths perforate a lone plane).
+  **THT pads are fully modeled**: their exact copper shapes (incl.
+  oblong pads, fetched from KiCad; the outer shape stands in for inner
+  rings) are stamped onto every included layer, and every **populated**
+  pad carries its full **soldered joint** on its SOLDER side (opposite
+  the component; the component-side pad face stays bare). The hole
+  holds the **component lead** (a cylinder of drill −
+  `THT_LEAD_CLEARANCE_MM`, resistivity `THT_LEAD_RHO_OHM_M`, copper by
+  default; raise it for brass/steel leads) **plus solder** in the
+  remaining annulus, both in parallel with the plating. The filled
+  hole also conducts **in-plane on every spanned layer** (component
+  side and inner layers included): the mouth keeps its copper and
+  additionally carries the plug — lead disc plus solder bore — as
+  conduction-equivalent copper of the **full hole depth** (the pin
+  continues beyond both mouths, so each layer sees the whole plug
+  cross-section). The joint is side-symmetric except for the solder:
+  coat and cone on the solder side only. On the raster map these
+  mouths render in a darker tin color. The pad face
+  gets the average-thickness solder coat (exact pad shape) and the
+  protruding-lead cone (see barrel contacts below; on oblong pads the
+  cone tapers to the pad's short dimension). **Slotted (oval) holes**
+  keep their true stadium shape: the barrel wall, drill mouth, contact
+  ring and lead cone all follow the slot (rotated with the pad), and
+  the barrel conducts over the slot's real perimeter/bore area — not a
+  circle of the slot's long dimension. Whether a hole is a via or
+  a THT pad, the owning footprint's side, and its **Do not populate**
+  flag are all read from KiCad. **DNP pads** get an **open hole** and
+  a plating-only barrel, no joint. At f > 0 the thickness scaling is
+  applied multiplicatively to the skin-corrected sheet conductance
   (approximation). Per layer a barrel attaches to
   the fill cell under it, or to the nearest copper cell within the pad
   footprint plus one grid cell — fills joined by **thermal-relief
   spokes** still connect; wider antipads do not, and the barrel bridges
   the layers above/below with the full barrel length. Barrels that reach
   fill on fewer than two layers carry no current and are reported.
+
+  ![How drilled holes are modeled — cross-section](docs/img/hole-model.png)
+  *The four hole types: capped small via, open large via, populated THT
+  pad with its full solder joint (lead ∥ solder ∥ plating in the hole,
+  one-sided pad coat, protruding-lead solder cone), and a DNP THT pad.*
+
 - The net's **traces** (straight and arc tracks, exact outline polygons
   incl. rounded ends) conduct together with the fills — dialog checkbox,
   on by default (`INCLUDE_TRACKS`). Traces narrower than
@@ -98,8 +173,13 @@ SWIG API. Requires KiCad **10.0.1+**.
   series resistance carries no discretization error and no cell-size
   tuning is needed for thin traces. 1D-modeled traces show potential,
   power density, and |J| (the true in-trace density from the link
-  currents, |ΔV|/(ρ·Δl)). Pad copper other than the selected
-  contacts is still **not** part of the conductor model.
+  currents, |ΔV|/(ρ·Δl)). Pad copper is part of the conductor: THT pad
+  shapes are stamped on every included layer (see above), **SMD** pad
+  shapes on their own layer (`INCLUDE_SMD_PADS`) — pads are the
+  junctions where traces and thermal-relief spokes actually meet, so
+  without them a multi-track junction necks down to the accidental
+  overlap of the track ends. Dead-end pads (component terminals) are
+  dropped with the other copper not connected to both contacts.
 - **Solder buildup on mask openings** (dialog checkbox, **off by
   default**; `INCLUDE_MASK_BUILDUP`): zones drawn on `F.Mask`/`B.Mask`
   are treated as mask openings that collect `SOLDER_THICKNESS_UM`
@@ -110,6 +190,29 @@ SWIG API. Requires KiCad **10.0.1+**.
   interface faces use harmonic-mean conductances. Buildup areas render
   tin-gray on the raster map; |J| in them is referenced to the
   conductance-equivalent copper thickness.
+- **Barrel contacts**: a selected **via or through-hole pad** injects at
+  the **drill-wall ring** on every layer the barrel spans — the current
+  physically enters through the lead/wire soldered into the hole, so
+  the spreading resistance across the pad and surrounding pour is part
+  of the result (both contact models; verified against
+  R = ρ/(π·t)·acosh(d/2a) for two circular contacts on a sheet).
+  Slotted holes inject along the stadium-shaped slot wall. A
+  soldered **THT joint** additionally assumes the **hole is filled with
+  solder** (core in parallel with the plating) and the **pad face on
+  the solder side carries an average-thickness solder coat**
+  (`SOLDER_THICKNESS_UM`, 50 µm) over the modeled copper under the pad
+  shape — the solder side is the side opposite the component (taken
+  from the owning footprint; assumed `B.Cu` if it cannot be found),
+  and the component-side pad face stays bare. There the **clipped
+  lead protrudes** `THT_LEAD_PROTRUSION_MM` (1.5 mm, 0 = off)
+  and a **solder cone** wraps
+  it: full protrusion height at the drill wall, tapering linearly to
+  zero at the pad edge, applied as extra conduction-equivalent copper
+  per cell. The tall solder column at the wall pulls the joint
+  vicinity to lead potential — equivalent to extending the barrel wall
+  vertically — while the taper carries the radial spreading. To model
+  a probe pressed onto the pad face instead, draw a marker rectangle
+  over the pad.
 - **Contact models** (dialog / `CONTACT_MODEL`): default **uniform
   injection** — a conductor pressed on top feeds the current orthogonally
   with uniform surface density, so |J| ramps across the contact area
@@ -120,13 +223,19 @@ SWIG API. Requires KiCad **10.0.1+**.
   terminals (e.g. planes joined only through the bolted lugs), only the
   equipotential model is well-defined; the uniform model stops with an
   error instead of prescribing an arbitrary split.
+
+  ![The two contact models bracket a real contact](docs/img/contact-models.png)
+  *|J| around the same 3×3 mm contact under both models: the ideal
+  bonded lug crowds the current at the contact edges (no in-sheet
+  current inside an equipotential region), the pressed conductor ramps
+  it across the contact area.*
 - Fields are reported at the dialog's test current; power scales with I².
 - **Skin effect (f > 0)**: per-layer effective sheet resistance from the
   exact 1D foil-diffusion solution `Zs = τρ·coth(τt)`, `τ = (1+j)/δ`
   (`SKIN_SIDES = 1` in config: plane facing a return plane; `2` =
   isolated foil), and the analogous correction for the 18 µm barrel wall.
   Enter one frequency per run (e.g. a switching harmonic, with its RMS
-  amplitude as the test current) — suffixes `k`/`M` accepted.
+  amplitude as the test current); suffixes `k`/`M` are accepted.
   **Caveat:** only through-thickness crowding is modeled. Lateral
   (proximity-effect) redistribution needs a magneto-quasistatic solver
   and is not captured — since the resistance-driven distribution is the
@@ -138,15 +247,15 @@ SWIG API. Requires KiCad **10.0.1+**.
   not the geometric foil thickness.
 - 5-point FDM per layer on an auto-sized shared grid (~2 M fine cells
   with the uniform grid; ~8 M with the adaptive grid, whose unknown
-  count no longer scales with them). Direct sparse solve up to 500 k
-  unknowns, AMG-preconditioned CG (pyamg) above — Jacobi-CG if pyamg is
-  missing. Discretization error typically ≲ 2 % at defaults — halve the
-  cell size and compare to judge convergence.
+  count no longer scales with the fine-cell count). Direct sparse solve
+  up to 500 k unknowns, AMG-preconditioned CG (pyamg) above (Jacobi-CG
+  if pyamg is missing). Discretization error typically ≲ 2 % at
+  defaults; halve the cell size and compare to judge convergence.
 - **Adaptive cells** (dialog checkbox, **on by default**;
   `ADAPTIVE_CELLS`):
   solves on a 2:1-balanced quadtree — fine cells at copper boundaries,
   electrodes, traces, via mouths and buildup, blocks up to
-  `ADAPTIVE_MAX_CELL_UM` (2 mm) in plane interiors (`ADAPTIVE_GUARD`
+  `ADAPTIVE_MAX_CELL_UM` (1 mm) in plane interiors (`ADAPTIVE_GUARD`
   sets the clearance a block needs to grow). The **minimum element size
   is the grid cell size itself** (auto / dialog / `CELL_UM_OVERRIDE`);
   the uniform limit reproduces the normal grid exactly. Large
@@ -158,25 +267,40 @@ SWIG API. Requires KiCad **10.0.1+**.
   uniform grid ≲ 0.03 %, with the power-balance identity intact. All
   fields are expanded back to the fine grid for the maps and reports.
 
+  ![Rasterized map with the adaptive mesh overlay](docs/img/demo-raster.png)
+  *The raster map of the demo net: quadtree leaves drawn on the copper
+  (fine at boundaries, electrodes, via mouths and pads; coarse blocks
+  in plane interiors), the tin-gray solder coat of the THT-pad contact
+  P1, and the via field with its pad copper.*
+
+**Measured vs. computed**: we tested the plugin on a few real boards
+against a UT3513+ micro-ohm meter; the measured resistances were within
+±20 % of the computed values. We attribute the deviation to
+imperfections of the testing setup (probe placement and probe contact
+resistance vs. the ideal modeled contacts) and to manufacturing
+inaccuracies — actual copper and plating thicknesses routinely deviate
+from nominal. Relative comparisons between layout variants are
+accordingly more trustworthy than absolute numbers.
+
 ## Offline / development
 
 Every run writes `geometry_dump.json`; re-solve without KiCad:
 
 ```powershell
-.venv\Scripts\python.exe -m fill_resistance.standalone dump.json `
+uv run python -m fill_resistance.standalone dump.json `
     [--current 40] [--cell-um 50] [--layers F.Cu,In1.Cu] [--no-show] `
     [--out DIR] [--force-iterative]
 ```
 
-Dev environment, tests, headless extraction (Windows shown; on
-Linux/macOS use `.venv/bin/python`):
+Dev environment, tests, headless extraction — [uv](https://docs.astral.sh/uv/)
+manages the venv from `pyproject.toml`/`uv.lock` (`requirements.txt`
+stays: KiCad builds the plugin's runtime venv from it):
 
 ```powershell
-uv venv --python 3.11 .venv
-uv pip install --python .venv\Scripts\python.exe kicad-python numpy scipy pyamg matplotlib pytest
-.venv\Scripts\python.exe -m pytest tests -q          # incl. exact analytic cases
-.venv\Scripts\python.exe tools\api_probe.py          # IPC API probe vs live KiCad
-.venv\Scripts\python.exe -m fill_resistance.board_io dump.json [NET]  # extract only
+uv sync                                              # one-time env setup
+uv run pytest -q                                     # incl. exact analytic cases
+uv run python tools/api_probe.py                     # IPC API probe vs live KiCad
+uv run python -m fill_resistance.board_io dump.json [NET]  # extract only
 ```
 
 ## Packaging / publishing
@@ -188,7 +312,10 @@ filled in. To publish: upload the zip to a release, set `download_url`
 (and the `homepage` resource in `metadata.json`), then submit the
 registry copy as `packages/th.co.b4l.fill-resistance/metadata.json` in a
 merge request to <https://gitlab.com/kicad/addons/metadata>. Icons are
-regenerated with `python tools/gen_icons.py`.
+regenerated with `python tools/gen_icons.py`; the README figures in
+`docs/img/` with `uv run python tools/gen_readme_figs.py`
+(real solver output on small synthetic boards, plus the hand-drawn
+hole cross-section).
 
 ## License
 
@@ -197,7 +324,9 @@ GPL-3.0-or-later — see [LICENSE](LICENSE).
 ## Troubleshooting
 
 - **No toolbar button**: venv still building (wait), or build failed →
-  *Recreate Plugin Environment*; check the interpreter path (setup 2).
+  *Recreate Plugin Environment* (right-click the plugin's row in
+  Preferences → *PCB Editor → Action Plugins*); check the interpreter
+  path (setup 2).
 - **"Could not connect to KiCad's IPC API"**: API server not enabled, or
   KiCad not running (no headless mode in KiCad 10).
 - **"KiCad is busy"**: a modal dialog is open in KiCad — close it, rerun.
@@ -205,3 +334,23 @@ GPL-3.0-or-later — see [LICENSE](LICENSE).
   best-effort); PNGs are always saved regardless.
 - **Result seems too low/high**: remember the model is fills + barrels
   only, with ideal contacts; measure electrode-to-electrode.
+
+## LLM disclaimer
+
+This plugin was developed with an LLM: Anthropic's **Claude** (Claude
+Code, model Claude Fable 5). The physics model, solver, tests, tooling
+and this documentation (including the figures; all but the hand-drawn
+hole cross-section are generated by the solver itself) were written by
+the model, feature by feature, under human direction and review
+(janik / B4L); most commits carry a `Co-Authored-By: Claude` trailer.
+
+What keeps this honest: the test suite pins the numerics to exact
+analytic references (strip and annulus resistances, the acosh spreading
+resistance of two circular contacts, skin-effect limits, power-balance
+identities) and to convergence/regression checks; run it with
+`uv run pytest`. Real boards were measured against a UT3513+ micro-ohm
+meter (see *Measured vs. computed* above). Nevertheless, an LLM wrote
+this: read *Model & limits*
+critically, treat surprising numbers with the usual engineering
+suspicion, and cross-check against a hand estimate before trusting a
+result with hardware. Bug reports are very welcome.

deploy.ps1

@@ -30,7 +30,8 @@ if ($Mode -eq 'Junction') {
     New-Item -ItemType Junction -Path $dst -Target $src | Out-Null
     Write-Host "junction created: $dst -> $src"
 } else {
-    $exclude = @('.venv', '.git', 'tests', 'tools', '__pycache__', '.pytest_cache')
+    $exclude = @('.venv', '.git', 'tests', 'tools', '__pycache__', '.pytest_cache',
+                 'pyproject.toml', 'uv.lock')
     New-Item -ItemType Directory -Force $dst | Out-Null
     Get-ChildItem $src -Force | Where-Object { $exclude -notcontains $_.Name } |
         ForEach-Object { Copy-Item $_.FullName -Destination $dst -Recurse -Force }

docs/RELEASING.md

@@ -0,0 +1,67 @@
+# Releasing a new version
+
+The PCM addon zip is built by CI (`.gitea/workflows/build-pcm.yml`).
+Every push to `main` builds it as a downloadable artifact; pushing a
+`v<version>` tag additionally creates a Gitea release with the zip
+attached. The release job checks that the tag matches `metadata.json`
+and fails on a mismatch.
+
+## Steps
+
+1. **Bump the version** in `metadata.json` — the single entry in
+   `versions` (plain `MAJOR.MINOR.PATCH`, no `v` prefix; the PCM schema
+   rejects anything else):
+
+   ```json
+   "versions": [
+       {
+           "version": "1.0.2",
+           "status": "stable",
+           "kicad_version": "10.0",
+           "runtime": "ipc"
+       }
+   ]
+   ```
+
+2. **Commit, tag, push** (tag = `v` + the manifest version):
+
+   ```powershell
+   git add metadata.json
+   git commit -m "Release 1.0.2"
+   git tag v1.0.2
+   git push
+   git push origin v1.0.2
+   ```
+
+3. **Verify**: the Actions run for the tag builds
+   `th.co.b4l.fill-resistance_<version>.zip` and publishes it at
+   <https://git.b4l.co.th/B4L/kicad-zone-resistance/releases>, together
+   with `metadata-registry.json`. The zip installs directly via
+   Plugin and Content Manager → *Install from File*.
+
+## Publishing to the official KiCad registry (optional)
+
+The attached `metadata-registry.json` already carries the release
+`download_url`, `download_sha256` and sizes. Submit it as
+`packages/th.co.b4l.fill-resistance/metadata.json` in a merge request
+to <https://gitlab.com/kicad/addons/metadata>. The registry keeps every
+published version: append the new entry to the `versions` array of the
+registry copy instead of replacing the previous one (the repo's own
+`metadata.json` only ever holds the current version —
+`tools/build_package.py` reads `versions[0]`).
+
+## Local build (no CI)
+
+```powershell
+python tools/build_package.py    # writes dist/<identifier>_<version>.zip
+```
+
+Pure stdlib — no venv needed. `dist/` is gitignored.
+
+## CI prerequisites (one-time, server side)
+
+- Actions enabled for the repo (Settings → Actions unit).
+- A runner registered with the `ubuntu-latest` label; the default
+  act_runner image works — the build needs only Python 3.
+- Workflow actions are pinned to commit SHAs; when bumping them, update
+  the SHA and the trailing version comment together.

fill_resistance/adaptive.py

@@ -94,6 +94,7 @@ def run_solve_adaptive(problem: Problem, stack: RasterStack,
 
     # --- leaves per layer -------------------------------------------------
     t0 = time.perf_counter()
+    links, dead_barrels = sv._barrel_links(stack, problem)
     keep = e1 | e2
     if stack.chain is not None:
         keep |= stack.chain
@@ -101,6 +102,13 @@ def run_solve_adaptive(problem: Problem, stack: RasterStack,
         keep |= stack.buildup
     if stack.thick_scale is not None:
         keep |= stack.thick_scale != 1.0
+    # pin every barrel attachment cell fine: a point-like barrel
+    # injection into a coarse leaf makes the whole leaf equipotential
+    # and deletes the local spreading resistance (via fields read up
+    # to ~13% low otherwise); the guard ring then grades around it
+    for _vi, la, ia_, ja_, lb, ib_, jb_, _r in links:
+        keep[la, ia_, ja_] = True
+        keep[lb, ib_, jb_] = True
     mb = _max_block(stack.h_nm)
     grids = [quadtree.build_leaves(stack.masks[li], keep_fine=keep[li],
                                    max_block=mb,
@@ -181,7 +189,6 @@ def run_solve_adaptive(problem: Problem, stack: RasterStack,
             xx.append(np.full(k, -1, dtype=np.int8))
             ee.append(np.full(k, -1, dtype=np.int16))
 
-    links, dead_barrels = sv._barrel_links(stack, problem)
     for vi, la, ia_, ja_, lb, ib_, jb_, r_dc in links:
         na = offs[la] + grids[la].id_grid[ia_, ja_]
         nb = offs[lb] + grids[lb].id_grid[ib_, jb_]

fill_resistance/board_io.py

@@ -6,12 +6,13 @@ KiCad to extract without the dialog (all layers of the net, defaults).
 """
 from __future__ import annotations
 
+import math
 from dataclasses import dataclass, field
 from pathlib import Path
 
 from kipy import KiCad
 from kipy.board import Board
-from kipy.board_types import ArcTrack, BoardRectangle, Pad
+from kipy.board_types import ArcTrack, BoardRectangle, Pad, Via
 from kipy.proto.board.board_pb2 import BoardStackupLayerType
 from kipy.proto.board.board_types_pb2 import ZoneType
 from kipy.util.board_layer import (canonical_name, is_copper_layer,
@@ -22,7 +23,9 @@ import numpy as np
 from . import config
 from .errors import ApiVersionError, CandidateError, SelectionError
 from .geometry import (Electrode, LayerFill, Polygon, Problem, Rect,
-                       SurfaceBuildup, TrackSeg, ViaLink, linearize_ring)
+                       SurfaceBuildup, TrackSeg, ViaLink,
+                       contact_solder_buildups, linearize_ring,
+                       tht_joint_buildups)
 
 MASK_TO_COPPER = {"F.Mask": "F.Cu", "B.Mask": "B.Cu"}
 
@@ -137,11 +140,33 @@ def _convert_poly(poly_with_holes) -> Polygon:
                    holes=[ring(h) for h in poly_with_holes.holes])
 
 
-def _pad_drill_nm(pad_or_via) -> int:
+def _drill_info(pad_or_via) -> tuple[int, int, int]:
+    """(width_nm, slot_dx_nm, slot_dy_nm) of a padstack drill. Round
+    holes: (diameter, 0, 0). Slotted (oblong) holes: width is the
+    NARROW dimension, (slot_dx, slot_dy) the board-frame offset from
+    the drill center to each end-cap center of the slot. The slot
+    follows the pad rotation (KiCad rotates CCW with y down:
+    x' = x cos + y sin, y' = y cos - x sin)."""
     try:
-        return int(pad_or_via.padstack.drill.diameter.x)
+        d = pad_or_via.padstack.drill.diameter
+        dx, dy = int(d.x), int(d.y)
     except Exception:
-        return 0
+        return 0, 0, 0
+    if dx <= 0 or dy <= 0 or dx == dy:
+        return max(dx, 0), 0, 0
+    half = (max(dx, dy) - min(dx, dy)) / 2.0
+    try:
+        th = math.radians(pad_or_via.padstack.angle.degrees)
+    except Exception:
+        th = 0.0
+    ux, uy = (1.0, 0.0) if dx > dy else (0.0, 1.0)
+    return (min(dx, dy),
+            int(round(half * (ux * math.cos(th) + uy * math.sin(th)))),
+            int(round(half * (uy * math.cos(th) - ux * math.sin(th)))))
+
+
+def _pad_drill_nm(pad_or_via) -> int:
+    return _drill_info(pad_or_via)[0]
 
 
 def _pad_default_contact(pad: Pad) -> str:
@@ -179,7 +204,42 @@ def _pad_polygons(board: Board, pad: Pad, contact: str) -> list[Polygon] | None:
     return None
 
 
-def _to_electrode(board: Board, item) -> Electrode:
+def _footprint_pad_map(footprints) -> dict:
+    """(x, y, number) -> owning FootprintInstance. Footprint pads are
+    stored with absolute positions, so the lookup is exact."""
+    out = {}
+    for fp in footprints or []:
+        try:
+            for fpad in fp.definition.pads:
+                out[(fpad.position.x, fpad.position.y, fpad.number)] = fp
+        except Exception:
+            continue
+    return out
+
+
+def _pad_owner(pad: Pad, pad_map: dict):
+    return pad_map.get((pad.position.x, pad.position.y, pad.number))
+
+
+def _tht_protrusion_side(pad: Pad, pad_map: dict, quiet: bool = False) -> str:
+    """Outer layer where the clipped THT lead protrudes (tent + solder
+    cone): the side OPPOSITE the component. Unknown owner -> assume the
+    component sits on F.Cu (lead tents on B.Cu)."""
+    fp = _pad_owner(pad, pad_map)
+    if fp is not None:
+        try:
+            side = canonical_name(fp.layer)
+            return "F.Cu" if side == "B.Cu" else "B.Cu"
+        except Exception:
+            pass
+    if not quiet:
+        print(f"note: no footprint found for pad {pad.number} - assuming "
+              f"its lead protrudes on B.Cu")
+    return "B.Cu"
+
+
+def _to_electrode(board: Board, item, stackup: StackupInfo | None = None,
+                  pad_map: dict | None = None) -> Electrode:
     if isinstance(item, BoardRectangle):
         tl, br = item.top_left, item.bottom_right
         rect = Rect.normalized(tl.x, tl.y, br.x, br.y,
@@ -188,6 +248,22 @@ def _to_electrode(board: Board, item) -> Electrode:
         cy = (rect.y0 + rect.y1) / 2e6
         return Electrode(rect=rect, contact="all",
                          label=f"rect({cx:.1f},{cy:.1f})")
+    if isinstance(item, Via):
+        via: Via = item
+        x, y = via.position.x, via.position.y
+        drill = int(via.drill_diameter or 0) or _pad_drill_nm(via)
+        if drill <= 0:
+            raise SelectionError(
+                f"Selected via at ({x / 1e6:.2f}, {y / 1e6:.2f}) mm has no "
+                f"drill diameter - cannot use it as a contact.")
+        pad_nm = _padstack_pad_nm(via)
+        r = max(pad_nm, drill) // 2
+        rect = Rect.normalized(x - r, y - r, x + r, y + r, "via")
+        return Electrode(
+            rect=rect, contact="all", label=f"via({x / 1e6:.1f},{y / 1e6:.1f})",
+            drill_nm=drill, pad_nm=pad_nm, center=(x, y),
+            barrel_z=(_padstack_span(via.padstack, stackup)
+                      if stackup is not None else None))
     # Pad
     pad: Pad = item
     contact = _pad_default_contact(pad)
@@ -197,37 +273,57 @@ def _to_electrode(board: Board, item) -> Electrode:
     if box is None:
         raise SelectionError(f"Could not get the bounding box of {label}.")
     rect = _box2_to_rect(box, "pad")
+    drill, slot_dx, slot_dy = _drill_info(pad)
     return Electrode(rect=rect, contact=contact,
-                     polygons=_pad_polygons(board, pad, contact), label=label)
+                     polygons=_pad_polygons(board, pad, contact), label=label,
+                     # through-hole pad: current enters at the soldered
+                     # barrel; the joint is solder-filled + pad-coated,
+                     # with a solder cone around the protruding lead
+                     drill_nm=drill, pad_nm=_padstack_pad_nm(pad),
+                     pad_min_nm=_padstack_pad_min_nm(pad),
+                     slot_dx_nm=slot_dx, slot_dy_nm=slot_dy,
+                     center=(pad.position.x, pad.position.y),
+                     solder=drill > 0,
+                     protrusion_side=(_tht_protrusion_side(pad, pad_map or {})
+                                      if drill > 0 else None))
 
 
-def _net_hint_of(pads: list[Pad]) -> str | None:
-    for pad in pads:
-        if pad.net is not None:
-            return pad.net.name
+def _net_hint_of(items: list) -> str | None:
+    for item in items:
+        if item.net is not None:
+            return item.net.name
     return None
 
 
-def get_electrodes(board: Board
+def get_electrodes(board: Board, stackup: StackupInfo | None = None
                    ) -> tuple[list[Electrode], list[Electrode], str | None]:
     """Terminals from the selection. Each terminal may have MULTIPLE parts
     (all merged into one externally-bonded contact):
 
     - rectangles on ELECTRODE_POS_LAYER -> V+ parts, on ELECTRODE_NEG_LAYER
-      -> V- parts; selected pads fill a side that has no rectangles;
+      -> V- parts; selected pads/vias fill a side that has no rectangles;
     - no marker rectangles selected: legacy mode, exactly 2 items
-      (rects/pads, any layer) -> one part each;
+      (rects/pads/vias, any layer) -> one part each;
     - empty selection: board-wide scan of both marker layers.
+
+    Selected vias and through-hole pads become BARREL contacts: current
+    enters at the drill-wall ring (the soldered lead/wire), not the pad
+    face. Draw a marker rectangle over the pad instead to model a probe
+    pressed onto the pad face.
     """
     pos_l = config.ELECTRODE_POS_LAYER
     neg_l = config.ELECTRODE_NEG_LAYER
     scheme = (f"Draw V+ rectangle(s) on {pos_l} and V- rectangle(s) on "
-              f"{neg_l} (axis-aligned), and/or select pads for a side "
+              f"{neg_l} (axis-aligned), and/or select pads/vias for a side "
               f"without rectangles.")
 
     selection = list(board.get_selection())
     rects = [s for s in selection if isinstance(s, BoardRectangle)]
-    pads = [s for s in selection if isinstance(s, Pad)]
+    pads = [s for s in selection if isinstance(s, (Pad, Via))]
+    # protrusion-side lookup needs the owning footprints (THT pads only)
+    pad_map = (_footprint_pad_map(board.get_footprints())
+               if any(isinstance(s, Pad) and _pad_drill_nm(s) > 0
+                      for s in pads) else {})
 
     if not selection:
         allr = [s for s in board.get_shapes() if isinstance(s, BoardRectangle)]
@@ -258,12 +354,13 @@ def get_electrodes(board: Board
         es2 = [_to_electrode(board, r) for r in neg]
         if pads and es1 and es2:
             raise SelectionError(
-                f"Cannot assign the {len(pads)} selected pad(s): both marker "
-                f"layers already provide rectangles. Use pads only for a "
-                f"side that has none."
+                f"Cannot assign the {len(pads)} selected pad(s)/via(s): both "
+                f"marker layers already provide rectangles. Use pads/vias "
+                f"only for a side that has none."
             )
         if pads:
-            pad_parts = [_to_electrode(board, p) for p in pads]
+            pad_parts = [_to_electrode(board, p, stackup, pad_map)
+                         for p in pads]
             if not es1:
                 es1 = pad_parts
             else:
@@ -277,12 +374,13 @@ def get_electrodes(board: Board
 
     items = rects + pads
     if len(items) == 2:
-        return ([_to_electrode(board, items[0])],
-                [_to_electrode(board, items[1])], _net_hint_of(pads))
+        return ([_to_electrode(board, items[0], stackup, pad_map)],
+                [_to_electrode(board, items[1], stackup, pad_map)],
+                _net_hint_of(pads))
     raise SelectionError(
         f"The selection has {len(rects)} rectangle(s) (none on the marker "
-        f"layers) and {len(pads)} pad(s); without marker layers exactly 2 "
-        f"contacts are needed.\n{scheme}"
+        f"layers) and {len(pads)} pad(s)/via(s); without marker layers "
+        f"exactly 2 contacts are needed.\n{scheme}"
     )
 
 
@@ -416,6 +514,18 @@ def _padstack_pad_nm(item) -> int:
         return 0
 
 
+def _padstack_pad_min_nm(item) -> int:
+    """Smallest dimension of the (largest) copper pad of a padstack; 0
+    if unknown. Bounds the lead-cone taper on oblong pads: the cone
+    stays within the inscribed circle."""
+    try:
+        sizes = [min(int(l.size.x), int(l.size.y))
+                 for l in item.padstack.copper_layers]
+        return max(sizes) if sizes else 0
+    except Exception:
+        return 0
+
+
 def _padstack_span(padstack, stackup: StackupInfo) -> tuple[int, int]:
     """(z_top, z_bot) of the barrel; falls back to the full stack."""
     try:
@@ -444,19 +554,154 @@ def gather_barrels(board: Board, net_name: str,
                                z_bot_nm=z_bot, kind="via",
                                pad_nm=_padstack_pad_nm(via)))
     if config.INCLUDE_TH_PADS:
-        for pad in board.get_pads():
-            if pad.net is None or pad.net.name != net_name:
-                continue
-            drill = _pad_drill_nm(pad)
-            if drill <= 0:
-                continue
-            barrels.append(ViaLink(x=pad.position.x, y=pad.position.y,
-                                   drill_nm=drill, z_top_nm=-1,
-                                   z_bot_nm=stackup.z_bot_nm + 1, kind="pad",
-                                   pad_nm=_padstack_pad_nm(pad)))
+        net_pads = [pad for pad in board.get_pads()
+                    if pad.net is not None and pad.net.name == net_name
+                    and _pad_drill_nm(pad) > 0]
+        # populated (non-DNP) THT pads carry a soldered joint: filled
+        # hole + coat + lead cone on the side opposite the component
+        pad_map = (_footprint_pad_map(board.get_footprints())
+                   if net_pads else {})
+        unknown = 0
+        for pad in net_pads:
+            fp = _pad_owner(pad, pad_map)
+            unknown += fp is None
+            populated = True
+            if fp is not None:
+                try:
+                    populated = not fp.attributes.do_not_populate
+                except Exception:
+                    pass
+            drill, slot_dx, slot_dy = _drill_info(pad)
+            barrels.append(ViaLink(
+                x=pad.position.x, y=pad.position.y,
+                drill_nm=drill, z_top_nm=-1,
+                z_bot_nm=stackup.z_bot_nm + 1, kind="pad",
+                pad_nm=_padstack_pad_nm(pad),
+                pad_min_nm=_padstack_pad_min_nm(pad),
+                slot_dx_nm=slot_dx, slot_dy_nm=slot_dy,
+                solder_filled=populated,
+                protrusion_side=(_tht_protrusion_side(pad, pad_map,
+                                                      quiet=True)
+                                 if populated else None)))
+        if unknown:
+            print(f"note: {unknown} THT pad(s) without an identifiable "
+                  f"footprint - assumed populated, leads on B.Cu")
     return barrels
 
 
+def gather_smd_pad_copper(board: Board, net_name: str
+                          ) -> dict[str, list[Polygon]]:
+    """layer name -> exact copper shape(s) of every SMD (undrilled) pad
+    on the net. Pads are junctions: traces and thermal-relief spokes
+    meet ON the pad copper, and without it the junction necks down to
+    the accidental overlap of the track ends - or is severed outright.
+    Dead-end pads (component terminals) become floating islands that
+    the solver's connectivity restriction drops. One API call per pad;
+    pads whose copper layer cannot be determined are skipped."""
+    shapes: dict[str, list[Polygon]] = {}
+    for pad in board.get_pads():
+        if pad.net is None or pad.net.name != net_name \
+                or _pad_drill_nm(pad) > 0:
+            continue
+        layer = _pad_default_contact(pad)      # SMD: its own copper layer
+        if layer == "all":
+            continue
+        polys = _pad_polygons(board, pad, layer)
+        if polys:
+            shapes.setdefault(layer, []).extend(polys)
+    return shapes
+
+
+def gather_tht_pad_copper(board: Board, net_name: str
+                          ) -> dict[tuple[int, int], list[Polygon]]:
+    """(x, y) -> exact copper shape(s) of every drilled (THT) pad on the
+    net. The annular-ring copper conducts on every layer the barrel
+    spans, so build_problem stamps these onto each included layer. One
+    API call per pad; the outer-layer shape stands in for the inner
+    rings (approximation - inner rings are usually the same or
+    smaller)."""
+    shapes: dict[tuple[int, int], list[Polygon]] = {}
+    for pad in board.get_pads():
+        if pad.net is None or pad.net.name != net_name \
+                or _pad_drill_nm(pad) <= 0:
+            continue
+        polys = _pad_polygons(board, pad, "all")
+        if polys:
+            shapes[(pad.position.x, pad.position.y)] = polys
+    return shapes
+
+
+# --- in-KiCad result overlays (EXPERIMENTAL) ---------------------------------
+
+# KiCad sizes reference images as pixels * (1 inch / PPI) * image_scale
+# and assumes 300 PPI for PNGs without a density chunk (BITMAP_BASE)
+OVERLAY_PIX_NM = 25.4e6 / 300
+
+
+def _create_reference_image(board: Board, ref) -> None:
+    """create_items with the per-item status surfaced (kipy <= 0.7.1
+    swallows it and returns an empty wrapper on failure)."""
+    from kipy.proto.common.commands.editor_commands_pb2 import (
+        CreateItems, CreateItemsResponse)
+    from kipy.util import pack_any
+
+    cmd = CreateItems()
+    cmd.header.document.CopyFrom(board._doc)
+    cmd.items.append(pack_any(ref.proto))
+    result = board._kicad.send(cmd, CreateItemsResponse).created_items[0]
+    if result.status.code != 1:                 # 1 = ISC_OK
+        raise RuntimeError(
+            f"KiCad rejected the image (status {result.status.code}) "
+            f"{result.status.error_message or ''} - is the layer enabled "
+            f"in Board Setup? (KiCad >= 10.0.1 required)")
+
+
+def remove_overlays(board: Board, layer) -> int:
+    """Remove every reference image on the given layer; returns count."""
+    ours = [r for r in board.get_reference_images() if r.layer == layer]
+    if ours:
+        board.remove_items(ours)
+    return len(ours)
+
+
+def push_result_overlays(board: Board, stack, result,
+                         lock: bool = False) -> None:
+    """EXPERIMENTAL: the solved |J| of every included copper layer as an
+    unlocked reference image on config.OVERLAY_LAYERS (stackup order,
+    top first; existing images there are replaced). Editor-only -
+    reference images never plot. Per-layer failures are reported and
+    skipped, never fatal to the run."""
+    from kipy.board_types import ReferenceImage
+    from kipy.geometry import Vector2
+
+    from .overlay import heatmap_png
+
+    names = stack.layer_names
+    pairs = list(zip(names, config.OVERLAY_LAYERS))
+    if len(names) > len(config.OVERLAY_LAYERS):
+        print(f"overlays: more copper layers than slots - "
+              f"{', '.join(names[len(config.OVERLAY_LAYERS):])} skipped")
+    ny, nx = stack.shape2d
+    w_nm, h_nm = nx * stack.h_nm, ny * stack.h_nm
+    for src, dest_name in pairs:
+        try:
+            dest = layer_from_canonical_name(dest_name)
+            png = heatmap_png(result.Jmag * 1e-6, names.index(src))
+            remove_overlays(board, dest)
+            ref = ReferenceImage()
+            ref.layer = dest
+            ref.position = Vector2.from_xy(round(stack.x0_nm + w_nm / 2),
+                                           round(stack.y0_nm + h_nm / 2))
+            ref.image_scale = w_nm / (nx * OVERLAY_PIX_NM)
+            ref.image_data = png
+            ref.locked = lock
+            _create_reference_image(board, ref)
+            print(f"overlay: |J| of {src} -> {dest_name} "
+                  f"({len(png) / 1024:.0f} kB)")
+        except Exception as e:
+            print(f"overlay: {src} -> {dest_name} failed: {e}")
+
+
 # --- top level ----------------------------------------------------------------
 
 def build_problem(board: Board, net: str, layer_names: list[str],
@@ -465,7 +710,8 @@ def build_problem(board: Board, net: str, layer_names: list[str],
                   buildups: dict[str, list[Polygon]] | None = None,
                   extra_cu_um: float | None = None,
                   tracks: dict | None = None,
-                  vias_capped: bool | None = None) -> Problem:
+                  vias_capped: bool | None = None,
+                  cap_max_drill_mm: float | None = None) -> Problem:
     per_layer = fills.get(net, {})
     per_layer_tracks = (tracks or {}).get(net, {})
     layers = []
@@ -490,7 +736,32 @@ def build_problem(board: Board, net: str, layer_names: list[str],
             f"Net {net} has no fill on any of the selected layers "
             f"({', '.join(layer_names)})."
         )
-    vias = gather_barrels(board, net, stackup) if len(layers) > 1 else []
+    # barrels matter on a single layer too: via rings + drill mouths
+    # perforate the plane, THT joints locally stiffen it
+    vias = gather_barrels(board, net, stackup)
+    # THT pad copper is part of the conductor: stamp the exact pad
+    # shapes onto every included layer (the barrel spans the stack)
+    pad_shapes = (gather_tht_pad_copper(board, net)
+                  if any(v.kind == "pad" for v in vias) else {})
+    if pad_shapes:
+        extra = [poly for polys in pad_shapes.values() for poly in polys]
+        for layer in layers:
+            layer.polygons = list(layer.polygons) + extra
+        print(f"{len(pad_shapes)} THT pad shape(s) stamped on every "
+              f"included layer")
+    # SMD pad copper too: pads are the junctions where traces/spokes
+    # meet (also gives selected SMD-pad contacts their real copper)
+    smd_shapes = (gather_smd_pad_copper(board, net)
+                  if config.INCLUDE_SMD_PADS else {})
+    if smd_shapes:
+        n = 0
+        for layer in layers:
+            polys = smd_shapes.get(layer.layer_name, [])
+            if polys:
+                layer.polygons = list(layer.polygons) + polys
+                n += len(polys)
+        if n:
+            print(f"{n} SMD pad shape(s) stamped on their layers")
     included = {l.layer_name for l in layers}
     buildup_list = [
         SurfaceBuildup(layer_name=name, polygons=polys)
@@ -502,7 +773,7 @@ def build_problem(board: Board, net: str, layer_names: list[str],
           + (f", solder buildup on "
              f"{', '.join(b.layer_name for b in buildup_list)}"
              if buildup_list else ""))
-    return Problem(
+    problem = Problem(
         board_path=board.name or "",
         net_name=net,
         rho_ohm_m=config.RHO_CU_OHM_M,
@@ -522,7 +793,34 @@ def build_problem(board: Board, net: str, layer_names: list[str],
         vias_capped=(vias_capped if vias_capped is not None
                      else config.VIAS_CAPPED),
         cap_plating_nm=int(config.CAP_PLATING_UM * 1000),
+        cap_max_drill_nm=int((cap_max_drill_mm if cap_max_drill_mm is not None
+                              else config.CAP_MAX_DRILL_MM) * 1e6),
+        tht_protrusion_nm=int(config.THT_LEAD_PROTRUSION_MM * 1e6),
+        tht_lead_clearance_nm=int(config.THT_LEAD_CLEARANCE_MM * 1e6),
+        tht_lead_rho_ohm_m=config.THT_LEAD_RHO_OHM_M,
     )
+    solder_layers = contact_solder_buildups(problem)
+    if solder_layers:
+        sides = sorted({e.protrusion_side
+                        for e in problem.electrodes1 + problem.electrodes2
+                        if e.solder and e.protrusion_side})
+        cone = (f", {config.THT_LEAD_PROTRUSION_MM:g} mm lead + solder cone "
+                f"on {', '.join(sides)}"
+                if sides and problem.tht_protrusion_nm > 0 else "")
+        print(f"THT contact(s): solder-filled hole + "
+              f"{config.SOLDER_THICKNESS_UM:g} um average solder coat on the "
+              f"pad face ({', '.join(solder_layers)}){cone}")
+    tht_joint_buildups(problem, pad_shapes)
+    n_joint = sum(1 for v in problem.vias
+                  if v.kind == "pad" and v.solder_filled)
+    n_dnp = sum(1 for v in problem.vias
+                if v.kind == "pad" and not v.solder_filled)
+    if n_joint or n_dnp:
+        print(f"{n_joint} populated THT pad joint(s): lead + solder in the "
+              f"hole, coat + cone on the solder side"
+              + (f"; {n_dnp} DNP pad(s): open hole, plating-only"
+                 if n_dnp else ""))
+    return problem
 
 
 if __name__ == "__main__":
@@ -533,7 +831,7 @@ if __name__ == "__main__":
     out = Path(sys.argv[1]) if len(sys.argv) > 1 else Path("geometry_dump.json")
     _, board = connect()
     stackup = get_stackup_info(board)
-    es1, es2, net_hint = get_electrodes(board)
+    es1, es2, net_hint = get_electrodes(board, stackup)
     if any_zone_unfilled(board):
         refill(board)
     fills = gather_net_fills(board)

fill_resistance/config.py

@@ -28,7 +28,30 @@ VIAS_CAPPED = True              # filled + capped vias (dialog checkbox):
                                 # Ring/pad copper of vias is modeled either
                                 # way; THT-pad copper/drills are not.
 CAP_PLATING_UM = 15.0           # cap plating thickness (fab spec)
-INCLUDE_TH_PADS = True          # plated through-hole pads stitch layers too
+CAP_MAX_DRILL_MM = 0.5          # fab caps only small vias: drills above this
+                                # stay open even with VIAS_CAPPED
+                                # (dialog-settable)
+INCLUDE_SMD_PADS = True         # the net's SMD pad copper conducts too (exact
+                                # shapes on the pad's layer): pads are the
+                                # junctions where traces/spokes meet, and
+                                # selected pad contacts get their real copper.
+                                # Dead-end pads are dropped as floating islands
+INCLUDE_TH_PADS = True          # plated through-hole pads stitch layers too;
+                                # their holes are modeled solder-filled (a
+                                # soldered component lead), so the solder core
+                                # conducts in parallel with the plating
+THT_LEAD_PROTRUSION_MM = 1.5    # clipped THT lead protrusion on the side
+                                # opposite the component: a solder cone of
+                                # this height at the drill wall (tapering to
+                                # zero at the pad edge) wraps the lead of
+                                # every populated THT pad. 0 = no cones
+THT_LEAD_CLEARANCE_MM = 0.25    # hole diameter minus lead diameter (fab
+                                # rule): a lead cylinder of drill - this
+                                # conducts inside every solder-filled hole
+THT_LEAD_RHO_OHM_M = 1.68e-8    # lead material resistivity: copper leads/
+                                # wires; brass ~6.4e-8, phosphor bronze
+                                # ~1.1e-7, copper-clad steel higher - raise
+                                # this if your components use such leads
 SKIN_SIDES = 1                  # skin-effect field config: 1 = plane facing a
                                 # return plane (conservative), 2 = isolated foil
 
@@ -55,6 +78,20 @@ ELECTRODE_POS_LAYER = "User.1"  # rectangles on this layer mark V+ contact parts
 ELECTRODE_NEG_LAYER = "User.2"  # rectangles on this layer mark V- contact parts
 ALWAYS_REFILL = False           # refill zones even if KiCad says they are filled
 
+# --- In-KiCad result overlays (EXPERIMENTAL) ---
+PUSH_OVERLAYS = False           # after solving, push the per-layer |J|
+                                # heatmaps into the open board as unlocked
+                                # reference images (editor-only, never
+                                # plotted); dialog-toggleable
+OVERLAY_LAYERS = ("User.9", "User.10", "User.11", "User.12")
+                                # copper layers map here in stackup order
+                                # (top first); existing reference images on
+                                # these layers are REPLACED on every push;
+                                # each must be enabled in Board Setup
+OVERLAY_ALPHA = 255             # overlay opacity over copper (0-255);
+                                # translucency washes out over bright
+                                # copper - toggle the User layer instead
+
 # --- Adaptive grid ---
 ADAPTIVE_CELLS = True           # solve on a 2:1-balanced quadtree: fine at
                                 # copper boundaries/electrodes/features,

fill_resistance/dialog.py

@@ -37,7 +37,9 @@ class Selection:
     extra_cu_um: float = 0.0
     include_tracks: bool = True
     vias_capped: bool = True
-    adaptive: bool = False
+    cap_max_drill_mm: float = 0.5
+    adaptive: bool = True
+    push_overlays: bool = False   # EXPERIMENTAL in-KiCad |J| overlays
 
 
 class _Dialog(QDialog):
@@ -73,9 +75,14 @@ class _Dialog(QDialog):
         self.capped_check.setChecked(config.VIAS_CAPPED)
         form.addRow("Vias:", self.capped_check)
 
+        self.cap_drill_edit = QLineEdit(f"{config.CAP_MAX_DRILL_MM:g}")
+        self.cap_drill_edit.setEnabled(config.VIAS_CAPPED)
+        self.capped_check.toggled.connect(self.cap_drill_edit.setEnabled)
+        form.addRow("Capped up to drill [mm]:", self.cap_drill_edit)
+
         self.adaptive_check = QCheckBox(
             "adaptive cells (coarsen plane interiors; faster on large "
-            "boards, corrected to ≲0.1 % of the uniform grid)")
+            "boards, corrected to ≲0.03 % of the uniform grid)")
         self.adaptive_check.setChecked(config.ADAPTIVE_CELLS)
         form.addRow("Grid:", self.adaptive_check)
 
@@ -115,6 +122,14 @@ class _Dialog(QDialog):
         self.extracu_edit.setEnabled(bool(buildup_layers))
         form.addRow("Extra Cu in openings [µm]:", self.extracu_edit)
 
+        first, last = config.OVERLAY_LAYERS[0], config.OVERLAY_LAYERS[-1]
+        self.overlay_check = QCheckBox(
+            f"experimental: push per-layer |J| heatmaps into the board as "
+            f"reference images on {first}..{last} (replaces images there; "
+            f"layers must be enabled in Board Setup)")
+        self.overlay_check.setChecked(config.PUSH_OVERLAYS)
+        form.addRow("Overlays:", self.overlay_check)
+
         buttons = QDialogButtonBox(QDialogButtonBox.Ok | QDialogButtonBox.Cancel)
         buttons.accepted.connect(self._try_accept)
         buttons.rejected.connect(self.reject)
@@ -181,10 +196,13 @@ class _Dialog(QDialog):
             raise ValueError("Check at least one layer.")
 
         def number(edit: QLineEdit, name: str) -> float:
+            text = edit.text().strip()
             try:
-                return float(edit.text().strip().replace(",", "."))
-            except ValueError:
-                raise ValueError(f"{name}: '{edit.text()}' is not a number.")
+                return float(skin.normalize_decimal(text))
+            except ValueError as exc:
+                if "separator" in str(exc):
+                    raise ValueError(f"{name}: {exc}")
+                raise ValueError(f"{name}: '{text}' is not a number.")
 
         current = number(self.current_edit, "Test current")
         if current <= 0:
@@ -205,6 +223,11 @@ class _Dialog(QDialog):
             extra_cu = number(self.extracu_edit, "Extra Cu")
             if extra_cu < 0:
                 raise ValueError("Extra Cu must be ≥ 0 µm.")
+        cap_max_drill = config.CAP_MAX_DRILL_MM
+        if self.capped_check.isChecked():
+            cap_max_drill = number(self.cap_drill_edit, "Capped up to drill")
+            if cap_max_drill <= 0:
+                raise ValueError("Capped-up-to drill must be > 0 mm.")
 
         def contact(box: QComboBox) -> str:
             t = box.currentText()
@@ -222,7 +245,9 @@ class _Dialog(QDialog):
                          extra_cu_um=extra_cu,
                          include_tracks=self.tracks_check.isChecked(),
                          vias_capped=self.capped_check.isChecked(),
-                         adaptive=self.adaptive_check.isChecked())
+                         cap_max_drill_mm=cap_max_drill,
+                         adaptive=self.adaptive_check.isChecked(),
+                         push_overlays=self.overlay_check.isChecked())
 
     def _try_accept(self) -> None:
         try:

fill_resistance/geometry.py

@@ -17,7 +17,7 @@ from pathlib import Path
 
 import numpy as np
 
-JSON_SCHEMA_VERSION = 5
+JSON_SCHEMA_VERSION = 6
 
 
 @dataclass(frozen=True)
@@ -93,16 +93,44 @@ class TrackSeg:
 
 @dataclass
 class Electrode:
-    """One PART of a current-injection terminal: a drawn rectangle or a
-    selected pad. A terminal (V+ or V-) is a LIST of parts, all merged
-    into one equipotential contact (externally bonded). `polygons`
-    (board nm) is the exact copper shape when known (pads); None means
-    the rectangle itself is the shape. `contact` = 'all' or a layer
-    name: which included layers this part touches."""
+    """One PART of a current-injection terminal: a drawn rectangle, a
+    selected pad, or a selected via. A terminal (V+ or V-) is a LIST of
+    parts, all merged into one equipotential contact (externally
+    bonded). `polygons` (board nm) is the exact copper shape when known
+    (pads); None means the rectangle itself is the shape. `contact` =
+    'all' or a layer name: which included layers this part touches.
+
+    drill_nm > 0 marks a BARREL contact (selected via or through-hole
+    pad): the current physically enters through the plated barrel (the
+    lead/wire soldered into the hole), so the contact cells are the
+    copper ring at the drill wall, not the whole pad face. `solder`
+    additionally models a soldered THT joint: the hole is filled with
+    solder and the pad face on the SOLDER side (protrusion_side,
+    opposite the component) carries an average-thickness solder coat
+    (Problem.solder_thickness_nm over `polygons`) plus the
+    protruding-lead cone."""
     rect: Rect                                # bounding box (labels/summary)
     contact: str = "all"
     polygons: list[Polygon] | None = None
     label: str = "rect"
+    drill_nm: int = 0                         # >0: barrel contact (slotted
+                                              # holes: the slot WIDTH)
+    pad_nm: int = 0                           # pad diameter (search bound;
+                                              # largest dimension if oblong)
+    pad_min_nm: int = 0                       # smallest pad dimension (cone
+                                              # taper bound); 0 = pad_nm
+    slot_dx_nm: int = 0                       # slotted (oblong) hole: offset
+    slot_dy_nm: int = 0                       # from `center` to each end-cap
+                                              # center of the slot, board
+                                              # frame; (0, 0) = round drill
+    center: tuple[int, int] | None = None     # drill center; None = rect center
+    barrel_z: tuple[int, int] | None = None   # (z_top, z_bot); None = full stack
+    solder: bool = False                      # soldered THT joint (see above)
+    protrusion_side: str | None = None        # outer layer where the clipped
+                                              # lead protrudes (opposite the
+                                              # component): a solder cone
+                                              # wraps it there, see
+                                              # Problem.tht_protrusion_nm
 
 
 @dataclass
@@ -111,21 +139,55 @@ class ViaLink:
     layers whose z lies within [z_top_nm, z_bot_nm]."""
     x: int
     y: int
-    drill_nm: int
+    drill_nm: int                             # slotted holes: the slot WIDTH
     z_top_nm: int
     z_bot_nm: int
     kind: str = "via"                         # "via" | "pad"
     pad_nm: int = 0                           # pad/annular diameter; 0 = unknown
+                                              # (oblong pads: LARGEST dimension,
+                                              # used as a search bound)
+    pad_min_nm: int = 0                       # smallest pad dimension (bounds
+                                              # the lead-cone taper on oblong
+                                              # pads); 0 = same as pad_nm
+    slot_dx_nm: int = 0                       # slotted (oblong) hole: offset
+    slot_dy_nm: int = 0                       # from (x, y) to each end-cap
+                                              # center of the slot, board
+                                              # frame; (0, 0) = round drill
+    solder_filled: bool = False               # populated THT pad: the hole
+                                              # holds lead + solder (in parallel
+                                              # with the plating); False for
+                                              # vias and DNP footprints
+    protrusion_side: str | None = None        # populated THT pad: outer layer
+                                              # where the clipped lead tents
+                                              # (solder cone), opposite the
+                                              # component side
 
     def spans(self, z_nm: int) -> bool:
         return self.z_top_nm - 1 <= z_nm <= self.z_bot_nm + 1
 
     def barrel_resistance(self, length_nm: int, rho_ohm_m: float,
-                          plating_nm: int) -> float:
+                          plating_nm: int,
+                          solder_rho_ohm_m: float | None = None,
+                          lead_nm: float = 0,
+                          lead_rho_ohm_m: float | None = None) -> float:
         """Barrel segment resistance over length_nm: thin-wall annulus of
-        plating around the drill."""
-        area_m2 = math.pi * (self.drill_nm * 1e-9) * (plating_nm * 1e-9)
-        return rho_ohm_m * (length_nm * 1e-9) / area_m2
+        plating around the drill (slotted holes: thin wall around the
+        stadium-shaped slot). With solder_rho_ohm_m the hole holds a
+        soldered THT joint: the component lead (a cylinder of lead_nm
+        diameter, resistivity lead_rho_ohm_m) and the solder filling the
+        remaining bore conduct in parallel with the plating."""
+        ext = 2.0 * math.hypot(self.slot_dx_nm, self.slot_dy_nm) * 1e-9
+        wall = math.pi * (self.drill_nm * 1e-9) + 2.0 * ext
+        ga = wall * (plating_nm * 1e-9) / rho_ohm_m
+                                              # conductance-area [m^2/ohm-m]
+        if solder_rho_ohm_m is not None:
+            r_core = max(self.drill_nm / 2.0 - plating_nm, 0.0) * 1e-9
+            r_lead = min(lead_nm * 1e-9 / 2.0, r_core)
+            if lead_rho_ohm_m is not None and r_lead > 0:
+                ga += math.pi * r_lead * r_lead / lead_rho_ohm_m
+            ga += (math.pi * r_core * r_core + 2.0 * r_core * ext
+                   - math.pi * r_lead * r_lead) / solder_rho_ohm_m
+        return (length_nm * 1e-9) / ga
 
 
 @dataclass
@@ -147,6 +209,23 @@ class Problem:
     vias_capped: bool = True                  # filled+capped vias: thin cap
     cap_plating_nm: int = 15_000              # over outer-layer mouths;
                                               # False = open mouths
+    cap_max_drill_nm: int = 500_000           # fab caps only small vias:
+                                              # drills above this stay open
+                                              # even with vias_capped
+    tht_protrusion_nm: int = 1_500_000        # clipped THT lead protrusion:
+                                              # a solder cone of this height
+                                              # at the drill wall (tapering
+                                              # to zero at the pad edge)
+                                              # wraps the lead on each solder
+                                              # contact's protrusion_side;
+                                              # 0 disables the cones
+    tht_lead_clearance_nm: int = 250_000      # hole minus lead diameter (fab
+                                              # rule): the lead cylinder of
+                                              # drill - this conducts inside
+                                              # every solder-filled hole
+    tht_lead_rho_ohm_m: float = 1.68e-8       # lead material resistivity
+                                              # (copper; brass ~6.4e-8,
+                                              # copper-clad steel higher)
 
     @property
     def layer_names(self) -> list[str]:
@@ -173,6 +252,101 @@ class Problem:
         return int(x.min()), int(y.min()), int(x.max()), int(y.max())
 
 
+def contact_solder_buildups(problem: Problem) -> list[str]:
+    """Soldered THT-joint contacts: the pad face on the SOLDER side (the
+    protrusion side, opposite the component - the component-side face
+    stays bare) is covered in solder of average thickness
+    solder_thickness_nm. Adds one SurfaceBuildup there for every
+    `solder` electrode's pad shape (the buildup machinery intersects
+    with actual copper at raster time). Returns the affected layer
+    names. Called once when the problem is built."""
+    included = {l.layer_name for l in problem.layers}
+    touched = []
+    for e in problem.electrodes1 + problem.electrodes2:
+        if not e.solder or not e.polygons \
+                or e.protrusion_side not in included:
+            continue
+        problem.buildups.append(
+            SurfaceBuildup(layer_name=e.protrusion_side,
+                           polygons=list(e.polygons)))
+        touched.append(e.protrusion_side)
+    return sorted(set(touched))
+
+
+def slot_distance(xg, yg, dx_nm: int, dy_nm: int):
+    """Distance from points (xg, yg) (numpy-broadcastable, coordinates
+    RELATIVE to the hole center) to a slotted hole's axis - the segment
+    (-dx, -dy)..(+dx, +dy) between the end-cap centers. The slot wall
+    sits at distance width/2. Round drills (dx = dy = 0) reduce to the
+    plain radius, so callers need no special case."""
+    if dx_nm == 0 and dy_nm == 0:
+        return np.hypot(xg, yg)
+    l2 = float(dx_nm) * dx_nm + float(dy_nm) * dy_nm
+    t = np.clip((xg * dx_nm + yg * dy_nm) / l2, -1.0, 1.0)
+    return np.hypot(xg - t * dx_nm, yg - t * dy_nm)
+
+
+def _disc_polygon(x_nm: float, y_nm: float, r_nm: float,
+                  n: int = 32) -> Polygon:
+    th = np.linspace(0.0, 2.0 * math.pi, n, endpoint=False)
+    return Polygon(outline=np.round(np.stack(
+        [x_nm + r_nm * np.cos(th), y_nm + r_nm * np.sin(th)],
+        axis=1)).astype(np.int64))
+
+
+def _capsule_polygon(x_nm: float, y_nm: float, dx_nm: float, dy_nm: float,
+                     r_nm: float, n: int = 16) -> Polygon:
+    """Stadium: two half-circle caps of radius r_nm centered at
+    (x +- dx, y +- dy), joined by straight flanks."""
+    a0 = math.atan2(dy_nm, dx_nm)
+    th = np.linspace(-0.5 * math.pi, 0.5 * math.pi, n) + a0
+    cap1 = np.stack([x_nm + dx_nm + r_nm * np.cos(th),
+                     y_nm + dy_nm + r_nm * np.sin(th)], axis=1)
+    cap2 = np.stack([x_nm - dx_nm + r_nm * np.cos(th + math.pi),
+                     y_nm - dy_nm + r_nm * np.sin(th + math.pi)], axis=1)
+    return Polygon(outline=np.round(np.vstack([cap1, cap2])).astype(np.int64))
+
+
+def tht_joint_buildups(problem: Problem,
+                       shapes: dict | None = None) -> list[str]:
+    """Solder coat of the net's populated STITCHING through-hole pads
+    (ViaLink kind 'pad' with solder_filled), on the pad's SOLDER side
+    (the protrusion side, opposite the component; the component-side
+    face stays bare). `shapes` maps (x, y) to the exact pad polygons
+    (fetched from KiCad); pads without one fall back to a pad-diameter
+    disc. Contact pads are skipped: contact_solder_buildups already
+    coats them with the exact pad shape. Returns the affected layer
+    names."""
+    included = {l.layer_name for l in problem.layers}
+    contacts = {e.center for e in problem.electrodes1 + problem.electrodes2
+                if e.drill_nm > 0 and e.center is not None}
+    touched = []
+    for v in problem.vias:
+        if v.kind != "pad" or not v.solder_filled \
+                or (v.x, v.y) in contacts \
+                or v.protrusion_side not in included:
+            continue
+        polys = (shapes or {}).get((v.x, v.y))
+        if polys is None:
+            if v.pad_nm <= v.drill_nm:
+                continue
+            # oblong pads: never coat past the pad - a capsule along the
+            # slot axis, or the inscribed disc when the axis is unknown
+            w = v.pad_min_nm or v.pad_nm
+            hl = math.hypot(v.slot_dx_nm, v.slot_dy_nm)
+            if hl > 0.0 and v.pad_nm > w:
+                s = (v.pad_nm - w) / 2.0 / hl
+                polys = [_capsule_polygon(v.x, v.y, v.slot_dx_nm * s,
+                                          v.slot_dy_nm * s, w / 2.0)]
+            else:
+                polys = [_disc_polygon(v.x, v.y, w / 2.0)]
+        problem.buildups.append(
+            SurfaceBuildup(layer_name=v.protrusion_side,
+                           polygons=list(polys)))
+        touched.append(v.protrusion_side)
+    return sorted(set(touched))
+
+
 def _arc_params(start, mid, end) -> tuple[float, float, float, float, float] | None:
     """Circle through three points: (cx, cy, r, a0, sweep) with a0 the
     start angle and sweep signed; None if the points are collinear."""
@@ -321,6 +495,15 @@ def _electrode_to_json(e: Electrode) -> dict:
         "label": e.label,
         "polygons": (None if e.polygons is None
                      else [_poly_to_json(poly) for poly in e.polygons]),
+        "drill_nm": e.drill_nm,
+        "pad_nm": e.pad_nm,
+        "pad_min_nm": e.pad_min_nm,
+        "slot_dx_nm": e.slot_dx_nm,
+        "slot_dy_nm": e.slot_dy_nm,
+        "center": (None if e.center is None else list(e.center)),
+        "barrel_z": (None if e.barrel_z is None else list(e.barrel_z)),
+        "solder": e.solder,
+        "protrusion_side": e.protrusion_side,
     }
 
 
@@ -331,6 +514,17 @@ def _electrode_from_json(d: dict) -> Electrode:
         label=d.get("label", "rect"),
         polygons=(None if d.get("polygons") is None
                   else [_poly_from_json(pd) for pd in d["polygons"]]),
+        drill_nm=int(d.get("drill_nm", 0)),
+        pad_nm=int(d.get("pad_nm", 0)),
+        pad_min_nm=int(d.get("pad_min_nm", 0)),
+        slot_dx_nm=int(d.get("slot_dx_nm", 0)),
+        slot_dy_nm=int(d.get("slot_dy_nm", 0)),
+        center=(None if d.get("center") is None
+                else (int(d["center"][0]), int(d["center"][1]))),
+        barrel_z=(None if d.get("barrel_z") is None
+                  else (int(d["barrel_z"][0]), int(d["barrel_z"][1]))),
+        solder=bool(d.get("solder", False)),
+        protrusion_side=d.get("protrusion_side"),
     )
 
 
@@ -369,6 +563,10 @@ def problem_to_json(p: Problem) -> dict:
         "extra_cu_nm": p.extra_cu_nm,
         "vias_capped": p.vias_capped,
         "cap_plating_nm": p.cap_plating_nm,
+        "cap_max_drill_nm": p.cap_max_drill_nm,
+        "tht_protrusion_nm": p.tht_protrusion_nm,
+        "tht_lead_clearance_nm": p.tht_lead_clearance_nm,
+        "tht_lead_rho_ohm_m": p.tht_lead_rho_ohm_m,
     }
 
 
@@ -415,7 +613,14 @@ def problem_from_json(d: dict) -> Problem:
             ViaLink(x=int(vd["x"]), y=int(vd["y"]), drill_nm=int(vd["drill_nm"]),
                     z_top_nm=int(vd["z_top_nm"]), z_bot_nm=int(vd["z_bot_nm"]),
                     kind=vd.get("kind", "via"),
-                    pad_nm=int(vd.get("pad_nm", 0)))
+                    pad_nm=int(vd.get("pad_nm", 0)),
+                    pad_min_nm=int(vd.get("pad_min_nm", 0)),
+                    slot_dx_nm=int(vd.get("slot_dx_nm", 0)),
+                    slot_dy_nm=int(vd.get("slot_dy_nm", 0)),
+                    # older dumps: every THT pad counted as solder-filled
+                    solder_filled=bool(vd.get(
+                        "solder_filled", vd.get("kind", "via") == "pad")),
+                    protrusion_side=vd.get("protrusion_side"))
             for vd in d["vias"]
         ],
         electrodes1=(
@@ -442,6 +647,10 @@ def problem_from_json(d: dict) -> Problem:
         ],
         vias_capped=bool(d.get("vias_capped", True)),
         cap_plating_nm=int(d.get("cap_plating_nm", 15_000)),
+        cap_max_drill_nm=int(d.get("cap_max_drill_nm", 500_000)),
+        tht_protrusion_nm=int(d.get("tht_protrusion_nm", 1_500_000)),
+        tht_lead_clearance_nm=int(d.get("tht_lead_clearance_nm", 250_000)),
+        tht_lead_rho_ohm_m=float(d.get("tht_lead_rho_ohm_m", 1.68e-8)),
     )
 
 

fill_resistance/main.py

@@ -33,7 +33,7 @@ def main() -> None:
         try:
             kicad, board = board_io.connect()
             stackup = board_io.get_stackup_info(board)
-            es1, es2, net_hint = board_io.get_electrodes(board)
+            es1, es2, net_hint = board_io.get_electrodes(board, stackup)
             if board_io.any_zone_unfilled(board) or config.ALWAYS_REFILL:
                 board_io.refill(board)
             fills = board_io.gather_net_fills(board)
@@ -95,15 +95,21 @@ def main() -> None:
                 buildups=(buildups if selection.include_buildup else None),
                 extra_cu_um=selection.extra_cu_um,
                 tracks=(tracks if selection.include_tracks else None),
-                vias_capped=selection.vias_capped)
+                vias_capped=selection.vias_capped,
+                cap_max_drill_mm=selection.cap_max_drill_mm)
             outdir = report.make_output_dir(board_io.board_dir(board))
         except ApiError as e:
             raise UserFacingError(f"KiCad API error: {e}")
 
         report.write_geometry_dump(outdir, problem)
+        overlay_cb = None
+        if selection.push_overlays:
+            def overlay_cb(stack, result):
+                board_io.push_result_overlays(board, stack, result)
         pipeline.run(problem, outdir, show=True, i_test=selection.current_a,
                      freq_hz=selection.freq_hz,
-                     contact_model=selection.contact_model)
+                     contact_model=selection.contact_model,
+                     overlay=overlay_cb)
     except UserFacingError as e:
         _fail(str(e), outdir)
     except Exception:

fill_resistance/overlay.py

@@ -0,0 +1,64 @@
+"""Rendering for the experimental in-KiCad result overlays: a solved
+field (|J|) as an RGBA PNG, one pixel per grid cell, transparent where
+there is no copper. The pushing side (ReferenceImages via the IPC API)
+lives in board_io; this module stays KiCad-free so it is testable
+headless.
+"""
+from __future__ import annotations
+
+import io
+
+import numpy as np
+
+from . import config
+
+# the colormap's near-black bottom must stay distinguishable from
+# KiCad's dark canvas (matplotlib figures sit on a light background
+# instead), so the log scale starts this far up the colormap
+FLOOR = 0.18
+
+
+def heatmap_png(data3: np.ndarray, li: int, alpha: int | None = None,
+                bleed: bool = True) -> bytes:
+    """One layer of a field (e.g. |J|, NaN = no copper) as opaque-over-
+    copper RGBA PNG bytes. Color scale matches the plugin's log figure
+    (global vmax across layers). `bleed` extends the edge color one
+    pixel outward at half opacity: the raster mask covers cells whose
+    CENTER is inside the copper, so without it the overlay stops half a
+    cell short of the outline KiCad draws."""
+    import matplotlib
+    from PIL import Image
+    from scipy import ndimage
+
+    if alpha is None:
+        alpha = config.OVERLAY_ALPHA
+    if not np.isfinite(data3).any():
+        raise ValueError("field is empty - nothing to overlay")
+    vmax = float(np.nanmax(data3))
+    if vmax <= 0:
+        raise ValueError("field is empty - nothing to overlay")
+    vmin = vmax / config.CURRENT_DYNAMIC_RANGE
+    d = np.clip(data3[li], vmin, vmax)
+    if config.LOG_CURRENT_SCALE:
+        u = (np.log(d) - np.log(vmin)) / (np.log(vmax) - np.log(vmin))
+    else:
+        u = d / vmax
+    u = FLOOR + (1.0 - FLOOR) * u
+    cmap = matplotlib.colormaps[config.CMAP_CURRENT]
+    rgba = (cmap(np.nan_to_num(u)) * 255).astype(np.uint8)
+    copper = ~np.isnan(data3[li])
+    rgba[..., 3] = np.where(copper, alpha, 0)
+
+    if bleed and copper.any() and not copper.all():
+        ring = ndimage.binary_dilation(
+            copper, structure=np.ones((3, 3), dtype=bool)) & ~copper
+        iy, ix = ndimage.distance_transform_edt(
+            ~copper, return_distances=False, return_indices=True)
+        rgba[ring, :3] = rgba[iy[ring], ix[ring], :3]
+        rgba[ring, 3] = alpha // 2
+
+    buf = io.BytesIO()
+    # no dpi metadata: KiCad assumes its 300 PPI default, which the
+    # pusher's scale computation relies on
+    Image.fromarray(rgba, "RGBA").save(buf, format="PNG")
+    return buf.getvalue()

fill_resistance/pipeline.py

@@ -12,7 +12,9 @@ from .solver import Result
 
 def run(problem: Problem, outdir: Path | None, show: bool = True,
         i_test: float | None = None, freq_hz: float = 0.0,
-        contact_model: str | None = None) -> Result:
+        contact_model: str | None = None, overlay=None) -> Result:
+    """overlay: optional callback(stack, result) run after the solve
+    (EXPERIMENTAL in-KiCad overlays); its failures are non-fatal."""
     if i_test is None:
         i_test = config.TEST_CURRENT_A
     if i_test <= 0:
@@ -43,6 +45,12 @@ def run(problem: Problem, outdir: Path | None, show: bool = True,
         report.write_summary(outdir, problem, stack, result)
     print(report.result_line(result, problem, stack))
 
+    if overlay is not None:
+        try:
+            overlay(stack, result)
+        except Exception as e:
+            print(f"overlay push failed: {e}")
+
     figs = [
         (plots.fig_raster(stack, e1, e2, problem, result), "1_raster_map"),
         (plots.fig_potential(result, stack, e1, e2, problem), "2_potential"),

fill_resistance/plots.py

@@ -50,7 +50,9 @@ _COPPER = "#c98b4e"
 _E1_COLOR = "#c8385a"
 _E2_COLOR = "#2f6fb0"
 _VIA_COLOR = "#2d6b45"
+_PAD_COLOR = "#5b4a8a"   # THT pad barrels (kind='pad'), violet-ink
 _SOLDER = "#9aa3ad"      # tin-gray: solder buildup areas
+_PLUG = "#6e7885"        # darker tin: solder-filled THT holes (lead + plug)
 _MESH = "#a56c33"        # darker copper: adaptive leaf boundaries
 _INK = "#3a3a3a"
 _GRID_INK = "#b8b4ae"
@@ -187,10 +189,14 @@ def _electrode_labels(ax, stack, e1_l, e2_l):
 
 
 def _via_markers(ax, problem, layer):
-    xs = [v.x * 1e-6 for v in problem.vias if v.spans(layer.z_nm)]
-    ys = [v.y * 1e-6 for v in problem.vias if v.spans(layer.z_nm)]
-    if xs:
-        ax.plot(xs, ys, ".", ms=2.5, color=_VIA_COLOR, alpha=0.7)
+    """One dot per barrel spanning the layer: vias green, THT pad
+    barrels violet (same joint markers, different physics)."""
+    for kind, color in (("via", _VIA_COLOR), ("pad", _PAD_COLOR)):
+        pts = [(v.x * 1e-6, v.y * 1e-6) for v in problem.vias
+               if v.kind == kind and v.spans(layer.z_nm)]
+        if pts:
+            xs, ys = zip(*pts)
+            ax.plot(xs, ys, ".", ms=2.5, color=color, alpha=0.7)
 
 
 def area_tag(sign: str, index: int) -> str:
@@ -219,8 +225,9 @@ def _injection_area_labels(ax, li, layer_name, problem, result):
 
 def fig_raster(stack, e1, e2, problem, result=None):
     cmap = ListedColormap([_BG, _COPPER, _E1_COLOR, _E2_COLOR, _SOLDER,
-                           _MESH])
+                           _MESH, _PLUG])
     has_buildup = stack.buildup is not None and stack.buildup.any()
+    has_plug = stack.plug is not None and stack.plug.any()
     has_mesh = stack.mesh is not None and stack.mesh.any()
 
     def paint(ax, li):
@@ -228,11 +235,13 @@ def fig_raster(stack, e1, e2, problem, result=None):
         codes[stack.masks[li]] = 1
         if has_buildup:
             codes[stack.buildup[li]] = 4
+        if has_plug:
+            codes[stack.plug[li]] = 6
         if has_mesh:
             codes[stack.mesh[li]] = 5
         codes[e1[li]] = 2
         codes[e2[li]] = 3
-        ax.imshow(codes, cmap=cmap, vmin=0, vmax=5, origin="upper",
+        ax.imshow(codes, cmap=cmap, vmin=0, vmax=6, origin="upper",
                   extent=stack.extent_mm(), interpolation="nearest")
         _via_markers(ax, problem, problem.layers[li])
         if result is not None and (result.part_currents1
@@ -243,8 +252,12 @@ def fig_raster(stack, e1, e2, problem, result=None):
             _electrode_labels(ax, stack, e1[li], e2[li])
 
     def finalize(fig, rows):
-        handles = [Patch(fc=_COPPER, label="copper"),
-                   Patch(fc=_VIA_COLOR, label="vias")]
+        kinds = {v.kind for v in problem.vias}
+        handles = [Patch(fc=_COPPER, label="copper")]
+        if "via" in kinds or not kinds:
+            handles.append(Patch(fc=_VIA_COLOR, label="vias"))
+        if "pad" in kinds:
+            handles.append(Patch(fc=_PAD_COLOR, label="THT pad barrels"))
         if has_mesh:
             handles.append(Patch(fc=_MESH,
                                  label="adaptive mesh (coarse leaves)"))
@@ -255,6 +268,9 @@ def fig_raster(stack, e1, e2, problem, result=None):
                       f"({problem.solder_thickness_nm / 1000:.0f} µm"
                       + (f" + {problem.extra_cu_nm / 1000:.0f} µm Cu"
                          if problem.extra_cu_nm else "") + ")"))
+        if has_plug:
+            handles.append(Patch(
+                fc=_PLUG, label="solder-filled THT hole (lead + solder)"))
         if result is not None and (result.part_currents1
                                    or result.part_currents2):
             entries = ([("+", _E1_COLOR, i, amps)
@@ -416,7 +432,7 @@ def fig_power(result, stack, e1, e2, problem):
 def fig_error(message: str):
     fig, ax = plt.subplots(figsize=(9, 4.5), layout="constrained")
     ax.axis("off")
-    ax.set_title("Fill Resistance — ERROR", color="#b02a2a",
+    ax.set_title("Fill Resistance - ERROR", color="#b02a2a",
                  fontsize=14, fontweight="bold", loc="left")
     wrapped = "\n".join(
         textwrap.fill(line, width=90) for line in message.splitlines()

fill_resistance/raster.py

@@ -23,7 +23,7 @@ from scipy import ndimage
 
 from . import config
 from .errors import ElectrodeError, GridSizeError
-from .geometry import Electrode, Problem, Rect
+from .geometry import Electrode, Problem, Rect, slot_distance
 
 # 4-connectivity: matches the in-plane 5-point stencil of the solver
 _STRUCT4 = ndimage.generate_binary_structure(2, 1)
@@ -46,7 +46,16 @@ class RasterStack:
     thick_scale: np.ndarray | None = None   # float (L, ny, nx): per-cell
                                             # copper-thickness factor (via
                                             # mouths: cap-thin or partially
-                                            # drilled cells); None = all 1
+                                            # drilled cells; folded-in cone
+                                            # and plug extras); None = all 1
+    t_extra_nm: np.ndarray | None = None    # float (L, ny, nx): additive
+                                            # conduction-equivalent copper
+                                            # (lead cones + hole plugs),
+                                            # folded into thick_scale at the
+                                            # end of rasterize_stack
+    plug: np.ndarray | None = None      # bool (L, ny, nx): solder-filled THT
+                                        # hole mouths (lead + solder plug,
+                                        # drawn on the raster map)
     mesh: np.ndarray | None = None      # bool (L, ny, nx): adaptive leaf
                                         # boundaries (drawn on the raster map)
 
@@ -232,9 +241,89 @@ def rasterize_stack(problem: Problem, h_nm: float) -> RasterStack:
                     _paint_ring(stack, hole, False, pmask)
                 stack.buildup[li] |= pmask
         stack.buildup &= stack.masks        # solder wets exposed copper only
+
+    _paint_lead_fillets(stack, problem)
+
+    if stack.t_extra_nm is not None:
+        # cones + plugs are ADDITIVE conduction-equivalent copper; fold
+        # them into the multiplicative per-cell scale once (multiplying
+        # per contribution would overstate cells carrying both)
+        if stack.thick_scale is None:
+            stack.thick_scale = np.ones(stack.masks.shape)
+        for li, layer in enumerate(problem.layers):
+            stack.thick_scale[li] *= np.where(
+                stack.masks[li],
+                1.0 + stack.t_extra_nm[li] / layer.thickness_nm, 1.0)
     return stack
 
 
+def _paint_lead_fillets(stack: RasterStack, problem: Problem) -> None:
+    """Protruding THT leads (barrel contacts AND the net's populated
+    stitching through-hole pads): the clipped lead sticks
+    tht_protrusion_nm out of the hole on the side opposite the
+    component, wrapped by a solder cone - full protrusion height at
+    the drill wall, tapering linearly to zero at the pad edge. Modeled
+    as extra conduction-equivalent copper (stack.t_extra_nm, ADDITIVE
+    with the hole plug, folded into thick_scale by rasterize_stack): the
+    tall solder column next to the wall pulls those cells to lead
+    potential (equivalent to extending the barrel wall vertically), the
+    taper carries the radial spreading. At f > 0 the factor multiplies
+    the skin-corrected sheet conductance, like the via mouths
+    (approximation)."""
+    H = problem.tht_protrusion_nm
+    if H <= 0:
+        return
+    ny, nx = stack.shape2d
+    h = stack.h_nm
+    index = {name: li for li, name in enumerate(stack.layer_names)}
+
+    # one cone per joint: contact electrodes first (exact data), then the
+    # net's populated stitching THT pads, skipping the contacts' barrels
+    jobs = []
+    seen = set()
+    for e in problem.electrodes1 + problem.electrodes2:
+        if e.drill_nm <= 0:
+            continue
+        if e.center is not None:
+            x, y = e.center
+        else:
+            x = (e.rect.x0 + e.rect.x1) / 2.0
+            y = (e.rect.y0 + e.rect.y1) / 2.0
+        seen.add((int(x), int(y)))
+        if e.solder and e.protrusion_side:
+            # oblong pads: taper from the (slot) wall to the inscribed
+            # dimension (conservative)
+            jobs.append((x, y, e.drill_nm, e.pad_min_nm or e.pad_nm,
+                         e.protrusion_side, e.slot_dx_nm, e.slot_dy_nm))
+    for v in problem.vias:
+        if v.kind == "pad" and v.solder_filled and v.protrusion_side \
+                and (v.x, v.y) not in seen:
+            jobs.append((v.x, v.y, v.drill_nm, v.pad_min_nm or v.pad_nm,
+                         v.protrusion_side, v.slot_dx_nm, v.slot_dy_nm))
+
+    for x, y, drill_nm, pad_nm, side, sdx, sdy in jobs:
+        li = index.get(side)
+        if li is None or pad_nm <= drill_nm:
+            continue
+        ra, rb = drill_nm / 2.0, pad_nm / 2.0
+        ex, ey = rb + abs(sdx), rb + abs(sdy)
+        j0 = max(0, math.floor((x - ex - stack.x0_nm) / h))
+        j1 = min(nx, math.floor((x + ex - stack.x0_nm) / h) + 1)
+        i0 = max(0, math.floor((y - ey - stack.y0_nm) / h))
+        i1 = min(ny, math.floor((y + ey - stack.y0_nm) / h) + 1)
+        if i0 >= i1 or j0 >= j1:
+            continue
+        xs = stack.x0_nm + (np.arange(j0, j1) + 0.5) * h - x
+        ys = stack.y0_nm + (np.arange(i0, i1) + 0.5) * h - y
+        r = slot_distance(xs[None, :], ys[:, None], sdx, sdy)
+        t_sn = H * np.clip((rb - r) / (rb - ra), 0.0, 1.0)
+        t_eq = t_sn * (problem.rho_ohm_m / problem.solder_rho_ohm_m)
+        if stack.t_extra_nm is None:
+            stack.t_extra_nm = np.zeros(stack.masks.shape)
+        m = stack.masks[li, i0:i1, j0:j1]
+        stack.t_extra_nm[li, i0:i1, j0:j1] += np.where(m, t_eq, 0.0)
+
+
 def _via_span(problem: Problem, via) -> list[int]:
     return [li for li, layer in enumerate(problem.layers)
             if via.spans(layer.z_nm)]
@@ -268,39 +357,79 @@ def _apply_via_mouths(stack: RasterStack, problem: Problem) -> None:
     """Drill-mouth treatment, area-weighted per cell (4x4 supersampling):
     capped vias carry a cap_plating-thin copper cap over the mouth on the
     OUTER layers, uncapped vias (and inner layers either way) get an open
-    hole. Fully swallowed cells leave the mask; partially covered cells
-    keep a thickness-scaled sheet conductance via stack.thick_scale."""
+    hole. The fab caps only small vias: drills above cap_max_drill_nm
+    stay open even with vias_capped. THT pad mouths: populated pads are
+    solder-filled - the mouth keeps its copper and additionally carries
+    the PLUG (the component lead plus the solder filling the bore) as
+    in-plane conduction-equivalent copper of the FULL hole depth on
+    EVERY spanned layer (the pin continues beyond both mouths, so each
+    layer sees the whole plug cross-section); the joint is then
+    side-symmetric except for the solder: the solder-side coat and cone
+    come on top, additively (see _paint_lead_fillets).
+    DNP pad holes are cut open on every layer. Fully swallowed cells
+    leave the mask; partially covered cells keep a thickness-scaled
+    sheet conductance via stack.thick_scale."""
     ny, nx = stack.shape2d
     h = stack.h_nm
     outer = {li for li, n in enumerate(stack.layer_names)
              if n in ("F.Cu", "B.Cu")}
     sub = (np.arange(4) + 0.5) / 4.0
     for via in problem.vias:
-        if via.kind != "via" or via.drill_nm <= 0:
+        if via.drill_nm <= 0:
             continue
+        plugged = via.kind == "pad" and via.solder_filled
         r = via.drill_nm / 2.0
-        j0 = max(0, math.floor((via.x - r - stack.x0_nm) / h))
-        j1 = min(nx, math.floor((via.x + r - stack.x0_nm) / h) + 1)
-        i0 = max(0, math.floor((via.y - r - stack.y0_nm) / h))
-        i1 = min(ny, math.floor((via.y + r - stack.y0_nm) / h) + 1)
+        ex, ey = r + abs(via.slot_dx_nm), r + abs(via.slot_dy_nm)
+        j0 = max(0, math.floor((via.x - ex - stack.x0_nm) / h))
+        j1 = min(nx, math.floor((via.x + ex - stack.x0_nm) / h) + 1)
+        i0 = max(0, math.floor((via.y - ey - stack.y0_nm) / h))
+        i1 = min(ny, math.floor((via.y + ey - stack.y0_nm) / h) + 1)
         if i0 >= i1 or j0 >= j1:
             continue
         xs = stack.x0_nm + (np.arange(j0, j1)[:, None] + sub[None, :]) * h \
             - via.x
         ys = stack.y0_nm + (np.arange(i0, i1)[:, None] + sub[None, :]) * h \
             - via.y
-        cov = ((ys[:, None, :, None] ** 2 + xs[None, :, None, :] ** 2)
-               <= r * r).mean(axis=(2, 3))
+        cov = (slot_distance(xs[None, :, None, :], ys[:, None, :, None],
+                             via.slot_dx_nm, via.slot_dy_nm)
+               <= r).mean(axis=(2, 3))
         if not (cov > 0).any():
             continue                            # mouth far smaller than h
+        span = _via_span(problem, via)
+
+        if plugged:
+            # lead cylinder + solder bore: the pin continues beyond BOTH
+            # mouths (component body / clipped stickout), so every
+            # spanned layer sees the FULL plug depth for lateral
+            # spreading - no per-layer split
+            r_lead = max(via.drill_nm - problem.tht_lead_clearance_nm,
+                         0) / 2.0
+            cov_lead = (np.hypot(xs[None, :, None, :],
+                                 ys[:, None, :, None])
+                        <= r_lead).mean(axis=(2, 3))
+            t_sn = problem.rho_ohm_m / problem.solder_rho_ohm_m
+            t_pb = problem.rho_ohm_m / problem.tht_lead_rho_ohm_m
+            depth = max(float(via.z_bot_nm - via.z_top_nm), 0.0)
+            t_eq = depth * (cov_lead * t_pb + (cov - cov_lead) * t_sn)
+            if stack.t_extra_nm is None:
+                stack.t_extra_nm = np.zeros(stack.masks.shape)
+            if stack.plug is None:
+                stack.plug = np.zeros_like(stack.masks)
+            for li in span:
+                m = stack.masks[li, i0:i1, j0:j1]
+                stack.t_extra_nm[li, i0:i1, j0:j1] += np.where(m, t_eq, 0.0)
+                stack.plug[li, i0:i1, j0:j1] |= m & (cov > 0.5)
+            continue
+
         if stack.thick_scale is None:
             stack.thick_scale = np.ones(stack.masks.shape)
-        for li in _via_span(problem, via):
-            if problem.vias_capped and li in outer:
+        for li in span:
+            if via.kind == "via" and problem.vias_capped and li in outer \
+                    and via.drill_nm <= problem.cap_max_drill_nm:
                 ratio = min(problem.cap_plating_nm
                             / problem.layers[li].thickness_nm, 1.0)
             else:
-                ratio = 0.0
+                ratio = 0.0                 # open hole (also DNP THT holes)
             s = 1.0 - cov * (1.0 - ratio)
             gone = s <= 1e-9
             stack.masks[li, i0:i1, j0:j1] &= ~gone
@@ -399,26 +528,90 @@ def _electrode_cells2d(stack: RasterStack, e: Electrode) -> np.ndarray:
     return _rect_cells(stack, e.rect)
 
 
+def _barrel_ring2d(stack: RasterStack, e: Electrode,
+                   mask2d: np.ndarray) -> np.ndarray:
+    """Contact cells of a barrel electrode on one layer: the copper ring
+    at the drill wall (cell centers within one cell of radius drill/2;
+    slotted holes: within one cell of the stadium-shaped slot wall),
+    where the lead/wire soldered into the hole actually meets the layer.
+    If rasterization or an antipad leaves no copper there, fall back to
+    the nearest copper ring within the pad footprint (+1 cell of slop) -
+    the same search bound as the solver's barrel attachment."""
+    ny, nx = stack.shape2d
+    h = stack.h_nm
+    if e.center is not None:
+        x, y = e.center
+    else:
+        x = (e.rect.x0 + e.rect.x1) / 2.0
+        y = (e.rect.y0 + e.rect.y1) / 2.0
+    r = e.drill_nm / 2.0
+    rw = max(e.pad_nm, e.drill_nm + 300_000) / 2.0 + h
+    ex, ey = rw + abs(e.slot_dx_nm), rw + abs(e.slot_dy_nm)
+    out = np.zeros((ny, nx), dtype=bool)
+    j0 = max(0, math.floor((x - ex - stack.x0_nm) / h))
+    j1 = min(nx, math.floor((x + ex - stack.x0_nm) / h) + 1)
+    i0 = max(0, math.floor((y - ey - stack.y0_nm) / h))
+    i1 = min(ny, math.floor((y + ey - stack.y0_nm) / h) + 1)
+    if i0 >= i1 or j0 >= j1:
+        return out
+    xs = stack.x0_nm + (np.arange(j0, j1) + 0.5) * h - x
+    ys = stack.y0_nm + (np.arange(i0, i1) + 0.5) * h - y
+    d = slot_distance(xs[None, :], ys[:, None], e.slot_dx_nm, e.slot_dy_nm)
+    m = mask2d[i0:i1, j0:j1]
+    ring = m & (np.abs(d - r) <= h)
+    if not ring.any():
+        dc = np.where(m & (d <= rw), d, np.inf)
+        dmin = dc.min()
+        if np.isfinite(dmin):
+            ring = dc <= dmin + h              # e.g. thermal-spoke tips
+    out[i0:i1, j0:j1] = ring
+    return out
+
+
+def _part_mask3d(stack: RasterStack, problem: Problem,
+                 el: Electrode) -> np.ndarray:
+    """(L, ny, nx) contact cells of one electrode part: the barrel-wall
+    ring on every spanned layer for via/THT-pad contacts, else the
+    part's shape ∩ copper on its contact layer(s)."""
+    part = np.zeros_like(stack.masks)
+    if el.drill_nm > 0:
+        for li, name in enumerate(stack.layer_names):
+            if el.contact not in ("all", name):
+                continue
+            if el.barrel_z is not None:
+                z = problem.layers[li].z_nm
+                if not (el.barrel_z[0] - 1 <= z <= el.barrel_z[1] + 1):
+                    continue
+            part[li] = _barrel_ring2d(stack, el, stack.masks[li])
+        return part
+    cells2d = _electrode_cells2d(stack, el)
+    for li, name in enumerate(stack.layer_names):
+        if el.contact in ("all", name):
+            part[li] = cells2d & stack.masks[li]
+    return part
+
+
 def electrode_masks(stack: RasterStack, problem: Problem
                     ) -> tuple[np.ndarray, np.ndarray]:
     """Terminal mask = OR over its parts; part = shape ∩ copper on the
-    part's contact layer(s). contact 'all' = every included layer (bolted
-    lug / through pad); a layer name = that layer only. Every part must
-    individually land on copper (clear feedback). V+/V- must not overlap;
-    touching is checked later, only for the equipotential contact model."""
+    part's contact layer(s), or the barrel-wall ring for via/THT-pad
+    contacts (current enters through the soldered barrel, not the pad
+    face). contact 'all' = every included layer (bolted lug / through
+    pad); a layer name = that layer only. Every part must individually
+    land on copper (clear feedback). V+/V- must not overlap; touching is
+    checked later, only for the equipotential contact model."""
     def build(parts: list[Electrode], which: str) -> np.ndarray:
         e = np.zeros_like(stack.masks)
         for el in parts:
-            cells2d = _electrode_cells2d(stack, el)
-            part = np.zeros_like(stack.masks)
-            for li, name in enumerate(stack.layer_names):
-                if el.contact == "all" or el.contact == name:
-                    part[li] = cells2d & stack.masks[li]
+            part = _part_mask3d(stack, problem, el)
             if not part.any():
+                where = ("near its barrel (drill-wall ring / pad footprint)"
+                         if el.drill_nm > 0 else
+                         "(or is smaller than one grid cell)")
                 raise ElectrodeError(
                     f"A {which} contact part ({el.label}) does not overlap "
                     f"any copper of the selected fill on contact layer(s) "
-                    f"'{el.contact}' (or is smaller than one grid cell)."
+                    f"'{el.contact}' {where}."
                 )
             e |= part
         if not e.any():
@@ -446,11 +639,7 @@ def electrode_partition(stack: RasterStack, problem: Problem
         out = []
         claimed = np.zeros_like(stack.masks)
         for el in parts:
-            cells2d = _electrode_cells2d(stack, el)
-            m = np.zeros_like(stack.masks)
-            for li, name in enumerate(stack.layer_names):
-                if el.contact == "all" or el.contact == name:
-                    m[li] = cells2d & stack.masks[li]
+            m = _part_mask3d(stack, problem, el)
             m &= ~claimed
             claimed |= m
             out.append((el.label, m))

fill_resistance/skin.py

@@ -21,6 +21,7 @@ from __future__ import annotations
 
 import cmath
 import math
+import re
 
 MU0 = 4e-7 * math.pi
 
@@ -58,11 +59,25 @@ def resistance_factor(thickness_m: float, freq_hz: float,
             / (rho_ohm_m / thickness_m))
 
 
+def normalize_decimal(text: str) -> str:
+    """Accept a European decimal comma ('1,5' -> '1.5'); reject
+    thousands-separator commas ('1,500' would silently become 1.5,
+    a 1000x error that propagates unnoticed into the result)."""
+    if "," in text:
+        if "." in text or text.count(",") > 1 \
+                or re.search(r",\d{3}(?=\D|$)", text):
+            raise ValueError(
+                f"ambiguous comma in '{text}': use '.' as the decimal "
+                "separator and no thousands separators")
+        text = text.replace(",", ".")
+    return text
+
+
 def parse_frequency(text: str) -> float:
     """'0', '100k', '1.5M', '142500' -> Hz; empty -> 0 (DC).
-    Raises ValueError on unparseable or negative input (a typo silently
-    becoming DC would mislabel the result)."""
-    t = text.strip().lower().replace(",", ".").removesuffix("hz").strip()
+    Raises ValueError on unparseable, ambiguous or negative input (a
+    typo silently becoming DC would mislabel the result)."""
+    t = normalize_decimal(text.strip().lower()).removesuffix("hz").strip()
     if not t:
         return 0.0
     mult = 1.0

fill_resistance/solver.py

@@ -47,7 +47,7 @@ from scipy.sparse import linalg as sla
 
 from . import config, skin
 from .errors import ConnectivityError, ElectrodeError, SolverError
-from .geometry import Problem
+from .geometry import Problem, slot_distance
 from .raster import RasterStack, electrodes_touch
 
 
@@ -156,12 +156,14 @@ def _barrel_links(stack: RasterStack, problem: Problem
         span = [li for li, layer in enumerate(problem.layers)
                 if via.spans(layer.z_nm)]
         r_nm = max(via.pad_nm, via.drill_nm + 300_000) / 2.0 + h
-        win = int(r_nm // h) + 1
-        i0, i1 = max(0, i - win), min(ny, i + win + 1)
-        j0, j1 = max(0, j - win), min(nx, j + win + 1)
+        win_j = int((r_nm + abs(via.slot_dx_nm)) // h) + 1
+        win_i = int((r_nm + abs(via.slot_dy_nm)) // h) + 1
+        i0, i1 = max(0, i - win_i), min(ny, i + win_i + 1)
+        j0, j1 = max(0, j - win_j), min(nx, j + win_j + 1)
         xs = stack.x0_nm + (np.arange(j0, j1) + 0.5) * h - via.x
         ys = stack.y0_nm + (np.arange(i0, i1) + 0.5) * h - via.y
-        d2 = ys[:, None] ** 2 + xs[None, :] ** 2
+        d2 = slot_distance(xs[None, :], ys[:, None],
+                           via.slot_dx_nm, via.slot_dy_nm) ** 2
         d2 = np.where(d2 <= r_nm * r_nm, d2, np.inf)
         present = []                            # (layer, i, j) per layer
         for li in span:
@@ -178,8 +180,16 @@ def _barrel_links(stack: RasterStack, problem: Problem
             length = problem.layers[lb].z_nm - problem.layers[la].z_nm
             if length <= 0:
                 continue
-            r_dc = via.barrel_resistance(length, problem.rho_ohm_m,
-                                         problem.plating_nm)
+            # populated THT pads carry a soldered component lead: lead
+            # cylinder (drill minus the fab clearance) + solder annulus
+            # in parallel with the plating (DNP pads and vias stay
+            # plating-only)
+            r_dc = via.barrel_resistance(
+                length, problem.rho_ohm_m, problem.plating_nm,
+                solder_rho_ohm_m=(problem.solder_rho_ohm_m
+                                  if via.solder_filled else None),
+                lead_nm=max(via.drill_nm - problem.tht_lead_clearance_nm, 0),
+                lead_rho_ohm_m=problem.tht_lead_rho_ohm_m)
             links.append((vi, la, ia, ja, lb, ib, jb, r_dc))
     return links, dead
 

fill_resistance/standalone.py

@@ -42,6 +42,10 @@ def main(argv=None) -> int:
     ap.add_argument("--uncapped", action="store_true",
                     help="treat vias as uncapped (open drill mouths on "
                          "all layers)")
+    ap.add_argument("--cap-max-drill", type=float, default=None,
+                    metavar="MM",
+                    help="cap only vias with drill <= this [mm]; larger "
+                         "drills stay open (default: from the dump)")
     ap.add_argument("--extra-cu-um", type=float, default=None,
                     help="override the added copper in mask openings [um]")
     ap.add_argument("--force-iterative", action="store_true",
@@ -68,6 +72,8 @@ def main(argv=None) -> int:
         problem.buildups = []
     if args.uncapped:
         problem.vias_capped = False
+    if args.cap_max_drill is not None:
+        problem.cap_max_drill_nm = int(args.cap_max_drill * 1e6)
     if args.extra_cu_um is not None:
         problem.extra_cu_nm = int(args.extra_cu_um * 1000)
     if args.layers:

metadata.json

@@ -2,7 +2,7 @@
     "$schema": "https://go.kicad.org/pcm/schemas/v2",
     "name": "Fill Resistance",
     "description": "DC/AC resistance of copper zone fills and traces between two contacts, single- or multi-layer with via coupling; current and power density maps.",
-    "description_full": "Computes the DC or AC resistance of copper zone fills and traces between two contacts (marker rectangles on User.1/User.2 and/or selected pads), single- or multi-layer: the chosen net's fills and tracks are solved as coupled finite-difference sheets linked by the net's via and through-hole-pad barrels; traces narrower than the grid become exact 1D resistor chains, and an adaptive multi-resolution grid (fine at features, coarse plane interiors, deferred-corrected) keeps large boards fast.\n\nShows per-layer rasterized maps, potential, current density and power density, reports per-via currents (via ampacity) and total dissipation at a selectable test current. At a user-set frequency the exact 1D foil/barrel skin-effect correction is applied (AC results are a rigorous lower bound). PNGs, a text summary and a re-solvable geometry dump are saved per run.\n\nNote: the first load builds the plugin's Python environment (numpy, scipy, pyamg, matplotlib, PySide6) and can take several minutes.",
+    "description_full": "Computes the DC or AC resistance of copper zone fills and traces between two contacts (marker rectangles on User.1/User.2 and/or selected pads/vias), single- or multi-layer: the chosen net's fills and tracks are solved as coupled finite-difference sheets linked by the net's via and through-hole-pad barrels. Selected vias/THT pads inject at the drill-wall barrel, and every populated THT hole carries its full solder joint (component lead, solder fill, one-sided pad coat and protruding-lead cone) with exact pad shapes and do-not-populate flags read from KiCad; traces narrower than the grid become exact 1D resistor chains, and an adaptive multi-resolution grid (fine at features, coarse plane interiors, deferred-corrected) keeps large boards fast.\n\nShows per-layer rasterized maps, potential, current density and power density, reports per-via currents (via ampacity) and total dissipation at a selectable test current. At a user-set frequency the exact 1D foil/barrel skin-effect correction is applied (AC results are a rigorous lower bound). PNGs, a text summary and a re-solvable geometry dump are saved per run.\n\nNote: the first load builds the plugin's Python environment (numpy, scipy, pyamg, matplotlib, PySide6) and can take several minutes.",
     "identifier": "th.co.b4l.fill-resistance",
     "type": "plugin",
     "author": {
@@ -17,7 +17,7 @@
     },
     "versions": [
         {
-            "version": "1.0.1",
+            "version": "1.1.0",
             "status": "stable",
             "kicad_version": "10.0",
             "runtime": "ipc"

pyproject.toml

@@ -0,0 +1,28 @@
+# Development environment only (uv sync / uv run). The KiCad plugin
+# manager builds the runtime venv itself from requirements.txt — keep
+# the dependency list there in sync with [project.dependencies].
+[project]
+name = "fill-resistance"
+version = "1.1.0"
+description = "DC/AC resistance of copper zone fills and traces between two contacts (KiCad 10 plugin)"
+license = "GPL-3.0-or-later"
+requires-python = ">=3.11"
+dependencies = [
+    "kicad-python>=0.7.0",
+    "numpy",
+    "scipy",
+    "pyamg",
+    "matplotlib",
+    "PySide6",
+]
+
+[dependency-groups]
+dev = [
+    "pytest",
+]
+
+[tool.uv]
+package = false
+
+[tool.pytest.ini_options]
+testpaths = ["tests"]

tests/test_adaptive.py

@@ -166,6 +166,32 @@ def test_part_currents_and_ac(monkeypatch):
     assert ada.rs_ratios == ref.rs_ratios
 
 
+def test_stitching_pad_mid_plane_close(monkeypatch):
+    """A solder-filled THT stitching pad mid-pour leaves no keep-fine
+    marker of its own (mouth not cut, thick_scale untouched, no copper
+    boundary nearby): without the barrel-attachment pinning its links
+    landed in a coarse equipotential leaf and the local spreading
+    resistance vanished - R read ~20% low on this exact case."""
+    sq = [(0, 0), (40, 0), (40, 40), (0, 40)]
+
+    def prob():
+        p = make_multilayer(
+            [[(sq, [])], [(sq, [])]],
+            rect1_mm=(0, 15, 2, 25), rect2_mm=(38, 15, 40, 25),
+            contact1="L0", contact2="L1",
+            vias_mm=[(20, 20)], gap_mm=1.6, drill_mm=1.0)
+        v = p.vias[0]
+        v.kind = "pad"
+        v.pad_nm = int(1.8 * NM)
+        v.solder_filled = True
+        return p
+
+    ref = _run(prob(), 0.15, adaptive=False, monkeypatch=monkeypatch)
+    ada = _run(prob(), 0.15, adaptive=True, monkeypatch=monkeypatch)
+    assert ada.n_free < 0.4 * ref.n_free       # pour still coarsens
+    assert ada.R_ohm == pytest.approx(ref.R_ohm, rel=2e-3)
+
+
 def test_auto_cell_size_finer_with_adaptive(monkeypatch):
     """The auto sizer affords a larger fine-cell budget (finer h) when
     the adaptive grid is on."""

tests/test_barrel_contacts.py

@@ -0,0 +1,542 @@
+"""Barrel (via / through-hole pad) contact tests: current enters at the
+drill-wall ring, not the pad face, and soldered THT joints carry a
+solder-filled hole plus an average-thickness solder coat on the pad."""
+import math
+
+import numpy as np
+import pytest
+
+from fill_resistance import raster, solver
+from fill_resistance.geometry import (Electrode, Polygon, ViaLink,
+                                      contact_solder_buildups, load_problem,
+                                      problem_from_json, problem_to_json,
+                                      save_problem, tht_joint_buildups)
+from tests.util import NM, make_multilayer, make_problem, rect_mm, ring_mm
+
+PLATE20 = [(0, 0), (20, 0), (20, 20), (0, 20)]
+
+
+def _barrel(x_mm, y_mm, drill_mm, pad_mm=0.0, solder=False, polygons=None):
+    r = max(pad_mm, drill_mm) / 2
+    return Electrode(
+        rect=rect_mm((x_mm - r, y_mm - r, x_mm + r, y_mm + r)),
+        contact="all", label=f"via({x_mm},{y_mm})",
+        drill_nm=int(drill_mm * NM), pad_nm=int(pad_mm * NM),
+        center=(int(x_mm * NM), int(y_mm * NM)), solder=solder,
+        polygons=polygons)
+
+
+def _disc(x_mm, y_mm, r_mm, n=64) -> Polygon:
+    ang = np.linspace(0, 2 * np.pi, n, endpoint=False)
+    return Polygon(outline=ring_mm(
+        [(x_mm + r_mm * np.cos(a), y_mm + r_mm * np.sin(a)) for a in ang]))
+
+
+def _solve(p, h_mm, model="equipotential"):
+    stack = raster.rasterize_stack(p, h_mm * NM)
+    e1, e2 = raster.electrode_masks(stack, p)
+    return solver.run_solve(p, stack, e1, e2, 1.0, contact_model=model), stack
+
+
+def test_ring_cells_at_drill_wall():
+    """The contact cells of a barrel electrode form a ring at the drill
+    wall (one-cell tolerance), not the pad face."""
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    p.electrodes1 = [_barrel(10, 10, drill_mm=1.0, pad_mm=1.6)]
+    stack = raster.rasterize_stack(p, 0.1 * NM)
+    e1, _ = raster.electrode_masks(stack, p)
+    ii, jj = np.nonzero(e1[0])
+    xs = stack.x0_nm + (jj + 0.5) * stack.h_nm - 10 * NM
+    ys = stack.y0_nm + (ii + 0.5) * stack.h_nm - 10 * NM
+    d = np.hypot(xs, ys)
+    assert len(ii) >= 8
+    assert (np.abs(d - 0.5 * NM) <= stack.h_nm + 1).all()
+    # far fewer cells than the full 1.6 mm pad disc
+    assert len(ii) < 0.5 * math.pi * (0.8 * NM / stack.h_nm) ** 2
+
+
+def test_two_barrel_contacts_match_acosh():
+    """Two equipotential circular contacts of radius a, centers d apart,
+    on a large sheet: R = rho/(pi t) * acosh(d / 2a). The barrel-ring
+    contact must reproduce the analytic spreading resistance."""
+    t_um, rho = 70.0, 1.68e-8
+    plate = [(0, 0), (80, 0), (80, 60), (0, 60)]
+    p = make_problem([(plate, [])], rect1_mm=(0, 0, 1, 1),
+                     rect2_mm=(79, 59, 80, 60), t_um=t_um, rho=rho)
+    p.electrodes1 = [_barrel(30, 30, drill_mm=2.0)]
+    p.electrodes2 = [_barrel(50, 30, drill_mm=2.0)]
+    res, _ = _solve(p, 0.15)
+    r_ref = rho / (math.pi * t_um * 1e-6) * math.acosh(20e-3 / (2 * 1e-3))
+    assert res.R_ohm == pytest.approx(r_ref, rel=0.08)
+
+
+def test_barrel_includes_pad_spreading_resistance():
+    """Injecting at the barrel wall (0.5 mm ring) sees the spreading
+    resistance the whole-pad-face contact (2.4 mm equipotential disc)
+    short-circuits: R_barrel > R_pad_face."""
+    p1 = make_problem([(PLATE20, [])],
+                      rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    p1.electrodes1 = [_barrel(10, 10, drill_mm=1.0, pad_mm=2.4)]
+    r_barrel, _ = _solve(p1, 0.1)
+
+    p2 = make_problem([(PLATE20, [])],
+                      rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    p2.electrodes1 = [Electrode(rect=rect_mm((8.8, 8.8, 11.2, 11.2)),
+                                contact="all", label="pad face",
+                                polygons=[_disc(10, 10, 1.2)])]
+    r_face, _ = _solve(p2, 0.1)
+    assert r_barrel.R_ohm > r_face.R_ohm * 1.05
+
+
+def test_ring_fallback_nearest_copper():
+    """Antipad bigger than the drill: no copper at the wall ring, the
+    contact falls back to the nearest copper ring inside the pad
+    footprint (e.g. thermal-spoke tips / hole edge)."""
+    hole = [(10 + 1.2 * np.cos(a), 10 + 1.2 * np.sin(a))
+            for a in np.linspace(0, 2 * np.pi, 64, endpoint=False)]
+    p = make_problem([(PLATE20, [hole])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    p.electrodes1 = [_barrel(10, 10, drill_mm=0.6, pad_mm=4.0)]
+    res, stack = _solve(p, 0.1)
+    e1, _ = raster.electrode_masks(stack, p)
+    ii, jj = np.nonzero(e1[0])
+    d = np.hypot(stack.x0_nm + (jj + 0.5) * stack.h_nm - 10 * NM,
+                 stack.y0_nm + (ii + 0.5) * stack.h_nm - 10 * NM)
+    assert len(ii) >= 8
+    assert (d >= 1.2 * NM - stack.h_nm).all()
+    assert (d <= 1.2 * NM + 2.5 * stack.h_nm).all()
+    assert np.isfinite(res.R_ohm) and res.R_ohm > 0
+
+
+def test_solder_filled_barrel_resistance():
+    """THT joints: the solder core conducts in parallel with the plating.
+    Exact parallel-area formula, and a sanity ratio for a 1 mm drill."""
+    v = ViaLink(x=0, y=0, drill_nm=1_000_000, z_top_nm=-1, z_bot_nm=1)
+    rho, sn = 1.68e-8, 1.32e-7
+    r_plain = v.barrel_resistance(1_600_000, rho, 18_000)
+    r_fill = v.barrel_resistance(1_600_000, rho, 18_000,
+                                 solder_rho_ohm_m=sn)
+    ga = math.pi * 1e-3 * 18e-6 / rho
+    ga += math.pi * (0.5e-3 - 18e-6) ** 2 / sn
+    assert r_fill == pytest.approx(1.6e-3 / ga, rel=1e-12)
+    assert 1.5 < r_plain / r_fill < 4.0
+
+
+def test_contact_solder_coat():
+    """A soldered THT contact adds an average-thickness solder buildup
+    over the pad face on its SOLDER side only (opposite the component),
+    lowering the spreading resistance vs the bare barrel contact."""
+    def prob():
+        p = make_problem([(PLATE20, [])],
+                         rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+        p.electrodes1 = [_barrel(10, 10, drill_mm=1.0, pad_mm=2.4,
+                                 solder=True, polygons=[_disc(10, 10, 1.2)])]
+        p.electrodes1[0].protrusion_side = "F.Cu"
+        return p
+
+    # solder side not among the included layers -> no coat there
+    q = prob()
+    q.electrodes1[0].protrusion_side = "B.Cu"
+    assert contact_solder_buildups(q) == []
+
+    p = prob()
+    assert contact_solder_buildups(p) == ["F.Cu"]
+    assert len(p.buildups) == 1 and p.buildups[0].layer_name == "F.Cu"
+    r_coat, stack = _solve(p, 0.1)
+    assert stack.buildup is not None and stack.buildup.any()
+
+    r_bare, _ = _solve(prob(), 0.1)          # helper not called: no coat
+    assert r_coat.R_ohm < r_bare.R_ohm
+
+
+def test_lead_fillet_profile():
+    """The protruding-lead solder cone paints thick_scale with the exact
+    per-cell formula: 1 + H*clip((rb-r)/(rb-ra), 0, 1)*(rho_cu/rho_sn)/t
+    on copper of the protrusion side; nothing elsewhere."""
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    p.electrodes1 = [_barrel(10, 10, drill_mm=1.0, pad_mm=2.4, solder=True)]
+    p.electrodes1[0].protrusion_side = "F.Cu"
+    stack = raster.rasterize_stack(p, 0.1 * NM)
+    assert stack.thick_scale is not None
+    ny, nx = stack.shape2d
+    jj, ii = np.meshgrid(np.arange(nx), np.arange(ny))
+    r = np.hypot(stack.x0_nm + (jj + 0.5) * stack.h_nm - 10 * NM,
+                 stack.y0_nm + (ii + 0.5) * stack.h_nm - 10 * NM)
+    ra, rb, H = 0.5 * NM, 1.2 * NM, p.tht_protrusion_nm
+    t_eq = H * np.clip((rb - r) / (rb - ra), 0, 1) \
+        * (p.rho_ohm_m / p.solder_rho_ohm_m)
+    expect = np.where(stack.masks[0],
+                      1.0 + t_eq / p.layers[0].thickness_nm, 1.0)
+    assert np.allclose(stack.thick_scale[0], expect, rtol=1e-12)
+    # 1.5 mm of solder at the wall ~ 191 um copper: factor ~ 3.7 on 70 um
+    assert stack.thick_scale[0].max() > 3.0
+
+    p.electrodes1[0].protrusion_side = None    # e.g. via contact: no cone
+    s2 = raster.rasterize_stack(p, 0.1 * NM)
+    assert s2.thick_scale is None
+
+
+def test_lead_fillet_lowers_resistance(monkeypatch):
+    """The cone shorts the joint vicinity: R(with cone) < R(coat-less
+    bare barrel); the adaptive grid pins the cone cells fine and
+    matches the uniform grid."""
+    def prob(protrude=True):
+        p = make_problem([(PLATE20, [])],
+                         rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+        p.electrodes1 = [_barrel(10, 10, drill_mm=1.0, pad_mm=2.4,
+                                 solder=True)]
+        p.electrodes1[0].protrusion_side = "F.Cu"
+        if not protrude:
+            p.tht_protrusion_nm = 0
+        return p
+
+    r_cone, _ = _solve(prob(), 0.1)
+    r_bare, _ = _solve(prob(protrude=False), 0.1)
+    assert r_cone.R_ohm < r_bare.R_ohm
+
+    from fill_resistance import config
+    monkeypatch.setattr(config, "ADAPTIVE_CELLS", True)
+    r_ada, _ = _solve(prob(), 0.1)
+    assert r_ada.R_ohm == pytest.approx(r_cone.R_ohm, rel=2e-3)
+
+
+def _pad_link(populated=True):
+    return ViaLink(x=10 * NM, y=10 * NM, drill_nm=1_000_000, z_top_nm=-1,
+                   z_bot_nm=1, kind="pad", pad_nm=2_400_000,
+                   solder_filled=populated,
+                   protrusion_side="F.Cu" if populated else None)
+
+
+def test_stitching_pad_joint():
+    """A populated THT pad on the net (not a contact) gets the full
+    joint: solder-side coat, cone, and a conducting (plugged) mouth;
+    a DNP pad gets an open hole and nothing else."""
+    def prob(populated=True):
+        p = make_problem([(PLATE20, [])],
+                         rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+        p.vias = [_pad_link(populated)]
+        return p
+
+    p = prob()
+    assert tht_joint_buildups(p) == ["F.Cu"]
+    assert len(p.buildups) == 1
+    r_joint, stack = _solve(p, 0.1)
+    assert stack.thick_scale is not None and stack.thick_scale.max() > 3.0
+    assert stack.buildup is not None and stack.buildup.any()
+    assert stack.masks[0][stack.cell_of(10 * NM, 10 * NM)]   # plugged mouth
+
+    q = prob(populated=False)
+    assert tht_joint_buildups(q) == []
+    r_bare, s2 = _solve(q, 0.1)
+    assert s2.buildup is None
+    assert not s2.masks[0][s2.cell_of(10 * NM, 10 * NM)]     # DNP: open hole
+    assert r_joint.R_ohm < r_bare.R_ohm
+
+
+def test_cone_not_doubled_at_contact():
+    """A contact THT pad also appears in the net's pad list (ViaLink):
+    the cone and coat must be applied once, not squared/stacked. The
+    hole plug (this synthetic barrel spans z = -1..1, so 2 nm of lead)
+    ADDS to the cone at the mouth instead of multiplying it."""
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    p.electrodes1 = [_barrel(10, 10, drill_mm=1.0, pad_mm=2.4, solder=True,
+                             polygons=[_disc(10, 10, 1.2)])]
+    p.electrodes1[0].protrusion_side = "F.Cu"
+    p.vias = [_pad_link()]
+    assert contact_solder_buildups(p) == ["F.Cu"]
+    assert tht_joint_buildups(p) == []       # contact center is skipped
+    stack = raster.rasterize_stack(p, 0.1 * NM)
+    t_cone = p.tht_protrusion_nm * (p.rho_ohm_m / p.solder_rho_ohm_m)
+    t_plug = 2.0 * (p.rho_ohm_m / p.tht_lead_rho_ohm_m)
+    wall = 1.0 + (t_cone + t_plug) / p.layers[0].thickness_nm
+    assert stack.thick_scale.max() == pytest.approx(wall, rel=1e-12)
+
+
+def test_lead_in_barrel_resistance():
+    """Populated hole: plating || lead cylinder || solder annulus, with
+    the lead clipped to the plating bore."""
+    v = ViaLink(x=0, y=0, drill_nm=1_000_000, z_top_nm=-1, z_bot_nm=1)
+    rho, sn = 1.68e-8, 1.32e-7
+    r_solder = v.barrel_resistance(1_600_000, rho, 18_000,
+                                   solder_rho_ohm_m=sn)
+    r_lead = v.barrel_resistance(1_600_000, rho, 18_000,
+                                 solder_rho_ohm_m=sn,
+                                 lead_nm=750_000, lead_rho_ohm_m=rho)
+    rl, rc = 0.375e-3, 0.5e-3 - 18e-6
+    ga = math.pi * 1e-3 * 18e-6 / rho
+    ga += math.pi * rl ** 2 / rho + math.pi * (rc ** 2 - rl ** 2) / sn
+    assert r_lead == pytest.approx(1.6e-3 / ga, rel=1e-12)
+    assert r_lead < r_solder
+    # a lead wider than the bore is clipped to it
+    r_big = v.barrel_resistance(1_600_000, rho, 18_000,
+                                solder_rho_ohm_m=sn,
+                                lead_nm=2_000_000, lead_rho_ohm_m=rho)
+    ga2 = math.pi * 1e-3 * 18e-6 / rho + math.pi * rc ** 2 / rho
+    assert r_big == pytest.approx(1.6e-3 / ga2, rel=1e-12)
+
+
+def test_oblong_pad_cone_uses_inscribed_dim():
+    """Oblong pads: the cone tapers to the inscribed circle (pad_min),
+    never past it, so the long pad axis is not overstated sideways."""
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    p.vias = [_pad_link()]
+    p.vias[0].pad_min_nm = 1_600_000          # 2.4 mm max, 1.6 mm min
+    stack = raster.rasterize_stack(p, 0.1 * NM)
+    ii, jj = np.nonzero(stack.thick_scale[0] != 1.0)
+    d = np.hypot(stack.x0_nm + (jj + 0.5) * stack.h_nm - 10 * NM,
+                 stack.y0_nm + (ii + 0.5) * stack.h_nm - 10 * NM)
+    assert len(d) and d.max() < 0.8 * NM
+
+
+def test_stitching_coat_exact_shape():
+    """When KiCad supplies the exact pad polygon, the coat uses it
+    instead of the pad-diameter disc (oblong pads stay honest)."""
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    p.vias = [_pad_link()]
+    shape = _disc(10, 10, 0.9)
+    assert tht_joint_buildups(p, {(10 * NM, 10 * NM): [shape]}) == ["F.Cu"]
+    assert p.buildups[0].polygons[0] is shape
+
+
+def test_vialink_solder_json():
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    p.vias = [_pad_link()]
+    d = problem_to_json(p)
+    q = problem_from_json(d)
+    assert q.vias[0].solder_filled is True
+    assert q.vias[0].protrusion_side == "F.Cu"
+    # legacy dumps without the flag: THT pads counted as solder-filled,
+    # vias as plating-only
+    del d["vias"][0]["solder_filled"], d["vias"][0]["protrusion_side"]
+    q = problem_from_json(d)
+    assert q.vias[0].solder_filled is True
+    assert q.vias[0].protrusion_side is None
+    d["vias"][0]["kind"] = "via"
+    assert problem_from_json(d).vias[0].solder_filled is False
+
+
+# --- slotted (oblong) holes --------------------------------------------------
+# The lead/barrel of a slotted hole is a stadium, not a circle: modeling
+# it as a circle of the slot's LONG dimension painted contact rings,
+# mouths and cones bigger than the oblong pad itself.
+
+def _slot_dist_mm(stack, ii, jj, x_mm, y_mm, dx_nm):
+    """Distance of cells (ii, jj) to a slot axis (+-dx_nm along x)."""
+    xs = stack.x0_nm + (jj + 0.5) * stack.h_nm - x_mm * NM
+    ys = stack.y0_nm + (ii + 0.5) * stack.h_nm - y_mm * NM
+    t = np.clip(xs / dx_nm, -1.0, 1.0)
+    return np.hypot(xs - t * dx_nm, ys), xs, ys
+
+
+def test_slot_ring_hugs_slot_wall():
+    """The contact ring of a slotted THT pad follows the stadium-shaped
+    slot wall: it reaches around the end caps but never pokes past the
+    oblong pad's short side (the old circular model of the slot's long
+    dimension put cells at radius 1.5 mm straight above/below)."""
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    e = _barrel(10, 10, drill_mm=1.0, pad_mm=3.6)   # slot 3.0 x 1.0 mm
+    e.pad_min_nm = int(1.6 * NM)                    # pad 3.6 x 1.6 mm
+    e.slot_dx_nm = 1 * NM
+    p.electrodes1 = [e]
+    stack = raster.rasterize_stack(p, 0.1 * NM)
+    e1, _ = raster.electrode_masks(stack, p)
+    ii, jj = np.nonzero(e1[0])
+    d, xs, ys = _slot_dist_mm(stack, ii, jj, 10, 10, 1 * NM)
+    assert len(ii) >= 16
+    assert (np.abs(d - 0.5 * NM) <= stack.h_nm + 1).all()
+    assert xs.max() > 1.2 * NM and xs.min() < -1.2 * NM   # rings the caps
+    assert np.abs(ys).max() < 0.8 * NM      # stays inside the 1.6 mm side
+
+
+def test_slot_mouth_is_stadium():
+    """A DNP slotted pad cuts a stadium-shaped hole: open along the whole
+    slot, copper kept just past the slot width and the end caps."""
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    v = _pad_link(populated=False)
+    v.slot_dx_nm = 1 * NM                   # slot 3.0 x 1.0 mm along x
+    p.vias = [v]
+    stack = raster.rasterize_stack(p, 0.1 * NM)
+    m = stack.masks[0]
+    assert not m[stack.cell_of(10 * NM, 10 * NM)]
+    assert not m[stack.cell_of(int(10.9 * NM), 10 * NM)]  # slot end: open
+    assert not m[stack.cell_of(int(9.1 * NM), 10 * NM)]
+    assert m[stack.cell_of(10 * NM, int(10.8 * NM))]  # past the width: copper
+    assert m[stack.cell_of(10 * NM, int(9.2 * NM))]
+    assert m[stack.cell_of(int(11.8 * NM), 10 * NM)]  # past the cap: copper
+
+
+def test_slot_cone_follows_slot():
+    """The lead cone of a slotted oblong pad tapers from the slot WALL
+    to the pad's short dimension. The old circular-drill model (diameter
+    = the slot's long dimension) skipped the cone entirely
+    (pad_min <= drill) and, for the mouth, ate the pad's short side."""
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    e = _barrel(10, 10, drill_mm=1.0, pad_mm=3.6, solder=True)
+    e.pad_min_nm = int(1.6 * NM)
+    e.slot_dx_nm = 1 * NM
+    e.protrusion_side = "F.Cu"
+    p.electrodes1 = [e]
+    stack = raster.rasterize_stack(p, 0.1 * NM)
+    assert stack.thick_scale is not None
+    ny, nx = stack.shape2d
+    jj, ii = np.meshgrid(np.arange(nx), np.arange(ny))
+    r, _, _ = _slot_dist_mm(stack, ii, jj, 10, 10, 1 * NM)
+    ra, rb, H = 0.5 * NM, 0.8 * NM, p.tht_protrusion_nm
+    t_eq = H * np.clip((rb - r) / (rb - ra), 0, 1) \
+        * (p.rho_ohm_m / p.solder_rho_ohm_m)
+    expect = np.where(stack.masks[0],
+                      1.0 + t_eq / p.layers[0].thickness_nm, 1.0)
+    assert np.allclose(stack.thick_scale[0], expect, rtol=1e-12)
+    assert stack.thick_scale[0].max() > 3.0
+
+
+def test_plug_conducts_on_component_side():
+    """A populated THT pad's filled hole (lead + solder plug) conducts
+    IN-PLANE across the mouth on EVERY spanned layer - the component
+    side is not bare foil. Each layer carries the FULL hole depth (the
+    pin continues beyond both mouths, so the whole plug cross-section
+    spreads current at every layer; side-to-side the only difference
+    is the solder coat + cone), converted to conduction-equivalent
+    copper: lead disc at lead resistivity, solder bore around it."""
+    p = make_multilayer([[(PLATE20, [])], [(PLATE20, [])]],
+                        rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    p.vias = [ViaLink(x=10 * NM, y=10 * NM, drill_nm=1_000_000, z_top_nm=-1,
+                      z_bot_nm=1 * NM + 1, kind="pad", pad_nm=2_400_000,
+                      solder_filled=True, protrusion_side="L0")]
+    stack = raster.rasterize_stack(p, 0.1 * NM)
+    c = stack.cell_of(10 * NM, 10 * NM)
+    assert stack.masks[0][c] and stack.masks[1][c]     # plugged, not open
+    assert stack.plug[0][c] and stack.plug[1][c]       # drawn on both sides
+    t = p.layers[0].thickness_nm
+    # full hole depth z = -1 .. 1 mm + 1 on both layers; the mouth
+    # center lies inside the 0.75 mm lead (copper resistivity)
+    depth = 1 * NM + 2.0
+    t_cone = p.tht_protrusion_nm * (p.rho_ohm_m / p.solder_rho_ohm_m)
+    assert stack.thick_scale[0][c] == pytest.approx(
+        1.0 + (t_cone + depth * (p.rho_ohm_m / p.tht_lead_rho_ohm_m)) / t,
+        rel=1e-9)                                      # solder side: + cone
+    assert stack.thick_scale[1][c] == pytest.approx(
+        1.0 + depth * (p.rho_ohm_m / p.tht_lead_rho_ohm_m) / t, rel=1e-9)
+    # far from the joint: untouched foil
+    assert stack.thick_scale[1][stack.cell_of(14 * NM, 10 * NM)] == 1.0
+
+    # clearance swallowing the bore -> no lead, solder-only plug
+    p.tht_lead_clearance_nm = 1_000_000
+    s_sn = raster.rasterize_stack(p, 0.1 * NM)
+    assert s_sn.thick_scale[1][c] == pytest.approx(
+        1.0 + depth * (p.rho_ohm_m / p.solder_rho_ohm_m) / t, rel=1e-9)
+    p.tht_lead_clearance_nm = 250_000
+
+    # a DNP pad still cuts an open hole and gets no plug
+    p.vias[0].solder_filled = False
+    p.vias[0].protrusion_side = None
+    s2 = raster.rasterize_stack(p, 0.1 * NM)
+    assert not s2.masks[0][c] and not s2.masks[1][c]
+    assert s2.plug is None
+
+
+def test_slot_barrel_resistance():
+    """Slotted barrel: plating wall = stadium perimeter, solder core =
+    stadium bore area (both reduce to the circle for dx = dy = 0)."""
+    v = ViaLink(x=0, y=0, drill_nm=1_000_000, z_top_nm=-1, z_bot_nm=1,
+                slot_dx_nm=800_000, slot_dy_nm=600_000)   # ext = 2 mm
+    rho, sn = 1.68e-8, 1.32e-7
+    ga = (math.pi * 1e-3 + 2 * 2e-3) * 18e-6 / rho
+    r_plain = v.barrel_resistance(1_600_000, rho, 18_000)
+    assert r_plain == pytest.approx(1.6e-3 / ga, rel=1e-12)
+    rc = 0.5e-3 - 18e-6
+    ga += (math.pi * rc * rc + 2 * rc * 2e-3) / sn
+    r_fill = v.barrel_resistance(1_600_000, rho, 18_000, solder_rho_ohm_m=sn)
+    assert r_fill == pytest.approx(1.6e-3 / ga, rel=1e-12)
+
+
+def test_slot_coat_fallback_within_pad():
+    """Without an exact pad shape the stitching coat falls back to a
+    capsule along the slot (width = pad_min), not the old pad_nm disc
+    that stuck out past an oblong pad's short side."""
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    v = _pad_link()                          # pad_nm = 2.4 mm
+    v.pad_min_nm = 1_600_000
+    v.slot_dx_nm = 1 * NM
+    p.vias = [v]
+    assert tht_joint_buildups(p) == ["F.Cu"]
+    pts = p.buildups[0].polygons[0].outline.astype(float)
+    xs, ys = pts[:, 0] - 10 * NM, pts[:, 1] - 10 * NM
+    t = np.clip(xs / (0.4 * NM), -1.0, 1.0)  # caps at +-(2.4-1.6)/2 mm
+    d = np.hypot(xs - t * 0.4 * NM, ys)
+    assert np.allclose(d, 0.8 * NM, atol=2)
+    assert np.abs(xs).max() <= 1.2 * NM + 2  # never past pad_nm / 2
+    assert np.abs(ys).max() <= 0.8 * NM + 2  # never past pad_min / 2
+
+
+def test_slot_json_roundtrip():
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    e = _barrel(10, 10, drill_mm=1.0, pad_mm=3.6)
+    e.slot_dx_nm, e.slot_dy_nm = 700_000, -700_000
+    p.electrodes1 = [e]
+    v = _pad_link()
+    v.slot_dx_nm = 1 * NM
+    p.vias = [v]
+    q = problem_from_json(problem_to_json(p))
+    assert (q.electrodes1[0].slot_dx_nm, q.electrodes1[0].slot_dy_nm) \
+        == (700_000, -700_000)
+    assert (q.vias[0].slot_dx_nm, q.vias[0].slot_dy_nm) == (1 * NM, 0)
+    # legacy dumps: round drills
+    d = problem_to_json(p)
+    for vd in d["vias"]:
+        del vd["slot_dx_nm"], vd["slot_dy_nm"]
+    assert problem_from_json(d).vias[0].slot_dx_nm == 0
+
+
+def test_drill_info_slot_rotation():
+    """_drill_info: slot axis from the drill x/y sizes, rotated with the
+    pad (KiCad angles are CCW with y down: 90 deg sends +x to -y)."""
+    from types import SimpleNamespace as NS
+
+    from fill_resistance.board_io import _drill_info
+
+    def pad(dx_mm, dy_mm, angle_deg):
+        return NS(padstack=NS(
+            drill=NS(diameter=NS(x=int(dx_mm * NM), y=int(dy_mm * NM))),
+            angle=NS(degrees=angle_deg)))
+
+    assert _drill_info(pad(1.0, 1.0, 0.0)) == (1 * NM, 0, 0)    # round
+    assert _drill_info(pad(3.0, 1.0, 0.0)) == (1 * NM, 1 * NM, 0)
+    assert _drill_info(pad(1.0, 3.0, 0.0)) == (1 * NM, 0, 1 * NM)
+    w, dx, dy = _drill_info(pad(3.0, 1.0, 90.0))
+    assert (w, dx, dy) == (1 * NM, 0, -1 * NM)
+    w, dx, dy = _drill_info(pad(3.0, 1.0, 45.0))
+    assert w == 1 * NM
+    assert dx == pytest.approx(1 * NM / math.sqrt(2), abs=2)
+    assert dy == pytest.approx(-1 * NM / math.sqrt(2), abs=2)
+
+
+def test_barrel_electrode_json_roundtrip(tmp_path):
+    p = make_problem([(PLATE20, [])],
+                     rect1_mm=(0, 0, 1, 20), rect2_mm=(19, 0, 20, 20))
+    p.electrodes1 = [_barrel(10, 10, drill_mm=0.6, pad_mm=1.2, solder=True,
+                             polygons=[_disc(10, 10, 0.6)])]
+    p.electrodes1[0].barrel_z = (-1, 1_600_001)
+    p.electrodes1[0].protrusion_side = "B.Cu"
+    p.tht_protrusion_nm = 1_200_000
+    f = tmp_path / "d.json"
+    save_problem(p, f)
+    q = load_problem(f)
+    e = q.electrodes1[0]
+    assert e.drill_nm == 600_000 and e.pad_nm == 1_200_000
+    assert e.center == (10 * NM, 10 * NM)
+    assert e.barrel_z == (-1, 1_600_001)
+    assert e.solder is True and len(e.polygons) == 1
+    assert e.protrusion_side == "B.Cu"
+    assert q.tht_protrusion_nm == 1_200_000

tests/test_capping.py

@@ -9,10 +9,13 @@ from tests.util import NM, make_multilayer
 
 
 def _two_layer(width_mm=5.0, drill_mm=0.3, pad_mm=0.6, kind="via",
-               capped=True, cap_um=15.0, hole_mm=None):
+               capped=True, cap_um=15.0, hole_mm=None,
+               cap_max_drill_mm=10.0):
     """10 x width strip on both (outer-named) layers, e1 left on F.Cu,
     e2 right on B.Cu, one via mid-strip. Optionally a circular hole in
-    the F.Cu fill around the via (ring-bridging scenario)."""
+    the F.Cu fill around the via (ring-bridging scenario). The cap-drill
+    threshold defaults to 10 mm here (= every drill capped) so the tests
+    exercise the mouth treatment itself; the threshold has its own test."""
     y = width_mm / 2
     strip = [(0, 0), (10, 0), (10, width_mm), (0, width_mm)]
     holes = []
@@ -31,6 +34,7 @@ def _two_layer(width_mm=5.0, drill_mm=0.3, pad_mm=0.6, kind="via",
     p.vias[0].pad_nm = int(pad_mm * NM)
     p.vias_capped = capped
     p.cap_plating_nm = int(cap_um * 1000)
+    p.cap_max_drill_nm = int(cap_max_drill_mm * NM)
     return p
 
 
@@ -46,7 +50,13 @@ def test_cap_at_foil_thickness_is_identity():
     so the result equals the feature-off reference (a 'pad'-kind barrel,
     which skips rings and mouths) with the mouth fully inside copper."""
     r_cap, _ = _solve(_two_layer(capped=True, cap_um=70.0), 0.1)
-    r_ref, _ = _solve(_two_layer(kind="pad"), 0.1)
+    ref = _two_layer(kind="pad")
+    # populated pads skip rings and mouths; kill the lead + solder core
+    # so the reference barrel matches the via's plating-only resistance
+    ref.vias[0].solder_filled = True
+    ref.solder_rho_ohm_m = 1e30
+    ref.tht_lead_clearance_nm = 10 ** 9
+    r_ref, _ = _solve(ref, 0.1)
     assert r_cap.R_ohm == pytest.approx(r_ref.R_ohm, rel=1e-9)
 
 
@@ -90,14 +100,31 @@ def test_subcell_mouth_perturbs_gently():
     assert r_solid.R_ohm <= r_open.R_ohm <= 1.05 * r_solid.R_ohm
 
 
+def test_cap_drill_threshold():
+    """Drills above cap_max_drill_nm stay open even with capping on: a
+    2 mm drill over a 0.5 mm threshold behaves exactly like uncapped,
+    while a threshold above the drill restores the cap."""
+    kw = dict(drill_mm=2.0, pad_mm=2.6)
+    r_big, s_big = _solve(_two_layer(capped=True, cap_max_drill_mm=0.5,
+                                     **kw), 0.25)
+    r_open, s_open = _solve(_two_layer(capped=False, **kw), 0.25)
+    assert r_big.R_ohm == pytest.approx(r_open.R_ohm, rel=1e-12)
+    assert int(s_big.masks.sum()) == int(s_open.masks.sum())
+
+    r_cap, _ = _solve(_two_layer(capped=True, cap_max_drill_mm=2.1,
+                                 **kw), 0.25)
+    assert r_cap.R_ohm < r_open.R_ohm
+
+
 def test_capping_json_roundtrip(tmp_path):
     from fill_resistance.geometry import load_problem, save_problem
-    p = _two_layer(capped=False, cap_um=12.0)
+    p = _two_layer(capped=False, cap_um=12.0, cap_max_drill_mm=0.8)
     f = tmp_path / "d.json"
     save_problem(p, f)
     q = load_problem(f)
     assert q.vias_capped is False
     assert q.cap_plating_nm == 12_000
+    assert q.cap_max_drill_nm == 800_000
     r_p, _ = _solve(p, 0.25)
     r_q, _ = _solve(q, 0.25)
     assert r_q.R_ohm == pytest.approx(r_p.R_ohm, rel=1e-12)

tests/test_overlay.py

@@ -0,0 +1,60 @@
+"""In-KiCad overlay rendering (fill_resistance.overlay): copper-shaped
+RGBA heatmaps with a visibility floor and a soft edge bleed. The kipy
+pushing side is exercised only against a live KiCad (tools/)."""
+import io
+
+import numpy as np
+import pytest
+from PIL import Image
+
+from fill_resistance import config, overlay
+
+
+def _field(ny=20, nx=30):
+    """Two-layer |J| field: copper disc on layer 0, NaN elsewhere."""
+    data = np.full((2, ny, nx), np.nan)
+    yy, xx = np.mgrid[:ny, :nx]
+    disc = (yy - ny / 2) ** 2 + (xx - nx / 2) ** 2 <= 8 ** 2
+    data[0][disc] = 1.0 + xx[disc]              # spans the log range
+    data[1][disc] = 1e-12                       # below the global log floor
+    return data, disc
+
+
+def test_heatmap_png_shape_and_alpha():
+    data, disc = _field()
+    img = Image.open(io.BytesIO(overlay.heatmap_png(data, 0, bleed=False)))
+    assert img.size == (30, 20)
+    rgba = np.asarray(img)
+    assert (rgba[..., 3][disc] == config.OVERLAY_ALPHA).all()
+    assert (rgba[..., 3][~disc] == 0).all()
+
+
+def test_heatmap_floor_not_black():
+    """The coldest copper must stay distinguishable from a dark canvas:
+    the colormap starts FLOOR up, never at its near-black bottom."""
+    data, disc = _field()
+    rgba = np.asarray(Image.open(io.BytesIO(
+        overlay.heatmap_png(data, 1, bleed=False))))   # layer 1: all-cold
+    floor = np.array(__import__("matplotlib").colormaps[
+        config.CMAP_CURRENT](overlay.FLOOR)[:3]) * 255
+    assert np.abs(rgba[..., :3][disc] - floor).max() <= 1
+    assert rgba[..., :3][disc].sum(axis=-1).min() > 30   # not near-black
+
+
+def test_heatmap_bleed_ring():
+    """bleed=True: one pixel of half-alpha edge color outside the copper
+    (the mask stops half a cell short of the drawn outline)."""
+    from scipy import ndimage
+    data, disc = _field()
+    rgba = np.asarray(Image.open(io.BytesIO(overlay.heatmap_png(data, 0))))
+    ring = ndimage.binary_dilation(
+        disc, structure=np.ones((3, 3), dtype=bool)) & ~disc
+    assert (rgba[..., 3][ring] == config.OVERLAY_ALPHA // 2).all()
+    outside = ~disc & ~ring
+    assert (rgba[..., 3][outside] == 0).all()
+    assert (rgba[..., 3][disc] == config.OVERLAY_ALPHA).all()
+
+
+def test_heatmap_empty_field():
+    with pytest.raises(ValueError):
+        overlay.heatmap_png(np.full((1, 4, 4), np.nan), 0)

tests/test_skin.py

@@ -64,6 +64,24 @@ def test_parse_frequency():
         skin.parse_frequency("-5k")
 
 
+def test_normalize_decimal():
+    """European decimal commas parse; thousands-separator patterns are
+    rejected ('1,500' silently becoming 1.5 was a 1000x input error)."""
+    assert skin.normalize_decimal("1,5") == "1.5"
+    assert skin.normalize_decimal("0,25") == "0.25"
+    assert skin.normalize_decimal("1,5000") == "1.5000"   # 4 digits: decimal
+    assert skin.normalize_decimal("2.5") == "2.5"
+    for bad in ("1,500", "1.500,5", "1,000,000", "12,345"):
+        with pytest.raises(ValueError, match="separator"):
+            skin.normalize_decimal(bad)
+
+
+def test_parse_frequency_decimal_comma():
+    assert skin.parse_frequency("1,5k") == 1500.0
+    with pytest.raises(ValueError):
+        skin.parse_frequency("1,500")     # ambiguous, not 1.5 Hz
+
+
 def test_single_layer_ac_scales_exactly():
     """Uniform conductance scaling leaves the field shape unchanged:
     R_AC = R_DC * factor to solver precision."""

tests/test_tracks.py

@@ -259,3 +259,25 @@ def test_track_unions_with_fill():
     assert int(s_both.masks.sum()) > int(s_plate.masks.sum())
     assert r_both.R_ohm < 0.75 * r_plate.R_ohm    # bridge shortens the detour
     assert r_both.power_balance_rel < 1e-9
+
+
+def test_pad_copper_bridges_track_junction():
+    """Two traces meet ON an SMD pad, their rounded ends 0.5 mm apart:
+    the junction only exists through the pad copper (board_io stamps
+    the net's pad shapes onto their layers). Without the pad the net
+    is severed - at both track models (rasterized and 1D chain)."""
+    from fill_resistance.errors import ConnectivityError
+    tabs = [[(0, 4.5), (1, 4.5), (1, 5.5), (0, 5.5)],
+            [(19, 4.5), (20, 4.5), (20, 5.5), (19, 5.5)]]
+    pad = [(9.25, 4.4), (10.75, 4.4), (10.75, 5.6), (9.25, 5.6)]
+    segs = [_seg([(0.5, 5), (9.5, 5)], 0.5),
+            _seg([(10.5, 5), (19.5, 5)], 0.5)]
+    r1, r2 = (0, 4.5, 1, 5.5), (19, 4.5, 20, 5.5)
+
+    for h in (0.1, 0.25):          # 5 cells: outlines; 2 cells: 1D chains
+        res, _ = _solve(_seg_problem(segs, r1, r2, fills_mm=tabs + [pad]), h)
+        # ~36 squares of 0.5 mm trace + tabs/pad: sanity-band the value
+        assert 0.007 < res.R_ohm < 0.011
+
+        with pytest.raises(ConnectivityError):
+            _solve(_seg_problem(segs, r1, r2, fills_mm=tabs), h)

tools/deploy.py

@@ -18,7 +18,7 @@ from pathlib import Path
 ROOT = Path(__file__).resolve().parent.parent
 COPY_EXCLUDE = {".venv", ".git", "tests", "tools", "dist", "resources",
                 "__pycache__", ".pytest_cache", "conftest.py", "deploy.ps1",
-                ".gitignore", "metadata.json"}
+                ".gitignore", "metadata.json", "pyproject.toml", "uv.lock"}
 
 
 def plugins_dir(kicad_version: str) -> Path:

tools/gen_readme_figs.py

@@ -0,0 +1,317 @@
+"""Generate the README model figures in docs/img/.
+
+    .venv\\Scripts\\python.exe tools\\gen_readme_figs.py
+
+Real solver output wherever possible: the demo-board maps (raster map
+with the adaptive mesh, current density) and the contact-model
+comparison come straight from the plugin's own pipeline on small
+synthetic boards; only the hole-anatomy cross-section is drawn by hand.
+"""
+import sys
+from pathlib import Path
+
+import numpy as np
+
+ROOT = Path(__file__).resolve().parents[1]
+sys.path.insert(0, str(ROOT))
+
+from fill_resistance import config  # noqa: E402
+
+config.INTERACTIVE = False
+
+import matplotlib  # noqa: E402
+
+matplotlib.use("Agg", force=True)
+
+from fill_resistance import plots, raster, solver  # noqa: E402
+from fill_resistance.geometry import (Electrode, LayerFill, Polygon,  # noqa: E402
+                                      Problem, Rect, ViaLink,
+                                      contact_solder_buildups)
+
+import matplotlib.pyplot as plt  # noqa: E402
+
+plt.switch_backend("Agg")
+plots.INTERACTIVE_BACKEND = None      # save-only: no window panels
+
+NM = 1_000_000
+OUT = ROOT / "docs" / "img"
+
+COPPER = plots._COPPER
+SOLDER = plots._SOLDER
+LEAD = "#5c6570"
+CORE = "#ccd6b3"                      # FR-4
+FILLER = "#eae5dc"                    # non-conductive via fill
+INK = plots._INK
+
+
+def _poly(pts_mm, holes_mm=()) -> Polygon:
+    ring = lambda pts: np.asarray(  # noqa: E731
+        [[int(x * NM), int(y * NM)] for x, y in pts], dtype=np.int64)
+    return Polygon(outline=ring(pts_mm),
+                   holes=[ring(h) for h in holes_mm])
+
+
+def _rect(x0, y0, x1, y1, layer="User.1") -> Rect:
+    return Rect.normalized(int(x0 * NM), int(y0 * NM),
+                           int(x1 * NM), int(y1 * NM), layer)
+
+
+def _disc(x_mm, y_mm, r_mm, n=64) -> Polygon:
+    ang = np.linspace(0, 2 * np.pi, n, endpoint=False)
+    return _poly([(x_mm + r_mm * np.cos(a), y_mm + r_mm * np.sin(a))
+                  for a in ang])
+
+
+def _solve(problem, h_mm, i_test=10.0, model=None):
+    stack = raster.rasterize_stack(problem, int(h_mm * NM))
+    e1, e2 = raster.electrode_masks(stack, problem)
+    p1, p2 = raster.electrode_partition(stack, problem)
+    res = solver.run_solve(problem, stack, e1, e2, i_test,
+                           contact_model=model, parts1=p1, parts2=p2)
+    return res, stack, e1, e2
+
+
+# --- demo board: 2 layers, notched F.Cu pour, half-size B.Cu pour, ----
+# --- a soldered THT-pad contact and a stitching-via field -------------
+
+def demo_problem() -> Problem:
+    # F.Cu: full pour with a notch from the top edge down to y=6 - the
+    # current from the left contact must squeeze through the channel
+    top = _poly([(0, 0), (40, 0), (40, 22), (13, 22), (13, 6),
+                 (10, 6), (10, 22), (0, 22)])
+    # B.Cu: pour on the right half only - vias must carry the transfer
+    bot = _poly([(18, 0), (40, 0), (40, 22), (18, 22)])
+    z_bot = int(1.6 * NM)
+    vias = [ViaLink(x=int(x * NM), y=int(y * NM), drill_nm=300_000,
+                    z_top_nm=-1, z_bot_nm=z_bot + 1, pad_nm=600_000)
+            for x in (20.5, 24.5, 28.5, 32.5, 36) for y in (3, 7, 11, 15, 19)]
+    tht = Electrode(
+        rect=_rect(3.4, 9.9, 5.6, 12.1), contact="F.Cu",
+        label="THT pad (4.5, 11)", drill_nm=1_000_000,
+        pad_nm=int(2.2 * NM), pad_min_nm=int(2.2 * NM),
+        center=(int(4.5 * NM), int(11 * NM)), solder=True,
+        protrusion_side="F.Cu", polygons=[_disc(4.5, 11, 1.1)])
+    lug = Electrode(rect=_rect(37.5, 3, 39.5, 19, "User.2"),
+                    contact="B.Cu", label="lug")
+    p = Problem(
+        board_path="synthetic", net_name="DEMO",
+        rho_ohm_m=1.68e-8, plating_nm=18_000,
+        layers=[LayerFill("F.Cu", 70_000, 0, [top]),
+                LayerFill("B.Cu", 70_000, z_bot, [bot])],
+        vias=vias, electrodes1=[tht], electrodes2=[lug],
+        thickness_source="override")
+    contact_solder_buildups(p)
+    return p
+
+
+def gen_demo_maps():
+    p = demo_problem()
+    res, stack, e1, e2 = _solve(p, 0.05)
+    figs = [
+        (plots.fig_raster(stack, e1, e2, p, res), "demo-raster"),
+        (plots.fig_current(res, stack, e1, e2, p), "demo-current"),
+        (plots.fig_potential(res, stack, e1, e2, p), "demo-potential"),
+    ]
+    plots.save_and_show(figs, OUT, show=False)
+
+
+# --- contact models: equipotential vs uniform injection ---------------
+
+def gen_contact_models():
+    plate = [(0, 0), (24, 0), (24, 18), (0, 18)]
+    r_ohm, zooms = {}, {}
+    # uniform grid: the coarse adaptive leaves would pixelate the |J| zoom
+    config.ADAPTIVE_CELLS = False
+    for model in ("equipotential", "uniform"):
+        p = Problem(
+            board_path="synthetic", net_name="DEMO",
+            rho_ohm_m=1.68e-8, plating_nm=18_000,
+            layers=[LayerFill("F.Cu", 70_000, 0, [_poly(plate)])],
+            vias=[],
+            electrodes1=[Electrode(rect=_rect(4.5, 7.5, 7.5, 10.5))],
+            electrodes2=[Electrode(rect=_rect(22, 1, 23.5, 17, "User.2"))],
+            thickness_source="override")
+        res, stack, _, _ = _solve(p, 0.05, model=model)
+        h_mm = stack.h_nm / NM
+        j = res.Jmag[0] * 1e-6                       # A/mm^2
+        x0, y0 = stack.x0_nm / NM, stack.y0_nm / NM
+        c0, c1 = int((2 - x0) / h_mm), int((13 - x0) / h_mm)
+        r0, r1 = int((3 - y0) / h_mm), int((15 - y0) / h_mm)
+        zooms[model] = (j[r0:r1, c0:c1],
+                        (x0 + c0 * h_mm, x0 + c1 * h_mm,
+                         y0 + r1 * h_mm, y0 + r0 * h_mm))
+        r_ohm[model] = res.R_ohm
+    config.ADAPTIVE_CELLS = True
+
+    vmax = float(np.percentile(
+        zooms["equipotential"][0][np.isfinite(zooms["equipotential"][0])],
+        99.0))
+    fig, axes = plt.subplots(1, 2, figsize=(9.5, 4.2), sharey=True,
+                             layout="constrained")
+    titles = {"equipotential": "equipotential (ideal bonded lug):\n"
+                               "|J| crowds at the contact edges",
+              "uniform": "uniform injection (pressed conductor):\n"
+                         "|J| ramps across the contact"}
+    for ax, model in zip(axes, ("equipotential", "uniform")):
+        data, extent = zooms[model]
+        cmap = matplotlib.colormaps[config.CMAP_CURRENT].copy()
+        cmap.set_bad(plots._BG)
+        im = ax.imshow(data, cmap=cmap, vmin=0, vmax=vmax, origin="upper",
+                       extent=extent, interpolation="nearest")
+        ax.add_patch(plt.Rectangle((4.5, 7.5), 3, 3, fill=False,
+                                   ec="white", ls="--", lw=1.0))
+        ax.set_title(f"{titles[model]}\nR = {r_ohm[model] * 1e3:.3f} mΩ",
+                     fontsize=9, color=INK)
+        ax.set_xlabel("x [mm]", fontsize=8)
+        ax.tick_params(labelsize=8, colors=INK)
+    axes[0].set_ylabel("y [mm]", fontsize=8)
+    cb = fig.colorbar(im, ax=axes, shrink=0.85)
+    cb.set_label("|J| [A/mm²] @ 10 A", fontsize=9)
+    fig.suptitle("The two contact models bracket a real contact:  "
+                 "R$_{equipotential}$ ≤ R$_{real}$ ≤ R$_{uniform}$",
+                 fontsize=10, color=INK)
+    fig.savefig(OUT / "contact-models.png", dpi=config.DPI,
+                facecolor="white", bbox_inches="tight")
+    plt.close(fig)
+    print(f"saved {OUT / 'contact-models.png'}")
+
+
+# --- hole anatomy: hand-drawn cross-section of the four hole types ----
+
+CORE_T = 1.6          # substrate thickness [drawing units ~ mm]
+FOIL_T = 0.18         # foil thickness, exaggerated
+PLATE_W = 0.12        # barrel plating, exaggerated
+CAP_T = 0.07          # via cap
+COAT_T = 0.1          # pad-face solder coat
+Y_TOP = CORE_T + FOIL_T
+
+
+def _board_segment(ax, x0, x1):
+    ax.add_patch(plt.Rectangle((x0, 0), x1 - x0, CORE_T, fc=CORE, ec="none"))
+    for y in (CORE_T, -FOIL_T):
+        ax.add_patch(plt.Rectangle((x0, y), x1 - x0, FOIL_T,
+                                   fc=COPPER, ec="none"))
+
+
+def _barrel(ax, xc, drill):
+    for s in (-1, 1):
+        x = xc + s * drill / 2 - (PLATE_W if s > 0 else 0)
+        ax.add_patch(plt.Rectangle((x, -FOIL_T), PLATE_W,
+                                   CORE_T + 2 * FOIL_T, fc=COPPER,
+                                   ec="none"))
+
+
+def _label(ax, text, xy, xytext, ha="left"):
+    ax.annotate(text, xy, xytext=xytext, fontsize=7.5, color=INK, ha=ha,
+                va="center",
+                arrowprops=dict(arrowstyle="-", color=INK, lw=0.7,
+                                shrinkA=2, shrinkB=1))
+
+
+def gen_hole_anatomy():
+    fig, ax = plt.subplots(figsize=(12.5, 5.2), layout="constrained")
+    holes = [(3.0, 0.7), (10.0, 1.6), (17.5, 1.6), (25.0, 1.6)]
+    edges = [0.0]
+    for xc, d in holes:
+        edges += [xc - d / 2, xc + d / 2]
+    edges.append(28.5)
+    for x0, x1 in zip(edges[::2], edges[1::2]):
+        _board_segment(ax, x0, x1)
+    for xc, d in holes:
+        _barrel(ax, xc, d)
+
+    # 1: small via, filled + capped
+    xc, d = holes[0]
+    ax.add_patch(plt.Rectangle((xc - d / 2 + PLATE_W, -FOIL_T),
+                               d - 2 * PLATE_W, CORE_T + 2 * FOIL_T,
+                               fc=FILLER, ec="none"))
+    for y in (Y_TOP, -FOIL_T - CAP_T):
+        ax.add_patch(plt.Rectangle((xc - d / 2 - 0.12, y), d + 0.24, CAP_T,
+                                   fc=COPPER, ec="none"))
+    _label(ax, "cap, CAP_PLATING_UM (15 µm)\non both outer mouths",
+           (xc, Y_TOP + CAP_T), (xc, 3.3), ha="center")
+    _label(ax, "non-conductive fill", (xc, 0.8), (5.6, -0.9))
+
+    # 2: big via, mouth open
+    xc, d = holes[1]
+    _label(ax, "open mouth: covered cells\nremoved, sub-cell mouths\n"
+               "scale the sheet conductance",
+           (xc, Y_TOP - FOIL_T / 2), (xc, 3.2), ha="center")
+
+    # 3: populated THT pad - full solder joint
+    xc, d = holes[2]
+    lead_w = d - 0.5                      # drill - clearance, exaggerated
+    pad_r = 1.7
+    prot = 1.5
+    sn_edge = "#7d8791"                   # delineate solder sub-shapes
+    # solder fill between plating and lead
+    for s in (-1, 1):
+        x0 = xc + s * lead_w / 2 if s > 0 else xc - d / 2 + PLATE_W
+        ax.add_patch(plt.Rectangle((x0, -FOIL_T),
+                                   d / 2 - PLATE_W - lead_w / 2,
+                                   CORE_T + 2 * FOIL_T, fc=SOLDER,
+                                   ec="none"))
+    # pad-face coat, solder side only
+    ax.add_patch(plt.Rectangle((xc - pad_r, Y_TOP), 2 * pad_r, COAT_T,
+                               fc=SOLDER, ec=sn_edge, lw=0.5))
+    # solder cone: protrusion height at the wall -> 0 at the pad edge
+    for s in (-1, 1):
+        wall = xc + s * lead_w / 2
+        ax.add_patch(plt.Polygon(
+            [(wall, Y_TOP + prot), (wall, Y_TOP + COAT_T),
+             (xc + s * pad_r, Y_TOP + COAT_T)],
+            closed=True, fc=SOLDER, ec=sn_edge, lw=0.5))
+    # lead: through the hole, protruding on top, component below
+    ax.add_patch(plt.Rectangle((xc - lead_w / 2, -2.05), lead_w,
+                               2.05 + Y_TOP + prot, fc=LEAD, ec="none"))
+    ax.add_patch(plt.Rectangle((xc - 1.5, -2.75), 3.0, 0.7,
+                               fc="#8a8f96", ec="none"))
+    ax.text(xc, -2.4, "component", fontsize=7.5, color="white",
+            ha="center", va="center")
+    _label(ax, "clipped lead protrudes\nTHT_LEAD_PROTRUSION_MM (1.5 mm)",
+           (xc + lead_w / 2, Y_TOP + prot - 0.2), (xc + 3.4, 4.15))
+    _label(ax, "solder cone: full height at the\nwall, tapers to 0 at the "
+               "pad edge", (xc - (lead_w / 2 + pad_r) / 2, Y_TOP + 0.7),
+           (13.6, 4.2), ha="center")
+    _label(ax, "pad-face solder coat (50 µm),\nSOLDER side only",
+           (xc - pad_r + 0.2, Y_TOP + COAT_T / 2), (12.9, 2.35),
+           ha="center")
+    _label(ax, "solder-filled hole: lead ∥ solder ∥ plating\n"
+               "lead ⌀ = drill − THT_LEAD_CLEARANCE_MM",
+           (xc - d / 2 + PLATE_W + 0.07, 0.5), (12.3, -1.5), ha="center")
+    _label(ax, "component side:\npad face stays bare",
+           (xc + pad_r - 0.3, -FOIL_T), (xc + 4.0, -1.05))
+
+    # 4: DNP THT pad
+    xc, d = holes[3]
+    _label(ax, "open hole on every layer,\nplating-only barrel, no joint",
+           (xc, 0.8), (xc + 1.3, -2.45), ha="center")
+
+    for (xc, _), title in zip(holes, (
+            "via ≤ cap-drill\n(capped)", "via > cap-drill\n(open)",
+            "THT pad, populated\n(read from KiCad)", "THT pad, DNP")):
+        ax.text(xc, 5.6, title, fontsize=9, color=INK, ha="center",
+                va="top", fontweight="bold")
+
+    handles = [plt.Rectangle((0, 0), 1, 1, fc=c) for c in
+               (COPPER, SOLDER, LEAD, CORE, FILLER)]
+    fig.legend(handles, ("copper (foil / plating / pad)", "solder",
+                         "component lead", "FR-4", "non-conductive fill"),
+               loc="outside right center", fontsize=8, framealpha=0.95)
+    ax.set_title("How drilled holes are modeled — cross-section "
+                 "(vertical scale exaggerated)", fontsize=11, color=INK)
+    ax.set_xlim(-0.3, 29.0)
+    ax.set_ylim(-3.1, 5.7)
+    ax.set_aspect("equal")
+    ax.axis("off")
+    fig.savefig(OUT / "hole-model.png", dpi=config.DPI,
+                facecolor="white", bbox_inches="tight")
+    plt.close(fig)
+    print(f"saved {OUT / 'hole-model.png'}")
+
+
+if __name__ == "__main__":
+    OUT.mkdir(parents=True, exist_ok=True)
+    gen_hole_anatomy()
+    gen_contact_models()
+    gen_demo_maps()

tools/kicad_heatmap_overlay.py

@@ -0,0 +1,141 @@
+"""Standalone runner for the EXPERIMENTAL in-KiCad result overlays
+(also available as a dialog checkbox in the plugin): solve the open
+board headlessly and push per-layer |J| heatmaps as unlocked
+ReferenceImages, transparent outside copper.
+
+    python tools/kicad_heatmap_overlay.py --net VOUT+ --amps 45
+        -> all included copper layers onto config.OVERLAY_LAYERS
+           (User.9..User.12, stackup order, top first)
+    python tools/kicad_heatmap_overlay.py --net X --source B.Cu --dest Eco1.User
+        -> a single layer wherever you want
+
+Needs KiCad >= 10.0.1 with the board open, electrode markers or a
+selection as in a normal plugin run, and the destination layers enabled
+in Board Setup. Re-running replaces the previous overlays. Remove with
+tools/kicad_overlay_test.py --remove --layer <dest>.
+"""
+import argparse
+import sys
+from pathlib import Path
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from kipy.board_types import ReferenceImage
+from kipy.geometry import Vector2
+from kipy.util.board_layer import layer_from_canonical_name
+
+from fill_resistance import config, raster, solver
+from fill_resistance import board_io as bio
+from fill_resistance.overlay import heatmap_png
+
+
+def extract_problem(board, net_arg=None):
+    """Same flow as `python -m fill_resistance.board_io` (dump path)."""
+    # a clicked overlay must not switch the electrode scan into
+    # selection mode - reference images can never be contacts
+    sel = list(board.get_selection())
+    if sel and all(isinstance(s, ReferenceImage) for s in sel):
+        board.clear_selection()
+    stackup = bio.get_stackup_info(board)
+    es1, es2, net_hint = bio.get_electrodes(board, stackup)
+    if bio.any_zone_unfilled(board):
+        bio.refill(board)
+    fills = bio.gather_net_fills(board)
+    tracks = bio.gather_net_tracks(board) if config.INCLUDE_TRACKS else {}
+    copper = bio.merge_copper(fills, bio.tracks_as_polygons(tracks))
+    nets = bio.nets_overlapping(copper, es1, es2)
+    if net_arg:
+        net = net_arg
+    elif net_hint in nets:
+        net = net_hint
+    elif len(nets) == 1:
+        net = nets[0]
+    else:
+        raise SystemExit(f"candidate nets: {nets}; pass one with --net")
+    # marker rectangles may exist for SEVERAL nets (board-wide scan):
+    # keep only the parts overlapping the chosen net's copper
+    per_layer = copper.get(net, {})
+    def on_net(e):
+        return any(bio._rect_overlaps(e.rect, polys)
+                   for polys in per_layer.values())
+    es1, es2 = [e for e in es1 if on_net(e)], [e for e in es2 if on_net(e)]
+    if not es1 or not es2:
+        raise SystemExit(f"no V+/V- marker overlaps {net} copper")
+    print(f"{len(es1)} V+ / {len(es2)} V- marker(s) on {net}")
+    return bio.build_problem(board, net, list(per_layer), es1, es2,
+                             stackup, fills, tracks=tracks)
+
+
+def main() -> None:
+    ap = argparse.ArgumentParser(description=__doc__)
+    ap.add_argument("--source", default=None,
+                    help="single copper layer to overlay (default: ALL "
+                         "included layers onto config.OVERLAY_LAYERS)")
+    ap.add_argument("--dest", default=None,
+                    help="destination layer for --source (default User.9; "
+                         "must be enabled in Board Setup)")
+    ap.add_argument("--net", default=None, help="net name (default: auto)")
+    ap.add_argument("--amps", type=float, default=None,
+                    help="test current [A] (default: config)")
+    ap.add_argument("--lock", action="store_true",
+                    help="lock the overlays (default unlocked: easier to "
+                         "delete; reruns replace them either way)")
+    ap.add_argument("--alpha", type=int, default=None,
+                    help="overlay opacity over copper, 0-255 (default "
+                         "config.OVERLAY_ALPHA)")
+    args = ap.parse_args()
+    if args.alpha is not None:
+        config.OVERLAY_ALPHA = args.alpha
+
+    _, board = bio.connect()
+    problem = extract_problem(board, args.net)
+
+    h = raster.choose_cell_size(problem.copper_bbox(), len(problem.layers))
+    print(f"rasterizing at {h / 1000:.1f} um ...")
+    stack = raster.rasterize_stack(problem, h)
+    # board-wide marker scan: drop parts that land on no copper of THIS
+    # net (markers belonging to other nets' analyses)
+    for name in ("electrodes1", "electrodes2"):
+        parts = getattr(problem, name)
+        keep = [e for e in parts
+                if raster._part_mask3d(stack, problem, e).any()]
+        if len(keep) != len(parts):
+            print(f"ignoring {len(parts) - len(keep)} marker(s) off-net "
+                  f"({name[-1] == '1' and 'V+' or 'V-'})")
+        if not keep:
+            raise SystemExit(f"no {name} marker lands on this net's copper")
+        setattr(problem, name, keep)
+    e1, e2 = raster.electrode_masks(stack, problem)
+    i_test = args.amps if args.amps is not None else config.TEST_CURRENT_A
+    print(f"solving @ {i_test:g} A DC ...")
+    result = solver.run_solve(problem, stack, e1, e2, i_test)
+    print(f"R = {result.R_ohm * 1e3:.4f} mOhm, P = {result.P_total:.3f} W "
+          f"@ {i_test:g} A")
+
+    if args.source is None:
+        bio.push_result_overlays(board, stack, result, lock=args.lock)
+        return
+
+    names = stack.layer_names
+    if args.source not in names:
+        raise SystemExit(f"layer {args.source} not in solve ({names})")
+    png = heatmap_png(result.Jmag * 1e-6, names.index(args.source))
+    ny, nx = stack.shape2d
+    w_nm, h_nm = nx * stack.h_nm, ny * stack.h_nm
+    dest_name = args.dest or "User.9"
+    dest = layer_from_canonical_name(dest_name)
+    n = bio.remove_overlays(board, dest)
+    ref = ReferenceImage()
+    ref.layer = dest
+    ref.position = Vector2.from_xy(round(stack.x0_nm + w_nm / 2),
+                                   round(stack.y0_nm + h_nm / 2))
+    ref.image_scale = w_nm / (nx * bio.OVERLAY_PIX_NM)
+    ref.image_data = png
+    ref.locked = args.lock
+    bio._create_reference_image(board, ref)
+    print(f"{args.source} -> {dest_name} ({nx}x{ny} px, "
+          f"{len(png) / 1024:.0f} kB" + (f", replaced {n}" if n else "") + ")")
+
+
+if __name__ == "__main__":
+    main()

tools/kicad_overlay_test.py

@@ -0,0 +1,156 @@
+"""Route-A experiment: push a bitmap overlay into the open KiCad board as
+a locked ReferenceImage on a User layer via the IPC API.
+
+Pushes a fiducial test pattern (corner + center crosshairs, 10 mm grid,
+translucent gradient) sized to the board outline so alignment and scale
+can be verified by eye in the editor. Re-running replaces the previous
+overlay. Requires KiCad >= 10.0.1 (ReferenceImage over the API).
+
+    python tools/kicad_overlay_test.py [--layer Cmts.User] [--remove]
+    python tools/kicad_overlay_test.py --image heat.png --bbox x0,y0,x1,y1
+                                       (mm; push an arbitrary PNG instead)
+
+The overlay is editor-only: reference images never plot to gerbers.
+Delete it any time by selecting it in KiCad (it sits on the chosen
+layer) or with --remove. The layer must be enabled in Board Setup:
+User.1..User.45 usually are NOT (KiCad refuses the item with 'no
+overlapping layers with the board'); Cmts.User/Eco1.User always exist.
+"""
+import argparse
+import io
+import sys
+from pathlib import Path
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1]))
+
+from kipy.board_types import ReferenceImage
+from kipy.geometry import Vector2
+from kipy.util.board_layer import canonical_name, layer_from_canonical_name
+
+from fill_resistance.board_io import (OVERLAY_PIX_NM as PIX_NM,
+                                      _create_reference_image, connect,
+                                      remove_overlays)
+
+NM = 1_000_000
+
+
+def board_bbox_nm(board):
+    """Union bbox of the Edge.Cuts shapes (fallback: all pads)."""
+    items = [s for s in board.get_shapes()
+             if canonical_name(s.layer) == "Edge.Cuts"]
+    if not items:
+        items = list(board.get_pads())
+    if not items:
+        raise SystemExit("board has no Edge.Cuts shapes and no pads")
+    x0 = y0 = None
+    x1 = y1 = None
+    for it in items:
+        box = board.get_item_bounding_box(it)
+        if box is None:
+            continue
+        lo_x, lo_y = box.pos.x, box.pos.y
+        hi_x, hi_y = lo_x + box.size.x, lo_y + box.size.y
+        x0 = lo_x if x0 is None else min(x0, lo_x)
+        y0 = lo_y if y0 is None else min(y0, lo_y)
+        x1 = hi_x if x1 is None else max(x1, hi_x)
+        y1 = hi_y if y1 is None else max(y1, hi_y)
+    return x0, y0, x1, y1
+
+
+def fiducial_png(w_nm: float, h_nm: float, px_per_mm: float = 16.0):
+    """RGBA test pattern: translucent gradient, 10 mm grid, opaque
+    crosshairs at the four corners and the center."""
+    import numpy as np
+    from PIL import Image
+
+    w_px = max(2, round(w_nm / NM * px_per_mm))
+    h_px = max(2, round(h_nm / NM * px_per_mm))
+    xx = np.linspace(0.0, 1.0, w_px)[None, :]
+    yy = np.linspace(0.0, 1.0, h_px)[:, None]
+    rgba = np.zeros((h_px, w_px, 4), dtype=np.uint8)
+    rgba[..., 0] = (255 * xx).astype(np.uint8)          # red ramp ->
+    rgba[..., 2] = (255 * yy).astype(np.uint8)          # blue ramp v
+    rgba[..., 1] = 60
+    rgba[..., 3] = 70                                   # mostly see-through
+
+    step = round(10.0 * px_per_mm)                      # 10 mm grid
+    for x in range(0, w_px, step):
+        rgba[:, x:x + 2, :3] = 255
+        rgba[:, x:x + 2, 3] = 150
+    for y in range(0, h_px, step):
+        rgba[y:y + 2, :, :3] = 255
+        rgba[y:y + 2, :, 3] = 150
+
+    def cross(cx, cy, arm=round(3 * px_per_mm)):
+        x_lo, x_hi = max(0, cx - arm), min(w_px, cx + arm + 1)
+        y_lo, y_hi = max(0, cy - arm), min(h_px, cy + arm + 1)
+        cy2 = np.clip(cy, 0, h_px - 2)
+        cx2 = np.clip(cx, 0, w_px - 2)
+        rgba[cy2:cy2 + 2, x_lo:x_hi] = (255, 0, 0, 255)
+        rgba[y_lo:y_hi, cx2:cx2 + 2] = (255, 0, 0, 255)
+
+    for cx in (0, w_px - 1):
+        for cy in (0, h_px - 1):
+            cross(cx, cy)
+    cross(w_px // 2, h_px // 2)
+
+    buf = io.BytesIO()
+    # no dpi= : without a density chunk KiCad assumes the 300 PPI default
+    Image.fromarray(rgba, "RGBA").save(buf, format="PNG")
+    return buf.getvalue(), w_px, h_px
+
+
+def main() -> None:
+    ap = argparse.ArgumentParser(description=__doc__)
+    ap.add_argument("--layer", default="Cmts.User",
+                    help="destination layer (default Cmts.User; must be "
+                         "enabled in Board Setup)")
+    ap.add_argument("--remove", action="store_true",
+                    help="only remove existing overlays on the layer")
+    ap.add_argument("--image", help="push this PNG instead of the pattern")
+    ap.add_argument("--bbox", help="x0,y0,x1,y1 [mm] for --image")
+    args = ap.parse_args()
+
+    _, board = connect()
+    layer = layer_from_canonical_name(args.layer)
+
+    n = remove_overlays(board, layer)
+    if n:
+        print(f"removed {n} previous overlay(s) on {args.layer}")
+    if args.remove:
+        return
+
+    if args.image:
+        if not args.bbox:
+            raise SystemExit("--image needs --bbox x0,y0,x1,y1 [mm]")
+        x0, y0, x1, y1 = (float(v) * NM for v in args.bbox.split(","))
+        png = Path(args.image).read_bytes()
+        from PIL import Image
+        w_px, h_px = Image.open(io.BytesIO(png)).size
+    else:
+        x0, y0, x1, y1 = board_bbox_nm(board)
+        png, w_px, h_px = fiducial_png(x1 - x0, y1 - y0)
+
+    scale = (x1 - x0) / (w_px * PIX_NM)
+
+    ref = ReferenceImage()
+    ref.layer = layer
+    ref.position = Vector2.from_xy(round((x0 + x1) / 2), round((y0 + y1) / 2))
+    ref.image_scale = scale
+    ref.image_data = png
+    ref.locked = False          # unlocked: easy to delete; reruns replace
+    _create_reference_image(board, ref)
+
+    got = [r for r in board.get_reference_images() if r.layer == layer]
+    print(f"pushed {len(png) / 1024:.0f} kB PNG ({w_px}x{h_px} px) onto "
+          f"{args.layer}: {(x1 - x0) / NM:.2f} x {(y1 - y0) / NM:.2f} mm at "
+          f"({x0 / NM:.2f}, {y0 / NM:.2f}) mm, scale {scale:.4f}")
+    for r in got:
+        print(f"readback: {r!r}")
+    print(f"-> enable layer '{args.layer}' in the Appearance panel; the "
+          f"red crosshairs must sit on the board bbox corners/center and "
+          f"the white grid must be 10 mm. Remove with --remove.")
+
+
+if __name__ == "__main__":
+    main()