fix: correct Zeo++ engine arguments and output handling

docs/unified-api.md

@@ -94,6 +94,15 @@ Paths in CLI specifications are relative to the specification file. Model names
 remain identifiers unless they resolve to local checkpoint files. Python paths
 are relative to the caller's working directory.
 
+For Zeo++ surface area, volume, and PSD, `channel_radius` controls accessibility
+and `probe_radius` controls sampling. The probe radius must not exceed the
+channel radius; equal radii remain the default. Both APIs pass these values in
+Zeo++'s documented channel/probe order. All requested outputs must be present,
+complete, and finite. Collection accepts multiline channel output and the
+default `.psd_histo` histogram filename, retaining explicitly named legacy
+`.psd` files and existing result keys. See the
+[Zeo++ command and output examples](https://www.zeoplusplus.org/examples.html).
+
 ```bash
 matkit pores --spec pores.json --outdir runs/pores
 matkit prepare --spec pores.json --outdir runs/prepared

src/matkit/api/adapters.py

@@ -289,47 +289,30 @@ def external_command(root, execution, engine, arguments):
 
 
 def run_pores(root, request, execution):
-    args = ["-ha"] if request.high_accuracy else []
-    args.extend(["-r", request.radii_file])
-    for analysis in request.analyses:
-        args.append(f"-{analysis}")
-        if analysis in {"sa", "vol", "psd"}:
-            args.extend(
-                [
-                    str(request.probe_radius),
-                    str(request.channel_radius),
-                    str(request.num_samples),
-                ]
-            )
-        elif analysis == "chan":
-            args.append(str(request.probe_radius))
-    args.append("structure.cif")
+    from matkit.zeopp.zeopp import _analysis_arguments
+
+    args = _analysis_arguments(
+        request.analyses,
+        request.probe_radius,
+        request.channel_radius,
+        request.num_samples,
+        request.high_accuracy,
+        request.radii_file,
+        "structure.cif",
+    )
     external_command(root, execution, "zeopp", args)
     return parse_pores(root, request)
 
 
 def parse_pores(root, request):
-    from matkit.zeopp.zeopp import _PARSERS
-
-    required = {
-        "res": {"Di", "Df", "Dif"},
-        "sa": {"ASA", "NASA", "density", "unitcell_volume"},
-        "vol": {"AV", "NAV", "density", "unitcell_volume"},
-        "psd": {"bin_lower", "counts"},
-        "chan": {"num_channels", "dimensionalities"},
+    from matkit.zeopp.zeopp import _output_path, _parse_output
+
+    results = {
+        analysis: _parse_output(
+            _output_path(root / "work", "structure", analysis), analysis
+        )
+        for analysis in request.analyses
     }
-    results = {}
-    for analysis in request.analyses:
-        data = _PARSERS[analysis](root / "work" / f"structure.{analysis}")
-        if not required[analysis] <= data.keys():
-            raise ValueError(f"Incomplete requested Zeo++ analysis: {analysis}")
-        if analysis == "psd" and not data["counts"]:
-            raise ValueError("Empty pore-size distribution")
-        if analysis == "chan" and data["num_channels"] != len(
-            data["dimensionalities"]
-        ):
-            raise ValueError("Channel count does not match dimensionalities")
-        results[analysis] = data
     return (
         PorePayload(results=results),
         [

src/matkit/api/models.py

@@ -231,6 +231,10 @@ class PoreRequest(RequestBase):
     def unique_analyses(self):
         if len(set(self.analyses)) != len(self.analyses):
             raise ValueError("Analyses must be unique")
+        if {"sa", "vol", "psd"}.intersection(
+            self.analyses
+        ) and self.probe_radius > self.channel_radius:
+            raise ValueError("probe_radius must not exceed channel_radius")
         return self
 
 

src/matkit/zeopp/zeopp.py

@@ -3,6 +3,8 @@
 import logging
 from concurrent.futures import ThreadPoolExecutor
 import json
+import math
+from numbers import Integral, Real
 from pathlib import Path
 import shutil
 import subprocess
@@ -14,6 +16,61 @@
 
 
 VALID_ANALYSES = {"res", "sa", "vol", "psd", "chan"}
+_OUTPUT_SUFFIXES = {analysis: (f".{analysis}",) for analysis in VALID_ANALYSES}
+# Zeo++ uses this suffix by default; retain explicitly named legacy files.
+_OUTPUT_SUFFIXES["psd"] = (".psd_histo", ".psd")
+
+
+def _analysis_arguments(
+    analyses,
+    probe_radius,
+    channel_radius,
+    num_samples,
+    high_accuracy,
+    radii_file,
+    structure_file,
+):
+    """Build the engine syntax shared by legacy and unified execution."""
+    if not analyses or set(analyses) - VALID_ANALYSES:
+        raise ValueError(f"Invalid analysis types: {analyses}")
+    if len(set(analyses)) != len(analyses):
+        raise ValueError("Analyses must be unique")
+    for name, value in (
+        ("probe_radius", probe_radius),
+        ("channel_radius", channel_radius),
+    ):
+        if (
+            isinstance(value, bool)
+            or not isinstance(value, Real)
+            or not math.isfinite(value)
+            or value <= 0
+        ):
+            raise ValueError(f"{name} must be positive and finite")
+    if (
+        isinstance(num_samples, bool)
+        or not isinstance(num_samples, Integral)
+        or num_samples < 1
+    ):
+        raise ValueError("num_samples must be a positive integer")
+    sampled = {"sa", "vol", "psd"}
+    if sampled.intersection(analyses) and probe_radius > channel_radius:
+        raise ValueError("probe_radius must not exceed channel_radius")
+    args = ["-ha"] if high_accuracy else []
+    args.extend(["-r", str(radii_file)])
+    for analysis in analyses:
+        args.append(f"-{analysis}")
+        if analysis in sampled:
+            args.extend(map(str, (channel_radius, probe_radius, num_samples)))
+        elif analysis == "chan":
+            args.append(str(probe_radius))
+    return [*args, str(structure_file)]
+
+
+def _output_path(directory: Path, stem: str, analysis: str) -> Path:
+    candidates = [
+        directory / f"{stem}{suffix}" for suffix in _OUTPUT_SUFFIXES[analysis]
+    ]
+    return next((path for path in candidates if path.is_file()), candidates[0])
 
 
 def _find_network_binary(network_path: str | None = None) -> str:
@@ -181,16 +238,34 @@ def _parse_psd(filepath: Path) -> dict:
         counts = []
         for line in lines:
             parts = line.strip().split()
-            if len(parts) >= 2:
-                try:
-                    low = float(parts[0])
-                    count = float(parts[1])
-                    bin_lower.append(low)
-                    counts.append(count)
-                except ValueError:
-                    continue
+            if not parts:
+                continue
+            try:
+                low = float(parts[0])
+            except ValueError:
+                continue  # Zeo++ histogram headers are not numeric rows.
+            if len(parts) < 2:
+                raise ValueError("Incomplete PSD histogram row")
+            values = [float(value) for value in parts]
+            if not all(math.isfinite(value) for value in values):
+                raise ValueError("PSD histogram values must be finite")
+            count = values[1]
+            if low < 0 or count < 0:
+                raise ValueError("PSD bins and counts must be nonnegative")
+            if bin_lower and low <= bin_lower[-1]:
+                raise ValueError("PSD bins must be strictly increasing")
+            bin_lower.append(low)
+            counts.append(count)
+
+        if not counts:
+            raise ValueError("Empty pore-size distribution")
 
         bin_size = bin_lower[1] - bin_lower[0] if len(bin_lower) > 1 else 0.0
+        if any(
+            not math.isclose(b - a, bin_size, rel_tol=1e-6, abs_tol=1e-9)
+            for a, b in zip(bin_lower, bin_lower[1:])
+        ):
+            raise ValueError("PSD histogram requires uniformly spaced bins")
         bin_upper = [b + bin_size for b in bin_lower]
 
         return {
@@ -215,12 +290,41 @@ def _parse_chan(filepath: Path) -> dict:
         Dict with num_channels and dimensionalities list.
     """
     try:
-        text = filepath.read_text().strip()
-        parts = text.split()
+        lines = [
+            line.strip()
+            for line in filepath.read_text().splitlines()
+            if line.strip() and not line.lstrip().startswith("#")
+        ]
+        parts = lines[0].split()
         channels_idx = parts.index("channels")
         num_channels = int(parts[channels_idx - 1])
         dim_idx = parts.index("dimensionality")
         dimensionalities = [int(d) for d in parts[dim_idx + 1 :]]
+        if (
+            num_channels < 0
+            or num_channels != len(dimensionalities)
+            or any(d not in {1, 2, 3} for d in dimensionalities)
+            or parts[channels_idx + 1 : dim_idx] != ["identified", "of"]
+        ):
+            raise ValueError("Invalid channel count or dimensionalities")
+        channels = set()
+        for line in lines[1:]:
+            row = line.split()
+            if row[0] != "Channel":
+                continue  # Optional summary rows do not list dimensionalities.
+            if len(row) != 5:
+                raise ValueError("Incomplete channel detail row")
+            index = int(row[1])
+            diameters = [float(value) for value in row[2:]]
+            if index in channels or not 0 <= index < num_channels:
+                raise ValueError("Invalid or duplicate channel index")
+            if any(not math.isfinite(d) or d < 0 for d in diameters):
+                raise ValueError(
+                    "Channel diameters must be finite and nonnegative"
+                )
+            channels.add(index)
+        if channels and len(channels) != num_channels:
+            raise ValueError("Incomplete channel detail rows")
         return {
             "num_channels": num_channels,
             "dimensionalities": dimensionalities,
@@ -238,6 +342,25 @@ def _parse_chan(filepath: Path) -> dict:
 }
 
 
+def _parse_output(path: Path, analysis: str) -> dict:
+    """Require complete, finite results in every public Zeo++ interface."""
+    required = {
+        "res": {"Di", "Df", "Dif"},
+        "sa": {"ASA", "NASA", "density", "unitcell_volume"},
+        "vol": {"AV", "NAV", "density", "unitcell_volume"},
+        "psd": {"bin_lower", "counts"},
+        "chan": {"num_channels", "dimensionalities"},
+    }
+    data = _PARSERS[analysis](path)
+    if not required[analysis] <= data.keys():
+        raise ValueError(f"Incomplete requested Zeo++ analysis: {analysis}")
+    try:
+        json.dumps(data, allow_nan=False)
+    except ValueError as exc:
+        raise ValueError(f"Nonfinite Zeo++ analysis: {analysis}") from exc
+    return data
+
+
 def get_output_data(
     output_path: str,
     analyses: list[str] | None = None,
@@ -247,7 +370,7 @@ def get_output_data(
     Args:
         output_path: Path to directory containing Zeo++ output files.
         analyses: Which analyses to parse. If None, auto-detects from
-            available files (.res, .sa, .vol, .psd, .chan).
+            available files (.res, .sa, .vol, .psd_histo, .psd, .chan).
 
     Returns:
         Dict with 'success' key and per-analysis result sub-dicts.
@@ -271,17 +394,27 @@ def get_output_data(
 
     if outdir.is_file():
         ext = outdir.suffix.lstrip(".")
+        if ext == "psd_histo":
+            ext = "psd"
+        if analyses is not None and set(analyses) != {ext}:
+            raise ValueError("File does not contain all requested analyses")
         if ext in _PARSERS:
-            results[ext] = _PARSERS[ext](outdir)
+            results[ext] = _parse_output(outdir, ext)
             results["success"] = True
         return results
 
     # Directory mode: find output files
     detect = analyses if analyses is not None else list(VALID_ANALYSES)
     for analysis in detect:
-        matches = list(outdir.glob(f"*.{analysis}"))
+        matches = [
+            path
+            for suffix in _OUTPUT_SUFFIXES[analysis]
+            for path in sorted(outdir.glob(f"*{suffix}"))
+        ]
         if matches:
-            results[analysis] = _PARSERS[analysis](matches[0])
+            results[analysis] = _parse_output(matches[0], analysis)
+        elif analyses is not None:
+            raise ValueError(f"Missing requested Zeo++ analysis: {analysis}")
 
     if len(results) > 1:  # has at least one analysis result beyond 'success'
         results["success"] = True
@@ -367,47 +500,18 @@ def run_zeopp(
         rad_dest = workdir / radii_path.name
         shutil.copy(radii_path, rad_dest)
 
-        # Build command
-        cmd = [binary]
-
-        if ha:
-            cmd.append("-ha")
-
-        cmd.extend(["-r", str(rad_dest)])
-
-        for analysis in analyses:
-            if analysis == "res":
-                cmd.extend(["-res"])
-            elif analysis == "sa":
-                cmd.extend(
-                    [
-                        "-sa",
-                        str(probe_radius),
-                        str(chan_radius),
-                        str(num_samples),
-                    ]
-                )
-            elif analysis == "vol":
-                cmd.extend(
-                    [
-                        "-vol",
-                        str(probe_radius),
-                        str(chan_radius),
-                        str(num_samples),
-                    ]
-                )
-            elif analysis == "psd":
-                cmd.extend(
-                    [
-                        "-psd",
-                        str(probe_radius),
-                        str(chan_radius),
-                        str(num_samples),
-                    ]
-                )
-            elif analysis == "chan":
-                cmd.extend(["-chan", str(probe_radius)])
-        cmd.append(str(cif_dest))
+        cmd = [
+            binary,
+            *_analysis_arguments(
+                analyses,
+                probe_radius,
+                chan_radius,
+                num_samples,
+                ha,
+                rad_dest,
+                cif_dest,
+            ),
+        ]
 
         proc = subprocess.run(cmd, capture_output=True, text=True)
         if proc.returncode != 0:
@@ -420,9 +524,12 @@ def run_zeopp(
         result = {"success": False, "results": {}, "error": None}
         stem = cif_dest.stem
         for analysis in analyses:
-            out_file = workdir / f"{stem}.{analysis}"
-            if out_file.exists():
-                result["results"][analysis] = _PARSERS[analysis](out_file)
+            out_file = _output_path(workdir, stem, analysis)
+            if not out_file.is_file():
+                raise ValueError(
+                    f"Missing requested Zeo++ analysis: {analysis}"
+                )
+            result["results"][analysis] = _parse_output(out_file, analysis)
 
         result["success"] = len(result["results"]) > 0
         return result

tests/data/zeopp/test_structure.chan

@@ -1 +1,3 @@
 test_structure 2 channels identified of dimensionality 3 3
+Channel 0 18.569 8.023 10.928
+Channel 1 12.0 6.0 8.0