2.0.0b2: robust bilayer fitting, thickness quality flag, and config defaults

pyproject.toml

@@ -4,7 +4,7 @@ build-backend = "setuptools.build_meta"
 
 [project]
 name = "surface-morphometrics"
-version = "2.0.0b1"
+version = "2.0.0b2"
 description = "Quantification of membrane surfaces segmented from cryo-ET or other volumetric imaging."
 readme = "README.md"
 requires-python = ">=3.11"
@@ -57,4 +57,4 @@ morphometrics = "surface_morphometrics.cli:main"
 include = ["surface_morphometrics*"]
 
 [tool.setuptools.package-data]
-surface_morphometrics = ["config_template.yml"]
+surface_morphometrics = ["config_template.yml", "config_template_simple.yml"]

surface_morphometrics/accept_refinement.py

@@ -25,12 +25,30 @@
 __license__ = "GPLv3"
 
 import os
+import re
 from glob import glob
 from pathlib import Path
 
 import click
 import yaml
 
+from .config_utils import load_config
+
+
+def _iterations_by_surface(work_dir):
+    """Map each surface basename to its sorted available refinement iterations.
+
+    Refinement stops early once a surface converges, so different surfaces may
+    have different numbers of `_refined_iter*` files.
+    """
+    pattern = re.compile(r"^(.*)_refined_iter(\d+)\.surface\.vtp$")
+    out = {}
+    for path in glob(f"{work_dir}*_refined_iter*.surface.vtp"):
+        match = pattern.match(os.path.basename(path))
+        if match:
+            out.setdefault(match.group(1), []).append(int(match.group(2)))
+    return {base: sorted(iters) for base, iters in out.items()}
+
 
 # The canonical surface files the downstream pipeline consumes.  Only these are
 # backed up (from the original) and promoted (from the accepted iteration).
@@ -185,8 +203,7 @@ def accept_refinement_cli(configfile, step, component_name, tomogram, dry_run):
     graph, the command warns that pycurv must be re-run on the promoted surface
     before any downstream analysis.
     """
-    with open(configfile) as f:
-        config = yaml.safe_load(f)
+    config = load_config(configfile, require=("work_dir",))
 
     work_dir = config.get("work_dir", config.get("seg_dir", "./"))
     if not work_dir.endswith("/"):
@@ -200,11 +217,10 @@ def accept_refinement_cli(configfile, step, component_name, tomogram, dry_run):
         print("DRY RUN - no files will be changed")
     print("-" * 60)
 
-    # Discover surfaces that have the requested iteration by finding their
-    # accepted-step surface files, then strip the tag to recover the basename.
-    iter_tag = f"_refined_iter{step}.surface.vtp"
-    candidates = sorted(glob(f"{work_dir}*{iter_tag}"))
-    basenames = [os.path.basename(p)[:-len(iter_tag)] for p in candidates]
+    # Discover every refined surface and the iterations it actually has (a surface
+    # that converged early stops before `step`).
+    iters_by_surface = _iterations_by_surface(work_dir)
+    basenames = sorted(iters_by_surface)
 
     # Apply optional component / tomogram filters on the basename.
     if tomogram:
@@ -213,14 +229,23 @@ def accept_refinement_cli(configfile, step, component_name, tomogram, dry_run):
         basenames = [b for b in basenames if b.endswith(f"_{component_name}")]
 
     if not basenames:
-        print(f"No surfaces with a refined iteration {step} found in {work_dir}"
-              + (f" matching the given filters." if (component_name or tomogram) else "."))
+        print(f"No refined surfaces found in {work_dir}"
+              + (" matching the given filters." if (component_name or tomogram) else "."))
         return
 
     accepted = 0
     needs_pycurv = []
     for basename in basenames:
-        result = accept_one(work_dir, basename, step, radius_hit, dry_run)
+        available = iters_by_surface[basename]
+        # Use the requested iteration, or the last available one if refinement
+        # converged before it (largest available iteration <= step).
+        usable = [i for i in available if i <= step] or available
+        effective_step = usable[-1]
+        if effective_step != step:
+            print(f"  WARNING: {basename}: iteration {step} not available - refinement "
+                  f"converged early at iteration {effective_step}; using iteration "
+                  f"{effective_step}.")
+        result = accept_one(work_dir, basename, effective_step, radius_hit, dry_run)
         if result is not None:
             accepted += 1
             if result is False:
@@ -229,7 +254,8 @@ def accept_refinement_cli(configfile, step, component_name, tomogram, dry_run):
 
     print("-" * 60)
     verb = "Would accept" if dry_run else "Accepted"
-    print(f"{verb} iteration {step} for {accepted} surface(s).")
+    print(f"{verb} refinement for {accepted} surface(s) (requested iteration {step}; "
+          "early-converged surfaces used their last iteration).")
     if accepted:
         print("Originals backed up with a .orig.bak suffix; only canonical surfaces "
               "and refinement summaries were kept.")
@@ -239,7 +265,7 @@ def accept_refinement_cli(configfile, step, component_name, tomogram, dry_run):
         print("iteration and have no curvature graph. Re-run pycurv before any downstream")
         print("step (distances, density sampling, thickness):")
         for basename in needs_pycurv:
-            print(f"  python run_pycurv.py {configfile} {basename}.surface.vtp")
+            print(f"  morphometrics pycurv {configfile} {basename}.surface.vtp")
 
 
 if __name__ == "__main__":

surface_morphometrics/cli.py

@@ -173,20 +173,34 @@ def cli():
 @click.option("-o", "--output", "output_path", default="config.yml",
               type=click.Path(), show_default=True,
               help="Where to write the config file.")
+@click.option("--simple/--verbose", "simple", default=False,
+              help="Write a minimal config with only the commonly-adjusted settings "
+                   "(--simple), or the fully annotated template exposing every option "
+                   "(--verbose, the default).")
 @click.option("--force", is_flag=True, default=False,
               help="Overwrite the output file if it already exists.")
-def new_config(output_path, force):
-    """Write a template config.yml into the current directory."""
+def new_config(output_path, simple, force):
+    """Write a template config.yml into the current directory.
+
+    Omitted settings fall back to documented defaults, so a --simple config is fully
+    runnable; use --verbose to see and tune every available option.
+    """
     from importlib.resources import files
 
     if os.path.exists(output_path) and not force:
         raise click.ClickException(
             f"{output_path} already exists; use --force to overwrite."
         )
-    template = files("surface_morphometrics").joinpath("config_template.yml").read_text()
+    template_name = "config_template_simple.yml" if simple else "config_template.yml"
+    template = files("surface_morphometrics").joinpath(template_name).read_text()
     with open(output_path, "w") as handle:
         handle.write(template)
-    click.echo(f"Wrote {output_path}. Edit it to point at your data and set parameters.")
+    kind = "minimal" if simple else "full"
+    click.echo(f"Wrote {kind} config to {output_path}. Edit it to point at your data "
+               "and set parameters.")
+    if simple:
+        click.echo("Omitted settings use sensible defaults; run "
+                   "`morphometrics new_config --verbose` to see every option.")
 
 
 cli.add_command(new_config)

surface_morphometrics/config_template_simple.yml

@@ -0,0 +1,36 @@
+# Minimal surface_morphometrics configuration.
+#
+# Only the most commonly-adjusted settings are listed here; every other option uses
+# a sensible default. Run `morphometrics new_config --verbose` for the fully
+# annotated template that exposes every setting.
+
+# --- Directories (required) ---
+seg_dir: "./segmentations/"   # Directory of segmentation MRC files
+tomo_dir: "./tomograms/"      # Directory of tomogram MRC files (density sampling, thickness, refinement)
+work_dir: "./morphometrics/"  # Directory for output files
+
+# --- Features (required): map each feature name to its segmentation label value ---
+segmentation_values:
+  OMM: 1
+  IMM: 2
+  ER: 3
+
+cores: 6  # Number of parallel processes
+
+# --- Mesh generation ---
+surface_generation:
+  angstroms: false                 # false rescales surfaces to nm; true keeps angstrom scale
+  isotropic_remesh: true           # Remesh to near-equilateral triangles of ~target_area
+  target_area: 1.0                 # Target triangle area in nm^2 (when isotropic_remesh is true)
+  simplify: false                  # Instead of remeshing, decimate to simplify_max_triangles
+  simplify_max_triangles: 300000   # Only used when simplify is true
+
+# --- Thickness measurement ---
+thickness_measurements:
+  average_radius: 12   # Radius in nm for local averaging of density profiles
+
+# --- Mesh refinement (optional step) ---
+mesh_refinement:
+  average_radius: 25             # Radius in nm for local averaging of density profiles
+  iterations: 6                  # Number of refinement iterations
+  xcorr_iterations: [1, 2, 3]    # Iterations that cross-correlate to sharpen the bilayer before fitting