Update reg presets: mode knob + FOV evaluator + dimension-aware patch; models.json fetched from HF (unbundled)

CBCT_CT_HeadNeck/Model.py

@@ -32,6 +32,7 @@ NOTE: do NOT add ``from __future__ import annotations`` here — KonfAI's config
 runtime-evaluated annotations (``get_origin``); PEP 563 stringized annotations break arg resolution.
 """
 
+import json
 import os
 import re
 import shutil
@@ -52,6 +53,212 @@ from konfai.utils.dataset import Attribute, data_to_image, image_to_data
 # Set KONFAI_ELASTIX_DIR to point at an existing install and skip the download.
 ELASTIX_CACHE = Path.home() / ".cache" / "konfai" / "elastix-impact"
 
+# ---------------------------------------------------------------------------------------------------
+# Per-resolution model matrix (the config is the source of truth) -> generated IMPACT parameter map.
+# IMPACT field docs: https://github.com/vboussot/ImpactLoss/tree/main/ParameterMaps
+# The forced per-model props (dimension/channels/FOV formula) live in a registry (models.json on
+# VBoussot/impact-torchscript-models); the config carries the FREE knobs (which models per resolution,
+# feature voxel size, iterations, per-model layer weights/mask/subset/pca/distance) and the global
+# ``mode``. PatchSize follows ImpactMode: Static -> "0 0 0" (whole image); Jacobian -> the model FOV
+# evaluated from the registry formula (MIND 2*r*d+1, TS/MRSeg 2^l+3, SAM 29, DINOv2 14) as a cube.
+# ---------------------------------------------------------------------------------------------------
+
+_IMPACT_MODELS_REGISTRY = "VBoussot/impact-torchscript-models:models.json"
+
+# ``2^l+3`` grows with depth but the segmenters' receptive field plateaus: layers 7-8 share layer 6's
+# FOV (the "ramp max"). A config that deep should really run in Static (whole image) anyway; in Jacobian
+# we clamp ``l`` to this plateau so the patch stays finite and matches the real FOV.
+_FOV_RAMP_MAX_LAYER = 6
+
+
+def _num(x: object) -> str:
+    """Format a number the elastix way: integers without a trailing '.0' (6.0 -> '6', 0.2 -> '0.2')."""
+    return "%g" % float(x)
+
+
+class ModelSpec:
+    """One feature model at one resolution, with its OWN config (several models may share a resolution).
+
+    ``ref`` selects the model; ``voxel_size`` / ``layers_weight`` / ``subset_features`` / ``pca`` /
+    ``distance`` are its free per-(resolution, model) tuning knobs (the doc's per-model *tuning* fields).
+    The intrinsic per-model props — dimension, channels, ``layers_mask``, patch-size (FOV) — come from the
+    registry (read-only); ``layers_mask`` / ``distance`` left empty fall back to the registry default.
+    """
+
+    def __init__(
+        self,
+        ref: str,
+        voxel_size: list[float] = [],
+        layers_weight: list[float] = [1.0],
+        subset_features: int = 0,
+        pca: int = 0,
+        distance: str = "",
+        layers_mask: str = "",
+    ) -> None:
+        self.ref = ref
+        self.voxel_size = voxel_size
+        self.layers_weight = layers_weight
+        self.subset_features = subset_features
+        self.pca = pca
+        self.distance = distance
+        self.layers_mask = layers_mask
+
+
+class ResolutionSpec:
+    """One elastix resolution level: its iteration budget and the models compared there (each self-configured)."""
+
+    def __init__(self, max_iterations: int, models: dict[str, ModelSpec]) -> None:
+        self.max_iterations = max_iterations
+        self.models = models
+
+
+def _sorted_specs(mapping: dict) -> list:
+    """dict keyed by string indices ('0','1',...) -> values in numeric order (well-defined res/model order)."""
+    return [mapping[k] for k in sorted(mapping, key=lambda key: int(key))]
+
+
+def load_models_registry(ref: str = _IMPACT_MODELS_REGISTRY) -> dict:
+    """Load models.json (forced params per model) from the model repo on Hugging Face.
+
+    The registry is NOT bundled with the preset — it lives on the models repo and is fetched from there.
+    Resolution: the ``KONFAI_IMPACT_MODELS_REGISTRY`` env path wins (dev/offline); otherwise ``ref`` must be
+    a ``repo:file`` Hugging Face reference.
+    """
+    local = os.environ.get("KONFAI_IMPACT_MODELS_REGISTRY", "")
+    if local:
+        path = Path(local)
+    elif ":" in ref:
+        repo, filename = ref.split(":", 1)
+        path = Path(hf_hub_download(repo_id=repo, filename=filename, repo_type="model"))  # nosec B615
+    else:
+        raise ValueError(
+            f"models_registry '{ref}' must be a 'repo:file' Hugging Face reference (the registry is fetched "
+            f"from HF, not bundled) — or set KONFAI_IMPACT_MODELS_REGISTRY to a local file for offline use."
+        )
+    return json.loads(path.read_text(encoding="utf-8"))
+
+
+def _model_key(ref: str) -> str:
+    """Registry key / staged relative path = the model file within the models repo (strip a 'repo:' prefix)."""
+    return ref.split(":", 1)[1] if ":" in ref else ref
+
+
+def _deepest_active_layer(layers_mask: str) -> int:
+    """Deepest (largest-FOV) layer selected by ``layers_mask``, as a 0-based index read left-to-right.
+
+    A model returns its feature layers shallow->deep (``[layer_0, layer_1, ...]``, see the model repo's
+    build scripts); ``layers_mask`` has one char per returned layer, position ``i`` == ``layer_i``, ``'1'``
+    = selected. In Jacobian the patch must cover the receptive field of the DEEPEST selected layer, so the
+    FOV is governed by the rightmost ``'1'``.
+    """
+    mask = layers_mask.strip().strip('"')
+    active = [i for i, char in enumerate(mask) if char == "1"]
+    if not active:
+        raise ValueError(f"LayersMask '{layers_mask}' selects no layer; cannot derive the model FOV.")
+    return max(active)
+
+
+def _fov_value(fov: dict, layers_mask: str) -> int:
+    """Evaluate a model's field-of-view (in voxels) from its registry ``fov`` spec.
+
+    Supported formulas (from the model repo, https://huggingface.co/VBoussot/impact-torchscript-models):
+      ``2*r*d+1``  MIND, from the handcrafted radius ``r`` / dilation ``d`` (e.g. R1D2 -> 5);
+      ``2^l+3``    TotalSegmentator / MRSegmentator, ``l`` = the deepest layer picked by ``layers_mask``,
+                   clamped to the receptive-field plateau ``_FOV_RAMP_MAX_LAYER`` (layers 7-8 -> layer 6);
+      a bare int  a fixed FOV (SAM2.1 -> 29, DINOv2 -> 14);
+      ``Global``   Anatomix — whole-image only (Static); has no finite Jacobian patch -> error.
+    An explicit ``value`` in the spec is honoured as a precomputed shortcut when the formula needs none.
+    """
+    formula = str(fov.get("formula", "")).strip()
+    key = re.sub(r"\s+", "", formula).lower()
+    if key.isdigit():
+        return int(key)
+    if key == "2*r*d+1":
+        return 2 * int(fov["r"]) * int(fov["d"]) + 1
+    if key == "2^l+3":
+        return 2 ** min(_deepest_active_layer(layers_mask), _FOV_RAMP_MAX_LAYER) + 3
+    if "global" in key:
+        raise ValueError(f"model FOV '{formula}' is whole-image only (Static); it has no Jacobian patch size.")
+    if fov.get("value") is not None:
+        return int(fov["value"])
+    raise ValueError(f"cannot evaluate model FOV formula '{formula}'.")
+
+
+def _patch_size(mode: str, entry: dict, layers_mask: str) -> str:
+    """PatchSize from the model FOV, one token per model axis (2D model -> 2 tokens, 3D -> 3): Static ->
+    whole image (all zeros); Jacobian -> the evaluated FOV repeated over the axes. A 2D model mixed with a
+    3D one at a resolution concatenates as e.g. ``29 29 11 11 11`` (SAM 2D + TS 3D), matching IMPACT."""
+    dim = int(entry.get("dimension", 3))
+    if mode.strip().strip('"').lower() != "jacobian":
+        return " ".join(["0"] * dim)
+    fov = _fov_value(entry.get("fov", {}), layers_mask)
+    return " ".join([str(fov)] * dim)
+
+
+def generate_impact_parameter_map(
+    template_text: str, resolutions: dict, registry: dict, mode: str = "Static"
+) -> str:
+    """Rewrite the resolution-dependent lines of ``template_text`` from the model matrix ``resolutions``.
+
+    Regenerated: MaximumNumberOfIterations, NumberOfResolutions, Fixed/MovingImagePyramidRescaleSchedule,
+    ImpactMode (from the config ``mode``), and the whole ImpactXxxK block; every other template line is
+    kept verbatim (optimizer, transform, metric weights, components...). N (number of resolutions) is
+    deduced from the config. ``mode`` drives PatchSize: Static -> ``0 0 0`` (whole image); Jacobian -> the
+    per-model FOV evaluated from the registry formula and the cell's ``layers_mask``.
+    """
+    res = _sorted_specs(resolutions)
+    n = len(res)
+    mode_clean = mode.strip().strip('"') or "Static"
+
+    impact: list[str] = []
+    for k, r in enumerate(res):
+        models = _sorted_specs(r.models)
+        entries = [registry[_model_key(m.ref)] for m in models]
+
+        def row(stem: str, values: list[str]) -> None:
+            impact.append(f"(Impact{stem}{k} " + " ".join(values) + ")")
+
+        # From the registry (models.json on the model repo) ONLY the 3 truly model-fixed props:
+        # Dimension, NumberOfChannels, PatchSize (the model FOV). Everything else is a per-model tuning knob
+        # taken straight from the cell: VoxelSize / LayersMask / SubsetFeatures / PCA / Distance / LayersWeight.
+        row("ModelsPath", [f'"{_model_key(m.ref)}"' for m in models])
+        row("Dimension", [e["dimension"] for e in entries])
+        row("NumberOfChannels", [e["numberofchannels"] for e in entries])
+        row("PatchSize", [_patch_size(mode_clean, e, m.layers_mask) for e, m in zip(entries, models)])
+        row("VoxelSize", [" ".join(_num(v) for v in m.voxel_size) for m in models])
+        row("LayersMask", [f'"{m.layers_mask}"' for m in models])
+        row("SubsetFeatures", [str(m.subset_features) for m in models])
+        row("PCA", [str(m.pca) for m in models])
+        row("Distance", [f'"{m.distance}"' for m in models])
+        row("LayersWeight", [" ".join(_num(w) for w in m.layers_weight) for m in models])
+        impact.append("")  # blank line between resolutions, mirroring the reference maps
+
+    # The per-resolution block is the contiguous span from the first to the last ``Impact<name><k>`` line
+    # (the blank lines the reference maps put BETWEEN resolutions fall inside that span). Replace the whole
+    # span in one shot with the generated block, so the reference blanks are not kept on top of ours.
+    lines = template_text.splitlines()
+    indexed = [(re.match(r"^\s*\((\S+?)\s+(.*?)\)\s*$", ln), ln) for ln in lines]
+    block_rows = [i for i, (m, _) in enumerate(indexed) if m and re.match(r"^Impact[A-Za-z]+\d+$", m.group(1))]
+    block_lo, block_hi = (block_rows[0], block_rows[-1]) if block_rows else (-1, -2)
+
+    out: list[str] = []
+    for i, (m, line) in enumerate(indexed):
+        key = m.group(1) if m else None
+        if block_lo <= i <= block_hi:
+            if i == block_lo:  # replace the whole span at its first line, drop the rest (incl. inner blanks)
+                out.extend(impact[:-1])
+        elif key == "MaximumNumberOfIterations":
+            out.append("(MaximumNumberOfIterations " + " ".join(_num(r.max_iterations) for r in res) + ")")
+        elif key == "NumberOfResolutions":
+            out.append(f"(NumberOfResolutions {n})")
+        elif key in ("FixedImagePyramidRescaleSchedule", "MovingImagePyramidRescaleSchedule"):
+            out.append(f"({key} " + " ".join(["1"] * 3 * n) + ")")
+        elif key == "ImpactMode":
+            out.append(f'(ImpactMode "{mode_clean}")')
+        else:
+            out.append(line)
+    return "\n".join(out)
+
 
 class ElastixEngine:
     """Run the elastix-IMPACT binary on a fixed/moving pair; return (moved, dvf) on the fixed grid.
@@ -60,14 +267,57 @@ class ElastixEngine:
     does NOT ship it), so registration is a subprocess call, not ``sitk.ElastixImageFilter``.
     """
 
-    def __init__(self, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         self._bundle_dir = Path(__file__).resolve().parent
         self._parameter_maps = [self._bundle_dir / Path(p).name for p in parameter_maps]
+        self._max_iterations = max_iterations
+        self._final_grid_spacing = final_grid_spacing
+        self._subset_features = subset_features
+        self._spatial_samples = spatial_samples
+        self._parameter_overrides = list(parameter_overrides)
+        # ImpactMode: Static computes features once per level (PatchSize 0 0 0 = whole image); Jacobian
+        # samples random patches sized to the model FOV each iteration. Global knob: one mode per preset.
+        self._mode = mode
+        # Matrix mode: when `resolutions` is given the parameter map is GENERATED from it (the config is the
+        # source of truth). An empty `resolutions` = an intensity preset (no IMPACT feature models): the fixed
+        # parameter maps are staged with only the global knob overrides.
+        self._resolutions = resolutions
+        self._registry = load_models_registry(models_registry) if resolutions else {}
+        # The feature models are DERIVED — the unique refs across the matrix cells (no flat `models` param).
+        models: list[str] = []
+        for res in _sorted_specs(resolutions):
+            for model in _sorted_specs(res.models):
+                if model.ref not in models:
+                    models.append(model.ref)
         self._models = models
-        self._iterations = iterations
+        # `iterations` (the progress-bar total) is NOT a config parameter — it is DERIVED: the sum of the
+        # per-resolution iteration budgets, read from the matrix (matrix mode) or the maps (legacy).
+        self._iterations = self._total_iterations()
         self._elastix_bin = self._ensure_binary()
         self._local_models = self._download_models()
 
+    def _total_iterations(self) -> int:
+        """Total iterations across all resolutions — the progress-bar budget, derived from the config."""
+        if self._resolutions:
+            return sum(int(res.max_iterations) for res in _sorted_specs(self._resolutions))
+        total = 0
+        for src in self._parameter_maps:
+            match = re.search(r"\(MaximumNumberOfIterations\s+([^)]*)\)", src.read_text(encoding="utf-8"))
+            if match:
+                total += sum(int(token) for token in match.group(1).split())
+        return total
+
     def _ensure_binary(self) -> Path:
         # Optional override: point at an existing elastix-IMPACT install (skips the download).
         override = os.environ.get("KONFAI_ELASTIX_DIR", "")
@@ -91,17 +341,90 @@ class ElastixEngine:
             models.append((filename, local))
         return models
 
+    def _parameter_map_overrides(self, global_only: bool = False) -> tuple[dict[str, str], list[tuple[str, str]]]:
+        """The tuned knobs as parameter-map overrides: ``(per_token, exact)``.
+
+        ``per_token`` maps an elastix key (or the ``ImpactSubsetFeatures`` prefix) to a value that replaces
+        **each** existing token, so per-resolution / per-model multiplicity is preserved (e.g.
+        ``(MaximumNumberOfIterations 500 250)`` -> ``(MaximumNumberOfIterations 300 300)``). ``exact``
+        entries (from ``parameter_overrides``, ``Key=value text``) replace the whole value verbatim and win
+        over the named knobs. Overrides only REPLACE keys already present in a map — never inject new ones.
+        ``global_only`` (matrix mode) keeps just the map-wide knobs and drops ``max_iterations`` /
+        ``subset_features`` — the per-resolution matrix already sets those per cell.
+        """
+        per_token: dict[str, str] = {}
+        if not global_only and self._max_iterations > 0:
+            per_token["MaximumNumberOfIterations"] = str(int(self._max_iterations))
+        if self._final_grid_spacing > 0:
+            per_token["FinalGridSpacingInPhysicalUnits"] = str(float(self._final_grid_spacing))
+        if not global_only and self._subset_features > 0:
+            per_token["ImpactSubsetFeatures"] = str(int(self._subset_features))  # prefix: indexed per metric
+        if self._spatial_samples > 0:
+            per_token["NumberOfSpatialSamples"] = str(int(self._spatial_samples))
+        exact: list[tuple[str, str]] = []
+        for entry in self._parameter_overrides:
+            key, sep, value = entry.partition("=")
+            if not sep or not key.strip():
+                raise ValueError(f"Invalid parameter_overrides entry '{entry}': expected 'Key=value text'.")
+            exact.append((key.strip(), value.strip()))
+        return per_token, exact
+
+    @staticmethod
+    def _apply_map_overrides(
+        text: str, per_token: dict[str, str], exact: list[tuple[str, str]], device_index: int
+    ) -> str:
+        """Patch a parameter map's text: set ImpactGPU to the device, apply exact key overrides, replace each
+        token of a per-token knob (preserving multiplicity), and warn for a requested key absent from the map.
+        """
+        entry_pattern = re.compile(r"^(\s*)\((\S+)((?:\s+[^)]*)?)\)\s*$")
+        requested = set(per_token) | {key for key, _ in exact}
+        seen: set[str] = set()
+        lines = []
+        for line in text.splitlines():
+            match = entry_pattern.match(line)
+            if match:
+                indent, key, values = match.group(1), match.group(2), match.group(3)
+                if key == "ImpactGPU":
+                    line = f"{indent}(ImpactGPU {device_index})"
+                else:
+                    exact_value = next((value for k, value in exact if k == key), None)
+                    if exact_value is not None:
+                        seen.add(key)
+                        line = f"{indent}({key} {exact_value})"
+                    else:
+                        token_key = "ImpactSubsetFeatures" if key.startswith("ImpactSubsetFeatures") else key
+                        if token_key in per_token:
+                            seen.add(token_key)
+                            replaced = " ".join(per_token[token_key] for _ in values.split())
+                            line = f"{indent}({key} {replaced})"
+            lines.append(line)
+        # Overrides never inject keys, so a knob set for a key absent from every map would silently do
+        # nothing — surface it (e.g. final_grid_spacing on a rigid-only preset).
+        for key in sorted(requested - seen):
+            print(f"[ImpactReg] note: override '{key}' matched no entry in the preset's parameter maps.")
+        return "\n".join(lines)
+
     def _stage_parameter_maps(self, work: Path, device_index: int) -> list[Path]:
-        """Copy parameter maps into the work dir, patching the ImpactGPU device field when present."""
+        """Stage the parameter maps into the work dir.
+
+        Matrix mode GENERATES each map from ``resolutions`` + the registry, then applies only the map-wide
+        knobs (grid spacing, spatial samples, exact overrides) — the matrix already sets iterations and
+        features per cell. Legacy mode copies the preset's maps and applies every per-token / exact override.
+        Both set the ImpactGPU device.
+        """
         staged = []
         for src in self._parameter_maps:
+            if self._resolutions:
+                text = generate_impact_parameter_map(
+                    src.read_text(encoding="utf-8"), self._resolutions, self._registry, self._mode
+                )
+                per_token, exact = self._parameter_map_overrides(global_only=True)
+            else:
+                text = src.read_text(encoding="utf-8")
+                per_token, exact = self._parameter_map_overrides()
+            text = self._apply_map_overrides(text, per_token, exact, device_index)
             dst = work / src.name
-            lines = []
-            for line in src.read_text(encoding="utf-8").splitlines():
-                if line.strip().startswith("(ImpactGPU"):
-                    line = f"(ImpactGPU {device_index})"
-                lines.append(line)
-            dst.write_text("\n".join(lines) + "\n", encoding="utf-8")
+            dst.write_text(text if text.endswith("\n") else text + "\n", encoding="utf-8")
             staged.append(dst)
         return staged
 
@@ -161,8 +484,10 @@ class ElastixEngine:
             captured: list[str] = []
             iteration_line = re.compile(r"^\d+\s")
             # ``iterations`` is the total iteration budget declared for the preset (summed over the
-            # chained parameter maps), so the bar spans the whole chain of registration stages.
-            progress = tqdm.tqdm(total=self._iterations or None, desc="Registration", ncols=0, leave=True)
+            # chained parameter maps), so the bar spans the whole chain of registration stages. A tuned
+            # ``max_iterations`` makes that declared budget stale — fall back to an open-ended bar.
+            budget = None if self._max_iterations > 0 else (self._iterations or None)
+            progress = tqdm.tqdm(total=budget, desc="Registration", ncols=0, leave=True)
             assert proc.stdout is not None
             resolution = 0
             for line in proc.stdout:
@@ -226,11 +551,33 @@ class ElastixRegistration(torch.nn.Module):
 
     accepts_attributes = True
 
-    def __init__(self, engine: str, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        engine: str,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         super().__init__()
         if engine != "elastix":
             raise NotImplementedError(f"ElastixRegistration engine '{engine}' is not implemented yet.")
-        self._engine = ElastixEngine(parameter_maps, models, iterations)
+        self._engine = ElastixEngine(
+            parameter_maps,
+            max_iterations,
+            final_grid_spacing,
+            subset_features,
+            spatial_samples,
+            parameter_overrides,
+            resolutions,
+            models_registry,
+            mode,
+        )
 
     def forward(
         self,
@@ -290,9 +637,23 @@ class RegistrationNet(network.Network):
         outputs_criterions: dict[str, network.TargetCriterionsLoader] = {"default": network.TargetCriterionsLoader()},
         engine: str = "elastix",
         parameter_maps: list[str] = [],
-        models: list[str] = [],
-        iterations: int = 1000,
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict[str, ResolutionSpec] = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
     ) -> None:
+        # The registration is fully described by the per-resolution model matrix ``resolutions`` (config =
+        # source of truth): each resolution lists its models, each model self-configured (ref, voxel_size,
+        # layers_mask, layers_weight, subset_features, pca, distance); intrinsic per-model props come from
+        # ``models_registry``. The feature-model download list is DERIVED from the matrix (no flat ``models``).
+        # Global knobs override the generated map: final_grid_spacing -> FinalGridSpacingInPhysicalUnits (mm),
+        # spatial_samples -> NumberOfSpatialSamples, parameter_overrides ('Key=value') -> any other entry.
+        # An empty ``resolutions`` = an intensity-only preset (no IMPACT models): the fixed maps are staged
+        # with just the global overrides. The total iteration count is derived (sum of per-resolution budgets).
         super().__init__(
             in_channels=1,
             optimizer=optimizer,
@@ -302,7 +663,18 @@ class RegistrationNet(network.Network):
         )
         self.add_module(
             "Registration",
-            ElastixRegistration(engine, parameter_maps, models, iterations),
+            ElastixRegistration(
+                engine,
+                parameter_maps,
+                max_iterations,
+                final_grid_spacing,
+                subset_features,
+                spatial_samples,
+                parameter_overrides,
+                resolutions,
+                models_registry,
+                mode,
+            ),
             in_branch=[0, 1, 2, 3],
             out_branch=["registration"],
         )

CBCT_CT_HeadNeck/ParameterMap_CBCT_HN.txt

@@ -1,4 +1,4 @@
-(MaximumNumberOfIterations 300  300 200 200 150)
+(MaximumNumberOfIterations 300 300 200 200 150)
 (NumberOfSpatialSamples 2000)
 (Transform "RecursiveBSplineTransform")
 (NumberOfResolutions 5)
@@ -25,7 +25,7 @@
 (ImpactModelsPath1 "TS/M731.pt")
 (ImpactDimension1 3)
 (ImpactNumberOfChannels1 1)
-(ImpactPatchSize1 0 0 0)  
+(ImpactPatchSize1 0 0 0)
 (ImpactVoxelSize1 3 3 3)
 (ImpactLayersMask1 "01")
 (ImpactSubsetFeatures1 64)
@@ -149,4 +149,4 @@
 (ResultImageFormat "mha")
 
 (ITKTransformOutputFileNameExtension "itk.txt")
-(WriteITKCompositeTransform "true")
+(WriteITKCompositeTransform "true")

CBCT_CT_HeadNeck/Prediction.yml

@@ -5,12 +5,90 @@ Predictor:
       engine: elastix
       parameter_maps:
       - ParameterMap_CBCT_HN.txt
-      models:
-      - VBoussot/impact-torchscript-models:TS/M732.pt
-      - VBoussot/impact-torchscript-models:TS/M731.pt
-      - VBoussot/impact-torchscript-models:TS/M730.pt
-      iterations: 1150
       outputs_criterions: None
+      max_iterations: 0
+      final_grid_spacing: 0.0
+      subset_features: 0
+      spatial_samples: 0
+      parameter_overrides: []
+      resolutions:
+        '0':
+          max_iterations: 300
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:TS/M732.pt
+              voxel_size:
+              - 6.0
+              - 6.0
+              - 6.0
+              layers_mask: '01'
+              layers_weight:
+              - 1.0
+              subset_features: 64
+              pca: 0
+              distance: L1
+        '1':
+          max_iterations: 300
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:TS/M731.pt
+              voxel_size:
+              - 3.0
+              - 3.0
+              - 3.0
+              layers_mask: '01'
+              layers_weight:
+              - 1.0
+              subset_features: 64
+              pca: 0
+              distance: L1
+        '2':
+          max_iterations: 200
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:TS/M731.pt
+              voxel_size:
+              - 3.0
+              - 3.0
+              - 3.0
+              layers_mask: '01'
+              layers_weight:
+              - 1.0
+              subset_features: 64
+              pca: 0
+              distance: L1
+        '3':
+          max_iterations: 200
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:TS/M730.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 3.0
+              layers_mask: '01'
+              layers_weight:
+              - 1.0
+              subset_features: 64
+              pca: 0
+              distance: L1
+        '4':
+          max_iterations: 150
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:TS/M730.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 3.0
+              layers_mask: '01'
+              layers_weight:
+              - 1.0
+              subset_features: 64
+              pca: 0
+              distance: L1
+      models_registry: VBoussot/impact-torchscript-models:models.json
+      mode: Static
   Dataset:
     groups_src:
       Volume_0:

CBCT_CT_MRSeg/Model.py

@@ -32,6 +32,7 @@ NOTE: do NOT add ``from __future__ import annotations`` here — KonfAI's config
 runtime-evaluated annotations (``get_origin``); PEP 563 stringized annotations break arg resolution.
 """
 
+import json
 import os
 import re
 import shutil
@@ -52,6 +53,212 @@ from konfai.utils.dataset import Attribute, data_to_image, image_to_data
 # Set KONFAI_ELASTIX_DIR to point at an existing install and skip the download.
 ELASTIX_CACHE = Path.home() / ".cache" / "konfai" / "elastix-impact"
 
+# ---------------------------------------------------------------------------------------------------
+# Per-resolution model matrix (the config is the source of truth) -> generated IMPACT parameter map.
+# IMPACT field docs: https://github.com/vboussot/ImpactLoss/tree/main/ParameterMaps
+# The forced per-model props (dimension/channels/FOV formula) live in a registry (models.json on
+# VBoussot/impact-torchscript-models); the config carries the FREE knobs (which models per resolution,
+# feature voxel size, iterations, per-model layer weights/mask/subset/pca/distance) and the global
+# ``mode``. PatchSize follows ImpactMode: Static -> "0 0 0" (whole image); Jacobian -> the model FOV
+# evaluated from the registry formula (MIND 2*r*d+1, TS/MRSeg 2^l+3, SAM 29, DINOv2 14) as a cube.
+# ---------------------------------------------------------------------------------------------------
+
+_IMPACT_MODELS_REGISTRY = "VBoussot/impact-torchscript-models:models.json"
+
+# ``2^l+3`` grows with depth but the segmenters' receptive field plateaus: layers 7-8 share layer 6's
+# FOV (the "ramp max"). A config that deep should really run in Static (whole image) anyway; in Jacobian
+# we clamp ``l`` to this plateau so the patch stays finite and matches the real FOV.
+_FOV_RAMP_MAX_LAYER = 6
+
+
+def _num(x: object) -> str:
+    """Format a number the elastix way: integers without a trailing '.0' (6.0 -> '6', 0.2 -> '0.2')."""
+    return "%g" % float(x)
+
+
+class ModelSpec:
+    """One feature model at one resolution, with its OWN config (several models may share a resolution).
+
+    ``ref`` selects the model; ``voxel_size`` / ``layers_weight`` / ``subset_features`` / ``pca`` /
+    ``distance`` are its free per-(resolution, model) tuning knobs (the doc's per-model *tuning* fields).
+    The intrinsic per-model props — dimension, channels, ``layers_mask``, patch-size (FOV) — come from the
+    registry (read-only); ``layers_mask`` / ``distance`` left empty fall back to the registry default.
+    """
+
+    def __init__(
+        self,
+        ref: str,
+        voxel_size: list[float] = [],
+        layers_weight: list[float] = [1.0],
+        subset_features: int = 0,
+        pca: int = 0,
+        distance: str = "",
+        layers_mask: str = "",
+    ) -> None:
+        self.ref = ref
+        self.voxel_size = voxel_size
+        self.layers_weight = layers_weight
+        self.subset_features = subset_features
+        self.pca = pca
+        self.distance = distance
+        self.layers_mask = layers_mask
+
+
+class ResolutionSpec:
+    """One elastix resolution level: its iteration budget and the models compared there (each self-configured)."""
+
+    def __init__(self, max_iterations: int, models: dict[str, ModelSpec]) -> None:
+        self.max_iterations = max_iterations
+        self.models = models
+
+
+def _sorted_specs(mapping: dict) -> list:
+    """dict keyed by string indices ('0','1',...) -> values in numeric order (well-defined res/model order)."""
+    return [mapping[k] for k in sorted(mapping, key=lambda key: int(key))]
+
+
+def load_models_registry(ref: str = _IMPACT_MODELS_REGISTRY) -> dict:
+    """Load models.json (forced params per model) from the model repo on Hugging Face.
+
+    The registry is NOT bundled with the preset — it lives on the models repo and is fetched from there.
+    Resolution: the ``KONFAI_IMPACT_MODELS_REGISTRY`` env path wins (dev/offline); otherwise ``ref`` must be
+    a ``repo:file`` Hugging Face reference.
+    """
+    local = os.environ.get("KONFAI_IMPACT_MODELS_REGISTRY", "")
+    if local:
+        path = Path(local)
+    elif ":" in ref:
+        repo, filename = ref.split(":", 1)
+        path = Path(hf_hub_download(repo_id=repo, filename=filename, repo_type="model"))  # nosec B615
+    else:
+        raise ValueError(
+            f"models_registry '{ref}' must be a 'repo:file' Hugging Face reference (the registry is fetched "
+            f"from HF, not bundled) — or set KONFAI_IMPACT_MODELS_REGISTRY to a local file for offline use."
+        )
+    return json.loads(path.read_text(encoding="utf-8"))
+
+
+def _model_key(ref: str) -> str:
+    """Registry key / staged relative path = the model file within the models repo (strip a 'repo:' prefix)."""
+    return ref.split(":", 1)[1] if ":" in ref else ref
+
+
+def _deepest_active_layer(layers_mask: str) -> int:
+    """Deepest (largest-FOV) layer selected by ``layers_mask``, as a 0-based index read left-to-right.
+
+    A model returns its feature layers shallow->deep (``[layer_0, layer_1, ...]``, see the model repo's
+    build scripts); ``layers_mask`` has one char per returned layer, position ``i`` == ``layer_i``, ``'1'``
+    = selected. In Jacobian the patch must cover the receptive field of the DEEPEST selected layer, so the
+    FOV is governed by the rightmost ``'1'``.
+    """
+    mask = layers_mask.strip().strip('"')
+    active = [i for i, char in enumerate(mask) if char == "1"]
+    if not active:
+        raise ValueError(f"LayersMask '{layers_mask}' selects no layer; cannot derive the model FOV.")
+    return max(active)
+
+
+def _fov_value(fov: dict, layers_mask: str) -> int:
+    """Evaluate a model's field-of-view (in voxels) from its registry ``fov`` spec.
+
+    Supported formulas (from the model repo, https://huggingface.co/VBoussot/impact-torchscript-models):
+      ``2*r*d+1``  MIND, from the handcrafted radius ``r`` / dilation ``d`` (e.g. R1D2 -> 5);
+      ``2^l+3``    TotalSegmentator / MRSegmentator, ``l`` = the deepest layer picked by ``layers_mask``,
+                   clamped to the receptive-field plateau ``_FOV_RAMP_MAX_LAYER`` (layers 7-8 -> layer 6);
+      a bare int  a fixed FOV (SAM2.1 -> 29, DINOv2 -> 14);
+      ``Global``   Anatomix — whole-image only (Static); has no finite Jacobian patch -> error.
+    An explicit ``value`` in the spec is honoured as a precomputed shortcut when the formula needs none.
+    """
+    formula = str(fov.get("formula", "")).strip()
+    key = re.sub(r"\s+", "", formula).lower()
+    if key.isdigit():
+        return int(key)
+    if key == "2*r*d+1":
+        return 2 * int(fov["r"]) * int(fov["d"]) + 1
+    if key == "2^l+3":
+        return 2 ** min(_deepest_active_layer(layers_mask), _FOV_RAMP_MAX_LAYER) + 3
+    if "global" in key:
+        raise ValueError(f"model FOV '{formula}' is whole-image only (Static); it has no Jacobian patch size.")
+    if fov.get("value") is not None:
+        return int(fov["value"])
+    raise ValueError(f"cannot evaluate model FOV formula '{formula}'.")
+
+
+def _patch_size(mode: str, entry: dict, layers_mask: str) -> str:
+    """PatchSize from the model FOV, one token per model axis (2D model -> 2 tokens, 3D -> 3): Static ->
+    whole image (all zeros); Jacobian -> the evaluated FOV repeated over the axes. A 2D model mixed with a
+    3D one at a resolution concatenates as e.g. ``29 29 11 11 11`` (SAM 2D + TS 3D), matching IMPACT."""
+    dim = int(entry.get("dimension", 3))
+    if mode.strip().strip('"').lower() != "jacobian":
+        return " ".join(["0"] * dim)
+    fov = _fov_value(entry.get("fov", {}), layers_mask)
+    return " ".join([str(fov)] * dim)
+
+
+def generate_impact_parameter_map(
+    template_text: str, resolutions: dict, registry: dict, mode: str = "Static"
+) -> str:
+    """Rewrite the resolution-dependent lines of ``template_text`` from the model matrix ``resolutions``.
+
+    Regenerated: MaximumNumberOfIterations, NumberOfResolutions, Fixed/MovingImagePyramidRescaleSchedule,
+    ImpactMode (from the config ``mode``), and the whole ImpactXxxK block; every other template line is
+    kept verbatim (optimizer, transform, metric weights, components...). N (number of resolutions) is
+    deduced from the config. ``mode`` drives PatchSize: Static -> ``0 0 0`` (whole image); Jacobian -> the
+    per-model FOV evaluated from the registry formula and the cell's ``layers_mask``.
+    """
+    res = _sorted_specs(resolutions)
+    n = len(res)
+    mode_clean = mode.strip().strip('"') or "Static"
+
+    impact: list[str] = []
+    for k, r in enumerate(res):
+        models = _sorted_specs(r.models)
+        entries = [registry[_model_key(m.ref)] for m in models]
+
+        def row(stem: str, values: list[str]) -> None:
+            impact.append(f"(Impact{stem}{k} " + " ".join(values) + ")")
+
+        # From the registry (models.json on the model repo) ONLY the 3 truly model-fixed props:
+        # Dimension, NumberOfChannels, PatchSize (the model FOV). Everything else is a per-model tuning knob
+        # taken straight from the cell: VoxelSize / LayersMask / SubsetFeatures / PCA / Distance / LayersWeight.
+        row("ModelsPath", [f'"{_model_key(m.ref)}"' for m in models])
+        row("Dimension", [e["dimension"] for e in entries])
+        row("NumberOfChannels", [e["numberofchannels"] for e in entries])
+        row("PatchSize", [_patch_size(mode_clean, e, m.layers_mask) for e, m in zip(entries, models)])
+        row("VoxelSize", [" ".join(_num(v) for v in m.voxel_size) for m in models])
+        row("LayersMask", [f'"{m.layers_mask}"' for m in models])
+        row("SubsetFeatures", [str(m.subset_features) for m in models])
+        row("PCA", [str(m.pca) for m in models])
+        row("Distance", [f'"{m.distance}"' for m in models])
+        row("LayersWeight", [" ".join(_num(w) for w in m.layers_weight) for m in models])
+        impact.append("")  # blank line between resolutions, mirroring the reference maps
+
+    # The per-resolution block is the contiguous span from the first to the last ``Impact<name><k>`` line
+    # (the blank lines the reference maps put BETWEEN resolutions fall inside that span). Replace the whole
+    # span in one shot with the generated block, so the reference blanks are not kept on top of ours.
+    lines = template_text.splitlines()
+    indexed = [(re.match(r"^\s*\((\S+?)\s+(.*?)\)\s*$", ln), ln) for ln in lines]
+    block_rows = [i for i, (m, _) in enumerate(indexed) if m and re.match(r"^Impact[A-Za-z]+\d+$", m.group(1))]
+    block_lo, block_hi = (block_rows[0], block_rows[-1]) if block_rows else (-1, -2)
+
+    out: list[str] = []
+    for i, (m, line) in enumerate(indexed):
+        key = m.group(1) if m else None
+        if block_lo <= i <= block_hi:
+            if i == block_lo:  # replace the whole span at its first line, drop the rest (incl. inner blanks)
+                out.extend(impact[:-1])
+        elif key == "MaximumNumberOfIterations":
+            out.append("(MaximumNumberOfIterations " + " ".join(_num(r.max_iterations) for r in res) + ")")
+        elif key == "NumberOfResolutions":
+            out.append(f"(NumberOfResolutions {n})")
+        elif key in ("FixedImagePyramidRescaleSchedule", "MovingImagePyramidRescaleSchedule"):
+            out.append(f"({key} " + " ".join(["1"] * 3 * n) + ")")
+        elif key == "ImpactMode":
+            out.append(f'(ImpactMode "{mode_clean}")')
+        else:
+            out.append(line)
+    return "\n".join(out)
+
 
 class ElastixEngine:
     """Run the elastix-IMPACT binary on a fixed/moving pair; return (moved, dvf) on the fixed grid.
@@ -60,14 +267,57 @@ class ElastixEngine:
     does NOT ship it), so registration is a subprocess call, not ``sitk.ElastixImageFilter``.
     """
 
-    def __init__(self, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         self._bundle_dir = Path(__file__).resolve().parent
         self._parameter_maps = [self._bundle_dir / Path(p).name for p in parameter_maps]
+        self._max_iterations = max_iterations
+        self._final_grid_spacing = final_grid_spacing
+        self._subset_features = subset_features
+        self._spatial_samples = spatial_samples
+        self._parameter_overrides = list(parameter_overrides)
+        # ImpactMode: Static computes features once per level (PatchSize 0 0 0 = whole image); Jacobian
+        # samples random patches sized to the model FOV each iteration. Global knob: one mode per preset.
+        self._mode = mode
+        # Matrix mode: when `resolutions` is given the parameter map is GENERATED from it (the config is the
+        # source of truth). An empty `resolutions` = an intensity preset (no IMPACT feature models): the fixed
+        # parameter maps are staged with only the global knob overrides.
+        self._resolutions = resolutions
+        self._registry = load_models_registry(models_registry) if resolutions else {}
+        # The feature models are DERIVED — the unique refs across the matrix cells (no flat `models` param).
+        models: list[str] = []
+        for res in _sorted_specs(resolutions):
+            for model in _sorted_specs(res.models):
+                if model.ref not in models:
+                    models.append(model.ref)
         self._models = models
-        self._iterations = iterations
+        # `iterations` (the progress-bar total) is NOT a config parameter — it is DERIVED: the sum of the
+        # per-resolution iteration budgets, read from the matrix (matrix mode) or the maps (legacy).
+        self._iterations = self._total_iterations()
         self._elastix_bin = self._ensure_binary()
         self._local_models = self._download_models()
 
+    def _total_iterations(self) -> int:
+        """Total iterations across all resolutions — the progress-bar budget, derived from the config."""
+        if self._resolutions:
+            return sum(int(res.max_iterations) for res in _sorted_specs(self._resolutions))
+        total = 0
+        for src in self._parameter_maps:
+            match = re.search(r"\(MaximumNumberOfIterations\s+([^)]*)\)", src.read_text(encoding="utf-8"))
+            if match:
+                total += sum(int(token) for token in match.group(1).split())
+        return total
+
     def _ensure_binary(self) -> Path:
         # Optional override: point at an existing elastix-IMPACT install (skips the download).
         override = os.environ.get("KONFAI_ELASTIX_DIR", "")
@@ -91,17 +341,90 @@ class ElastixEngine:
             models.append((filename, local))
         return models
 
+    def _parameter_map_overrides(self, global_only: bool = False) -> tuple[dict[str, str], list[tuple[str, str]]]:
+        """The tuned knobs as parameter-map overrides: ``(per_token, exact)``.
+
+        ``per_token`` maps an elastix key (or the ``ImpactSubsetFeatures`` prefix) to a value that replaces
+        **each** existing token, so per-resolution / per-model multiplicity is preserved (e.g.
+        ``(MaximumNumberOfIterations 500 250)`` -> ``(MaximumNumberOfIterations 300 300)``). ``exact``
+        entries (from ``parameter_overrides``, ``Key=value text``) replace the whole value verbatim and win
+        over the named knobs. Overrides only REPLACE keys already present in a map — never inject new ones.
+        ``global_only`` (matrix mode) keeps just the map-wide knobs and drops ``max_iterations`` /
+        ``subset_features`` — the per-resolution matrix already sets those per cell.
+        """
+        per_token: dict[str, str] = {}
+        if not global_only and self._max_iterations > 0:
+            per_token["MaximumNumberOfIterations"] = str(int(self._max_iterations))
+        if self._final_grid_spacing > 0:
+            per_token["FinalGridSpacingInPhysicalUnits"] = str(float(self._final_grid_spacing))
+        if not global_only and self._subset_features > 0:
+            per_token["ImpactSubsetFeatures"] = str(int(self._subset_features))  # prefix: indexed per metric
+        if self._spatial_samples > 0:
+            per_token["NumberOfSpatialSamples"] = str(int(self._spatial_samples))
+        exact: list[tuple[str, str]] = []
+        for entry in self._parameter_overrides:
+            key, sep, value = entry.partition("=")
+            if not sep or not key.strip():
+                raise ValueError(f"Invalid parameter_overrides entry '{entry}': expected 'Key=value text'.")
+            exact.append((key.strip(), value.strip()))
+        return per_token, exact
+
+    @staticmethod
+    def _apply_map_overrides(
+        text: str, per_token: dict[str, str], exact: list[tuple[str, str]], device_index: int
+    ) -> str:
+        """Patch a parameter map's text: set ImpactGPU to the device, apply exact key overrides, replace each
+        token of a per-token knob (preserving multiplicity), and warn for a requested key absent from the map.
+        """
+        entry_pattern = re.compile(r"^(\s*)\((\S+)((?:\s+[^)]*)?)\)\s*$")
+        requested = set(per_token) | {key for key, _ in exact}
+        seen: set[str] = set()
+        lines = []
+        for line in text.splitlines():
+            match = entry_pattern.match(line)
+            if match:
+                indent, key, values = match.group(1), match.group(2), match.group(3)
+                if key == "ImpactGPU":
+                    line = f"{indent}(ImpactGPU {device_index})"
+                else:
+                    exact_value = next((value for k, value in exact if k == key), None)
+                    if exact_value is not None:
+                        seen.add(key)
+                        line = f"{indent}({key} {exact_value})"
+                    else:
+                        token_key = "ImpactSubsetFeatures" if key.startswith("ImpactSubsetFeatures") else key
+                        if token_key in per_token:
+                            seen.add(token_key)
+                            replaced = " ".join(per_token[token_key] for _ in values.split())
+                            line = f"{indent}({key} {replaced})"
+            lines.append(line)
+        # Overrides never inject keys, so a knob set for a key absent from every map would silently do
+        # nothing — surface it (e.g. final_grid_spacing on a rigid-only preset).
+        for key in sorted(requested - seen):
+            print(f"[ImpactReg] note: override '{key}' matched no entry in the preset's parameter maps.")
+        return "\n".join(lines)
+
     def _stage_parameter_maps(self, work: Path, device_index: int) -> list[Path]:
-        """Copy parameter maps into the work dir, patching the ImpactGPU device field when present."""
+        """Stage the parameter maps into the work dir.
+
+        Matrix mode GENERATES each map from ``resolutions`` + the registry, then applies only the map-wide
+        knobs (grid spacing, spatial samples, exact overrides) — the matrix already sets iterations and
+        features per cell. Legacy mode copies the preset's maps and applies every per-token / exact override.
+        Both set the ImpactGPU device.
+        """
         staged = []
         for src in self._parameter_maps:
+            if self._resolutions:
+                text = generate_impact_parameter_map(
+                    src.read_text(encoding="utf-8"), self._resolutions, self._registry, self._mode
+                )
+                per_token, exact = self._parameter_map_overrides(global_only=True)
+            else:
+                text = src.read_text(encoding="utf-8")
+                per_token, exact = self._parameter_map_overrides()
+            text = self._apply_map_overrides(text, per_token, exact, device_index)
             dst = work / src.name
-            lines = []
-            for line in src.read_text(encoding="utf-8").splitlines():
-                if line.strip().startswith("(ImpactGPU"):
-                    line = f"(ImpactGPU {device_index})"
-                lines.append(line)
-            dst.write_text("\n".join(lines) + "\n", encoding="utf-8")
+            dst.write_text(text if text.endswith("\n") else text + "\n", encoding="utf-8")
             staged.append(dst)
         return staged
 
@@ -161,8 +484,10 @@ class ElastixEngine:
             captured: list[str] = []
             iteration_line = re.compile(r"^\d+\s")
             # ``iterations`` is the total iteration budget declared for the preset (summed over the
-            # chained parameter maps), so the bar spans the whole chain of registration stages.
-            progress = tqdm.tqdm(total=self._iterations or None, desc="Registration", ncols=0, leave=True)
+            # chained parameter maps), so the bar spans the whole chain of registration stages. A tuned
+            # ``max_iterations`` makes that declared budget stale — fall back to an open-ended bar.
+            budget = None if self._max_iterations > 0 else (self._iterations or None)
+            progress = tqdm.tqdm(total=budget, desc="Registration", ncols=0, leave=True)
             assert proc.stdout is not None
             resolution = 0
             for line in proc.stdout:
@@ -226,11 +551,33 @@ class ElastixRegistration(torch.nn.Module):
 
     accepts_attributes = True
 
-    def __init__(self, engine: str, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        engine: str,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         super().__init__()
         if engine != "elastix":
             raise NotImplementedError(f"ElastixRegistration engine '{engine}' is not implemented yet.")
-        self._engine = ElastixEngine(parameter_maps, models, iterations)
+        self._engine = ElastixEngine(
+            parameter_maps,
+            max_iterations,
+            final_grid_spacing,
+            subset_features,
+            spatial_samples,
+            parameter_overrides,
+            resolutions,
+            models_registry,
+            mode,
+        )
 
     def forward(
         self,
@@ -290,9 +637,23 @@ class RegistrationNet(network.Network):
         outputs_criterions: dict[str, network.TargetCriterionsLoader] = {"default": network.TargetCriterionsLoader()},
         engine: str = "elastix",
         parameter_maps: list[str] = [],
-        models: list[str] = [],
-        iterations: int = 1000,
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict[str, ResolutionSpec] = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
     ) -> None:
+        # The registration is fully described by the per-resolution model matrix ``resolutions`` (config =
+        # source of truth): each resolution lists its models, each model self-configured (ref, voxel_size,
+        # layers_mask, layers_weight, subset_features, pca, distance); intrinsic per-model props come from
+        # ``models_registry``. The feature-model download list is DERIVED from the matrix (no flat ``models``).
+        # Global knobs override the generated map: final_grid_spacing -> FinalGridSpacingInPhysicalUnits (mm),
+        # spatial_samples -> NumberOfSpatialSamples, parameter_overrides ('Key=value') -> any other entry.
+        # An empty ``resolutions`` = an intensity-only preset (no IMPACT models): the fixed maps are staged
+        # with just the global overrides. The total iteration count is derived (sum of per-resolution budgets).
         super().__init__(
             in_channels=1,
             optimizer=optimizer,
@@ -302,7 +663,18 @@ class RegistrationNet(network.Network):
         )
         self.add_module(
             "Registration",
-            ElastixRegistration(engine, parameter_maps, models, iterations),
+            ElastixRegistration(
+                engine,
+                parameter_maps,
+                max_iterations,
+                final_grid_spacing,
+                subset_features,
+                spatial_samples,
+                parameter_overrides,
+                resolutions,
+                models_registry,
+                mode,
+            ),
             in_branch=[0, 1, 2, 3],
             out_branch=["registration"],
         )

CBCT_CT_MRSeg/ParameterMap_CBCT_generic_MRSeg.txt

@@ -23,7 +23,7 @@
 (ImpactModelsPath1 "MIND/R1D2_3D.pt" "MRSeg/MRSeg.pt")
 (ImpactDimension1 3 3)
 (ImpactNumberOfChannels1 1 1)
-(ImpactPatchSize1 0 0 0 0 0 0)  
+(ImpactPatchSize1 0 0 0 0 0 0)
 (ImpactVoxelSize1 3 3 3 3 3 3)
 (ImpactLayersMask1 "1" "1")
 (ImpactSubsetFeatures1 32 64)
@@ -36,7 +36,7 @@
 (ImpactNumberOfChannels2 1 1)
 (ImpactPatchSize2 0 0 0 0 0 0)
 (ImpactVoxelSize2 2 2 2 2 2 2)
-(ImpactLayersMask2 "1"  "1")
+(ImpactLayersMask2 "1" "1")
 (ImpactSubsetFeatures2 32 64)
 (ImpactPCA2 0 0)
 (ImpactDistance2 "L1" "Dice")
@@ -47,7 +47,7 @@
 (ImpactNumberOfChannels3 1 1)
 (ImpactPatchSize3 0 0 0 0 0 0)
 (ImpactVoxelSize3 2 2 2 2 2 2)
-(ImpactLayersMask3 "1"  "1")
+(ImpactLayersMask3 "1" "1")
 (ImpactSubsetFeatures3 32 64)
 (ImpactPCA3 0 0)
 (ImpactDistance3 "L1" "Dice")
@@ -134,4 +134,4 @@
 (ResultImageFormat "mha")
 
 (ITKTransformOutputFileNameExtension "itk.txt")
-(WriteITKCompositeTransform "true")
+(WriteITKCompositeTransform "true")

CBCT_CT_MRSeg/Prediction.yml

@@ -5,11 +5,123 @@ Predictor:
       engine: elastix
       parameter_maps:
       - ParameterMap_CBCT_generic_MRSeg.txt
-      models:
-      - VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
-      - VBoussot/impact-torchscript-models:MRSeg/MRSeg.pt
-      iterations: 1100
       outputs_criterions: None
+      max_iterations: 0
+      final_grid_spacing: 0.0
+      subset_features: 0
+      spatial_samples: 0
+      parameter_overrides: []
+      resolutions:
+        '0':
+          max_iterations: 400
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 6.0
+              - 6.0
+              - 6.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.2
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:MRSeg/MRSeg.pt
+              voxel_size:
+              - 6.0
+              - 6.0
+              - 6.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.8
+              subset_features: 64
+              pca: 0
+              distance: Dice
+        '1':
+          max_iterations: 300
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 3.0
+              - 3.0
+              - 3.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.3
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:MRSeg/MRSeg.pt
+              voxel_size:
+              - 3.0
+              - 3.0
+              - 3.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.7
+              subset_features: 64
+              pca: 0
+              distance: Dice
+        '2':
+          max_iterations: 200
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.6
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:MRSeg/MRSeg.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.4
+              subset_features: 64
+              pca: 0
+              distance: Dice
+        '3':
+          max_iterations: 200
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.7
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:MRSeg/MRSeg.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.3
+              subset_features: 64
+              pca: 0
+              distance: Dice
+      models_registry: VBoussot/impact-torchscript-models:models.json
+      mode: Static
   Dataset:
     groups_src:
       Volume_0:

CBCT_CT_TS/Model.py

@@ -32,6 +32,7 @@ NOTE: do NOT add ``from __future__ import annotations`` here — KonfAI's config
 runtime-evaluated annotations (``get_origin``); PEP 563 stringized annotations break arg resolution.
 """
 
+import json
 import os
 import re
 import shutil
@@ -52,6 +53,212 @@ from konfai.utils.dataset import Attribute, data_to_image, image_to_data
 # Set KONFAI_ELASTIX_DIR to point at an existing install and skip the download.
 ELASTIX_CACHE = Path.home() / ".cache" / "konfai" / "elastix-impact"
 
+# ---------------------------------------------------------------------------------------------------
+# Per-resolution model matrix (the config is the source of truth) -> generated IMPACT parameter map.
+# IMPACT field docs: https://github.com/vboussot/ImpactLoss/tree/main/ParameterMaps
+# The forced per-model props (dimension/channels/FOV formula) live in a registry (models.json on
+# VBoussot/impact-torchscript-models); the config carries the FREE knobs (which models per resolution,
+# feature voxel size, iterations, per-model layer weights/mask/subset/pca/distance) and the global
+# ``mode``. PatchSize follows ImpactMode: Static -> "0 0 0" (whole image); Jacobian -> the model FOV
+# evaluated from the registry formula (MIND 2*r*d+1, TS/MRSeg 2^l+3, SAM 29, DINOv2 14) as a cube.
+# ---------------------------------------------------------------------------------------------------
+
+_IMPACT_MODELS_REGISTRY = "VBoussot/impact-torchscript-models:models.json"
+
+# ``2^l+3`` grows with depth but the segmenters' receptive field plateaus: layers 7-8 share layer 6's
+# FOV (the "ramp max"). A config that deep should really run in Static (whole image) anyway; in Jacobian
+# we clamp ``l`` to this plateau so the patch stays finite and matches the real FOV.
+_FOV_RAMP_MAX_LAYER = 6
+
+
+def _num(x: object) -> str:
+    """Format a number the elastix way: integers without a trailing '.0' (6.0 -> '6', 0.2 -> '0.2')."""
+    return "%g" % float(x)
+
+
+class ModelSpec:
+    """One feature model at one resolution, with its OWN config (several models may share a resolution).
+
+    ``ref`` selects the model; ``voxel_size`` / ``layers_weight`` / ``subset_features`` / ``pca`` /
+    ``distance`` are its free per-(resolution, model) tuning knobs (the doc's per-model *tuning* fields).
+    The intrinsic per-model props — dimension, channels, ``layers_mask``, patch-size (FOV) — come from the
+    registry (read-only); ``layers_mask`` / ``distance`` left empty fall back to the registry default.
+    """
+
+    def __init__(
+        self,
+        ref: str,
+        voxel_size: list[float] = [],
+        layers_weight: list[float] = [1.0],
+        subset_features: int = 0,
+        pca: int = 0,
+        distance: str = "",
+        layers_mask: str = "",
+    ) -> None:
+        self.ref = ref
+        self.voxel_size = voxel_size
+        self.layers_weight = layers_weight
+        self.subset_features = subset_features
+        self.pca = pca
+        self.distance = distance
+        self.layers_mask = layers_mask
+
+
+class ResolutionSpec:
+    """One elastix resolution level: its iteration budget and the models compared there (each self-configured)."""
+
+    def __init__(self, max_iterations: int, models: dict[str, ModelSpec]) -> None:
+        self.max_iterations = max_iterations
+        self.models = models
+
+
+def _sorted_specs(mapping: dict) -> list:
+    """dict keyed by string indices ('0','1',...) -> values in numeric order (well-defined res/model order)."""
+    return [mapping[k] for k in sorted(mapping, key=lambda key: int(key))]
+
+
+def load_models_registry(ref: str = _IMPACT_MODELS_REGISTRY) -> dict:
+    """Load models.json (forced params per model) from the model repo on Hugging Face.
+
+    The registry is NOT bundled with the preset — it lives on the models repo and is fetched from there.
+    Resolution: the ``KONFAI_IMPACT_MODELS_REGISTRY`` env path wins (dev/offline); otherwise ``ref`` must be
+    a ``repo:file`` Hugging Face reference.
+    """
+    local = os.environ.get("KONFAI_IMPACT_MODELS_REGISTRY", "")
+    if local:
+        path = Path(local)
+    elif ":" in ref:
+        repo, filename = ref.split(":", 1)
+        path = Path(hf_hub_download(repo_id=repo, filename=filename, repo_type="model"))  # nosec B615
+    else:
+        raise ValueError(
+            f"models_registry '{ref}' must be a 'repo:file' Hugging Face reference (the registry is fetched "
+            f"from HF, not bundled) — or set KONFAI_IMPACT_MODELS_REGISTRY to a local file for offline use."
+        )
+    return json.loads(path.read_text(encoding="utf-8"))
+
+
+def _model_key(ref: str) -> str:
+    """Registry key / staged relative path = the model file within the models repo (strip a 'repo:' prefix)."""
+    return ref.split(":", 1)[1] if ":" in ref else ref
+
+
+def _deepest_active_layer(layers_mask: str) -> int:
+    """Deepest (largest-FOV) layer selected by ``layers_mask``, as a 0-based index read left-to-right.
+
+    A model returns its feature layers shallow->deep (``[layer_0, layer_1, ...]``, see the model repo's
+    build scripts); ``layers_mask`` has one char per returned layer, position ``i`` == ``layer_i``, ``'1'``
+    = selected. In Jacobian the patch must cover the receptive field of the DEEPEST selected layer, so the
+    FOV is governed by the rightmost ``'1'``.
+    """
+    mask = layers_mask.strip().strip('"')
+    active = [i for i, char in enumerate(mask) if char == "1"]
+    if not active:
+        raise ValueError(f"LayersMask '{layers_mask}' selects no layer; cannot derive the model FOV.")
+    return max(active)
+
+
+def _fov_value(fov: dict, layers_mask: str) -> int:
+    """Evaluate a model's field-of-view (in voxels) from its registry ``fov`` spec.
+
+    Supported formulas (from the model repo, https://huggingface.co/VBoussot/impact-torchscript-models):
+      ``2*r*d+1``  MIND, from the handcrafted radius ``r`` / dilation ``d`` (e.g. R1D2 -> 5);
+      ``2^l+3``    TotalSegmentator / MRSegmentator, ``l`` = the deepest layer picked by ``layers_mask``,
+                   clamped to the receptive-field plateau ``_FOV_RAMP_MAX_LAYER`` (layers 7-8 -> layer 6);
+      a bare int  a fixed FOV (SAM2.1 -> 29, DINOv2 -> 14);
+      ``Global``   Anatomix — whole-image only (Static); has no finite Jacobian patch -> error.
+    An explicit ``value`` in the spec is honoured as a precomputed shortcut when the formula needs none.
+    """
+    formula = str(fov.get("formula", "")).strip()
+    key = re.sub(r"\s+", "", formula).lower()
+    if key.isdigit():
+        return int(key)
+    if key == "2*r*d+1":
+        return 2 * int(fov["r"]) * int(fov["d"]) + 1
+    if key == "2^l+3":
+        return 2 ** min(_deepest_active_layer(layers_mask), _FOV_RAMP_MAX_LAYER) + 3
+    if "global" in key:
+        raise ValueError(f"model FOV '{formula}' is whole-image only (Static); it has no Jacobian patch size.")
+    if fov.get("value") is not None:
+        return int(fov["value"])
+    raise ValueError(f"cannot evaluate model FOV formula '{formula}'.")
+
+
+def _patch_size(mode: str, entry: dict, layers_mask: str) -> str:
+    """PatchSize from the model FOV, one token per model axis (2D model -> 2 tokens, 3D -> 3): Static ->
+    whole image (all zeros); Jacobian -> the evaluated FOV repeated over the axes. A 2D model mixed with a
+    3D one at a resolution concatenates as e.g. ``29 29 11 11 11`` (SAM 2D + TS 3D), matching IMPACT."""
+    dim = int(entry.get("dimension", 3))
+    if mode.strip().strip('"').lower() != "jacobian":
+        return " ".join(["0"] * dim)
+    fov = _fov_value(entry.get("fov", {}), layers_mask)
+    return " ".join([str(fov)] * dim)
+
+
+def generate_impact_parameter_map(
+    template_text: str, resolutions: dict, registry: dict, mode: str = "Static"
+) -> str:
+    """Rewrite the resolution-dependent lines of ``template_text`` from the model matrix ``resolutions``.
+
+    Regenerated: MaximumNumberOfIterations, NumberOfResolutions, Fixed/MovingImagePyramidRescaleSchedule,
+    ImpactMode (from the config ``mode``), and the whole ImpactXxxK block; every other template line is
+    kept verbatim (optimizer, transform, metric weights, components...). N (number of resolutions) is
+    deduced from the config. ``mode`` drives PatchSize: Static -> ``0 0 0`` (whole image); Jacobian -> the
+    per-model FOV evaluated from the registry formula and the cell's ``layers_mask``.
+    """
+    res = _sorted_specs(resolutions)
+    n = len(res)
+    mode_clean = mode.strip().strip('"') or "Static"
+
+    impact: list[str] = []
+    for k, r in enumerate(res):
+        models = _sorted_specs(r.models)
+        entries = [registry[_model_key(m.ref)] for m in models]
+
+        def row(stem: str, values: list[str]) -> None:
+            impact.append(f"(Impact{stem}{k} " + " ".join(values) + ")")
+
+        # From the registry (models.json on the model repo) ONLY the 3 truly model-fixed props:
+        # Dimension, NumberOfChannels, PatchSize (the model FOV). Everything else is a per-model tuning knob
+        # taken straight from the cell: VoxelSize / LayersMask / SubsetFeatures / PCA / Distance / LayersWeight.
+        row("ModelsPath", [f'"{_model_key(m.ref)}"' for m in models])
+        row("Dimension", [e["dimension"] for e in entries])
+        row("NumberOfChannels", [e["numberofchannels"] for e in entries])
+        row("PatchSize", [_patch_size(mode_clean, e, m.layers_mask) for e, m in zip(entries, models)])
+        row("VoxelSize", [" ".join(_num(v) for v in m.voxel_size) for m in models])
+        row("LayersMask", [f'"{m.layers_mask}"' for m in models])
+        row("SubsetFeatures", [str(m.subset_features) for m in models])
+        row("PCA", [str(m.pca) for m in models])
+        row("Distance", [f'"{m.distance}"' for m in models])
+        row("LayersWeight", [" ".join(_num(w) for w in m.layers_weight) for m in models])
+        impact.append("")  # blank line between resolutions, mirroring the reference maps
+
+    # The per-resolution block is the contiguous span from the first to the last ``Impact<name><k>`` line
+    # (the blank lines the reference maps put BETWEEN resolutions fall inside that span). Replace the whole
+    # span in one shot with the generated block, so the reference blanks are not kept on top of ours.
+    lines = template_text.splitlines()
+    indexed = [(re.match(r"^\s*\((\S+?)\s+(.*?)\)\s*$", ln), ln) for ln in lines]
+    block_rows = [i for i, (m, _) in enumerate(indexed) if m and re.match(r"^Impact[A-Za-z]+\d+$", m.group(1))]
+    block_lo, block_hi = (block_rows[0], block_rows[-1]) if block_rows else (-1, -2)
+
+    out: list[str] = []
+    for i, (m, line) in enumerate(indexed):
+        key = m.group(1) if m else None
+        if block_lo <= i <= block_hi:
+            if i == block_lo:  # replace the whole span at its first line, drop the rest (incl. inner blanks)
+                out.extend(impact[:-1])
+        elif key == "MaximumNumberOfIterations":
+            out.append("(MaximumNumberOfIterations " + " ".join(_num(r.max_iterations) for r in res) + ")")
+        elif key == "NumberOfResolutions":
+            out.append(f"(NumberOfResolutions {n})")
+        elif key in ("FixedImagePyramidRescaleSchedule", "MovingImagePyramidRescaleSchedule"):
+            out.append(f"({key} " + " ".join(["1"] * 3 * n) + ")")
+        elif key == "ImpactMode":
+            out.append(f'(ImpactMode "{mode_clean}")')
+        else:
+            out.append(line)
+    return "\n".join(out)
+
 
 class ElastixEngine:
     """Run the elastix-IMPACT binary on a fixed/moving pair; return (moved, dvf) on the fixed grid.
@@ -60,14 +267,57 @@ class ElastixEngine:
     does NOT ship it), so registration is a subprocess call, not ``sitk.ElastixImageFilter``.
     """
 
-    def __init__(self, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         self._bundle_dir = Path(__file__).resolve().parent
         self._parameter_maps = [self._bundle_dir / Path(p).name for p in parameter_maps]
+        self._max_iterations = max_iterations
+        self._final_grid_spacing = final_grid_spacing
+        self._subset_features = subset_features
+        self._spatial_samples = spatial_samples
+        self._parameter_overrides = list(parameter_overrides)
+        # ImpactMode: Static computes features once per level (PatchSize 0 0 0 = whole image); Jacobian
+        # samples random patches sized to the model FOV each iteration. Global knob: one mode per preset.
+        self._mode = mode
+        # Matrix mode: when `resolutions` is given the parameter map is GENERATED from it (the config is the
+        # source of truth). An empty `resolutions` = an intensity preset (no IMPACT feature models): the fixed
+        # parameter maps are staged with only the global knob overrides.
+        self._resolutions = resolutions
+        self._registry = load_models_registry(models_registry) if resolutions else {}
+        # The feature models are DERIVED — the unique refs across the matrix cells (no flat `models` param).
+        models: list[str] = []
+        for res in _sorted_specs(resolutions):
+            for model in _sorted_specs(res.models):
+                if model.ref not in models:
+                    models.append(model.ref)
         self._models = models
-        self._iterations = iterations
+        # `iterations` (the progress-bar total) is NOT a config parameter — it is DERIVED: the sum of the
+        # per-resolution iteration budgets, read from the matrix (matrix mode) or the maps (legacy).
+        self._iterations = self._total_iterations()
         self._elastix_bin = self._ensure_binary()
         self._local_models = self._download_models()
 
+    def _total_iterations(self) -> int:
+        """Total iterations across all resolutions — the progress-bar budget, derived from the config."""
+        if self._resolutions:
+            return sum(int(res.max_iterations) for res in _sorted_specs(self._resolutions))
+        total = 0
+        for src in self._parameter_maps:
+            match = re.search(r"\(MaximumNumberOfIterations\s+([^)]*)\)", src.read_text(encoding="utf-8"))
+            if match:
+                total += sum(int(token) for token in match.group(1).split())
+        return total
+
     def _ensure_binary(self) -> Path:
         # Optional override: point at an existing elastix-IMPACT install (skips the download).
         override = os.environ.get("KONFAI_ELASTIX_DIR", "")
@@ -91,17 +341,90 @@ class ElastixEngine:
             models.append((filename, local))
         return models
 
+    def _parameter_map_overrides(self, global_only: bool = False) -> tuple[dict[str, str], list[tuple[str, str]]]:
+        """The tuned knobs as parameter-map overrides: ``(per_token, exact)``.
+
+        ``per_token`` maps an elastix key (or the ``ImpactSubsetFeatures`` prefix) to a value that replaces
+        **each** existing token, so per-resolution / per-model multiplicity is preserved (e.g.
+        ``(MaximumNumberOfIterations 500 250)`` -> ``(MaximumNumberOfIterations 300 300)``). ``exact``
+        entries (from ``parameter_overrides``, ``Key=value text``) replace the whole value verbatim and win
+        over the named knobs. Overrides only REPLACE keys already present in a map — never inject new ones.
+        ``global_only`` (matrix mode) keeps just the map-wide knobs and drops ``max_iterations`` /
+        ``subset_features`` — the per-resolution matrix already sets those per cell.
+        """
+        per_token: dict[str, str] = {}
+        if not global_only and self._max_iterations > 0:
+            per_token["MaximumNumberOfIterations"] = str(int(self._max_iterations))
+        if self._final_grid_spacing > 0:
+            per_token["FinalGridSpacingInPhysicalUnits"] = str(float(self._final_grid_spacing))
+        if not global_only and self._subset_features > 0:
+            per_token["ImpactSubsetFeatures"] = str(int(self._subset_features))  # prefix: indexed per metric
+        if self._spatial_samples > 0:
+            per_token["NumberOfSpatialSamples"] = str(int(self._spatial_samples))
+        exact: list[tuple[str, str]] = []
+        for entry in self._parameter_overrides:
+            key, sep, value = entry.partition("=")
+            if not sep or not key.strip():
+                raise ValueError(f"Invalid parameter_overrides entry '{entry}': expected 'Key=value text'.")
+            exact.append((key.strip(), value.strip()))
+        return per_token, exact
+
+    @staticmethod
+    def _apply_map_overrides(
+        text: str, per_token: dict[str, str], exact: list[tuple[str, str]], device_index: int
+    ) -> str:
+        """Patch a parameter map's text: set ImpactGPU to the device, apply exact key overrides, replace each
+        token of a per-token knob (preserving multiplicity), and warn for a requested key absent from the map.
+        """
+        entry_pattern = re.compile(r"^(\s*)\((\S+)((?:\s+[^)]*)?)\)\s*$")
+        requested = set(per_token) | {key for key, _ in exact}
+        seen: set[str] = set()
+        lines = []
+        for line in text.splitlines():
+            match = entry_pattern.match(line)
+            if match:
+                indent, key, values = match.group(1), match.group(2), match.group(3)
+                if key == "ImpactGPU":
+                    line = f"{indent}(ImpactGPU {device_index})"
+                else:
+                    exact_value = next((value for k, value in exact if k == key), None)
+                    if exact_value is not None:
+                        seen.add(key)
+                        line = f"{indent}({key} {exact_value})"
+                    else:
+                        token_key = "ImpactSubsetFeatures" if key.startswith("ImpactSubsetFeatures") else key
+                        if token_key in per_token:
+                            seen.add(token_key)
+                            replaced = " ".join(per_token[token_key] for _ in values.split())
+                            line = f"{indent}({key} {replaced})"
+            lines.append(line)
+        # Overrides never inject keys, so a knob set for a key absent from every map would silently do
+        # nothing — surface it (e.g. final_grid_spacing on a rigid-only preset).
+        for key in sorted(requested - seen):
+            print(f"[ImpactReg] note: override '{key}' matched no entry in the preset's parameter maps.")
+        return "\n".join(lines)
+
     def _stage_parameter_maps(self, work: Path, device_index: int) -> list[Path]:
-        """Copy parameter maps into the work dir, patching the ImpactGPU device field when present."""
+        """Stage the parameter maps into the work dir.
+
+        Matrix mode GENERATES each map from ``resolutions`` + the registry, then applies only the map-wide
+        knobs (grid spacing, spatial samples, exact overrides) — the matrix already sets iterations and
+        features per cell. Legacy mode copies the preset's maps and applies every per-token / exact override.
+        Both set the ImpactGPU device.
+        """
         staged = []
         for src in self._parameter_maps:
+            if self._resolutions:
+                text = generate_impact_parameter_map(
+                    src.read_text(encoding="utf-8"), self._resolutions, self._registry, self._mode
+                )
+                per_token, exact = self._parameter_map_overrides(global_only=True)
+            else:
+                text = src.read_text(encoding="utf-8")
+                per_token, exact = self._parameter_map_overrides()
+            text = self._apply_map_overrides(text, per_token, exact, device_index)
             dst = work / src.name
-            lines = []
-            for line in src.read_text(encoding="utf-8").splitlines():
-                if line.strip().startswith("(ImpactGPU"):
-                    line = f"(ImpactGPU {device_index})"
-                lines.append(line)
-            dst.write_text("\n".join(lines) + "\n", encoding="utf-8")
+            dst.write_text(text if text.endswith("\n") else text + "\n", encoding="utf-8")
             staged.append(dst)
         return staged
 
@@ -161,8 +484,10 @@ class ElastixEngine:
             captured: list[str] = []
             iteration_line = re.compile(r"^\d+\s")
             # ``iterations`` is the total iteration budget declared for the preset (summed over the
-            # chained parameter maps), so the bar spans the whole chain of registration stages.
-            progress = tqdm.tqdm(total=self._iterations or None, desc="Registration", ncols=0, leave=True)
+            # chained parameter maps), so the bar spans the whole chain of registration stages. A tuned
+            # ``max_iterations`` makes that declared budget stale — fall back to an open-ended bar.
+            budget = None if self._max_iterations > 0 else (self._iterations or None)
+            progress = tqdm.tqdm(total=budget, desc="Registration", ncols=0, leave=True)
             assert proc.stdout is not None
             resolution = 0
             for line in proc.stdout:
@@ -226,11 +551,33 @@ class ElastixRegistration(torch.nn.Module):
 
     accepts_attributes = True
 
-    def __init__(self, engine: str, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        engine: str,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         super().__init__()
         if engine != "elastix":
             raise NotImplementedError(f"ElastixRegistration engine '{engine}' is not implemented yet.")
-        self._engine = ElastixEngine(parameter_maps, models, iterations)
+        self._engine = ElastixEngine(
+            parameter_maps,
+            max_iterations,
+            final_grid_spacing,
+            subset_features,
+            spatial_samples,
+            parameter_overrides,
+            resolutions,
+            models_registry,
+            mode,
+        )
 
     def forward(
         self,
@@ -290,9 +637,23 @@ class RegistrationNet(network.Network):
         outputs_criterions: dict[str, network.TargetCriterionsLoader] = {"default": network.TargetCriterionsLoader()},
         engine: str = "elastix",
         parameter_maps: list[str] = [],
-        models: list[str] = [],
-        iterations: int = 1000,
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict[str, ResolutionSpec] = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
     ) -> None:
+        # The registration is fully described by the per-resolution model matrix ``resolutions`` (config =
+        # source of truth): each resolution lists its models, each model self-configured (ref, voxel_size,
+        # layers_mask, layers_weight, subset_features, pca, distance); intrinsic per-model props come from
+        # ``models_registry``. The feature-model download list is DERIVED from the matrix (no flat ``models``).
+        # Global knobs override the generated map: final_grid_spacing -> FinalGridSpacingInPhysicalUnits (mm),
+        # spatial_samples -> NumberOfSpatialSamples, parameter_overrides ('Key=value') -> any other entry.
+        # An empty ``resolutions`` = an intensity-only preset (no IMPACT models): the fixed maps are staged
+        # with just the global overrides. The total iteration count is derived (sum of per-resolution budgets).
         super().__init__(
             in_channels=1,
             optimizer=optimizer,
@@ -302,7 +663,18 @@ class RegistrationNet(network.Network):
         )
         self.add_module(
             "Registration",
-            ElastixRegistration(engine, parameter_maps, models, iterations),
+            ElastixRegistration(
+                engine,
+                parameter_maps,
+                max_iterations,
+                final_grid_spacing,
+                subset_features,
+                spatial_samples,
+                parameter_overrides,
+                resolutions,
+                models_registry,
+                mode,
+            ),
             in_branch=[0, 1, 2, 3],
             out_branch=["registration"],
         )

CBCT_CT_TS/Prediction.yml

@@ -5,11 +5,12 @@ Predictor:
       engine: elastix
       parameter_maps:
       - ParameterMap_CBCT_generic_TS.txt
-      models:
-      - VBoussot/impact-torchscript-models:TS/M852.pt
-      - VBoussot/impact-torchscript-models:TS/M850.pt
-      iterations: 1050
       outputs_criterions: None
+      max_iterations: 0
+      final_grid_spacing: 0.0
+      subset_features: 0
+      spatial_samples: 0
+      parameter_overrides: []
   Dataset:
     groups_src:
       Volume_0:

ConvexAdam_Coarse/Model.py

@@ -43,6 +43,9 @@ from konfai.network import network
 from konfai.utils.dataset import Attribute, data_to_image, image_to_data
 
 DIM = 3
+# The feature model's input channel count is an intrinsic property of the pretrained model (grayscale
+# medical images), not a tunable — so it's fixed here, never a config/signature parameter.
+NUM_CHANNELS = 1
 _IMAGE_F = itk.Image[itk.F, DIM]
 
 
@@ -106,7 +109,6 @@ class ConvexAdamEngine:
         self,
         models: list[str],
         voxel_size: list[float],
-        num_channels: int,
         overlap: int,
         layers_mask: list[bool],
         mixed_precision: bool,
@@ -131,7 +133,6 @@ class ConvexAdamEngine:
         # from disk in C++, so build the list once and reuse it across both stages and every case.
         self._configurations: "list[itk.ModelConfiguration] | None" = None
         self._voxel_size = voxel_size
-        self._num_channels = num_channels
         self._overlap = overlap
         self._layers_mask = layers_mask
         self._mixed_precision = mixed_precision
@@ -171,7 +172,7 @@ class ConvexAdamEngine:
                 itk.ModelConfiguration(
                     path,
                     DIM,
-                    self._num_channels,
+                    NUM_CHANNELS,
                     [0, 0, 0],
                     [float(v) for v in self._voxel_size],
                     self._overlap,
@@ -416,7 +417,6 @@ class RegistrationNet(network.Network):
         outputs_criterions: dict[str, network.TargetCriterionsLoader] = {"default": network.TargetCriterionsLoader()},
         models: list[str] = [],
         voxel_size: list[float] = [3.0, 3.0, 3.0],
-        num_channels: int = 1,
         overlap: int = 2,
         layers_mask: list[bool] = [True],
         mixed_precision: bool = False,
@@ -428,7 +428,7 @@ class RegistrationNet(network.Network):
         grid_shrink: int = 4,
         distance: list[str] = ["L1"],
         layers_weight: list[float] = [1.0],
-        subset_features: list[int] = [32],
+        subset_features: list[int] = [],  # feature-channel indices to keep (empty = all); NOT a count
         pca: list[int] = [0],
         stages: list[str] = ["coarse", "fine"],
         linear: bool = True,
@@ -445,7 +445,6 @@ class RegistrationNet(network.Network):
         engine = ConvexAdamEngine(
             models,
             voxel_size,
-            num_channels,
             overlap,
             layers_mask,
             mixed_precision,

ConvexAdam_Coarse/Prediction.yml

@@ -8,7 +8,6 @@ Predictor:
       - 3.0
       - 3.0
       - 3.0
-      num_channels: 1
       overlap: 2
       layers_mask:
       - true
@@ -23,8 +22,7 @@ Predictor:
       - L1
       layers_weight:
       - 1.0
-      subset_features:
-      - 32
+      subset_features: []
       pca:
       - 0
       linear: true

ConvexAdam_Composite/Model.py

@@ -43,6 +43,9 @@ from konfai.network import network
 from konfai.utils.dataset import Attribute, data_to_image, image_to_data
 
 DIM = 3
+# The feature model's input channel count is an intrinsic property of the pretrained model (grayscale
+# medical images), not a tunable — so it's fixed here, never a config/signature parameter.
+NUM_CHANNELS = 1
 _IMAGE_F = itk.Image[itk.F, DIM]
 
 
@@ -106,7 +109,6 @@ class ConvexAdamEngine:
         self,
         models: list[str],
         voxel_size: list[float],
-        num_channels: int,
         overlap: int,
         layers_mask: list[bool],
         mixed_precision: bool,
@@ -131,7 +133,6 @@ class ConvexAdamEngine:
         # from disk in C++, so build the list once and reuse it across both stages and every case.
         self._configurations: "list[itk.ModelConfiguration] | None" = None
         self._voxel_size = voxel_size
-        self._num_channels = num_channels
         self._overlap = overlap
         self._layers_mask = layers_mask
         self._mixed_precision = mixed_precision
@@ -171,7 +172,7 @@ class ConvexAdamEngine:
                 itk.ModelConfiguration(
                     path,
                     DIM,
-                    self._num_channels,
+                    NUM_CHANNELS,
                     [0, 0, 0],
                     [float(v) for v in self._voxel_size],
                     self._overlap,
@@ -416,7 +417,6 @@ class RegistrationNet(network.Network):
         outputs_criterions: dict[str, network.TargetCriterionsLoader] = {"default": network.TargetCriterionsLoader()},
         models: list[str] = [],
         voxel_size: list[float] = [3.0, 3.0, 3.0],
-        num_channels: int = 1,
         overlap: int = 2,
         layers_mask: list[bool] = [True],
         mixed_precision: bool = False,
@@ -428,7 +428,7 @@ class RegistrationNet(network.Network):
         grid_shrink: int = 4,
         distance: list[str] = ["L1"],
         layers_weight: list[float] = [1.0],
-        subset_features: list[int] = [32],
+        subset_features: list[int] = [],  # feature-channel indices to keep (empty = all); NOT a count
         pca: list[int] = [0],
         stages: list[str] = ["coarse", "fine"],
         linear: bool = True,
@@ -445,7 +445,6 @@ class RegistrationNet(network.Network):
         engine = ConvexAdamEngine(
             models,
             voxel_size,
-            num_channels,
             overlap,
             layers_mask,
             mixed_precision,

ConvexAdam_Composite/Prediction.yml

@@ -8,7 +8,6 @@ Predictor:
       - 3.0
       - 3.0
       - 3.0
-      num_channels: 1
       overlap: 2
       layers_mask:
       - true
@@ -23,8 +22,7 @@ Predictor:
       - L1
       layers_weight:
       - 1.0
-      subset_features:
-      - 32
+      subset_features: []
       pca:
       - 0
       linear: true

ConvexAdam_Fine/Model.py

@@ -43,6 +43,9 @@ from konfai.network import network
 from konfai.utils.dataset import Attribute, data_to_image, image_to_data
 
 DIM = 3
+# The feature model's input channel count is an intrinsic property of the pretrained model (grayscale
+# medical images), not a tunable — so it's fixed here, never a config/signature parameter.
+NUM_CHANNELS = 1
 _IMAGE_F = itk.Image[itk.F, DIM]
 
 
@@ -106,7 +109,6 @@ class ConvexAdamEngine:
         self,
         models: list[str],
         voxel_size: list[float],
-        num_channels: int,
         overlap: int,
         layers_mask: list[bool],
         mixed_precision: bool,
@@ -131,7 +133,6 @@ class ConvexAdamEngine:
         # from disk in C++, so build the list once and reuse it across both stages and every case.
         self._configurations: "list[itk.ModelConfiguration] | None" = None
         self._voxel_size = voxel_size
-        self._num_channels = num_channels
         self._overlap = overlap
         self._layers_mask = layers_mask
         self._mixed_precision = mixed_precision
@@ -171,7 +172,7 @@ class ConvexAdamEngine:
                 itk.ModelConfiguration(
                     path,
                     DIM,
-                    self._num_channels,
+                    NUM_CHANNELS,
                     [0, 0, 0],
                     [float(v) for v in self._voxel_size],
                     self._overlap,
@@ -416,7 +417,6 @@ class RegistrationNet(network.Network):
         outputs_criterions: dict[str, network.TargetCriterionsLoader] = {"default": network.TargetCriterionsLoader()},
         models: list[str] = [],
         voxel_size: list[float] = [3.0, 3.0, 3.0],
-        num_channels: int = 1,
         overlap: int = 2,
         layers_mask: list[bool] = [True],
         mixed_precision: bool = False,
@@ -428,7 +428,7 @@ class RegistrationNet(network.Network):
         grid_shrink: int = 4,
         distance: list[str] = ["L1"],
         layers_weight: list[float] = [1.0],
-        subset_features: list[int] = [32],
+        subset_features: list[int] = [],  # feature-channel indices to keep (empty = all); NOT a count
         pca: list[int] = [0],
         stages: list[str] = ["coarse", "fine"],
         linear: bool = True,
@@ -445,7 +445,6 @@ class RegistrationNet(network.Network):
         engine = ConvexAdamEngine(
             models,
             voxel_size,
-            num_channels,
             overlap,
             layers_mask,
             mixed_precision,

ConvexAdam_Fine/Prediction.yml

@@ -8,7 +8,6 @@ Predictor:
       - 3.0
       - 3.0
       - 3.0
-      num_channels: 1
       overlap: 2
       layers_mask:
       - true
@@ -23,8 +22,7 @@ Predictor:
       - L1
       layers_weight:
       - 1.0
-      subset_features:
-      - 32
+      subset_features: []
       pca:
       - 0
       linear: false

Generic_Rigid/Model.py

@@ -32,6 +32,7 @@ NOTE: do NOT add ``from __future__ import annotations`` here — KonfAI's config
 runtime-evaluated annotations (``get_origin``); PEP 563 stringized annotations break arg resolution.
 """
 
+import json
 import os
 import re
 import shutil
@@ -52,6 +53,212 @@ from konfai.utils.dataset import Attribute, data_to_image, image_to_data
 # Set KONFAI_ELASTIX_DIR to point at an existing install and skip the download.
 ELASTIX_CACHE = Path.home() / ".cache" / "konfai" / "elastix-impact"
 
+# ---------------------------------------------------------------------------------------------------
+# Per-resolution model matrix (the config is the source of truth) -> generated IMPACT parameter map.
+# IMPACT field docs: https://github.com/vboussot/ImpactLoss/tree/main/ParameterMaps
+# The forced per-model props (dimension/channels/FOV formula) live in a registry (models.json on
+# VBoussot/impact-torchscript-models); the config carries the FREE knobs (which models per resolution,
+# feature voxel size, iterations, per-model layer weights/mask/subset/pca/distance) and the global
+# ``mode``. PatchSize follows ImpactMode: Static -> "0 0 0" (whole image); Jacobian -> the model FOV
+# evaluated from the registry formula (MIND 2*r*d+1, TS/MRSeg 2^l+3, SAM 29, DINOv2 14) as a cube.
+# ---------------------------------------------------------------------------------------------------
+
+_IMPACT_MODELS_REGISTRY = "VBoussot/impact-torchscript-models:models.json"
+
+# ``2^l+3`` grows with depth but the segmenters' receptive field plateaus: layers 7-8 share layer 6's
+# FOV (the "ramp max"). A config that deep should really run in Static (whole image) anyway; in Jacobian
+# we clamp ``l`` to this plateau so the patch stays finite and matches the real FOV.
+_FOV_RAMP_MAX_LAYER = 6
+
+
+def _num(x: object) -> str:
+    """Format a number the elastix way: integers without a trailing '.0' (6.0 -> '6', 0.2 -> '0.2')."""
+    return "%g" % float(x)
+
+
+class ModelSpec:
+    """One feature model at one resolution, with its OWN config (several models may share a resolution).
+
+    ``ref`` selects the model; ``voxel_size`` / ``layers_weight`` / ``subset_features`` / ``pca`` /
+    ``distance`` are its free per-(resolution, model) tuning knobs (the doc's per-model *tuning* fields).
+    The intrinsic per-model props — dimension, channels, ``layers_mask``, patch-size (FOV) — come from the
+    registry (read-only); ``layers_mask`` / ``distance`` left empty fall back to the registry default.
+    """
+
+    def __init__(
+        self,
+        ref: str,
+        voxel_size: list[float] = [],
+        layers_weight: list[float] = [1.0],
+        subset_features: int = 0,
+        pca: int = 0,
+        distance: str = "",
+        layers_mask: str = "",
+    ) -> None:
+        self.ref = ref
+        self.voxel_size = voxel_size
+        self.layers_weight = layers_weight
+        self.subset_features = subset_features
+        self.pca = pca
+        self.distance = distance
+        self.layers_mask = layers_mask
+
+
+class ResolutionSpec:
+    """One elastix resolution level: its iteration budget and the models compared there (each self-configured)."""
+
+    def __init__(self, max_iterations: int, models: dict[str, ModelSpec]) -> None:
+        self.max_iterations = max_iterations
+        self.models = models
+
+
+def _sorted_specs(mapping: dict) -> list:
+    """dict keyed by string indices ('0','1',...) -> values in numeric order (well-defined res/model order)."""
+    return [mapping[k] for k in sorted(mapping, key=lambda key: int(key))]
+
+
+def load_models_registry(ref: str = _IMPACT_MODELS_REGISTRY) -> dict:
+    """Load models.json (forced params per model) from the model repo on Hugging Face.
+
+    The registry is NOT bundled with the preset — it lives on the models repo and is fetched from there.
+    Resolution: the ``KONFAI_IMPACT_MODELS_REGISTRY`` env path wins (dev/offline); otherwise ``ref`` must be
+    a ``repo:file`` Hugging Face reference.
+    """
+    local = os.environ.get("KONFAI_IMPACT_MODELS_REGISTRY", "")
+    if local:
+        path = Path(local)
+    elif ":" in ref:
+        repo, filename = ref.split(":", 1)
+        path = Path(hf_hub_download(repo_id=repo, filename=filename, repo_type="model"))  # nosec B615
+    else:
+        raise ValueError(
+            f"models_registry '{ref}' must be a 'repo:file' Hugging Face reference (the registry is fetched "
+            f"from HF, not bundled) — or set KONFAI_IMPACT_MODELS_REGISTRY to a local file for offline use."
+        )
+    return json.loads(path.read_text(encoding="utf-8"))
+
+
+def _model_key(ref: str) -> str:
+    """Registry key / staged relative path = the model file within the models repo (strip a 'repo:' prefix)."""
+    return ref.split(":", 1)[1] if ":" in ref else ref
+
+
+def _deepest_active_layer(layers_mask: str) -> int:
+    """Deepest (largest-FOV) layer selected by ``layers_mask``, as a 0-based index read left-to-right.
+
+    A model returns its feature layers shallow->deep (``[layer_0, layer_1, ...]``, see the model repo's
+    build scripts); ``layers_mask`` has one char per returned layer, position ``i`` == ``layer_i``, ``'1'``
+    = selected. In Jacobian the patch must cover the receptive field of the DEEPEST selected layer, so the
+    FOV is governed by the rightmost ``'1'``.
+    """
+    mask = layers_mask.strip().strip('"')
+    active = [i for i, char in enumerate(mask) if char == "1"]
+    if not active:
+        raise ValueError(f"LayersMask '{layers_mask}' selects no layer; cannot derive the model FOV.")
+    return max(active)
+
+
+def _fov_value(fov: dict, layers_mask: str) -> int:
+    """Evaluate a model's field-of-view (in voxels) from its registry ``fov`` spec.
+
+    Supported formulas (from the model repo, https://huggingface.co/VBoussot/impact-torchscript-models):
+      ``2*r*d+1``  MIND, from the handcrafted radius ``r`` / dilation ``d`` (e.g. R1D2 -> 5);
+      ``2^l+3``    TotalSegmentator / MRSegmentator, ``l`` = the deepest layer picked by ``layers_mask``,
+                   clamped to the receptive-field plateau ``_FOV_RAMP_MAX_LAYER`` (layers 7-8 -> layer 6);
+      a bare int  a fixed FOV (SAM2.1 -> 29, DINOv2 -> 14);
+      ``Global``   Anatomix — whole-image only (Static); has no finite Jacobian patch -> error.
+    An explicit ``value`` in the spec is honoured as a precomputed shortcut when the formula needs none.
+    """
+    formula = str(fov.get("formula", "")).strip()
+    key = re.sub(r"\s+", "", formula).lower()
+    if key.isdigit():
+        return int(key)
+    if key == "2*r*d+1":
+        return 2 * int(fov["r"]) * int(fov["d"]) + 1
+    if key == "2^l+3":
+        return 2 ** min(_deepest_active_layer(layers_mask), _FOV_RAMP_MAX_LAYER) + 3
+    if "global" in key:
+        raise ValueError(f"model FOV '{formula}' is whole-image only (Static); it has no Jacobian patch size.")
+    if fov.get("value") is not None:
+        return int(fov["value"])
+    raise ValueError(f"cannot evaluate model FOV formula '{formula}'.")
+
+
+def _patch_size(mode: str, entry: dict, layers_mask: str) -> str:
+    """PatchSize from the model FOV, one token per model axis (2D model -> 2 tokens, 3D -> 3): Static ->
+    whole image (all zeros); Jacobian -> the evaluated FOV repeated over the axes. A 2D model mixed with a
+    3D one at a resolution concatenates as e.g. ``29 29 11 11 11`` (SAM 2D + TS 3D), matching IMPACT."""
+    dim = int(entry.get("dimension", 3))
+    if mode.strip().strip('"').lower() != "jacobian":
+        return " ".join(["0"] * dim)
+    fov = _fov_value(entry.get("fov", {}), layers_mask)
+    return " ".join([str(fov)] * dim)
+
+
+def generate_impact_parameter_map(
+    template_text: str, resolutions: dict, registry: dict, mode: str = "Static"
+) -> str:
+    """Rewrite the resolution-dependent lines of ``template_text`` from the model matrix ``resolutions``.
+
+    Regenerated: MaximumNumberOfIterations, NumberOfResolutions, Fixed/MovingImagePyramidRescaleSchedule,
+    ImpactMode (from the config ``mode``), and the whole ImpactXxxK block; every other template line is
+    kept verbatim (optimizer, transform, metric weights, components...). N (number of resolutions) is
+    deduced from the config. ``mode`` drives PatchSize: Static -> ``0 0 0`` (whole image); Jacobian -> the
+    per-model FOV evaluated from the registry formula and the cell's ``layers_mask``.
+    """
+    res = _sorted_specs(resolutions)
+    n = len(res)
+    mode_clean = mode.strip().strip('"') or "Static"
+
+    impact: list[str] = []
+    for k, r in enumerate(res):
+        models = _sorted_specs(r.models)
+        entries = [registry[_model_key(m.ref)] for m in models]
+
+        def row(stem: str, values: list[str]) -> None:
+            impact.append(f"(Impact{stem}{k} " + " ".join(values) + ")")
+
+        # From the registry (models.json on the model repo) ONLY the 3 truly model-fixed props:
+        # Dimension, NumberOfChannels, PatchSize (the model FOV). Everything else is a per-model tuning knob
+        # taken straight from the cell: VoxelSize / LayersMask / SubsetFeatures / PCA / Distance / LayersWeight.
+        row("ModelsPath", [f'"{_model_key(m.ref)}"' for m in models])
+        row("Dimension", [e["dimension"] for e in entries])
+        row("NumberOfChannels", [e["numberofchannels"] for e in entries])
+        row("PatchSize", [_patch_size(mode_clean, e, m.layers_mask) for e, m in zip(entries, models)])
+        row("VoxelSize", [" ".join(_num(v) for v in m.voxel_size) for m in models])
+        row("LayersMask", [f'"{m.layers_mask}"' for m in models])
+        row("SubsetFeatures", [str(m.subset_features) for m in models])
+        row("PCA", [str(m.pca) for m in models])
+        row("Distance", [f'"{m.distance}"' for m in models])
+        row("LayersWeight", [" ".join(_num(w) for w in m.layers_weight) for m in models])
+        impact.append("")  # blank line between resolutions, mirroring the reference maps
+
+    # The per-resolution block is the contiguous span from the first to the last ``Impact<name><k>`` line
+    # (the blank lines the reference maps put BETWEEN resolutions fall inside that span). Replace the whole
+    # span in one shot with the generated block, so the reference blanks are not kept on top of ours.
+    lines = template_text.splitlines()
+    indexed = [(re.match(r"^\s*\((\S+?)\s+(.*?)\)\s*$", ln), ln) for ln in lines]
+    block_rows = [i for i, (m, _) in enumerate(indexed) if m and re.match(r"^Impact[A-Za-z]+\d+$", m.group(1))]
+    block_lo, block_hi = (block_rows[0], block_rows[-1]) if block_rows else (-1, -2)
+
+    out: list[str] = []
+    for i, (m, line) in enumerate(indexed):
+        key = m.group(1) if m else None
+        if block_lo <= i <= block_hi:
+            if i == block_lo:  # replace the whole span at its first line, drop the rest (incl. inner blanks)
+                out.extend(impact[:-1])
+        elif key == "MaximumNumberOfIterations":
+            out.append("(MaximumNumberOfIterations " + " ".join(_num(r.max_iterations) for r in res) + ")")
+        elif key == "NumberOfResolutions":
+            out.append(f"(NumberOfResolutions {n})")
+        elif key in ("FixedImagePyramidRescaleSchedule", "MovingImagePyramidRescaleSchedule"):
+            out.append(f"({key} " + " ".join(["1"] * 3 * n) + ")")
+        elif key == "ImpactMode":
+            out.append(f'(ImpactMode "{mode_clean}")')
+        else:
+            out.append(line)
+    return "\n".join(out)
+
 
 class ElastixEngine:
     """Run the elastix-IMPACT binary on a fixed/moving pair; return (moved, dvf) on the fixed grid.
@@ -60,14 +267,57 @@ class ElastixEngine:
     does NOT ship it), so registration is a subprocess call, not ``sitk.ElastixImageFilter``.
     """
 
-    def __init__(self, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         self._bundle_dir = Path(__file__).resolve().parent
         self._parameter_maps = [self._bundle_dir / Path(p).name for p in parameter_maps]
+        self._max_iterations = max_iterations
+        self._final_grid_spacing = final_grid_spacing
+        self._subset_features = subset_features
+        self._spatial_samples = spatial_samples
+        self._parameter_overrides = list(parameter_overrides)
+        # ImpactMode: Static computes features once per level (PatchSize 0 0 0 = whole image); Jacobian
+        # samples random patches sized to the model FOV each iteration. Global knob: one mode per preset.
+        self._mode = mode
+        # Matrix mode: when `resolutions` is given the parameter map is GENERATED from it (the config is the
+        # source of truth). An empty `resolutions` = an intensity preset (no IMPACT feature models): the fixed
+        # parameter maps are staged with only the global knob overrides.
+        self._resolutions = resolutions
+        self._registry = load_models_registry(models_registry) if resolutions else {}
+        # The feature models are DERIVED — the unique refs across the matrix cells (no flat `models` param).
+        models: list[str] = []
+        for res in _sorted_specs(resolutions):
+            for model in _sorted_specs(res.models):
+                if model.ref not in models:
+                    models.append(model.ref)
         self._models = models
-        self._iterations = iterations
+        # `iterations` (the progress-bar total) is NOT a config parameter — it is DERIVED: the sum of the
+        # per-resolution iteration budgets, read from the matrix (matrix mode) or the maps (legacy).
+        self._iterations = self._total_iterations()
         self._elastix_bin = self._ensure_binary()
         self._local_models = self._download_models()
 
+    def _total_iterations(self) -> int:
+        """Total iterations across all resolutions — the progress-bar budget, derived from the config."""
+        if self._resolutions:
+            return sum(int(res.max_iterations) for res in _sorted_specs(self._resolutions))
+        total = 0
+        for src in self._parameter_maps:
+            match = re.search(r"\(MaximumNumberOfIterations\s+([^)]*)\)", src.read_text(encoding="utf-8"))
+            if match:
+                total += sum(int(token) for token in match.group(1).split())
+        return total
+
     def _ensure_binary(self) -> Path:
         # Optional override: point at an existing elastix-IMPACT install (skips the download).
         override = os.environ.get("KONFAI_ELASTIX_DIR", "")
@@ -91,17 +341,90 @@ class ElastixEngine:
             models.append((filename, local))
         return models
 
+    def _parameter_map_overrides(self, global_only: bool = False) -> tuple[dict[str, str], list[tuple[str, str]]]:
+        """The tuned knobs as parameter-map overrides: ``(per_token, exact)``.
+
+        ``per_token`` maps an elastix key (or the ``ImpactSubsetFeatures`` prefix) to a value that replaces
+        **each** existing token, so per-resolution / per-model multiplicity is preserved (e.g.
+        ``(MaximumNumberOfIterations 500 250)`` -> ``(MaximumNumberOfIterations 300 300)``). ``exact``
+        entries (from ``parameter_overrides``, ``Key=value text``) replace the whole value verbatim and win
+        over the named knobs. Overrides only REPLACE keys already present in a map — never inject new ones.
+        ``global_only`` (matrix mode) keeps just the map-wide knobs and drops ``max_iterations`` /
+        ``subset_features`` — the per-resolution matrix already sets those per cell.
+        """
+        per_token: dict[str, str] = {}
+        if not global_only and self._max_iterations > 0:
+            per_token["MaximumNumberOfIterations"] = str(int(self._max_iterations))
+        if self._final_grid_spacing > 0:
+            per_token["FinalGridSpacingInPhysicalUnits"] = str(float(self._final_grid_spacing))
+        if not global_only and self._subset_features > 0:
+            per_token["ImpactSubsetFeatures"] = str(int(self._subset_features))  # prefix: indexed per metric
+        if self._spatial_samples > 0:
+            per_token["NumberOfSpatialSamples"] = str(int(self._spatial_samples))
+        exact: list[tuple[str, str]] = []
+        for entry in self._parameter_overrides:
+            key, sep, value = entry.partition("=")
+            if not sep or not key.strip():
+                raise ValueError(f"Invalid parameter_overrides entry '{entry}': expected 'Key=value text'.")
+            exact.append((key.strip(), value.strip()))
+        return per_token, exact
+
+    @staticmethod
+    def _apply_map_overrides(
+        text: str, per_token: dict[str, str], exact: list[tuple[str, str]], device_index: int
+    ) -> str:
+        """Patch a parameter map's text: set ImpactGPU to the device, apply exact key overrides, replace each
+        token of a per-token knob (preserving multiplicity), and warn for a requested key absent from the map.
+        """
+        entry_pattern = re.compile(r"^(\s*)\((\S+)((?:\s+[^)]*)?)\)\s*$")
+        requested = set(per_token) | {key for key, _ in exact}
+        seen: set[str] = set()
+        lines = []
+        for line in text.splitlines():
+            match = entry_pattern.match(line)
+            if match:
+                indent, key, values = match.group(1), match.group(2), match.group(3)
+                if key == "ImpactGPU":
+                    line = f"{indent}(ImpactGPU {device_index})"
+                else:
+                    exact_value = next((value for k, value in exact if k == key), None)
+                    if exact_value is not None:
+                        seen.add(key)
+                        line = f"{indent}({key} {exact_value})"
+                    else:
+                        token_key = "ImpactSubsetFeatures" if key.startswith("ImpactSubsetFeatures") else key
+                        if token_key in per_token:
+                            seen.add(token_key)
+                            replaced = " ".join(per_token[token_key] for _ in values.split())
+                            line = f"{indent}({key} {replaced})"
+            lines.append(line)
+        # Overrides never inject keys, so a knob set for a key absent from every map would silently do
+        # nothing — surface it (e.g. final_grid_spacing on a rigid-only preset).
+        for key in sorted(requested - seen):
+            print(f"[ImpactReg] note: override '{key}' matched no entry in the preset's parameter maps.")
+        return "\n".join(lines)
+
     def _stage_parameter_maps(self, work: Path, device_index: int) -> list[Path]:
-        """Copy parameter maps into the work dir, patching the ImpactGPU device field when present."""
+        """Stage the parameter maps into the work dir.
+
+        Matrix mode GENERATES each map from ``resolutions`` + the registry, then applies only the map-wide
+        knobs (grid spacing, spatial samples, exact overrides) — the matrix already sets iterations and
+        features per cell. Legacy mode copies the preset's maps and applies every per-token / exact override.
+        Both set the ImpactGPU device.
+        """
         staged = []
         for src in self._parameter_maps:
+            if self._resolutions:
+                text = generate_impact_parameter_map(
+                    src.read_text(encoding="utf-8"), self._resolutions, self._registry, self._mode
+                )
+                per_token, exact = self._parameter_map_overrides(global_only=True)
+            else:
+                text = src.read_text(encoding="utf-8")
+                per_token, exact = self._parameter_map_overrides()
+            text = self._apply_map_overrides(text, per_token, exact, device_index)
             dst = work / src.name
-            lines = []
-            for line in src.read_text(encoding="utf-8").splitlines():
-                if line.strip().startswith("(ImpactGPU"):
-                    line = f"(ImpactGPU {device_index})"
-                lines.append(line)
-            dst.write_text("\n".join(lines) + "\n", encoding="utf-8")
+            dst.write_text(text if text.endswith("\n") else text + "\n", encoding="utf-8")
             staged.append(dst)
         return staged
 
@@ -161,8 +484,10 @@ class ElastixEngine:
             captured: list[str] = []
             iteration_line = re.compile(r"^\d+\s")
             # ``iterations`` is the total iteration budget declared for the preset (summed over the
-            # chained parameter maps), so the bar spans the whole chain of registration stages.
-            progress = tqdm.tqdm(total=self._iterations or None, desc="Registration", ncols=0, leave=True)
+            # chained parameter maps), so the bar spans the whole chain of registration stages. A tuned
+            # ``max_iterations`` makes that declared budget stale — fall back to an open-ended bar.
+            budget = None if self._max_iterations > 0 else (self._iterations or None)
+            progress = tqdm.tqdm(total=budget, desc="Registration", ncols=0, leave=True)
             assert proc.stdout is not None
             resolution = 0
             for line in proc.stdout:
@@ -226,11 +551,33 @@ class ElastixRegistration(torch.nn.Module):
 
     accepts_attributes = True
 
-    def __init__(self, engine: str, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        engine: str,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         super().__init__()
         if engine != "elastix":
             raise NotImplementedError(f"ElastixRegistration engine '{engine}' is not implemented yet.")
-        self._engine = ElastixEngine(parameter_maps, models, iterations)
+        self._engine = ElastixEngine(
+            parameter_maps,
+            max_iterations,
+            final_grid_spacing,
+            subset_features,
+            spatial_samples,
+            parameter_overrides,
+            resolutions,
+            models_registry,
+            mode,
+        )
 
     def forward(
         self,
@@ -290,9 +637,23 @@ class RegistrationNet(network.Network):
         outputs_criterions: dict[str, network.TargetCriterionsLoader] = {"default": network.TargetCriterionsLoader()},
         engine: str = "elastix",
         parameter_maps: list[str] = [],
-        models: list[str] = [],
-        iterations: int = 1000,
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict[str, ResolutionSpec] = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
     ) -> None:
+        # The registration is fully described by the per-resolution model matrix ``resolutions`` (config =
+        # source of truth): each resolution lists its models, each model self-configured (ref, voxel_size,
+        # layers_mask, layers_weight, subset_features, pca, distance); intrinsic per-model props come from
+        # ``models_registry``. The feature-model download list is DERIVED from the matrix (no flat ``models``).
+        # Global knobs override the generated map: final_grid_spacing -> FinalGridSpacingInPhysicalUnits (mm),
+        # spatial_samples -> NumberOfSpatialSamples, parameter_overrides ('Key=value') -> any other entry.
+        # An empty ``resolutions`` = an intensity-only preset (no IMPACT models): the fixed maps are staged
+        # with just the global overrides. The total iteration count is derived (sum of per-resolution budgets).
         super().__init__(
             in_channels=1,
             optimizer=optimizer,
@@ -302,7 +663,18 @@ class RegistrationNet(network.Network):
         )
         self.add_module(
             "Registration",
-            ElastixRegistration(engine, parameter_maps, models, iterations),
+            ElastixRegistration(
+                engine,
+                parameter_maps,
+                max_iterations,
+                final_grid_spacing,
+                subset_features,
+                spatial_samples,
+                parameter_overrides,
+                resolutions,
+                models_registry,
+                mode,
+            ),
             in_branch=[0, 1, 2, 3],
             out_branch=["registration"],
         )

Generic_Rigid/Prediction.yml

@@ -5,9 +5,10 @@ Predictor:
       engine: elastix
       parameter_maps:
       - Parameters_Rigid.txt
-      models: []
-      iterations: 250
       outputs_criterions: None
+      max_iterations: 0
+      spatial_samples: 0
+      parameter_overrides: []
   Dataset:
     groups_src:
       Volume_0:

Generic_Rigid_BSpline/Model.py

@@ -32,6 +32,7 @@ NOTE: do NOT add ``from __future__ import annotations`` here — KonfAI's config
 runtime-evaluated annotations (``get_origin``); PEP 563 stringized annotations break arg resolution.
 """
 
+import json
 import os
 import re
 import shutil
@@ -52,6 +53,212 @@ from konfai.utils.dataset import Attribute, data_to_image, image_to_data
 # Set KONFAI_ELASTIX_DIR to point at an existing install and skip the download.
 ELASTIX_CACHE = Path.home() / ".cache" / "konfai" / "elastix-impact"
 
+# ---------------------------------------------------------------------------------------------------
+# Per-resolution model matrix (the config is the source of truth) -> generated IMPACT parameter map.
+# IMPACT field docs: https://github.com/vboussot/ImpactLoss/tree/main/ParameterMaps
+# The forced per-model props (dimension/channels/FOV formula) live in a registry (models.json on
+# VBoussot/impact-torchscript-models); the config carries the FREE knobs (which models per resolution,
+# feature voxel size, iterations, per-model layer weights/mask/subset/pca/distance) and the global
+# ``mode``. PatchSize follows ImpactMode: Static -> "0 0 0" (whole image); Jacobian -> the model FOV
+# evaluated from the registry formula (MIND 2*r*d+1, TS/MRSeg 2^l+3, SAM 29, DINOv2 14) as a cube.
+# ---------------------------------------------------------------------------------------------------
+
+_IMPACT_MODELS_REGISTRY = "VBoussot/impact-torchscript-models:models.json"
+
+# ``2^l+3`` grows with depth but the segmenters' receptive field plateaus: layers 7-8 share layer 6's
+# FOV (the "ramp max"). A config that deep should really run in Static (whole image) anyway; in Jacobian
+# we clamp ``l`` to this plateau so the patch stays finite and matches the real FOV.
+_FOV_RAMP_MAX_LAYER = 6
+
+
+def _num(x: object) -> str:
+    """Format a number the elastix way: integers without a trailing '.0' (6.0 -> '6', 0.2 -> '0.2')."""
+    return "%g" % float(x)
+
+
+class ModelSpec:
+    """One feature model at one resolution, with its OWN config (several models may share a resolution).
+
+    ``ref`` selects the model; ``voxel_size`` / ``layers_weight`` / ``subset_features`` / ``pca`` /
+    ``distance`` are its free per-(resolution, model) tuning knobs (the doc's per-model *tuning* fields).
+    The intrinsic per-model props — dimension, channels, ``layers_mask``, patch-size (FOV) — come from the
+    registry (read-only); ``layers_mask`` / ``distance`` left empty fall back to the registry default.
+    """
+
+    def __init__(
+        self,
+        ref: str,
+        voxel_size: list[float] = [],
+        layers_weight: list[float] = [1.0],
+        subset_features: int = 0,
+        pca: int = 0,
+        distance: str = "",
+        layers_mask: str = "",
+    ) -> None:
+        self.ref = ref
+        self.voxel_size = voxel_size
+        self.layers_weight = layers_weight
+        self.subset_features = subset_features
+        self.pca = pca
+        self.distance = distance
+        self.layers_mask = layers_mask
+
+
+class ResolutionSpec:
+    """One elastix resolution level: its iteration budget and the models compared there (each self-configured)."""
+
+    def __init__(self, max_iterations: int, models: dict[str, ModelSpec]) -> None:
+        self.max_iterations = max_iterations
+        self.models = models
+
+
+def _sorted_specs(mapping: dict) -> list:
+    """dict keyed by string indices ('0','1',...) -> values in numeric order (well-defined res/model order)."""
+    return [mapping[k] for k in sorted(mapping, key=lambda key: int(key))]
+
+
+def load_models_registry(ref: str = _IMPACT_MODELS_REGISTRY) -> dict:
+    """Load models.json (forced params per model) from the model repo on Hugging Face.
+
+    The registry is NOT bundled with the preset — it lives on the models repo and is fetched from there.
+    Resolution: the ``KONFAI_IMPACT_MODELS_REGISTRY`` env path wins (dev/offline); otherwise ``ref`` must be
+    a ``repo:file`` Hugging Face reference.
+    """
+    local = os.environ.get("KONFAI_IMPACT_MODELS_REGISTRY", "")
+    if local:
+        path = Path(local)
+    elif ":" in ref:
+        repo, filename = ref.split(":", 1)
+        path = Path(hf_hub_download(repo_id=repo, filename=filename, repo_type="model"))  # nosec B615
+    else:
+        raise ValueError(
+            f"models_registry '{ref}' must be a 'repo:file' Hugging Face reference (the registry is fetched "
+            f"from HF, not bundled) — or set KONFAI_IMPACT_MODELS_REGISTRY to a local file for offline use."
+        )
+    return json.loads(path.read_text(encoding="utf-8"))
+
+
+def _model_key(ref: str) -> str:
+    """Registry key / staged relative path = the model file within the models repo (strip a 'repo:' prefix)."""
+    return ref.split(":", 1)[1] if ":" in ref else ref
+
+
+def _deepest_active_layer(layers_mask: str) -> int:
+    """Deepest (largest-FOV) layer selected by ``layers_mask``, as a 0-based index read left-to-right.
+
+    A model returns its feature layers shallow->deep (``[layer_0, layer_1, ...]``, see the model repo's
+    build scripts); ``layers_mask`` has one char per returned layer, position ``i`` == ``layer_i``, ``'1'``
+    = selected. In Jacobian the patch must cover the receptive field of the DEEPEST selected layer, so the
+    FOV is governed by the rightmost ``'1'``.
+    """
+    mask = layers_mask.strip().strip('"')
+    active = [i for i, char in enumerate(mask) if char == "1"]
+    if not active:
+        raise ValueError(f"LayersMask '{layers_mask}' selects no layer; cannot derive the model FOV.")
+    return max(active)
+
+
+def _fov_value(fov: dict, layers_mask: str) -> int:
+    """Evaluate a model's field-of-view (in voxels) from its registry ``fov`` spec.
+
+    Supported formulas (from the model repo, https://huggingface.co/VBoussot/impact-torchscript-models):
+      ``2*r*d+1``  MIND, from the handcrafted radius ``r`` / dilation ``d`` (e.g. R1D2 -> 5);
+      ``2^l+3``    TotalSegmentator / MRSegmentator, ``l`` = the deepest layer picked by ``layers_mask``,
+                   clamped to the receptive-field plateau ``_FOV_RAMP_MAX_LAYER`` (layers 7-8 -> layer 6);
+      a bare int  a fixed FOV (SAM2.1 -> 29, DINOv2 -> 14);
+      ``Global``   Anatomix — whole-image only (Static); has no finite Jacobian patch -> error.
+    An explicit ``value`` in the spec is honoured as a precomputed shortcut when the formula needs none.
+    """
+    formula = str(fov.get("formula", "")).strip()
+    key = re.sub(r"\s+", "", formula).lower()
+    if key.isdigit():
+        return int(key)
+    if key == "2*r*d+1":
+        return 2 * int(fov["r"]) * int(fov["d"]) + 1
+    if key == "2^l+3":
+        return 2 ** min(_deepest_active_layer(layers_mask), _FOV_RAMP_MAX_LAYER) + 3
+    if "global" in key:
+        raise ValueError(f"model FOV '{formula}' is whole-image only (Static); it has no Jacobian patch size.")
+    if fov.get("value") is not None:
+        return int(fov["value"])
+    raise ValueError(f"cannot evaluate model FOV formula '{formula}'.")
+
+
+def _patch_size(mode: str, entry: dict, layers_mask: str) -> str:
+    """PatchSize from the model FOV, one token per model axis (2D model -> 2 tokens, 3D -> 3): Static ->
+    whole image (all zeros); Jacobian -> the evaluated FOV repeated over the axes. A 2D model mixed with a
+    3D one at a resolution concatenates as e.g. ``29 29 11 11 11`` (SAM 2D + TS 3D), matching IMPACT."""
+    dim = int(entry.get("dimension", 3))
+    if mode.strip().strip('"').lower() != "jacobian":
+        return " ".join(["0"] * dim)
+    fov = _fov_value(entry.get("fov", {}), layers_mask)
+    return " ".join([str(fov)] * dim)
+
+
+def generate_impact_parameter_map(
+    template_text: str, resolutions: dict, registry: dict, mode: str = "Static"
+) -> str:
+    """Rewrite the resolution-dependent lines of ``template_text`` from the model matrix ``resolutions``.
+
+    Regenerated: MaximumNumberOfIterations, NumberOfResolutions, Fixed/MovingImagePyramidRescaleSchedule,
+    ImpactMode (from the config ``mode``), and the whole ImpactXxxK block; every other template line is
+    kept verbatim (optimizer, transform, metric weights, components...). N (number of resolutions) is
+    deduced from the config. ``mode`` drives PatchSize: Static -> ``0 0 0`` (whole image); Jacobian -> the
+    per-model FOV evaluated from the registry formula and the cell's ``layers_mask``.
+    """
+    res = _sorted_specs(resolutions)
+    n = len(res)
+    mode_clean = mode.strip().strip('"') or "Static"
+
+    impact: list[str] = []
+    for k, r in enumerate(res):
+        models = _sorted_specs(r.models)
+        entries = [registry[_model_key(m.ref)] for m in models]
+
+        def row(stem: str, values: list[str]) -> None:
+            impact.append(f"(Impact{stem}{k} " + " ".join(values) + ")")
+
+        # From the registry (models.json on the model repo) ONLY the 3 truly model-fixed props:
+        # Dimension, NumberOfChannels, PatchSize (the model FOV). Everything else is a per-model tuning knob
+        # taken straight from the cell: VoxelSize / LayersMask / SubsetFeatures / PCA / Distance / LayersWeight.
+        row("ModelsPath", [f'"{_model_key(m.ref)}"' for m in models])
+        row("Dimension", [e["dimension"] for e in entries])
+        row("NumberOfChannels", [e["numberofchannels"] for e in entries])
+        row("PatchSize", [_patch_size(mode_clean, e, m.layers_mask) for e, m in zip(entries, models)])
+        row("VoxelSize", [" ".join(_num(v) for v in m.voxel_size) for m in models])
+        row("LayersMask", [f'"{m.layers_mask}"' for m in models])
+        row("SubsetFeatures", [str(m.subset_features) for m in models])
+        row("PCA", [str(m.pca) for m in models])
+        row("Distance", [f'"{m.distance}"' for m in models])
+        row("LayersWeight", [" ".join(_num(w) for w in m.layers_weight) for m in models])
+        impact.append("")  # blank line between resolutions, mirroring the reference maps
+
+    # The per-resolution block is the contiguous span from the first to the last ``Impact<name><k>`` line
+    # (the blank lines the reference maps put BETWEEN resolutions fall inside that span). Replace the whole
+    # span in one shot with the generated block, so the reference blanks are not kept on top of ours.
+    lines = template_text.splitlines()
+    indexed = [(re.match(r"^\s*\((\S+?)\s+(.*?)\)\s*$", ln), ln) for ln in lines]
+    block_rows = [i for i, (m, _) in enumerate(indexed) if m and re.match(r"^Impact[A-Za-z]+\d+$", m.group(1))]
+    block_lo, block_hi = (block_rows[0], block_rows[-1]) if block_rows else (-1, -2)
+
+    out: list[str] = []
+    for i, (m, line) in enumerate(indexed):
+        key = m.group(1) if m else None
+        if block_lo <= i <= block_hi:
+            if i == block_lo:  # replace the whole span at its first line, drop the rest (incl. inner blanks)
+                out.extend(impact[:-1])
+        elif key == "MaximumNumberOfIterations":
+            out.append("(MaximumNumberOfIterations " + " ".join(_num(r.max_iterations) for r in res) + ")")
+        elif key == "NumberOfResolutions":
+            out.append(f"(NumberOfResolutions {n})")
+        elif key in ("FixedImagePyramidRescaleSchedule", "MovingImagePyramidRescaleSchedule"):
+            out.append(f"({key} " + " ".join(["1"] * 3 * n) + ")")
+        elif key == "ImpactMode":
+            out.append(f'(ImpactMode "{mode_clean}")')
+        else:
+            out.append(line)
+    return "\n".join(out)
+
 
 class ElastixEngine:
     """Run the elastix-IMPACT binary on a fixed/moving pair; return (moved, dvf) on the fixed grid.
@@ -60,14 +267,57 @@ class ElastixEngine:
     does NOT ship it), so registration is a subprocess call, not ``sitk.ElastixImageFilter``.
     """
 
-    def __init__(self, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         self._bundle_dir = Path(__file__).resolve().parent
         self._parameter_maps = [self._bundle_dir / Path(p).name for p in parameter_maps]
+        self._max_iterations = max_iterations
+        self._final_grid_spacing = final_grid_spacing
+        self._subset_features = subset_features
+        self._spatial_samples = spatial_samples
+        self._parameter_overrides = list(parameter_overrides)
+        # ImpactMode: Static computes features once per level (PatchSize 0 0 0 = whole image); Jacobian
+        # samples random patches sized to the model FOV each iteration. Global knob: one mode per preset.
+        self._mode = mode
+        # Matrix mode: when `resolutions` is given the parameter map is GENERATED from it (the config is the
+        # source of truth). An empty `resolutions` = an intensity preset (no IMPACT feature models): the fixed
+        # parameter maps are staged with only the global knob overrides.
+        self._resolutions = resolutions
+        self._registry = load_models_registry(models_registry) if resolutions else {}
+        # The feature models are DERIVED — the unique refs across the matrix cells (no flat `models` param).
+        models: list[str] = []
+        for res in _sorted_specs(resolutions):
+            for model in _sorted_specs(res.models):
+                if model.ref not in models:
+                    models.append(model.ref)
         self._models = models
-        self._iterations = iterations
+        # `iterations` (the progress-bar total) is NOT a config parameter — it is DERIVED: the sum of the
+        # per-resolution iteration budgets, read from the matrix (matrix mode) or the maps (legacy).
+        self._iterations = self._total_iterations()
         self._elastix_bin = self._ensure_binary()
         self._local_models = self._download_models()
 
+    def _total_iterations(self) -> int:
+        """Total iterations across all resolutions — the progress-bar budget, derived from the config."""
+        if self._resolutions:
+            return sum(int(res.max_iterations) for res in _sorted_specs(self._resolutions))
+        total = 0
+        for src in self._parameter_maps:
+            match = re.search(r"\(MaximumNumberOfIterations\s+([^)]*)\)", src.read_text(encoding="utf-8"))
+            if match:
+                total += sum(int(token) for token in match.group(1).split())
+        return total
+
     def _ensure_binary(self) -> Path:
         # Optional override: point at an existing elastix-IMPACT install (skips the download).
         override = os.environ.get("KONFAI_ELASTIX_DIR", "")
@@ -91,17 +341,90 @@ class ElastixEngine:
             models.append((filename, local))
         return models
 
+    def _parameter_map_overrides(self, global_only: bool = False) -> tuple[dict[str, str], list[tuple[str, str]]]:
+        """The tuned knobs as parameter-map overrides: ``(per_token, exact)``.
+
+        ``per_token`` maps an elastix key (or the ``ImpactSubsetFeatures`` prefix) to a value that replaces
+        **each** existing token, so per-resolution / per-model multiplicity is preserved (e.g.
+        ``(MaximumNumberOfIterations 500 250)`` -> ``(MaximumNumberOfIterations 300 300)``). ``exact``
+        entries (from ``parameter_overrides``, ``Key=value text``) replace the whole value verbatim and win
+        over the named knobs. Overrides only REPLACE keys already present in a map — never inject new ones.
+        ``global_only`` (matrix mode) keeps just the map-wide knobs and drops ``max_iterations`` /
+        ``subset_features`` — the per-resolution matrix already sets those per cell.
+        """
+        per_token: dict[str, str] = {}
+        if not global_only and self._max_iterations > 0:
+            per_token["MaximumNumberOfIterations"] = str(int(self._max_iterations))
+        if self._final_grid_spacing > 0:
+            per_token["FinalGridSpacingInPhysicalUnits"] = str(float(self._final_grid_spacing))
+        if not global_only and self._subset_features > 0:
+            per_token["ImpactSubsetFeatures"] = str(int(self._subset_features))  # prefix: indexed per metric
+        if self._spatial_samples > 0:
+            per_token["NumberOfSpatialSamples"] = str(int(self._spatial_samples))
+        exact: list[tuple[str, str]] = []
+        for entry in self._parameter_overrides:
+            key, sep, value = entry.partition("=")
+            if not sep or not key.strip():
+                raise ValueError(f"Invalid parameter_overrides entry '{entry}': expected 'Key=value text'.")
+            exact.append((key.strip(), value.strip()))
+        return per_token, exact
+
+    @staticmethod
+    def _apply_map_overrides(
+        text: str, per_token: dict[str, str], exact: list[tuple[str, str]], device_index: int
+    ) -> str:
+        """Patch a parameter map's text: set ImpactGPU to the device, apply exact key overrides, replace each
+        token of a per-token knob (preserving multiplicity), and warn for a requested key absent from the map.
+        """
+        entry_pattern = re.compile(r"^(\s*)\((\S+)((?:\s+[^)]*)?)\)\s*$")
+        requested = set(per_token) | {key for key, _ in exact}
+        seen: set[str] = set()
+        lines = []
+        for line in text.splitlines():
+            match = entry_pattern.match(line)
+            if match:
+                indent, key, values = match.group(1), match.group(2), match.group(3)
+                if key == "ImpactGPU":
+                    line = f"{indent}(ImpactGPU {device_index})"
+                else:
+                    exact_value = next((value for k, value in exact if k == key), None)
+                    if exact_value is not None:
+                        seen.add(key)
+                        line = f"{indent}({key} {exact_value})"
+                    else:
+                        token_key = "ImpactSubsetFeatures" if key.startswith("ImpactSubsetFeatures") else key
+                        if token_key in per_token:
+                            seen.add(token_key)
+                            replaced = " ".join(per_token[token_key] for _ in values.split())
+                            line = f"{indent}({key} {replaced})"
+            lines.append(line)
+        # Overrides never inject keys, so a knob set for a key absent from every map would silently do
+        # nothing — surface it (e.g. final_grid_spacing on a rigid-only preset).
+        for key in sorted(requested - seen):
+            print(f"[ImpactReg] note: override '{key}' matched no entry in the preset's parameter maps.")
+        return "\n".join(lines)
+
     def _stage_parameter_maps(self, work: Path, device_index: int) -> list[Path]:
-        """Copy parameter maps into the work dir, patching the ImpactGPU device field when present."""
+        """Stage the parameter maps into the work dir.
+
+        Matrix mode GENERATES each map from ``resolutions`` + the registry, then applies only the map-wide
+        knobs (grid spacing, spatial samples, exact overrides) — the matrix already sets iterations and
+        features per cell. Legacy mode copies the preset's maps and applies every per-token / exact override.
+        Both set the ImpactGPU device.
+        """
         staged = []
         for src in self._parameter_maps:
+            if self._resolutions:
+                text = generate_impact_parameter_map(
+                    src.read_text(encoding="utf-8"), self._resolutions, self._registry, self._mode
+                )
+                per_token, exact = self._parameter_map_overrides(global_only=True)
+            else:
+                text = src.read_text(encoding="utf-8")
+                per_token, exact = self._parameter_map_overrides()
+            text = self._apply_map_overrides(text, per_token, exact, device_index)
             dst = work / src.name
-            lines = []
-            for line in src.read_text(encoding="utf-8").splitlines():
-                if line.strip().startswith("(ImpactGPU"):
-                    line = f"(ImpactGPU {device_index})"
-                lines.append(line)
-            dst.write_text("\n".join(lines) + "\n", encoding="utf-8")
+            dst.write_text(text if text.endswith("\n") else text + "\n", encoding="utf-8")
             staged.append(dst)
         return staged
 
@@ -161,8 +484,10 @@ class ElastixEngine:
             captured: list[str] = []
             iteration_line = re.compile(r"^\d+\s")
             # ``iterations`` is the total iteration budget declared for the preset (summed over the
-            # chained parameter maps), so the bar spans the whole chain of registration stages.
-            progress = tqdm.tqdm(total=self._iterations or None, desc="Registration", ncols=0, leave=True)
+            # chained parameter maps), so the bar spans the whole chain of registration stages. A tuned
+            # ``max_iterations`` makes that declared budget stale — fall back to an open-ended bar.
+            budget = None if self._max_iterations > 0 else (self._iterations or None)
+            progress = tqdm.tqdm(total=budget, desc="Registration", ncols=0, leave=True)
             assert proc.stdout is not None
             resolution = 0
             for line in proc.stdout:
@@ -226,11 +551,33 @@ class ElastixRegistration(torch.nn.Module):
 
     accepts_attributes = True
 
-    def __init__(self, engine: str, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        engine: str,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         super().__init__()
         if engine != "elastix":
             raise NotImplementedError(f"ElastixRegistration engine '{engine}' is not implemented yet.")
-        self._engine = ElastixEngine(parameter_maps, models, iterations)
+        self._engine = ElastixEngine(
+            parameter_maps,
+            max_iterations,
+            final_grid_spacing,
+            subset_features,
+            spatial_samples,
+            parameter_overrides,
+            resolutions,
+            models_registry,
+            mode,
+        )
 
     def forward(
         self,
@@ -290,9 +637,23 @@ class RegistrationNet(network.Network):
         outputs_criterions: dict[str, network.TargetCriterionsLoader] = {"default": network.TargetCriterionsLoader()},
         engine: str = "elastix",
         parameter_maps: list[str] = [],
-        models: list[str] = [],
-        iterations: int = 1000,
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict[str, ResolutionSpec] = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
     ) -> None:
+        # The registration is fully described by the per-resolution model matrix ``resolutions`` (config =
+        # source of truth): each resolution lists its models, each model self-configured (ref, voxel_size,
+        # layers_mask, layers_weight, subset_features, pca, distance); intrinsic per-model props come from
+        # ``models_registry``. The feature-model download list is DERIVED from the matrix (no flat ``models``).
+        # Global knobs override the generated map: final_grid_spacing -> FinalGridSpacingInPhysicalUnits (mm),
+        # spatial_samples -> NumberOfSpatialSamples, parameter_overrides ('Key=value') -> any other entry.
+        # An empty ``resolutions`` = an intensity-only preset (no IMPACT models): the fixed maps are staged
+        # with just the global overrides. The total iteration count is derived (sum of per-resolution budgets).
         super().__init__(
             in_channels=1,
             optimizer=optimizer,
@@ -302,7 +663,18 @@ class RegistrationNet(network.Network):
         )
         self.add_module(
             "Registration",
-            ElastixRegistration(engine, parameter_maps, models, iterations),
+            ElastixRegistration(
+                engine,
+                parameter_maps,
+                max_iterations,
+                final_grid_spacing,
+                subset_features,
+                spatial_samples,
+                parameter_overrides,
+                resolutions,
+                models_registry,
+                mode,
+            ),
             in_branch=[0, 1, 2, 3],
             out_branch=["registration"],
         )

Generic_Rigid_BSpline/Prediction.yml

@@ -6,9 +6,11 @@ Predictor:
       parameter_maps:
       - Parameters_Rigid.txt
       - Parameters_BSpline.txt
-      models: []
-      iterations: 3000
       outputs_criterions: None
+      max_iterations: 0
+      final_grid_spacing: 0.0
+      spatial_samples: 0
+      parameter_overrides: []
   Dataset:
     groups_src:
       Volume_0:

MR_CT_HeadNeck/Model.py

@@ -32,6 +32,7 @@ NOTE: do NOT add ``from __future__ import annotations`` here — KonfAI's config
 runtime-evaluated annotations (``get_origin``); PEP 563 stringized annotations break arg resolution.
 """
 
+import json
 import os
 import re
 import shutil
@@ -52,6 +53,212 @@ from konfai.utils.dataset import Attribute, data_to_image, image_to_data
 # Set KONFAI_ELASTIX_DIR to point at an existing install and skip the download.
 ELASTIX_CACHE = Path.home() / ".cache" / "konfai" / "elastix-impact"
 
+# ---------------------------------------------------------------------------------------------------
+# Per-resolution model matrix (the config is the source of truth) -> generated IMPACT parameter map.
+# IMPACT field docs: https://github.com/vboussot/ImpactLoss/tree/main/ParameterMaps
+# The forced per-model props (dimension/channels/FOV formula) live in a registry (models.json on
+# VBoussot/impact-torchscript-models); the config carries the FREE knobs (which models per resolution,
+# feature voxel size, iterations, per-model layer weights/mask/subset/pca/distance) and the global
+# ``mode``. PatchSize follows ImpactMode: Static -> "0 0 0" (whole image); Jacobian -> the model FOV
+# evaluated from the registry formula (MIND 2*r*d+1, TS/MRSeg 2^l+3, SAM 29, DINOv2 14) as a cube.
+# ---------------------------------------------------------------------------------------------------
+
+_IMPACT_MODELS_REGISTRY = "VBoussot/impact-torchscript-models:models.json"
+
+# ``2^l+3`` grows with depth but the segmenters' receptive field plateaus: layers 7-8 share layer 6's
+# FOV (the "ramp max"). A config that deep should really run in Static (whole image) anyway; in Jacobian
+# we clamp ``l`` to this plateau so the patch stays finite and matches the real FOV.
+_FOV_RAMP_MAX_LAYER = 6
+
+
+def _num(x: object) -> str:
+    """Format a number the elastix way: integers without a trailing '.0' (6.0 -> '6', 0.2 -> '0.2')."""
+    return "%g" % float(x)
+
+
+class ModelSpec:
+    """One feature model at one resolution, with its OWN config (several models may share a resolution).
+
+    ``ref`` selects the model; ``voxel_size`` / ``layers_weight`` / ``subset_features`` / ``pca`` /
+    ``distance`` are its free per-(resolution, model) tuning knobs (the doc's per-model *tuning* fields).
+    The intrinsic per-model props — dimension, channels, ``layers_mask``, patch-size (FOV) — come from the
+    registry (read-only); ``layers_mask`` / ``distance`` left empty fall back to the registry default.
+    """
+
+    def __init__(
+        self,
+        ref: str,
+        voxel_size: list[float] = [],
+        layers_weight: list[float] = [1.0],
+        subset_features: int = 0,
+        pca: int = 0,
+        distance: str = "",
+        layers_mask: str = "",
+    ) -> None:
+        self.ref = ref
+        self.voxel_size = voxel_size
+        self.layers_weight = layers_weight
+        self.subset_features = subset_features
+        self.pca = pca
+        self.distance = distance
+        self.layers_mask = layers_mask
+
+
+class ResolutionSpec:
+    """One elastix resolution level: its iteration budget and the models compared there (each self-configured)."""
+
+    def __init__(self, max_iterations: int, models: dict[str, ModelSpec]) -> None:
+        self.max_iterations = max_iterations
+        self.models = models
+
+
+def _sorted_specs(mapping: dict) -> list:
+    """dict keyed by string indices ('0','1',...) -> values in numeric order (well-defined res/model order)."""
+    return [mapping[k] for k in sorted(mapping, key=lambda key: int(key))]
+
+
+def load_models_registry(ref: str = _IMPACT_MODELS_REGISTRY) -> dict:
+    """Load models.json (forced params per model) from the model repo on Hugging Face.
+
+    The registry is NOT bundled with the preset — it lives on the models repo and is fetched from there.
+    Resolution: the ``KONFAI_IMPACT_MODELS_REGISTRY`` env path wins (dev/offline); otherwise ``ref`` must be
+    a ``repo:file`` Hugging Face reference.
+    """
+    local = os.environ.get("KONFAI_IMPACT_MODELS_REGISTRY", "")
+    if local:
+        path = Path(local)
+    elif ":" in ref:
+        repo, filename = ref.split(":", 1)
+        path = Path(hf_hub_download(repo_id=repo, filename=filename, repo_type="model"))  # nosec B615
+    else:
+        raise ValueError(
+            f"models_registry '{ref}' must be a 'repo:file' Hugging Face reference (the registry is fetched "
+            f"from HF, not bundled) — or set KONFAI_IMPACT_MODELS_REGISTRY to a local file for offline use."
+        )
+    return json.loads(path.read_text(encoding="utf-8"))
+
+
+def _model_key(ref: str) -> str:
+    """Registry key / staged relative path = the model file within the models repo (strip a 'repo:' prefix)."""
+    return ref.split(":", 1)[1] if ":" in ref else ref
+
+
+def _deepest_active_layer(layers_mask: str) -> int:
+    """Deepest (largest-FOV) layer selected by ``layers_mask``, as a 0-based index read left-to-right.
+
+    A model returns its feature layers shallow->deep (``[layer_0, layer_1, ...]``, see the model repo's
+    build scripts); ``layers_mask`` has one char per returned layer, position ``i`` == ``layer_i``, ``'1'``
+    = selected. In Jacobian the patch must cover the receptive field of the DEEPEST selected layer, so the
+    FOV is governed by the rightmost ``'1'``.
+    """
+    mask = layers_mask.strip().strip('"')
+    active = [i for i, char in enumerate(mask) if char == "1"]
+    if not active:
+        raise ValueError(f"LayersMask '{layers_mask}' selects no layer; cannot derive the model FOV.")
+    return max(active)
+
+
+def _fov_value(fov: dict, layers_mask: str) -> int:
+    """Evaluate a model's field-of-view (in voxels) from its registry ``fov`` spec.
+
+    Supported formulas (from the model repo, https://huggingface.co/VBoussot/impact-torchscript-models):
+      ``2*r*d+1``  MIND, from the handcrafted radius ``r`` / dilation ``d`` (e.g. R1D2 -> 5);
+      ``2^l+3``    TotalSegmentator / MRSegmentator, ``l`` = the deepest layer picked by ``layers_mask``,
+                   clamped to the receptive-field plateau ``_FOV_RAMP_MAX_LAYER`` (layers 7-8 -> layer 6);
+      a bare int  a fixed FOV (SAM2.1 -> 29, DINOv2 -> 14);
+      ``Global``   Anatomix — whole-image only (Static); has no finite Jacobian patch -> error.
+    An explicit ``value`` in the spec is honoured as a precomputed shortcut when the formula needs none.
+    """
+    formula = str(fov.get("formula", "")).strip()
+    key = re.sub(r"\s+", "", formula).lower()
+    if key.isdigit():
+        return int(key)
+    if key == "2*r*d+1":
+        return 2 * int(fov["r"]) * int(fov["d"]) + 1
+    if key == "2^l+3":
+        return 2 ** min(_deepest_active_layer(layers_mask), _FOV_RAMP_MAX_LAYER) + 3
+    if "global" in key:
+        raise ValueError(f"model FOV '{formula}' is whole-image only (Static); it has no Jacobian patch size.")
+    if fov.get("value") is not None:
+        return int(fov["value"])
+    raise ValueError(f"cannot evaluate model FOV formula '{formula}'.")
+
+
+def _patch_size(mode: str, entry: dict, layers_mask: str) -> str:
+    """PatchSize from the model FOV, one token per model axis (2D model -> 2 tokens, 3D -> 3): Static ->
+    whole image (all zeros); Jacobian -> the evaluated FOV repeated over the axes. A 2D model mixed with a
+    3D one at a resolution concatenates as e.g. ``29 29 11 11 11`` (SAM 2D + TS 3D), matching IMPACT."""
+    dim = int(entry.get("dimension", 3))
+    if mode.strip().strip('"').lower() != "jacobian":
+        return " ".join(["0"] * dim)
+    fov = _fov_value(entry.get("fov", {}), layers_mask)
+    return " ".join([str(fov)] * dim)
+
+
+def generate_impact_parameter_map(
+    template_text: str, resolutions: dict, registry: dict, mode: str = "Static"
+) -> str:
+    """Rewrite the resolution-dependent lines of ``template_text`` from the model matrix ``resolutions``.
+
+    Regenerated: MaximumNumberOfIterations, NumberOfResolutions, Fixed/MovingImagePyramidRescaleSchedule,
+    ImpactMode (from the config ``mode``), and the whole ImpactXxxK block; every other template line is
+    kept verbatim (optimizer, transform, metric weights, components...). N (number of resolutions) is
+    deduced from the config. ``mode`` drives PatchSize: Static -> ``0 0 0`` (whole image); Jacobian -> the
+    per-model FOV evaluated from the registry formula and the cell's ``layers_mask``.
+    """
+    res = _sorted_specs(resolutions)
+    n = len(res)
+    mode_clean = mode.strip().strip('"') or "Static"
+
+    impact: list[str] = []
+    for k, r in enumerate(res):
+        models = _sorted_specs(r.models)
+        entries = [registry[_model_key(m.ref)] for m in models]
+
+        def row(stem: str, values: list[str]) -> None:
+            impact.append(f"(Impact{stem}{k} " + " ".join(values) + ")")
+
+        # From the registry (models.json on the model repo) ONLY the 3 truly model-fixed props:
+        # Dimension, NumberOfChannels, PatchSize (the model FOV). Everything else is a per-model tuning knob
+        # taken straight from the cell: VoxelSize / LayersMask / SubsetFeatures / PCA / Distance / LayersWeight.
+        row("ModelsPath", [f'"{_model_key(m.ref)}"' for m in models])
+        row("Dimension", [e["dimension"] for e in entries])
+        row("NumberOfChannels", [e["numberofchannels"] for e in entries])
+        row("PatchSize", [_patch_size(mode_clean, e, m.layers_mask) for e, m in zip(entries, models)])
+        row("VoxelSize", [" ".join(_num(v) for v in m.voxel_size) for m in models])
+        row("LayersMask", [f'"{m.layers_mask}"' for m in models])
+        row("SubsetFeatures", [str(m.subset_features) for m in models])
+        row("PCA", [str(m.pca) for m in models])
+        row("Distance", [f'"{m.distance}"' for m in models])
+        row("LayersWeight", [" ".join(_num(w) for w in m.layers_weight) for m in models])
+        impact.append("")  # blank line between resolutions, mirroring the reference maps
+
+    # The per-resolution block is the contiguous span from the first to the last ``Impact<name><k>`` line
+    # (the blank lines the reference maps put BETWEEN resolutions fall inside that span). Replace the whole
+    # span in one shot with the generated block, so the reference blanks are not kept on top of ours.
+    lines = template_text.splitlines()
+    indexed = [(re.match(r"^\s*\((\S+?)\s+(.*?)\)\s*$", ln), ln) for ln in lines]
+    block_rows = [i for i, (m, _) in enumerate(indexed) if m and re.match(r"^Impact[A-Za-z]+\d+$", m.group(1))]
+    block_lo, block_hi = (block_rows[0], block_rows[-1]) if block_rows else (-1, -2)
+
+    out: list[str] = []
+    for i, (m, line) in enumerate(indexed):
+        key = m.group(1) if m else None
+        if block_lo <= i <= block_hi:
+            if i == block_lo:  # replace the whole span at its first line, drop the rest (incl. inner blanks)
+                out.extend(impact[:-1])
+        elif key == "MaximumNumberOfIterations":
+            out.append("(MaximumNumberOfIterations " + " ".join(_num(r.max_iterations) for r in res) + ")")
+        elif key == "NumberOfResolutions":
+            out.append(f"(NumberOfResolutions {n})")
+        elif key in ("FixedImagePyramidRescaleSchedule", "MovingImagePyramidRescaleSchedule"):
+            out.append(f"({key} " + " ".join(["1"] * 3 * n) + ")")
+        elif key == "ImpactMode":
+            out.append(f'(ImpactMode "{mode_clean}")')
+        else:
+            out.append(line)
+    return "\n".join(out)
+
 
 class ElastixEngine:
     """Run the elastix-IMPACT binary on a fixed/moving pair; return (moved, dvf) on the fixed grid.
@@ -60,14 +267,57 @@ class ElastixEngine:
     does NOT ship it), so registration is a subprocess call, not ``sitk.ElastixImageFilter``.
     """
 
-    def __init__(self, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         self._bundle_dir = Path(__file__).resolve().parent
         self._parameter_maps = [self._bundle_dir / Path(p).name for p in parameter_maps]
+        self._max_iterations = max_iterations
+        self._final_grid_spacing = final_grid_spacing
+        self._subset_features = subset_features
+        self._spatial_samples = spatial_samples
+        self._parameter_overrides = list(parameter_overrides)
+        # ImpactMode: Static computes features once per level (PatchSize 0 0 0 = whole image); Jacobian
+        # samples random patches sized to the model FOV each iteration. Global knob: one mode per preset.
+        self._mode = mode
+        # Matrix mode: when `resolutions` is given the parameter map is GENERATED from it (the config is the
+        # source of truth). An empty `resolutions` = an intensity preset (no IMPACT feature models): the fixed
+        # parameter maps are staged with only the global knob overrides.
+        self._resolutions = resolutions
+        self._registry = load_models_registry(models_registry) if resolutions else {}
+        # The feature models are DERIVED — the unique refs across the matrix cells (no flat `models` param).
+        models: list[str] = []
+        for res in _sorted_specs(resolutions):
+            for model in _sorted_specs(res.models):
+                if model.ref not in models:
+                    models.append(model.ref)
         self._models = models
-        self._iterations = iterations
+        # `iterations` (the progress-bar total) is NOT a config parameter — it is DERIVED: the sum of the
+        # per-resolution iteration budgets, read from the matrix (matrix mode) or the maps (legacy).
+        self._iterations = self._total_iterations()
         self._elastix_bin = self._ensure_binary()
         self._local_models = self._download_models()
 
+    def _total_iterations(self) -> int:
+        """Total iterations across all resolutions — the progress-bar budget, derived from the config."""
+        if self._resolutions:
+            return sum(int(res.max_iterations) for res in _sorted_specs(self._resolutions))
+        total = 0
+        for src in self._parameter_maps:
+            match = re.search(r"\(MaximumNumberOfIterations\s+([^)]*)\)", src.read_text(encoding="utf-8"))
+            if match:
+                total += sum(int(token) for token in match.group(1).split())
+        return total
+
     def _ensure_binary(self) -> Path:
         # Optional override: point at an existing elastix-IMPACT install (skips the download).
         override = os.environ.get("KONFAI_ELASTIX_DIR", "")
@@ -91,17 +341,90 @@ class ElastixEngine:
             models.append((filename, local))
         return models
 
+    def _parameter_map_overrides(self, global_only: bool = False) -> tuple[dict[str, str], list[tuple[str, str]]]:
+        """The tuned knobs as parameter-map overrides: ``(per_token, exact)``.
+
+        ``per_token`` maps an elastix key (or the ``ImpactSubsetFeatures`` prefix) to a value that replaces
+        **each** existing token, so per-resolution / per-model multiplicity is preserved (e.g.
+        ``(MaximumNumberOfIterations 500 250)`` -> ``(MaximumNumberOfIterations 300 300)``). ``exact``
+        entries (from ``parameter_overrides``, ``Key=value text``) replace the whole value verbatim and win
+        over the named knobs. Overrides only REPLACE keys already present in a map — never inject new ones.
+        ``global_only`` (matrix mode) keeps just the map-wide knobs and drops ``max_iterations`` /
+        ``subset_features`` — the per-resolution matrix already sets those per cell.
+        """
+        per_token: dict[str, str] = {}
+        if not global_only and self._max_iterations > 0:
+            per_token["MaximumNumberOfIterations"] = str(int(self._max_iterations))
+        if self._final_grid_spacing > 0:
+            per_token["FinalGridSpacingInPhysicalUnits"] = str(float(self._final_grid_spacing))
+        if not global_only and self._subset_features > 0:
+            per_token["ImpactSubsetFeatures"] = str(int(self._subset_features))  # prefix: indexed per metric
+        if self._spatial_samples > 0:
+            per_token["NumberOfSpatialSamples"] = str(int(self._spatial_samples))
+        exact: list[tuple[str, str]] = []
+        for entry in self._parameter_overrides:
+            key, sep, value = entry.partition("=")
+            if not sep or not key.strip():
+                raise ValueError(f"Invalid parameter_overrides entry '{entry}': expected 'Key=value text'.")
+            exact.append((key.strip(), value.strip()))
+        return per_token, exact
+
+    @staticmethod
+    def _apply_map_overrides(
+        text: str, per_token: dict[str, str], exact: list[tuple[str, str]], device_index: int
+    ) -> str:
+        """Patch a parameter map's text: set ImpactGPU to the device, apply exact key overrides, replace each
+        token of a per-token knob (preserving multiplicity), and warn for a requested key absent from the map.
+        """
+        entry_pattern = re.compile(r"^(\s*)\((\S+)((?:\s+[^)]*)?)\)\s*$")
+        requested = set(per_token) | {key for key, _ in exact}
+        seen: set[str] = set()
+        lines = []
+        for line in text.splitlines():
+            match = entry_pattern.match(line)
+            if match:
+                indent, key, values = match.group(1), match.group(2), match.group(3)
+                if key == "ImpactGPU":
+                    line = f"{indent}(ImpactGPU {device_index})"
+                else:
+                    exact_value = next((value for k, value in exact if k == key), None)
+                    if exact_value is not None:
+                        seen.add(key)
+                        line = f"{indent}({key} {exact_value})"
+                    else:
+                        token_key = "ImpactSubsetFeatures" if key.startswith("ImpactSubsetFeatures") else key
+                        if token_key in per_token:
+                            seen.add(token_key)
+                            replaced = " ".join(per_token[token_key] for _ in values.split())
+                            line = f"{indent}({key} {replaced})"
+            lines.append(line)
+        # Overrides never inject keys, so a knob set for a key absent from every map would silently do
+        # nothing — surface it (e.g. final_grid_spacing on a rigid-only preset).
+        for key in sorted(requested - seen):
+            print(f"[ImpactReg] note: override '{key}' matched no entry in the preset's parameter maps.")
+        return "\n".join(lines)
+
     def _stage_parameter_maps(self, work: Path, device_index: int) -> list[Path]:
-        """Copy parameter maps into the work dir, patching the ImpactGPU device field when present."""
+        """Stage the parameter maps into the work dir.
+
+        Matrix mode GENERATES each map from ``resolutions`` + the registry, then applies only the map-wide
+        knobs (grid spacing, spatial samples, exact overrides) — the matrix already sets iterations and
+        features per cell. Legacy mode copies the preset's maps and applies every per-token / exact override.
+        Both set the ImpactGPU device.
+        """
         staged = []
         for src in self._parameter_maps:
+            if self._resolutions:
+                text = generate_impact_parameter_map(
+                    src.read_text(encoding="utf-8"), self._resolutions, self._registry, self._mode
+                )
+                per_token, exact = self._parameter_map_overrides(global_only=True)
+            else:
+                text = src.read_text(encoding="utf-8")
+                per_token, exact = self._parameter_map_overrides()
+            text = self._apply_map_overrides(text, per_token, exact, device_index)
             dst = work / src.name
-            lines = []
-            for line in src.read_text(encoding="utf-8").splitlines():
-                if line.strip().startswith("(ImpactGPU"):
-                    line = f"(ImpactGPU {device_index})"
-                lines.append(line)
-            dst.write_text("\n".join(lines) + "\n", encoding="utf-8")
+            dst.write_text(text if text.endswith("\n") else text + "\n", encoding="utf-8")
             staged.append(dst)
         return staged
 
@@ -161,8 +484,10 @@ class ElastixEngine:
             captured: list[str] = []
             iteration_line = re.compile(r"^\d+\s")
             # ``iterations`` is the total iteration budget declared for the preset (summed over the
-            # chained parameter maps), so the bar spans the whole chain of registration stages.
-            progress = tqdm.tqdm(total=self._iterations or None, desc="Registration", ncols=0, leave=True)
+            # chained parameter maps), so the bar spans the whole chain of registration stages. A tuned
+            # ``max_iterations`` makes that declared budget stale — fall back to an open-ended bar.
+            budget = None if self._max_iterations > 0 else (self._iterations or None)
+            progress = tqdm.tqdm(total=budget, desc="Registration", ncols=0, leave=True)
             assert proc.stdout is not None
             resolution = 0
             for line in proc.stdout:
@@ -226,11 +551,33 @@ class ElastixRegistration(torch.nn.Module):
 
     accepts_attributes = True
 
-    def __init__(self, engine: str, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        engine: str,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         super().__init__()
         if engine != "elastix":
             raise NotImplementedError(f"ElastixRegistration engine '{engine}' is not implemented yet.")
-        self._engine = ElastixEngine(parameter_maps, models, iterations)
+        self._engine = ElastixEngine(
+            parameter_maps,
+            max_iterations,
+            final_grid_spacing,
+            subset_features,
+            spatial_samples,
+            parameter_overrides,
+            resolutions,
+            models_registry,
+            mode,
+        )
 
     def forward(
         self,
@@ -290,9 +637,23 @@ class RegistrationNet(network.Network):
         outputs_criterions: dict[str, network.TargetCriterionsLoader] = {"default": network.TargetCriterionsLoader()},
         engine: str = "elastix",
         parameter_maps: list[str] = [],
-        models: list[str] = [],
-        iterations: int = 1000,
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict[str, ResolutionSpec] = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
     ) -> None:
+        # The registration is fully described by the per-resolution model matrix ``resolutions`` (config =
+        # source of truth): each resolution lists its models, each model self-configured (ref, voxel_size,
+        # layers_mask, layers_weight, subset_features, pca, distance); intrinsic per-model props come from
+        # ``models_registry``. The feature-model download list is DERIVED from the matrix (no flat ``models``).
+        # Global knobs override the generated map: final_grid_spacing -> FinalGridSpacingInPhysicalUnits (mm),
+        # spatial_samples -> NumberOfSpatialSamples, parameter_overrides ('Key=value') -> any other entry.
+        # An empty ``resolutions`` = an intensity-only preset (no IMPACT models): the fixed maps are staged
+        # with just the global overrides. The total iteration count is derived (sum of per-resolution budgets).
         super().__init__(
             in_channels=1,
             optimizer=optimizer,
@@ -302,7 +663,18 @@ class RegistrationNet(network.Network):
         )
         self.add_module(
             "Registration",
-            ElastixRegistration(engine, parameter_maps, models, iterations),
+            ElastixRegistration(
+                engine,
+                parameter_maps,
+                max_iterations,
+                final_grid_spacing,
+                subset_features,
+                spatial_samples,
+                parameter_overrides,
+                resolutions,
+                models_registry,
+                mode,
+            ),
             in_branch=[0, 1, 2, 3],
             out_branch=["registration"],
         )

MR_CT_HeadNeck/ParameterMap_MRI_HN.txt

@@ -11,7 +11,7 @@
 
 
 (ImpactModelsPath0 "MIND/R1D2_3D.pt")
-(ImpactDimension0 3 )
+(ImpactDimension0 3)
 (ImpactNumberOfChannels0 1)
 (ImpactPatchSize0 0 0 0)
 (ImpactVoxelSize0 6 6 6)
@@ -24,7 +24,7 @@
 (ImpactModelsPath1 "MIND/R1D2_3D.pt")
 (ImpactDimension1 3)
 (ImpactNumberOfChannels1 1)
-(ImpactPatchSize1 0 0 0)  
+(ImpactPatchSize1 0 0 0)
 (ImpactVoxelSize1 3 3 3)
 (ImpactLayersMask1 "1")
 (ImpactSubsetFeatures1 32)
@@ -136,4 +136,4 @@
 (ResultImageFormat "mha")
 
 (ITKTransformOutputFileNameExtension "itk.txt")
-(WriteITKCompositeTransform "true")
+(WriteITKCompositeTransform "true")

MR_CT_HeadNeck/Prediction.yml

@@ -5,10 +5,75 @@ Predictor:
       engine: elastix
       parameter_maps:
       - ParameterMap_MRI_HN.txt
-      models:
-      - VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
-      iterations: 1050
       outputs_criterions: None
+      max_iterations: 0
+      final_grid_spacing: 0.0
+      subset_features: 0
+      spatial_samples: 0
+      parameter_overrides: []
+      resolutions:
+        '0':
+          max_iterations: 300
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 6.0
+              - 6.0
+              - 6.0
+              layers_mask: '1'
+              layers_weight:
+              - 1.0
+              subset_features: 32
+              pca: 0
+              distance: L1
+        '1':
+          max_iterations: 300
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 3.0
+              - 3.0
+              - 3.0
+              layers_mask: '1'
+              layers_weight:
+              - 1.0
+              subset_features: 32
+              pca: 0
+              distance: L1
+        '2':
+          max_iterations: 250
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 1.0
+              subset_features: 32
+              pca: 0
+              distance: L1
+        '3':
+          max_iterations: 200
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 1.0
+              subset_features: 32
+              pca: 0
+              distance: L1
+      models_registry: VBoussot/impact-torchscript-models:models.json
+      mode: Static
   Dataset:
     groups_src:
       Volume_0:

MR_CT_MRSeg/Model.py

@@ -32,6 +32,7 @@ NOTE: do NOT add ``from __future__ import annotations`` here — KonfAI's config
 runtime-evaluated annotations (``get_origin``); PEP 563 stringized annotations break arg resolution.
 """
 
+import json
 import os
 import re
 import shutil
@@ -52,6 +53,212 @@ from konfai.utils.dataset import Attribute, data_to_image, image_to_data
 # Set KONFAI_ELASTIX_DIR to point at an existing install and skip the download.
 ELASTIX_CACHE = Path.home() / ".cache" / "konfai" / "elastix-impact"
 
+# ---------------------------------------------------------------------------------------------------
+# Per-resolution model matrix (the config is the source of truth) -> generated IMPACT parameter map.
+# IMPACT field docs: https://github.com/vboussot/ImpactLoss/tree/main/ParameterMaps
+# The forced per-model props (dimension/channels/FOV formula) live in a registry (models.json on
+# VBoussot/impact-torchscript-models); the config carries the FREE knobs (which models per resolution,
+# feature voxel size, iterations, per-model layer weights/mask/subset/pca/distance) and the global
+# ``mode``. PatchSize follows ImpactMode: Static -> "0 0 0" (whole image); Jacobian -> the model FOV
+# evaluated from the registry formula (MIND 2*r*d+1, TS/MRSeg 2^l+3, SAM 29, DINOv2 14) as a cube.
+# ---------------------------------------------------------------------------------------------------
+
+_IMPACT_MODELS_REGISTRY = "VBoussot/impact-torchscript-models:models.json"
+
+# ``2^l+3`` grows with depth but the segmenters' receptive field plateaus: layers 7-8 share layer 6's
+# FOV (the "ramp max"). A config that deep should really run in Static (whole image) anyway; in Jacobian
+# we clamp ``l`` to this plateau so the patch stays finite and matches the real FOV.
+_FOV_RAMP_MAX_LAYER = 6
+
+
+def _num(x: object) -> str:
+    """Format a number the elastix way: integers without a trailing '.0' (6.0 -> '6', 0.2 -> '0.2')."""
+    return "%g" % float(x)
+
+
+class ModelSpec:
+    """One feature model at one resolution, with its OWN config (several models may share a resolution).
+
+    ``ref`` selects the model; ``voxel_size`` / ``layers_weight`` / ``subset_features`` / ``pca`` /
+    ``distance`` are its free per-(resolution, model) tuning knobs (the doc's per-model *tuning* fields).
+    The intrinsic per-model props — dimension, channels, ``layers_mask``, patch-size (FOV) — come from the
+    registry (read-only); ``layers_mask`` / ``distance`` left empty fall back to the registry default.
+    """
+
+    def __init__(
+        self,
+        ref: str,
+        voxel_size: list[float] = [],
+        layers_weight: list[float] = [1.0],
+        subset_features: int = 0,
+        pca: int = 0,
+        distance: str = "",
+        layers_mask: str = "",
+    ) -> None:
+        self.ref = ref
+        self.voxel_size = voxel_size
+        self.layers_weight = layers_weight
+        self.subset_features = subset_features
+        self.pca = pca
+        self.distance = distance
+        self.layers_mask = layers_mask
+
+
+class ResolutionSpec:
+    """One elastix resolution level: its iteration budget and the models compared there (each self-configured)."""
+
+    def __init__(self, max_iterations: int, models: dict[str, ModelSpec]) -> None:
+        self.max_iterations = max_iterations
+        self.models = models
+
+
+def _sorted_specs(mapping: dict) -> list:
+    """dict keyed by string indices ('0','1',...) -> values in numeric order (well-defined res/model order)."""
+    return [mapping[k] for k in sorted(mapping, key=lambda key: int(key))]
+
+
+def load_models_registry(ref: str = _IMPACT_MODELS_REGISTRY) -> dict:
+    """Load models.json (forced params per model) from the model repo on Hugging Face.
+
+    The registry is NOT bundled with the preset — it lives on the models repo and is fetched from there.
+    Resolution: the ``KONFAI_IMPACT_MODELS_REGISTRY`` env path wins (dev/offline); otherwise ``ref`` must be
+    a ``repo:file`` Hugging Face reference.
+    """
+    local = os.environ.get("KONFAI_IMPACT_MODELS_REGISTRY", "")
+    if local:
+        path = Path(local)
+    elif ":" in ref:
+        repo, filename = ref.split(":", 1)
+        path = Path(hf_hub_download(repo_id=repo, filename=filename, repo_type="model"))  # nosec B615
+    else:
+        raise ValueError(
+            f"models_registry '{ref}' must be a 'repo:file' Hugging Face reference (the registry is fetched "
+            f"from HF, not bundled) — or set KONFAI_IMPACT_MODELS_REGISTRY to a local file for offline use."
+        )
+    return json.loads(path.read_text(encoding="utf-8"))
+
+
+def _model_key(ref: str) -> str:
+    """Registry key / staged relative path = the model file within the models repo (strip a 'repo:' prefix)."""
+    return ref.split(":", 1)[1] if ":" in ref else ref
+
+
+def _deepest_active_layer(layers_mask: str) -> int:
+    """Deepest (largest-FOV) layer selected by ``layers_mask``, as a 0-based index read left-to-right.
+
+    A model returns its feature layers shallow->deep (``[layer_0, layer_1, ...]``, see the model repo's
+    build scripts); ``layers_mask`` has one char per returned layer, position ``i`` == ``layer_i``, ``'1'``
+    = selected. In Jacobian the patch must cover the receptive field of the DEEPEST selected layer, so the
+    FOV is governed by the rightmost ``'1'``.
+    """
+    mask = layers_mask.strip().strip('"')
+    active = [i for i, char in enumerate(mask) if char == "1"]
+    if not active:
+        raise ValueError(f"LayersMask '{layers_mask}' selects no layer; cannot derive the model FOV.")
+    return max(active)
+
+
+def _fov_value(fov: dict, layers_mask: str) -> int:
+    """Evaluate a model's field-of-view (in voxels) from its registry ``fov`` spec.
+
+    Supported formulas (from the model repo, https://huggingface.co/VBoussot/impact-torchscript-models):
+      ``2*r*d+1``  MIND, from the handcrafted radius ``r`` / dilation ``d`` (e.g. R1D2 -> 5);
+      ``2^l+3``    TotalSegmentator / MRSegmentator, ``l`` = the deepest layer picked by ``layers_mask``,
+                   clamped to the receptive-field plateau ``_FOV_RAMP_MAX_LAYER`` (layers 7-8 -> layer 6);
+      a bare int  a fixed FOV (SAM2.1 -> 29, DINOv2 -> 14);
+      ``Global``   Anatomix — whole-image only (Static); has no finite Jacobian patch -> error.
+    An explicit ``value`` in the spec is honoured as a precomputed shortcut when the formula needs none.
+    """
+    formula = str(fov.get("formula", "")).strip()
+    key = re.sub(r"\s+", "", formula).lower()
+    if key.isdigit():
+        return int(key)
+    if key == "2*r*d+1":
+        return 2 * int(fov["r"]) * int(fov["d"]) + 1
+    if key == "2^l+3":
+        return 2 ** min(_deepest_active_layer(layers_mask), _FOV_RAMP_MAX_LAYER) + 3
+    if "global" in key:
+        raise ValueError(f"model FOV '{formula}' is whole-image only (Static); it has no Jacobian patch size.")
+    if fov.get("value") is not None:
+        return int(fov["value"])
+    raise ValueError(f"cannot evaluate model FOV formula '{formula}'.")
+
+
+def _patch_size(mode: str, entry: dict, layers_mask: str) -> str:
+    """PatchSize from the model FOV, one token per model axis (2D model -> 2 tokens, 3D -> 3): Static ->
+    whole image (all zeros); Jacobian -> the evaluated FOV repeated over the axes. A 2D model mixed with a
+    3D one at a resolution concatenates as e.g. ``29 29 11 11 11`` (SAM 2D + TS 3D), matching IMPACT."""
+    dim = int(entry.get("dimension", 3))
+    if mode.strip().strip('"').lower() != "jacobian":
+        return " ".join(["0"] * dim)
+    fov = _fov_value(entry.get("fov", {}), layers_mask)
+    return " ".join([str(fov)] * dim)
+
+
+def generate_impact_parameter_map(
+    template_text: str, resolutions: dict, registry: dict, mode: str = "Static"
+) -> str:
+    """Rewrite the resolution-dependent lines of ``template_text`` from the model matrix ``resolutions``.
+
+    Regenerated: MaximumNumberOfIterations, NumberOfResolutions, Fixed/MovingImagePyramidRescaleSchedule,
+    ImpactMode (from the config ``mode``), and the whole ImpactXxxK block; every other template line is
+    kept verbatim (optimizer, transform, metric weights, components...). N (number of resolutions) is
+    deduced from the config. ``mode`` drives PatchSize: Static -> ``0 0 0`` (whole image); Jacobian -> the
+    per-model FOV evaluated from the registry formula and the cell's ``layers_mask``.
+    """
+    res = _sorted_specs(resolutions)
+    n = len(res)
+    mode_clean = mode.strip().strip('"') or "Static"
+
+    impact: list[str] = []
+    for k, r in enumerate(res):
+        models = _sorted_specs(r.models)
+        entries = [registry[_model_key(m.ref)] for m in models]
+
+        def row(stem: str, values: list[str]) -> None:
+            impact.append(f"(Impact{stem}{k} " + " ".join(values) + ")")
+
+        # From the registry (models.json on the model repo) ONLY the 3 truly model-fixed props:
+        # Dimension, NumberOfChannels, PatchSize (the model FOV). Everything else is a per-model tuning knob
+        # taken straight from the cell: VoxelSize / LayersMask / SubsetFeatures / PCA / Distance / LayersWeight.
+        row("ModelsPath", [f'"{_model_key(m.ref)}"' for m in models])
+        row("Dimension", [e["dimension"] for e in entries])
+        row("NumberOfChannels", [e["numberofchannels"] for e in entries])
+        row("PatchSize", [_patch_size(mode_clean, e, m.layers_mask) for e, m in zip(entries, models)])
+        row("VoxelSize", [" ".join(_num(v) for v in m.voxel_size) for m in models])
+        row("LayersMask", [f'"{m.layers_mask}"' for m in models])
+        row("SubsetFeatures", [str(m.subset_features) for m in models])
+        row("PCA", [str(m.pca) for m in models])
+        row("Distance", [f'"{m.distance}"' for m in models])
+        row("LayersWeight", [" ".join(_num(w) for w in m.layers_weight) for m in models])
+        impact.append("")  # blank line between resolutions, mirroring the reference maps
+
+    # The per-resolution block is the contiguous span from the first to the last ``Impact<name><k>`` line
+    # (the blank lines the reference maps put BETWEEN resolutions fall inside that span). Replace the whole
+    # span in one shot with the generated block, so the reference blanks are not kept on top of ours.
+    lines = template_text.splitlines()
+    indexed = [(re.match(r"^\s*\((\S+?)\s+(.*?)\)\s*$", ln), ln) for ln in lines]
+    block_rows = [i for i, (m, _) in enumerate(indexed) if m and re.match(r"^Impact[A-Za-z]+\d+$", m.group(1))]
+    block_lo, block_hi = (block_rows[0], block_rows[-1]) if block_rows else (-1, -2)
+
+    out: list[str] = []
+    for i, (m, line) in enumerate(indexed):
+        key = m.group(1) if m else None
+        if block_lo <= i <= block_hi:
+            if i == block_lo:  # replace the whole span at its first line, drop the rest (incl. inner blanks)
+                out.extend(impact[:-1])
+        elif key == "MaximumNumberOfIterations":
+            out.append("(MaximumNumberOfIterations " + " ".join(_num(r.max_iterations) for r in res) + ")")
+        elif key == "NumberOfResolutions":
+            out.append(f"(NumberOfResolutions {n})")
+        elif key in ("FixedImagePyramidRescaleSchedule", "MovingImagePyramidRescaleSchedule"):
+            out.append(f"({key} " + " ".join(["1"] * 3 * n) + ")")
+        elif key == "ImpactMode":
+            out.append(f'(ImpactMode "{mode_clean}")')
+        else:
+            out.append(line)
+    return "\n".join(out)
+
 
 class ElastixEngine:
     """Run the elastix-IMPACT binary on a fixed/moving pair; return (moved, dvf) on the fixed grid.
@@ -60,14 +267,57 @@ class ElastixEngine:
     does NOT ship it), so registration is a subprocess call, not ``sitk.ElastixImageFilter``.
     """
 
-    def __init__(self, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         self._bundle_dir = Path(__file__).resolve().parent
         self._parameter_maps = [self._bundle_dir / Path(p).name for p in parameter_maps]
+        self._max_iterations = max_iterations
+        self._final_grid_spacing = final_grid_spacing
+        self._subset_features = subset_features
+        self._spatial_samples = spatial_samples
+        self._parameter_overrides = list(parameter_overrides)
+        # ImpactMode: Static computes features once per level (PatchSize 0 0 0 = whole image); Jacobian
+        # samples random patches sized to the model FOV each iteration. Global knob: one mode per preset.
+        self._mode = mode
+        # Matrix mode: when `resolutions` is given the parameter map is GENERATED from it (the config is the
+        # source of truth). An empty `resolutions` = an intensity preset (no IMPACT feature models): the fixed
+        # parameter maps are staged with only the global knob overrides.
+        self._resolutions = resolutions
+        self._registry = load_models_registry(models_registry) if resolutions else {}
+        # The feature models are DERIVED — the unique refs across the matrix cells (no flat `models` param).
+        models: list[str] = []
+        for res in _sorted_specs(resolutions):
+            for model in _sorted_specs(res.models):
+                if model.ref not in models:
+                    models.append(model.ref)
         self._models = models
-        self._iterations = iterations
+        # `iterations` (the progress-bar total) is NOT a config parameter — it is DERIVED: the sum of the
+        # per-resolution iteration budgets, read from the matrix (matrix mode) or the maps (legacy).
+        self._iterations = self._total_iterations()
         self._elastix_bin = self._ensure_binary()
         self._local_models = self._download_models()
 
+    def _total_iterations(self) -> int:
+        """Total iterations across all resolutions — the progress-bar budget, derived from the config."""
+        if self._resolutions:
+            return sum(int(res.max_iterations) for res in _sorted_specs(self._resolutions))
+        total = 0
+        for src in self._parameter_maps:
+            match = re.search(r"\(MaximumNumberOfIterations\s+([^)]*)\)", src.read_text(encoding="utf-8"))
+            if match:
+                total += sum(int(token) for token in match.group(1).split())
+        return total
+
     def _ensure_binary(self) -> Path:
         # Optional override: point at an existing elastix-IMPACT install (skips the download).
         override = os.environ.get("KONFAI_ELASTIX_DIR", "")
@@ -91,17 +341,90 @@ class ElastixEngine:
             models.append((filename, local))
         return models
 
+    def _parameter_map_overrides(self, global_only: bool = False) -> tuple[dict[str, str], list[tuple[str, str]]]:
+        """The tuned knobs as parameter-map overrides: ``(per_token, exact)``.
+
+        ``per_token`` maps an elastix key (or the ``ImpactSubsetFeatures`` prefix) to a value that replaces
+        **each** existing token, so per-resolution / per-model multiplicity is preserved (e.g.
+        ``(MaximumNumberOfIterations 500 250)`` -> ``(MaximumNumberOfIterations 300 300)``). ``exact``
+        entries (from ``parameter_overrides``, ``Key=value text``) replace the whole value verbatim and win
+        over the named knobs. Overrides only REPLACE keys already present in a map — never inject new ones.
+        ``global_only`` (matrix mode) keeps just the map-wide knobs and drops ``max_iterations`` /
+        ``subset_features`` — the per-resolution matrix already sets those per cell.
+        """
+        per_token: dict[str, str] = {}
+        if not global_only and self._max_iterations > 0:
+            per_token["MaximumNumberOfIterations"] = str(int(self._max_iterations))
+        if self._final_grid_spacing > 0:
+            per_token["FinalGridSpacingInPhysicalUnits"] = str(float(self._final_grid_spacing))
+        if not global_only and self._subset_features > 0:
+            per_token["ImpactSubsetFeatures"] = str(int(self._subset_features))  # prefix: indexed per metric
+        if self._spatial_samples > 0:
+            per_token["NumberOfSpatialSamples"] = str(int(self._spatial_samples))
+        exact: list[tuple[str, str]] = []
+        for entry in self._parameter_overrides:
+            key, sep, value = entry.partition("=")
+            if not sep or not key.strip():
+                raise ValueError(f"Invalid parameter_overrides entry '{entry}': expected 'Key=value text'.")
+            exact.append((key.strip(), value.strip()))
+        return per_token, exact
+
+    @staticmethod
+    def _apply_map_overrides(
+        text: str, per_token: dict[str, str], exact: list[tuple[str, str]], device_index: int
+    ) -> str:
+        """Patch a parameter map's text: set ImpactGPU to the device, apply exact key overrides, replace each
+        token of a per-token knob (preserving multiplicity), and warn for a requested key absent from the map.
+        """
+        entry_pattern = re.compile(r"^(\s*)\((\S+)((?:\s+[^)]*)?)\)\s*$")
+        requested = set(per_token) | {key for key, _ in exact}
+        seen: set[str] = set()
+        lines = []
+        for line in text.splitlines():
+            match = entry_pattern.match(line)
+            if match:
+                indent, key, values = match.group(1), match.group(2), match.group(3)
+                if key == "ImpactGPU":
+                    line = f"{indent}(ImpactGPU {device_index})"
+                else:
+                    exact_value = next((value for k, value in exact if k == key), None)
+                    if exact_value is not None:
+                        seen.add(key)
+                        line = f"{indent}({key} {exact_value})"
+                    else:
+                        token_key = "ImpactSubsetFeatures" if key.startswith("ImpactSubsetFeatures") else key
+                        if token_key in per_token:
+                            seen.add(token_key)
+                            replaced = " ".join(per_token[token_key] for _ in values.split())
+                            line = f"{indent}({key} {replaced})"
+            lines.append(line)
+        # Overrides never inject keys, so a knob set for a key absent from every map would silently do
+        # nothing — surface it (e.g. final_grid_spacing on a rigid-only preset).
+        for key in sorted(requested - seen):
+            print(f"[ImpactReg] note: override '{key}' matched no entry in the preset's parameter maps.")
+        return "\n".join(lines)
+
     def _stage_parameter_maps(self, work: Path, device_index: int) -> list[Path]:
-        """Copy parameter maps into the work dir, patching the ImpactGPU device field when present."""
+        """Stage the parameter maps into the work dir.
+
+        Matrix mode GENERATES each map from ``resolutions`` + the registry, then applies only the map-wide
+        knobs (grid spacing, spatial samples, exact overrides) — the matrix already sets iterations and
+        features per cell. Legacy mode copies the preset's maps and applies every per-token / exact override.
+        Both set the ImpactGPU device.
+        """
         staged = []
         for src in self._parameter_maps:
+            if self._resolutions:
+                text = generate_impact_parameter_map(
+                    src.read_text(encoding="utf-8"), self._resolutions, self._registry, self._mode
+                )
+                per_token, exact = self._parameter_map_overrides(global_only=True)
+            else:
+                text = src.read_text(encoding="utf-8")
+                per_token, exact = self._parameter_map_overrides()
+            text = self._apply_map_overrides(text, per_token, exact, device_index)
             dst = work / src.name
-            lines = []
-            for line in src.read_text(encoding="utf-8").splitlines():
-                if line.strip().startswith("(ImpactGPU"):
-                    line = f"(ImpactGPU {device_index})"
-                lines.append(line)
-            dst.write_text("\n".join(lines) + "\n", encoding="utf-8")
+            dst.write_text(text if text.endswith("\n") else text + "\n", encoding="utf-8")
             staged.append(dst)
         return staged
 
@@ -161,8 +484,10 @@ class ElastixEngine:
             captured: list[str] = []
             iteration_line = re.compile(r"^\d+\s")
             # ``iterations`` is the total iteration budget declared for the preset (summed over the
-            # chained parameter maps), so the bar spans the whole chain of registration stages.
-            progress = tqdm.tqdm(total=self._iterations or None, desc="Registration", ncols=0, leave=True)
+            # chained parameter maps), so the bar spans the whole chain of registration stages. A tuned
+            # ``max_iterations`` makes that declared budget stale — fall back to an open-ended bar.
+            budget = None if self._max_iterations > 0 else (self._iterations or None)
+            progress = tqdm.tqdm(total=budget, desc="Registration", ncols=0, leave=True)
             assert proc.stdout is not None
             resolution = 0
             for line in proc.stdout:
@@ -226,11 +551,33 @@ class ElastixRegistration(torch.nn.Module):
 
     accepts_attributes = True
 
-    def __init__(self, engine: str, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        engine: str,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         super().__init__()
         if engine != "elastix":
             raise NotImplementedError(f"ElastixRegistration engine '{engine}' is not implemented yet.")
-        self._engine = ElastixEngine(parameter_maps, models, iterations)
+        self._engine = ElastixEngine(
+            parameter_maps,
+            max_iterations,
+            final_grid_spacing,
+            subset_features,
+            spatial_samples,
+            parameter_overrides,
+            resolutions,
+            models_registry,
+            mode,
+        )
 
     def forward(
         self,
@@ -290,9 +637,23 @@ class RegistrationNet(network.Network):
         outputs_criterions: dict[str, network.TargetCriterionsLoader] = {"default": network.TargetCriterionsLoader()},
         engine: str = "elastix",
         parameter_maps: list[str] = [],
-        models: list[str] = [],
-        iterations: int = 1000,
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict[str, ResolutionSpec] = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
     ) -> None:
+        # The registration is fully described by the per-resolution model matrix ``resolutions`` (config =
+        # source of truth): each resolution lists its models, each model self-configured (ref, voxel_size,
+        # layers_mask, layers_weight, subset_features, pca, distance); intrinsic per-model props come from
+        # ``models_registry``. The feature-model download list is DERIVED from the matrix (no flat ``models``).
+        # Global knobs override the generated map: final_grid_spacing -> FinalGridSpacingInPhysicalUnits (mm),
+        # spatial_samples -> NumberOfSpatialSamples, parameter_overrides ('Key=value') -> any other entry.
+        # An empty ``resolutions`` = an intensity-only preset (no IMPACT models): the fixed maps are staged
+        # with just the global overrides. The total iteration count is derived (sum of per-resolution budgets).
         super().__init__(
             in_channels=1,
             optimizer=optimizer,
@@ -302,7 +663,18 @@ class RegistrationNet(network.Network):
         )
         self.add_module(
             "Registration",
-            ElastixRegistration(engine, parameter_maps, models, iterations),
+            ElastixRegistration(
+                engine,
+                parameter_maps,
+                max_iterations,
+                final_grid_spacing,
+                subset_features,
+                spatial_samples,
+                parameter_overrides,
+                resolutions,
+                models_registry,
+                mode,
+            ),
             in_branch=[0, 1, 2, 3],
             out_branch=["registration"],
         )

MR_CT_MRSeg/ParameterMap_MRI_MRSeg.txt

@@ -24,7 +24,7 @@
 (ImpactModelsPath1 "MIND/R1D2_3D.pt" "MRSeg/MRSeg.pt")
 (ImpactDimension1 3 3)
 (ImpactNumberOfChannels1 1 1)
-(ImpactPatchSize1 0 0 0 0 0 0)  
+(ImpactPatchSize1 0 0 0 0 0 0)
 (ImpactVoxelSize1 3 3 3 3 3 3)
 (ImpactLayersMask1 "1" "1")
 (ImpactSubsetFeatures1 32 64)
@@ -37,7 +37,7 @@
 (ImpactNumberOfChannels2 1 1)
 (ImpactPatchSize2 0 0 0 0 0 0)
 (ImpactVoxelSize2 2 2 2 2 2 2)
-(ImpactLayersMask2 "1"  "1")
+(ImpactLayersMask2 "1" "1")
 (ImpactSubsetFeatures2 32 64)
 (ImpactPCA2 0 0)
 (ImpactDistance2 "L1" "Dice")
@@ -48,7 +48,7 @@
 (ImpactNumberOfChannels3 1 1)
 (ImpactPatchSize3 0 0 0 0 0 0)
 (ImpactVoxelSize3 2 2 2 2 2 2)
-(ImpactLayersMask3 "1"  "1")
+(ImpactLayersMask3 "1" "1")
 (ImpactSubsetFeatures3 32 64)
 (ImpactPCA3 0 0)
 (ImpactDistance3 "L1" "Dice")

MR_CT_MRSeg/Prediction.yml

@@ -5,11 +5,123 @@ Predictor:
       engine: elastix
       parameter_maps:
       - ParameterMap_MRI_MRSeg.txt
-      models:
-      - VBoussot/impact-torchscript-models:MRSeg/MRSeg.pt
-      - VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
-      iterations: 1100
       outputs_criterions: None
+      max_iterations: 0
+      final_grid_spacing: 0.0
+      subset_features: 0
+      spatial_samples: 0
+      parameter_overrides: []
+      resolutions:
+        '0':
+          max_iterations: 400
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 6.0
+              - 6.0
+              - 6.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.2
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:MRSeg/MRSeg.pt
+              voxel_size:
+              - 6.0
+              - 6.0
+              - 6.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.8
+              subset_features: 64
+              pca: 0
+              distance: Dice
+        '1':
+          max_iterations: 300
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 3.0
+              - 3.0
+              - 3.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.3
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:MRSeg/MRSeg.pt
+              voxel_size:
+              - 3.0
+              - 3.0
+              - 3.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.7
+              subset_features: 64
+              pca: 0
+              distance: Dice
+        '2':
+          max_iterations: 200
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.6
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:MRSeg/MRSeg.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.4
+              subset_features: 64
+              pca: 0
+              distance: Dice
+        '3':
+          max_iterations: 200
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.7
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:MRSeg/MRSeg.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.3
+              subset_features: 64
+              pca: 0
+              distance: Dice
+      models_registry: VBoussot/impact-torchscript-models:models.json
+      mode: Static
   Dataset:
     groups_src:
       Volume_0:

MR_CT_TS/Model.py

@@ -32,6 +32,7 @@ NOTE: do NOT add ``from __future__ import annotations`` here — KonfAI's config
 runtime-evaluated annotations (``get_origin``); PEP 563 stringized annotations break arg resolution.
 """
 
+import json
 import os
 import re
 import shutil
@@ -52,6 +53,212 @@ from konfai.utils.dataset import Attribute, data_to_image, image_to_data
 # Set KONFAI_ELASTIX_DIR to point at an existing install and skip the download.
 ELASTIX_CACHE = Path.home() / ".cache" / "konfai" / "elastix-impact"
 
+# ---------------------------------------------------------------------------------------------------
+# Per-resolution model matrix (the config is the source of truth) -> generated IMPACT parameter map.
+# IMPACT field docs: https://github.com/vboussot/ImpactLoss/tree/main/ParameterMaps
+# The forced per-model props (dimension/channels/FOV formula) live in a registry (models.json on
+# VBoussot/impact-torchscript-models); the config carries the FREE knobs (which models per resolution,
+# feature voxel size, iterations, per-model layer weights/mask/subset/pca/distance) and the global
+# ``mode``. PatchSize follows ImpactMode: Static -> "0 0 0" (whole image); Jacobian -> the model FOV
+# evaluated from the registry formula (MIND 2*r*d+1, TS/MRSeg 2^l+3, SAM 29, DINOv2 14) as a cube.
+# ---------------------------------------------------------------------------------------------------
+
+_IMPACT_MODELS_REGISTRY = "VBoussot/impact-torchscript-models:models.json"
+
+# ``2^l+3`` grows with depth but the segmenters' receptive field plateaus: layers 7-8 share layer 6's
+# FOV (the "ramp max"). A config that deep should really run in Static (whole image) anyway; in Jacobian
+# we clamp ``l`` to this plateau so the patch stays finite and matches the real FOV.
+_FOV_RAMP_MAX_LAYER = 6
+
+
+def _num(x: object) -> str:
+    """Format a number the elastix way: integers without a trailing '.0' (6.0 -> '6', 0.2 -> '0.2')."""
+    return "%g" % float(x)
+
+
+class ModelSpec:
+    """One feature model at one resolution, with its OWN config (several models may share a resolution).
+
+    ``ref`` selects the model; ``voxel_size`` / ``layers_weight`` / ``subset_features`` / ``pca`` /
+    ``distance`` are its free per-(resolution, model) tuning knobs (the doc's per-model *tuning* fields).
+    The intrinsic per-model props — dimension, channels, ``layers_mask``, patch-size (FOV) — come from the
+    registry (read-only); ``layers_mask`` / ``distance`` left empty fall back to the registry default.
+    """
+
+    def __init__(
+        self,
+        ref: str,
+        voxel_size: list[float] = [],
+        layers_weight: list[float] = [1.0],
+        subset_features: int = 0,
+        pca: int = 0,
+        distance: str = "",
+        layers_mask: str = "",
+    ) -> None:
+        self.ref = ref
+        self.voxel_size = voxel_size
+        self.layers_weight = layers_weight
+        self.subset_features = subset_features
+        self.pca = pca
+        self.distance = distance
+        self.layers_mask = layers_mask
+
+
+class ResolutionSpec:
+    """One elastix resolution level: its iteration budget and the models compared there (each self-configured)."""
+
+    def __init__(self, max_iterations: int, models: dict[str, ModelSpec]) -> None:
+        self.max_iterations = max_iterations
+        self.models = models
+
+
+def _sorted_specs(mapping: dict) -> list:
+    """dict keyed by string indices ('0','1',...) -> values in numeric order (well-defined res/model order)."""
+    return [mapping[k] for k in sorted(mapping, key=lambda key: int(key))]
+
+
+def load_models_registry(ref: str = _IMPACT_MODELS_REGISTRY) -> dict:
+    """Load models.json (forced params per model) from the model repo on Hugging Face.
+
+    The registry is NOT bundled with the preset — it lives on the models repo and is fetched from there.
+    Resolution: the ``KONFAI_IMPACT_MODELS_REGISTRY`` env path wins (dev/offline); otherwise ``ref`` must be
+    a ``repo:file`` Hugging Face reference.
+    """
+    local = os.environ.get("KONFAI_IMPACT_MODELS_REGISTRY", "")
+    if local:
+        path = Path(local)
+    elif ":" in ref:
+        repo, filename = ref.split(":", 1)
+        path = Path(hf_hub_download(repo_id=repo, filename=filename, repo_type="model"))  # nosec B615
+    else:
+        raise ValueError(
+            f"models_registry '{ref}' must be a 'repo:file' Hugging Face reference (the registry is fetched "
+            f"from HF, not bundled) — or set KONFAI_IMPACT_MODELS_REGISTRY to a local file for offline use."
+        )
+    return json.loads(path.read_text(encoding="utf-8"))
+
+
+def _model_key(ref: str) -> str:
+    """Registry key / staged relative path = the model file within the models repo (strip a 'repo:' prefix)."""
+    return ref.split(":", 1)[1] if ":" in ref else ref
+
+
+def _deepest_active_layer(layers_mask: str) -> int:
+    """Deepest (largest-FOV) layer selected by ``layers_mask``, as a 0-based index read left-to-right.
+
+    A model returns its feature layers shallow->deep (``[layer_0, layer_1, ...]``, see the model repo's
+    build scripts); ``layers_mask`` has one char per returned layer, position ``i`` == ``layer_i``, ``'1'``
+    = selected. In Jacobian the patch must cover the receptive field of the DEEPEST selected layer, so the
+    FOV is governed by the rightmost ``'1'``.
+    """
+    mask = layers_mask.strip().strip('"')
+    active = [i for i, char in enumerate(mask) if char == "1"]
+    if not active:
+        raise ValueError(f"LayersMask '{layers_mask}' selects no layer; cannot derive the model FOV.")
+    return max(active)
+
+
+def _fov_value(fov: dict, layers_mask: str) -> int:
+    """Evaluate a model's field-of-view (in voxels) from its registry ``fov`` spec.
+
+    Supported formulas (from the model repo, https://huggingface.co/VBoussot/impact-torchscript-models):
+      ``2*r*d+1``  MIND, from the handcrafted radius ``r`` / dilation ``d`` (e.g. R1D2 -> 5);
+      ``2^l+3``    TotalSegmentator / MRSegmentator, ``l`` = the deepest layer picked by ``layers_mask``,
+                   clamped to the receptive-field plateau ``_FOV_RAMP_MAX_LAYER`` (layers 7-8 -> layer 6);
+      a bare int  a fixed FOV (SAM2.1 -> 29, DINOv2 -> 14);
+      ``Global``   Anatomix — whole-image only (Static); has no finite Jacobian patch -> error.
+    An explicit ``value`` in the spec is honoured as a precomputed shortcut when the formula needs none.
+    """
+    formula = str(fov.get("formula", "")).strip()
+    key = re.sub(r"\s+", "", formula).lower()
+    if key.isdigit():
+        return int(key)
+    if key == "2*r*d+1":
+        return 2 * int(fov["r"]) * int(fov["d"]) + 1
+    if key == "2^l+3":
+        return 2 ** min(_deepest_active_layer(layers_mask), _FOV_RAMP_MAX_LAYER) + 3
+    if "global" in key:
+        raise ValueError(f"model FOV '{formula}' is whole-image only (Static); it has no Jacobian patch size.")
+    if fov.get("value") is not None:
+        return int(fov["value"])
+    raise ValueError(f"cannot evaluate model FOV formula '{formula}'.")
+
+
+def _patch_size(mode: str, entry: dict, layers_mask: str) -> str:
+    """PatchSize from the model FOV, one token per model axis (2D model -> 2 tokens, 3D -> 3): Static ->
+    whole image (all zeros); Jacobian -> the evaluated FOV repeated over the axes. A 2D model mixed with a
+    3D one at a resolution concatenates as e.g. ``29 29 11 11 11`` (SAM 2D + TS 3D), matching IMPACT."""
+    dim = int(entry.get("dimension", 3))
+    if mode.strip().strip('"').lower() != "jacobian":
+        return " ".join(["0"] * dim)
+    fov = _fov_value(entry.get("fov", {}), layers_mask)
+    return " ".join([str(fov)] * dim)
+
+
+def generate_impact_parameter_map(
+    template_text: str, resolutions: dict, registry: dict, mode: str = "Static"
+) -> str:
+    """Rewrite the resolution-dependent lines of ``template_text`` from the model matrix ``resolutions``.
+
+    Regenerated: MaximumNumberOfIterations, NumberOfResolutions, Fixed/MovingImagePyramidRescaleSchedule,
+    ImpactMode (from the config ``mode``), and the whole ImpactXxxK block; every other template line is
+    kept verbatim (optimizer, transform, metric weights, components...). N (number of resolutions) is
+    deduced from the config. ``mode`` drives PatchSize: Static -> ``0 0 0`` (whole image); Jacobian -> the
+    per-model FOV evaluated from the registry formula and the cell's ``layers_mask``.
+    """
+    res = _sorted_specs(resolutions)
+    n = len(res)
+    mode_clean = mode.strip().strip('"') or "Static"
+
+    impact: list[str] = []
+    for k, r in enumerate(res):
+        models = _sorted_specs(r.models)
+        entries = [registry[_model_key(m.ref)] for m in models]
+
+        def row(stem: str, values: list[str]) -> None:
+            impact.append(f"(Impact{stem}{k} " + " ".join(values) + ")")
+
+        # From the registry (models.json on the model repo) ONLY the 3 truly model-fixed props:
+        # Dimension, NumberOfChannels, PatchSize (the model FOV). Everything else is a per-model tuning knob
+        # taken straight from the cell: VoxelSize / LayersMask / SubsetFeatures / PCA / Distance / LayersWeight.
+        row("ModelsPath", [f'"{_model_key(m.ref)}"' for m in models])
+        row("Dimension", [e["dimension"] for e in entries])
+        row("NumberOfChannels", [e["numberofchannels"] for e in entries])
+        row("PatchSize", [_patch_size(mode_clean, e, m.layers_mask) for e, m in zip(entries, models)])
+        row("VoxelSize", [" ".join(_num(v) for v in m.voxel_size) for m in models])
+        row("LayersMask", [f'"{m.layers_mask}"' for m in models])
+        row("SubsetFeatures", [str(m.subset_features) for m in models])
+        row("PCA", [str(m.pca) for m in models])
+        row("Distance", [f'"{m.distance}"' for m in models])
+        row("LayersWeight", [" ".join(_num(w) for w in m.layers_weight) for m in models])
+        impact.append("")  # blank line between resolutions, mirroring the reference maps
+
+    # The per-resolution block is the contiguous span from the first to the last ``Impact<name><k>`` line
+    # (the blank lines the reference maps put BETWEEN resolutions fall inside that span). Replace the whole
+    # span in one shot with the generated block, so the reference blanks are not kept on top of ours.
+    lines = template_text.splitlines()
+    indexed = [(re.match(r"^\s*\((\S+?)\s+(.*?)\)\s*$", ln), ln) for ln in lines]
+    block_rows = [i for i, (m, _) in enumerate(indexed) if m and re.match(r"^Impact[A-Za-z]+\d+$", m.group(1))]
+    block_lo, block_hi = (block_rows[0], block_rows[-1]) if block_rows else (-1, -2)
+
+    out: list[str] = []
+    for i, (m, line) in enumerate(indexed):
+        key = m.group(1) if m else None
+        if block_lo <= i <= block_hi:
+            if i == block_lo:  # replace the whole span at its first line, drop the rest (incl. inner blanks)
+                out.extend(impact[:-1])
+        elif key == "MaximumNumberOfIterations":
+            out.append("(MaximumNumberOfIterations " + " ".join(_num(r.max_iterations) for r in res) + ")")
+        elif key == "NumberOfResolutions":
+            out.append(f"(NumberOfResolutions {n})")
+        elif key in ("FixedImagePyramidRescaleSchedule", "MovingImagePyramidRescaleSchedule"):
+            out.append(f"({key} " + " ".join(["1"] * 3 * n) + ")")
+        elif key == "ImpactMode":
+            out.append(f'(ImpactMode "{mode_clean}")')
+        else:
+            out.append(line)
+    return "\n".join(out)
+
 
 class ElastixEngine:
     """Run the elastix-IMPACT binary on a fixed/moving pair; return (moved, dvf) on the fixed grid.
@@ -60,14 +267,57 @@ class ElastixEngine:
     does NOT ship it), so registration is a subprocess call, not ``sitk.ElastixImageFilter``.
     """
 
-    def __init__(self, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         self._bundle_dir = Path(__file__).resolve().parent
         self._parameter_maps = [self._bundle_dir / Path(p).name for p in parameter_maps]
+        self._max_iterations = max_iterations
+        self._final_grid_spacing = final_grid_spacing
+        self._subset_features = subset_features
+        self._spatial_samples = spatial_samples
+        self._parameter_overrides = list(parameter_overrides)
+        # ImpactMode: Static computes features once per level (PatchSize 0 0 0 = whole image); Jacobian
+        # samples random patches sized to the model FOV each iteration. Global knob: one mode per preset.
+        self._mode = mode
+        # Matrix mode: when `resolutions` is given the parameter map is GENERATED from it (the config is the
+        # source of truth). An empty `resolutions` = an intensity preset (no IMPACT feature models): the fixed
+        # parameter maps are staged with only the global knob overrides.
+        self._resolutions = resolutions
+        self._registry = load_models_registry(models_registry) if resolutions else {}
+        # The feature models are DERIVED — the unique refs across the matrix cells (no flat `models` param).
+        models: list[str] = []
+        for res in _sorted_specs(resolutions):
+            for model in _sorted_specs(res.models):
+                if model.ref not in models:
+                    models.append(model.ref)
         self._models = models
-        self._iterations = iterations
+        # `iterations` (the progress-bar total) is NOT a config parameter — it is DERIVED: the sum of the
+        # per-resolution iteration budgets, read from the matrix (matrix mode) or the maps (legacy).
+        self._iterations = self._total_iterations()
         self._elastix_bin = self._ensure_binary()
         self._local_models = self._download_models()
 
+    def _total_iterations(self) -> int:
+        """Total iterations across all resolutions — the progress-bar budget, derived from the config."""
+        if self._resolutions:
+            return sum(int(res.max_iterations) for res in _sorted_specs(self._resolutions))
+        total = 0
+        for src in self._parameter_maps:
+            match = re.search(r"\(MaximumNumberOfIterations\s+([^)]*)\)", src.read_text(encoding="utf-8"))
+            if match:
+                total += sum(int(token) for token in match.group(1).split())
+        return total
+
     def _ensure_binary(self) -> Path:
         # Optional override: point at an existing elastix-IMPACT install (skips the download).
         override = os.environ.get("KONFAI_ELASTIX_DIR", "")
@@ -91,17 +341,90 @@ class ElastixEngine:
             models.append((filename, local))
         return models
 
+    def _parameter_map_overrides(self, global_only: bool = False) -> tuple[dict[str, str], list[tuple[str, str]]]:
+        """The tuned knobs as parameter-map overrides: ``(per_token, exact)``.
+
+        ``per_token`` maps an elastix key (or the ``ImpactSubsetFeatures`` prefix) to a value that replaces
+        **each** existing token, so per-resolution / per-model multiplicity is preserved (e.g.
+        ``(MaximumNumberOfIterations 500 250)`` -> ``(MaximumNumberOfIterations 300 300)``). ``exact``
+        entries (from ``parameter_overrides``, ``Key=value text``) replace the whole value verbatim and win
+        over the named knobs. Overrides only REPLACE keys already present in a map — never inject new ones.
+        ``global_only`` (matrix mode) keeps just the map-wide knobs and drops ``max_iterations`` /
+        ``subset_features`` — the per-resolution matrix already sets those per cell.
+        """
+        per_token: dict[str, str] = {}
+        if not global_only and self._max_iterations > 0:
+            per_token["MaximumNumberOfIterations"] = str(int(self._max_iterations))
+        if self._final_grid_spacing > 0:
+            per_token["FinalGridSpacingInPhysicalUnits"] = str(float(self._final_grid_spacing))
+        if not global_only and self._subset_features > 0:
+            per_token["ImpactSubsetFeatures"] = str(int(self._subset_features))  # prefix: indexed per metric
+        if self._spatial_samples > 0:
+            per_token["NumberOfSpatialSamples"] = str(int(self._spatial_samples))
+        exact: list[tuple[str, str]] = []
+        for entry in self._parameter_overrides:
+            key, sep, value = entry.partition("=")
+            if not sep or not key.strip():
+                raise ValueError(f"Invalid parameter_overrides entry '{entry}': expected 'Key=value text'.")
+            exact.append((key.strip(), value.strip()))
+        return per_token, exact
+
+    @staticmethod
+    def _apply_map_overrides(
+        text: str, per_token: dict[str, str], exact: list[tuple[str, str]], device_index: int
+    ) -> str:
+        """Patch a parameter map's text: set ImpactGPU to the device, apply exact key overrides, replace each
+        token of a per-token knob (preserving multiplicity), and warn for a requested key absent from the map.
+        """
+        entry_pattern = re.compile(r"^(\s*)\((\S+)((?:\s+[^)]*)?)\)\s*$")
+        requested = set(per_token) | {key for key, _ in exact}
+        seen: set[str] = set()
+        lines = []
+        for line in text.splitlines():
+            match = entry_pattern.match(line)
+            if match:
+                indent, key, values = match.group(1), match.group(2), match.group(3)
+                if key == "ImpactGPU":
+                    line = f"{indent}(ImpactGPU {device_index})"
+                else:
+                    exact_value = next((value for k, value in exact if k == key), None)
+                    if exact_value is not None:
+                        seen.add(key)
+                        line = f"{indent}({key} {exact_value})"
+                    else:
+                        token_key = "ImpactSubsetFeatures" if key.startswith("ImpactSubsetFeatures") else key
+                        if token_key in per_token:
+                            seen.add(token_key)
+                            replaced = " ".join(per_token[token_key] for _ in values.split())
+                            line = f"{indent}({key} {replaced})"
+            lines.append(line)
+        # Overrides never inject keys, so a knob set for a key absent from every map would silently do
+        # nothing — surface it (e.g. final_grid_spacing on a rigid-only preset).
+        for key in sorted(requested - seen):
+            print(f"[ImpactReg] note: override '{key}' matched no entry in the preset's parameter maps.")
+        return "\n".join(lines)
+
     def _stage_parameter_maps(self, work: Path, device_index: int) -> list[Path]:
-        """Copy parameter maps into the work dir, patching the ImpactGPU device field when present."""
+        """Stage the parameter maps into the work dir.
+
+        Matrix mode GENERATES each map from ``resolutions`` + the registry, then applies only the map-wide
+        knobs (grid spacing, spatial samples, exact overrides) — the matrix already sets iterations and
+        features per cell. Legacy mode copies the preset's maps and applies every per-token / exact override.
+        Both set the ImpactGPU device.
+        """
         staged = []
         for src in self._parameter_maps:
+            if self._resolutions:
+                text = generate_impact_parameter_map(
+                    src.read_text(encoding="utf-8"), self._resolutions, self._registry, self._mode
+                )
+                per_token, exact = self._parameter_map_overrides(global_only=True)
+            else:
+                text = src.read_text(encoding="utf-8")
+                per_token, exact = self._parameter_map_overrides()
+            text = self._apply_map_overrides(text, per_token, exact, device_index)
             dst = work / src.name
-            lines = []
-            for line in src.read_text(encoding="utf-8").splitlines():
-                if line.strip().startswith("(ImpactGPU"):
-                    line = f"(ImpactGPU {device_index})"
-                lines.append(line)
-            dst.write_text("\n".join(lines) + "\n", encoding="utf-8")
+            dst.write_text(text if text.endswith("\n") else text + "\n", encoding="utf-8")
             staged.append(dst)
         return staged
 
@@ -161,8 +484,10 @@ class ElastixEngine:
             captured: list[str] = []
             iteration_line = re.compile(r"^\d+\s")
             # ``iterations`` is the total iteration budget declared for the preset (summed over the
-            # chained parameter maps), so the bar spans the whole chain of registration stages.
-            progress = tqdm.tqdm(total=self._iterations or None, desc="Registration", ncols=0, leave=True)
+            # chained parameter maps), so the bar spans the whole chain of registration stages. A tuned
+            # ``max_iterations`` makes that declared budget stale — fall back to an open-ended bar.
+            budget = None if self._max_iterations > 0 else (self._iterations or None)
+            progress = tqdm.tqdm(total=budget, desc="Registration", ncols=0, leave=True)
             assert proc.stdout is not None
             resolution = 0
             for line in proc.stdout:
@@ -226,11 +551,33 @@ class ElastixRegistration(torch.nn.Module):
 
     accepts_attributes = True
 
-    def __init__(self, engine: str, parameter_maps: list[str], models: list[str], iterations: int) -> None:
+    def __init__(
+        self,
+        engine: str,
+        parameter_maps: list[str],
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
+    ) -> None:
         super().__init__()
         if engine != "elastix":
             raise NotImplementedError(f"ElastixRegistration engine '{engine}' is not implemented yet.")
-        self._engine = ElastixEngine(parameter_maps, models, iterations)
+        self._engine = ElastixEngine(
+            parameter_maps,
+            max_iterations,
+            final_grid_spacing,
+            subset_features,
+            spatial_samples,
+            parameter_overrides,
+            resolutions,
+            models_registry,
+            mode,
+        )
 
     def forward(
         self,
@@ -290,9 +637,23 @@ class RegistrationNet(network.Network):
         outputs_criterions: dict[str, network.TargetCriterionsLoader] = {"default": network.TargetCriterionsLoader()},
         engine: str = "elastix",
         parameter_maps: list[str] = [],
-        models: list[str] = [],
-        iterations: int = 1000,
+        max_iterations: int = 0,
+        final_grid_spacing: float = 0.0,
+        subset_features: int = 0,
+        spatial_samples: int = 0,
+        parameter_overrides: list[str] = [],
+        resolutions: dict[str, ResolutionSpec] = {},
+        models_registry: str = _IMPACT_MODELS_REGISTRY,
+        mode: str = "Static",
     ) -> None:
+        # The registration is fully described by the per-resolution model matrix ``resolutions`` (config =
+        # source of truth): each resolution lists its models, each model self-configured (ref, voxel_size,
+        # layers_mask, layers_weight, subset_features, pca, distance); intrinsic per-model props come from
+        # ``models_registry``. The feature-model download list is DERIVED from the matrix (no flat ``models``).
+        # Global knobs override the generated map: final_grid_spacing -> FinalGridSpacingInPhysicalUnits (mm),
+        # spatial_samples -> NumberOfSpatialSamples, parameter_overrides ('Key=value') -> any other entry.
+        # An empty ``resolutions`` = an intensity-only preset (no IMPACT models): the fixed maps are staged
+        # with just the global overrides. The total iteration count is derived (sum of per-resolution budgets).
         super().__init__(
             in_channels=1,
             optimizer=optimizer,
@@ -302,7 +663,18 @@ class RegistrationNet(network.Network):
         )
         self.add_module(
             "Registration",
-            ElastixRegistration(engine, parameter_maps, models, iterations),
+            ElastixRegistration(
+                engine,
+                parameter_maps,
+                max_iterations,
+                final_grid_spacing,
+                subset_features,
+                spatial_samples,
+                parameter_overrides,
+                resolutions,
+                models_registry,
+                mode,
+            ),
             in_branch=[0, 1, 2, 3],
             out_branch=["registration"],
         )

MR_CT_TS/ParameterMap_MRI_TS.txt

@@ -24,7 +24,7 @@
 (ImpactModelsPath1 "MIND/R1D2_3D.pt" "TS/M852.pt")
 (ImpactDimension1 3 3)
 (ImpactNumberOfChannels1 1 1)
-(ImpactPatchSize1 0 0 0 0 0 0)  
+(ImpactPatchSize1 0 0 0 0 0 0)
 (ImpactVoxelSize1 3 3 3 3 3 3)
 (ImpactLayersMask1 "1" "0000001")
 (ImpactSubsetFeatures1 32 64)
@@ -37,7 +37,7 @@
 (ImpactNumberOfChannels2 1 1)
 (ImpactPatchSize2 0 0 0 0 0 0)
 (ImpactVoxelSize2 2 2 2 2 2 2)
-(ImpactLayersMask2 "1"  "00000001")
+(ImpactLayersMask2 "1" "00000001")
 (ImpactSubsetFeatures2 32 64)
 (ImpactPCA2 0 0)
 (ImpactDistance2 "L1" "Dice")
@@ -48,7 +48,7 @@
 (ImpactNumberOfChannels3 1 1)
 (ImpactPatchSize3 0 0 0 0 0 0)
 (ImpactVoxelSize3 2 2 2 2 2 2)
-(ImpactLayersMask3 "1"  "00000001")
+(ImpactLayersMask3 "1" "00000001")
 (ImpactSubsetFeatures3 32 64)
 (ImpactPCA3 0 0)
 (ImpactDistance3 "L1" "Dice")
@@ -135,4 +135,4 @@
 (ResultImageFormat "mha")
 
 (ITKTransformOutputFileNameExtension "itk.txt")
-(WriteITKCompositeTransform "true")
+(WriteITKCompositeTransform "true")

MR_CT_TS/Prediction.yml

@@ -5,12 +5,123 @@ Predictor:
       engine: elastix
       parameter_maps:
       - ParameterMap_MRI_TS.txt
-      models:
-      - VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
-      - VBoussot/impact-torchscript-models:TS/M852.pt
-      - VBoussot/impact-torchscript-models:TS/M850.pt
-      iterations: 1100
       outputs_criterions: None
+      max_iterations: 0
+      final_grid_spacing: 0.0
+      subset_features: 0
+      spatial_samples: 0
+      parameter_overrides: []
+      resolutions:
+        '0':
+          max_iterations: 400
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 6.0
+              - 6.0
+              - 6.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.2
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:TS/M852.pt
+              voxel_size:
+              - 6.0
+              - 6.0
+              - 6.0
+              layers_mask: '0000001'
+              layers_weight:
+              - 0.8
+              subset_features: 64
+              pca: 0
+              distance: Dice
+        '1':
+          max_iterations: 300
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 3.0
+              - 3.0
+              - 3.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.3
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:TS/M852.pt
+              voxel_size:
+              - 3.0
+              - 3.0
+              - 3.0
+              layers_mask: '0000001'
+              layers_weight:
+              - 0.7
+              subset_features: 64
+              pca: 0
+              distance: Dice
+        '2':
+          max_iterations: 200
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.6
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:TS/M850.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '00000001'
+              layers_weight:
+              - 0.4
+              subset_features: 64
+              pca: 0
+              distance: Dice
+        '3':
+          max_iterations: 200
+          models:
+            '0':
+              ref: VBoussot/impact-torchscript-models:MIND/R1D2_3D.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '1'
+              layers_weight:
+              - 0.7
+              subset_features: 32
+              pca: 0
+              distance: L1
+            '1':
+              ref: VBoussot/impact-torchscript-models:TS/M850.pt
+              voxel_size:
+              - 2.0
+              - 2.0
+              - 2.0
+              layers_mask: '00000001'
+              layers_weight:
+              - 0.3
+              subset_features: 64
+              pca: 0
+              distance: Dice
+      models_registry: VBoussot/impact-torchscript-models:models.json
+      mode: Static
   Dataset:
     groups_src:
       Volume_0:

README.md

@@ -0,0 +1,81 @@
+---
+license: apache-2.0
+library_name: konfai
+pipeline_tag: image-to-image
+tags:
+- medical-imaging
+- registration
+- deformable
+- multimodal
+- ct
+- mri
+- cbct
+- konfai
+---
+
+# IMPACT-Reg — Multimodal Medical Image Registration
+
+Robust **multimodal** (**MR / CT / CBCT**) deformable **registration** presets, built with
+[**KonfAI**](https://github.com/vboussot/KonfAI). Alignment is driven by the **IMPACT** semantic
+similarity metric — deep features from pretrained segmentation / foundation models (**MIND**,
+**TotalSegmentator**, **MRSegmentator**) — so cross-modality pairs align while the deformation
+stays smooth and physically plausible.
+
+Each preset is a **self-contained KonfAI app**: on the **fixed** grid it produces the moving image
+resampled onto the fixed image (**`MovedImage`**) and the **`DisplacementField`**. Presets can be
+**ensembled** (their displacement fields are averaged into one transform).
+
+## 🧩 Presets
+
+| Preset | Pair | Engine | Description |
+|:--|:--|:--|:--|
+| `Generic_Rigid` | any | elastix | Rigid alignment (mutual information, multi-resolution) |
+| `Generic_Rigid_BSpline` | any | elastix | Rigid, then B-spline deformable refinement |
+| `MR_CT_HeadNeck` | MR/CT | elastix + IMPACT | MR/CT head & neck preset |
+| `MR_CT_TS` | MR/CT | elastix + IMPACT | MR/CT with MIND + TotalSegmentator features |
+| `MR_CT_MRSeg` | MR/CT | elastix + IMPACT | MR/CT with MIND + MRSegmentator features |
+| `CBCT_CT_HeadNeck` | CBCT/CT | elastix + IMPACT | CBCT/CT head & neck preset |
+| `CBCT_CT_TS` | CBCT/CT | elastix + IMPACT | CBCT/CT with TotalSegmentator features |
+| `CBCT_CT_MRSeg` | CBCT/CT | elastix + IMPACT | CBCT/CT with MRSegmentator features |
+| `ConvexAdam_Coarse` | any | itk-impact (native) | Global coarse coupled-convex init (IMPACT/MIND) |
+| `ConvexAdam_Fine` | any | itk-impact (native) | Adam instance-optimisation (tileable; expects a pre-aligned start) |
+| `ConvexAdam_Composite` | any | itk-impact (native) | Coarse + fine ConvexAdam in one app (IMPACT/MIND) |
+
+Inputs: **Fixed**, **Moving**, and optional **FixedMask** / **MovingMask** (restrict the metric region).
+
+## 🚀 Usage
+
+```bash
+pip install impact-reg-konfai
+# Register a moving image onto a fixed image (ensemble several presets by listing them):
+impact-reg-konfai register ConvexAdam_Composite -f fixed.nii.gz -m moving.nii.gz -o ./Output --gpu 0
+```
+
+- **Generic runner (single preset):** `konfai-apps infer VBoussot/ImpactReg:ConvexAdam_Composite -i fixed.nii.gz -i moving.nii.gz -o output/`
+- **Interactive:** [**SlicerImpactReg**](https://github.com/vboussot/SlicerImpactReg) — a 3D Slicer extension driving these presets.
+
+> The `ConvexAdam_*` presets depend on [`itk-impact`](https://pypi.org/project/itk-impact/); resolving the app installs it automatically (it reuses your existing PyTorch, CPU or GPU).
+
+## ⚡ Performance & VRAM
+
+`ConvexAdam` presets (native, GPU) benchmarked on an **NVIDIA RTX PRO 5000 (24 GB)** with a real
+abdomen **MR→CT** pair, **222 × 226 × 124 @ 2 mm** (single pass, no TTA):
+
+| Preset | Stages | Time / case | Peak VRAM |
+|:--|:--|:--:|:--:|
+| `ConvexAdam_Fine` | fine (150 Adam iters) | **≈ 0.5 s** | ~2.1 GB |
+| `ConvexAdam_Coarse` | linear + coarse | **≈ 4.6 s** | ~2.1 GB |
+| `ConvexAdam_Composite` | linear + coarse + fine | **≈ 5.1 s** | ~2.1 GB |
+
+Per-stage breakdown: linear pre-align **≈ 4.2 s** (ITK affine, MI) · coarse **≈ 0.4 s** · fine **≈ 0.5 s**.
+One-time TorchScript feature-model load **≈ 7 s** (amortised across a batch). Times scale with case size;
+`--tta k` multiplies runtime. The `elastix + IMPACT` presets run through elastix and scale differently.
+
+## 🔗 Links & Citation
+
+- 🧠 **KonfAI:** [github.com/vboussot/KonfAI](https://github.com/vboussot/KonfAI)
+- 📦 **PyPI:** [impact_reg_konfai](https://pypi.org/project/impact_reg_konfai/)
+- 🩻 **Slicer:** [SlicerImpactReg](https://github.com/vboussot/SlicerImpactReg)
+- 📄 **Paper:** KonfAI — [arXiv:2508.09823](https://arxiv.org/abs/2508.09823)
+</content>
+</invoke>