v2.5.0 release

README.md

@@ -56,8 +56,9 @@ Option | Default | Description
 Option | Default | Description
 --|--|--
 `-m`, `--mode` | `angular` | Destriping approach: `angular` (recommended) or `linear` (deprecated). See [Destriping mode](#destriping-mode) below.
-`--angle` | *(auto-estimated)* | Curtaining orientation, degrees from horizontal. Omit to estimate automatically from the tilt series' central frame[^estimation]; pass a value to override. A diagnostic plot is saved alongside the output when auto-estimated[^diagnosticplot].
-`--reference-frame` | *(mid-stack index)* | Stack index used for angle estimation and the destriping preview. See [Choosing a reference frame](#choosing-a-reference-frame) for more information.
+`--angle` | *(auto-estimated)* | Curtaining orientation, degrees from horizontal. Omit to estimate automatically per frame[^estimation] (see [Per-frame curtain angle](#per-frame-curtain-angle) below); pass a value to override and skip estimation entirely. A diagnostic plot for the reference frame is saved alongside the output when auto-estimated[^diagnosticplot].
+`--reference-frame` | *(mid-stack index)* | Stack index used to seed the per-frame angle walk and as the destriping preview frame. See [Choosing a reference frame](#choosing-a-reference-frame) for more information.
+`--anchor-tilts` | `5` | Number of frames nearest the reference frame used to seed the per-frame consensus walk. See [Per-frame curtain angle](#per-frame-curtain-angle) below.
 `--angular-width` | `8.0` | Angular width of the destriping notch, in degrees. Only used when `--mode angular`. Narrower keeps more real structure sharing a nearby angle to the curtains, at the cost of weaker curtain removal.
 `--notch-fraction` | `0.03` | Width of the destriping notch, as a fraction of image width. Narrower removes less real signal running parallel to the curtains, but leaves more curtaining behind. Only used when `--mode linear`, which is deprecated.
 `--protect-fraction` | `0.01` | Fraction of image width around the zero-frequency (DC) origin exempted from destriping. See [Destriping mode](#destriping-mode) below for why this exists and its trade-off. Only used when `--mode linear`, which is deprecated.
@@ -66,7 +67,7 @@ Option | Default | Description
 Option | Default | Description
 --|--|--
 `--output-dir` | *(required for directory input)* | Output directory for batch mode, mirroring the input directory's structure.
-`--angle-outlier-threshold` | `5.0` | Warn if an individual series' own angle estimate differs from the batch consensus by more than this many degrees.
+`--angle-outlier-threshold` | `5.0` | Within a stack, the maximum degrees a frame's own angle estimate may deviate from its nearest resolved neighbor before falling back to that neighbor's angle (see [Per-frame curtain angle](#per-frame-curtain-angle)). In batch mode, also used to warn (only) when a series' reference-frame estimate differs from the cross-stack diagnostic consensus by more than this many degrees.
 `--workers` | `0` (all CPUs) | Number of parallel processes to use in batch mode.
 
 ##### Other options
@@ -110,19 +111,30 @@ pylisc stack raw_tilt_series/ --output-dir cleared_tilt_series/
 
 Files already carrying a `_PyLisC_` suffix (i.e. previous PyLisC output) are skipped, so re-running against the same directory won't reprocess its own results.
 
-#### Shared curtain angle
-Unless `--angle` is given explicitly, single-file mode estimates the curtaining angle from one frame. In batch mode, PyLisC instead estimates an angle **per series** and combines them into a single shared angle, which is then applied to every series in the batch rather than letting each one drift independently.
+#### Per-frame curtain angle
+As curtaining angles can drift slightly across a tilt series, PyLisC only applies a single angle when one is provided explicitly with `--angle`. Otherwise, for each tilt series PyLisC:
+1. Estimates each frame's curtaining angle.
+2. Determines a 'seed' angle from the `--anchor-tilts` frames nearest the [reference frame](#choosing-a-reference-frame).
+3. Walks outward from the seed frames, comparing each frame's estimate to the nearest already-resolved frame and either:
+      1. Accepts the estimate: if it's within `--angle-outlier-threshold` degrees of the nearest already-resolved frame (itself becoming the new "nearest resolved" frame).
+      2. Rejects the estimate: replacing with that nearest resolved angle and logged as a warning naming the frame. 
+      - This allows for slight angle drift over a series while preventing one bad frame from throwing off its neighbours.
+4. Destripes each frame using the resolved angle.
 
-The combination is a confidence-weighted circular mean: each series' angle is weighted by its own confidence ratio (see [Curtain angle diagnostic plot](#curtain-angle-diagnostic-plot) below), so a series with a clear, sharp peak counts for more than one with a flat, uncertain profile.
+The reference frame's diagnostic plot (see [Curtain angle diagnostic plot](#curtain-angle-diagnostic-plot) below) is saved alongside the output.
 
-#### Outlier detection
+#### Multi-series consensus check
 
-If any individual series' own angle estimate differs from the batch consensus by more than `--angle-outlier-threshold` (default `5.0` degrees), a warning is printed naming that series:
+In batch mode, PyLisC additionally combines each tilt series' reference-frame estimate into a consensus angle to use as a diagnostic value. This is a confidence-weighted circular mean: each series' angle is weighted by its own confidence ratio (see [Curtain angle diagnostic plot](#curtain-angle-diagnostic-plot) below), clipped to a per-run cap first so a single sharp FFT peak can't dominate the consensus (the cap is the confidence distribution's median plus 3x its scaled median absolute deviation, falling back to 5x the median when every value is identical).
+
+If any individual series' own reference-frame estimate differs from this consensus by more than `--angle-outlier-threshold` (default `5.0` degrees), a warning is printed naming that series:
 
 ```
 WARNING: sample_07.mrc angle (58.3 deg) deviates 41.2 deg from consensus (17.1 deg) -- check its diagnostic plot
 ```
 
+This does not affect the actual destriping, however, with each tilt series destriped using its own per-frame walk described [above](#per-frame-curtain-angle).
+
 ### `pylisc frames`
 
 ```sh
@@ -132,11 +144,10 @@ pylisc frames [OPTIONS] --output-dir OUTPUT_DIR --filename-template TEMPLATE INP
 This command is aimed at destriping tilt images that exist as individual 2D MRC frames, i.e. not yet assembled/aligned into a stack. `INPUT_DIR` is not recursed into; every `*.mrc` directly inside it (excluding PyLisC's own `_PyLisC_` output) is treated as one tilt image. `--output-dir` is required, and mirrors the input's flat structure: each frame is written back out individually with a `_PyLisC_{mode}` suffix, same as `--output-dir` does for [batch mode](#batch-mode).
 
 #### Options
-`pylisc frames` uses many of the same options as `pylisc stack`, see [above](#options) or run `pylisc frames -h` for further information. Frames mode does not take `--reference-frame` or `--preview-strengths`, but adds:
+`pylisc frames` uses many of the same options as `pylisc stack`, see [above](#options) or run `pylisc frames -h` for further information, including `--anchor-tilts` (here, the number of tilt buckets rather than frames used to seed the consensus walk — see [Per-tilt curtain angle](#per-tilt-curtain-angle) below). Frames mode does not take `--reference-frame` or `--preview-strengths`, but adds:
 
 Option | Default | Description
 --|--|--
-`--anchor-tilts` | `5` | Number of tilt buckets nearest each series' median tilt used to seed its consensus walk. See [Per-tilt curtain angle](#per-tilt-curtain-angle) below.
 `--print-angles` | off | Print a diagnostic table of per-file and per-bucket angle estimation, and write it to CSV files in `--output-dir`.
 
 #### Filename template
@@ -158,7 +169,7 @@ Field boundaries default to underscore only. `--filename-delimiters` sets which
 Curtaining orientation drifts slightly with tilt angle, so unless `--angle` is given explicitly, frames mode does **not** use one consensus angle for the whole directory. Instead:
 1. Every frame's own angle is estimated.
 2. Each frame's confidence ratio is clipped to a per-run cap before use, so a single sharp FFT peak can't dominate its bucket's consensus. The cap is the confidence distribution's median plus 3x its (scaled) median absolute deviation, falling back to 5x the median when every value is identical, (i.e. zero deviation).
-3. Frames are grouped by tilt angle, rounded to the nearest whole degree (so e.g. two positions' `-30.00°` and `-29.98°` tilts fall in the same bucket), and each bucket's estimates are combined into a per-tilt consensus (the same confidence-weighted circular mean as [batch mode](#shared-curtain-angle)), using the clipped confidences.
+3. Frames are grouped by tilt angle, rounded to the nearest whole degree (so e.g. two positions' `-30.00°` and `-29.98°` tilts fall in the same bucket), and each bucket's estimates are combined into a per-tilt consensus (the same confidence-weighted circular mean used in [stack mode's cross-stack consensus check](#multi-series-consensus-check)), using the clipped confidences.
 4. Tilt buckets are clustered into separate acquisition series (by matching frame count per bucket, then splitting wherever the tilt spacing breaks step), so a directory holding more than one tilt series is handled independently per series rather than as one pool.
 5. Within each series, the `--anchor-tilts` buckets around the series' median tilt are used as a trusted seed. PyLisC then walks outward from this seed, and compares each tilt bucket's own consensus to the nearest already-resolved bucket. A tilt bucket's consensus angle is accepted (becoming the new "nearest resolved" point) if it's within `--angle-outlier-threshold` degrees of it, otherwise it is replaced with that nearest resolved angle and logged as a warning naming both tilts. This approach prevents buckets with incorrect angles from throwing off subsequent tilts, while allowing for slight angle drift over a tilt series.
 6. Pass `--print-angles` to see the full per-file and per-bucket breakdown (angle, confidence, resolved angle, seed/accepted/rejected status) as tables and CSV files written to `--output-dir`.
@@ -193,7 +204,7 @@ Curtaining removal works by finding curtaining's signature in Fourier space and
 - **`linear` (deprecated).** Dims frequencies by their *distance* from the curtain line rather than their direction. Below a radius set by `--notch-fraction`, distance alone can no longer distinguish direction at all, so without `--protect-fraction` exempting a small disc around the origin, large-scale contrast gets suppressed at every angle near that radius, not just along the curtains. Protecting that disc, in turn, risks letting broad, low-frequency curtaining pass through unfiltered if the curtaining's own frequency sits close to the protected radius. `angular` avoids this trade-off entirely.
 
 ### Choosing a reference frame
-The reference frame should be the tilt with the least foreshortening and the best signal-to-noise, since that gives the most reliable curtain angle estimate and the clearest destriping preview. In practice this is the 0° tilt (or the pretilt used during lamella imaging).
+The reference frame should be the tilt with the least foreshortening and the best signal-to-noise, since that gives the most reliable curtain angle estimate, the clearest destriping preview, and the most trustworthy seed for the [per-frame consensus walk](#per-frame-curtain-angle). In practice this is the pre-tilt used during lamella imaging.
 
 - **Dose-symmetric schemes** (0° acquired first, then alternating ±): use `--reference-frame 0`.
 - **Continuous sweeps** (most-negative tilt acquired first): 0° sits in the middle of the stack, so use roughly `--reference-frame <n//2>` for an n-tilt series.

pyproject.toml

@@ -1,6 +1,6 @@
 [project]
 name = "pylisc"
-version = "2.4.0"
+version = "2.5.0"
 description = "Python implementation of LisC algorithm"
 readme = "README.md"
 authors = [

src/pylisc/cli.py

@@ -91,7 +91,7 @@ def version_callback(value: bool | None) -> None:
 ]
 AnchorTiltsOpt = Annotated[
     int,
-    typer.Option('--anchor-tilts', help='Number of tilt buckets nearest each series\' median tilt used to seed its consensus walk (frames mode only).', rich_help_panel='Batch options', min=1),
+    typer.Option('--anchor-tilts', help='Number of tilt buckets (frames mode) or frames (stack mode) nearest the median/reference used to seed the angle consensus walk.', rich_help_panel='Batch options', min=1),
 ]
 
 # Define callback for pylisc (to allow version option)
@@ -139,6 +139,7 @@ def stack(
     notch_frac: NotchFracOpt = 0.03,
     dc_protect_frac: DcProtectFracOpt = 0.01,
     angle_outlier_threshold: AngleOutlierThresholdOpt = 5.0,
+    anchor_tilts: AnchorTiltsOpt = 5,
     version: VersionOpt = None,
     workers: WorkersOpt = 0,
 ):
@@ -181,6 +182,7 @@ def stack(
         notch_frac=notch_frac,
         dc_protect_frac=dc_protect_frac,
         angle_outlier_threshold=angle_outlier_threshold,
+        anchor_tilts=anchor_tilts,
         preview_strengths=preview_strengths,
         force=force,
         dry_run=dry_run,

src/pylisc/estimate_angle.py

@@ -89,6 +89,54 @@ def plot_angular_energy(
     plt.close(fig)
     return output_path
 
+def clip_confidence_outliers(confidences: dict) -> dict:
+    '''
+    Cap confidences above a median+MAD threshold so a single spuriously sharp FFT peak can't dominate a consensus
+    '''
+    conf_values = np.array(list(confidences.values()))
+    if len(conf_values) == 0:
+        return dict(confidences)
+    median_conf = np.median(conf_values)
+    mad = np.median(np.abs(conf_values - median_conf))
+    confidence_cap = median_conf + 3 * 1.4826 * mad if mad > 0 else median_conf * 5
+    if confidence_cap <= 0:
+        return dict(confidences)
+    return {k: min(v, confidence_cap) for k, v in confidences.items()}
+
+def resolve_walk(keys_sorted: list, values: dict, angle_outlier_threshold: float, anchor_window: int) -> tuple[dict, dict]:
+    '''
+    Seed a window of keys around the median key, and walk outward, checking each key's own value against the nearest already-resolved (seed/accepted) value
+    '''
+    if len(keys_sorted) == 1:
+        k = keys_sorted[0]
+        return {k: values[k]}, {k: 'seed'}
+    median_key = np.median(keys_sorted)
+    center_idx = int(np.argmin([abs(k - median_key) for k in keys_sorted]))
+    window = max(1, anchor_window)
+    half = window // 2
+    start = max(0, center_idx - half)
+    end = min(len(keys_sorted), start + window)
+    start = max(0, end - window)
+    seed_keys = keys_sorted[start:end]
+    resolved = {k: values[k] for k in seed_keys}
+    status = {k: 'seed' for k in seed_keys}
+    for direction, idx, edge in ((-1, start - 1, start), (1, end, end - 1)):
+        nearest = resolved[keys_sorted[edge]]
+        i = idx
+        while 0 <= i < len(keys_sorted):
+            k = keys_sorted[i]
+            own_value = values[k]
+            deviation = min(abs(own_value - nearest), 180 - abs(own_value - nearest))
+            if deviation <= angle_outlier_threshold:
+                resolved[k] = own_value
+                status[k] = 'accepted'
+                nearest = own_value
+            else:
+                resolved[k] = nearest
+                status[k] = 'rejected'
+            i += direction
+    return resolved, status
+
 def combine_angles(angles_deg: list, confidences: list) -> tuple:
     '''
     Confidence-weighted circular mean of curtain angles

src/pylisc/frames.py

@@ -9,7 +9,7 @@
 from rich.table import Table
 
 # Import internal PyLisC modules
-from pylisc.estimate_angle import combine_angles, estimate_curtain_angle
+from pylisc.estimate_angle import clip_confidence_outliers, combine_angles, estimate_curtain_angle, resolve_walk
 from pylisc.io import find_input_files, readMrcFile, writeMrcFile
 from pylisc.lisc import lisc_clear_frame
 from pylisc.log import logger, per_file_log
@@ -178,42 +178,14 @@ def _resolve_series(series_buckets, bucket_consensus, angle_outlier_threshold, a
     Seed buckets around the median tilt, and walk out, checking each against the nearest already-resolved (seed/accepted) bucket
     '''
     buckets_sorted = sorted(series_buckets)
-    if len(buckets_sorted) == 1:
-        b = buckets_sorted[0]
-        return {b: bucket_consensus[b]}, {b: 'seed'}
-
-    median_tilt = np.median(buckets_sorted)
-    center_idx = int(np.argmin([abs(b - median_tilt) for b in buckets_sorted]))
-    window = max(1, anchor_tilts)
-    half = window // 2
-    start = max(0, center_idx - half)
-    end = min(len(buckets_sorted), start + window)
-    start = max(0, end - window)
-
-    seed_buckets = buckets_sorted[start:end]
-    resolved = {b: bucket_consensus[b] for b in seed_buckets}
-    status = {b: 'seed' for b in seed_buckets}
-
-    for direction, idx, edge in ((-1, start - 1, start), (1, end, end - 1)):
-        nearest = resolved[buckets_sorted[edge]]
-        i = idx
-        while 0 <= i < len(buckets_sorted):
-            b = buckets_sorted[i]
-            own_angle = bucket_consensus[b]
-            deviation = min(abs(own_angle - nearest), 180 - abs(own_angle - nearest))
-            if deviation <= angle_outlier_threshold:
-                resolved[b] = own_angle
-                status[b] = 'accepted'
-                nearest = own_angle
-            else:
-                resolved[b] = nearest
-                status[b] = 'rejected'
-                logger.warning(
-                    "tilt {}° consensus angle ({}°) deviates {}° from nearest resolved angle ({}°) - using that angle instead",
-                    b, f'{own_angle:.1f}', f'{deviation:.1f}', f'{nearest:.1f}',
-                )
-            i += direction
-
+    resolved, status = resolve_walk(buckets_sorted, bucket_consensus, angle_outlier_threshold, anchor_tilts)
+    for b in buckets_sorted:
+        if status[b] == 'rejected':
+            deviation = min(abs(bucket_consensus[b] - resolved[b]), 180 - abs(bucket_consensus[b] - resolved[b]))
+            logger.warning(
+                "tilt {}° consensus angle ({}°) deviates {}° from nearest resolved angle ({}°) - using that angle instead",
+                b, f'{bucket_consensus[b]:.1f}', f'{deviation:.1f}', f'{resolved[b]:.1f}',
+            )
     return resolved, status
 
 def _estimate_per_tilt_angles(paths, tilt_of, angle_outlier_threshold, anchor_tilts=5, print_angles=False, output_dir=None):
@@ -230,18 +202,10 @@ def _estimate_per_tilt_angles(paths, tilt_of, angle_outlier_threshold, anchor_ti
 
     # A single spuriously sharp FFT peak can otherwise dominate its bucket's consensus
     raw_confidences = dict(confidences)
-    conf_values = np.array(list(confidences.values()))
-    if len(conf_values):
-        median_conf = np.median(conf_values)
-        mad = np.median(np.abs(conf_values - median_conf))
-        confidence_cap = median_conf + 3 * 1.4826 * mad if mad > 0 else median_conf * 5
-    else:
-        confidence_cap = 0.0
-    if confidence_cap > 0:
-        n_clipped = sum(1 for c in confidences.values() if c > confidence_cap)
-        if n_clipped:
-            logger.debug('clipping {} frame(s) with confidence above {}', n_clipped, f'{confidence_cap:.2f}')
-        confidences = {p: min(c, confidence_cap) for p, c in confidences.items()}
+    confidences = clip_confidence_outliers(confidences)
+    n_clipped = sum(1 for p in confidences if confidences[p] < raw_confidences[p])
+    if n_clipped:
+        logger.debug('clipping {} frame(s) with confidence above cap', n_clipped)
 
     tilt_buckets = {}
     for path in paths: