feat(ccd): publish writing-ready corrected analysis

Consolidate the validated Karman CCD chain into a bounded, navigable package with reproducible diagnostics and lightweight plotting contracts while keeping dense evidence external. Co-authored-by: Cursor <cursoragent@cursor.com>
2026-08-09 21:25:18 +08:00
parent a451dca961
commit 7eab04869c
32 changed files with 3407 additions and 579 deletions
@@ -141,3 +141,6 @@ src/drl_pinball/train/output/*/*
 src/drl_pinball/train/output/*/models/*
 !src/drl_pinball/train/output/*/models/best_model.zip
 !src/drl_pinball/train/calibrations/*/target.npy
+# CCD generated evidence/archive payloads remain external immutable data.
+src/CCD_analysis/evidence/*
+src/CCD_analysis/archive/*
@@ -0,0 +1,26 @@
+# CCD claim-status ledger
+
+## Current / supported
+
+- Mean-control contribution dominates the tested zero→DRL ROI improvement: `93.42%` versus `6.58%` dynamic increment.
+- The corrected DRL-minus-constant-template residual has a stable rank-3 phase-coherent CCD subspace under the declared cycle, bin, origin, and template checks.
+- In the pinball-inclusive domain, raw full-fit weighted POD has nearly identical rank-3 field subspace: principal cosines `[0.9999616609, 0.9998971005, 0.9976439804]`.
+- CCD supplies action-coordinate association (`front`, `rear-symmetric`, `rear-antisymmetric`) for descriptive interpretation of that residual.
+
+## Contextual / diagnostic
+
+- Modes 1/2 have predominantly antisymmetric field parity; mode 3 is predominantly symmetric.
+- Signed velocity maps and the body-local rotation proxy suggest near-body rotation combinations. They are descriptive reader diagnostics, not wall-vorticity measurements.
+- SR provides bounded contextual control-law term ranking: rear constant dominant tested element, rear-lift feedback secondary, tested front feedback weak over its deletion window.
+
+## Historical / superseded
+
+- Prior temporal CCD, Fourier-resampled phase CCD, wake-only interpretation siblings, and archived pre-reset routes are retained for provenance only. Their paths must not be presented as current authority.
+
+## Unsupported / prohibited wording
+
+CCD does not establish a paired counterfactual, causality, actuator response, mechanism, response time, independent-realization uncertainty, explained variance, global stability, CCD superiority over POD, or universal Reynolds-number law. Do not call individual SVD basis vectors uniquely stable physical objects; only the rank-3 subspace is stability-supported. Do not call the constant-template subtraction a causal DRL contribution.
+
+## Status vocabulary
+
+`Current` means supported within the exact artifact contract. `Contextual` means useful cross-analysis but not an independent claim. `Diagnostic` means descriptive evidence. `Historical` means retained provenance. `Unsupported` means do not use without new evidence.
@@ -0,0 +1,69 @@
+# CCD figures and redraw contract
+
+## Preferred reader publication
+
+Use `data/karman-dynamic/canonical/publication-pinball-domain-20260809T1925+0800/` for repository-visible figures. Its external machine-bound parent is `/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/publication-pinball-domain-20260809T1925+0800/`. The four PNG/PDF groups are derived presentation assets, not new evidence.
+
+| File | Purpose | Main source fields |
+|---|---|---|
+| `01_domain_estimand` | Frozen domain, bodies, partitions, and residual definition | `x_D`, `y_D`, `domain_mask` |
+| `02_ccd_mode_identity` | CCD `u_x`, `u_y`, action coordinates, body-local rotation proxy | `ccd_mode_velocity`, `ccd_mode_rotation_proxy`, `left_action_physical_coordinates` |
+| `03_phase_reconstruction` | Raw, CCD rank-2/rank-3, and POD rank-3 phase-bin means | `representative_*`, `representative_phase_bins` |
+| `04_raw_pod_vs_ccd` | Raw POD versus CCD modes and phase residual curves | `ccd_mode_velocity`, `pod_mode_velocity`, `reconstruction`, principal cosines |
+
+PNG/PDF files are reader views. For numeric redraw use `plot_data.npz`; do not digitize pixels.
+
+## `plot_data.npz` fields
+
+The attached `plot_data.json` is the machine-readable companion and records the source result manifest hash, geometry, domain, regions, and field semantics.
+
+- `x_D`: `(1280,)`, nondimensional x coordinates, float32, axis 0.
+- `y_D`: `(512,)`, nondimensional y coordinates, float32, axis 1.
+- `domain_mask`: `(1280,512)`, boolean common-fluid selected domain.
+- `point_weights`: `(99264,)`, positive coordinate quadrature weights in selected-point order.
+- `ccd_mode_velocity`: `(3,2,1280,512)`, zero-filled grid; component index 0=`u_x`, 1=`u_y`.
+- `ccd_mode_rotation_proxy`: `(3,1280,512)`, signed Gaussian body-local angular-momentum-like visualization projection, not wall vorticity.
+- `ccd_modes`: `(198528,3)`, component-major physical vectors: selected `u_x` then selected `u_y`.
+- `ccd_left_functions`: `(3,3)`, observable-side SVD vectors in native `front,upper,lower` order.
+- `left_action_physical_coordinates`: `(3,3)`, rows `front`, `rear-symmetric`, `rear-antisymmetric`; columns mode 1–3.
+- `representative_phase_bins`: `(4,)`, registered bins closest to `0, π/2, π, 3π/2`.
+
+The compact result additionally contains coefficients, raw/POD modes, all summary diagnostics, hashes, configuration, and manifest. It is the scientific result artifact; `plot_data.npz` is intentionally smaller and plotting-oriented.
+
+## Coordinate, mask, and geometry rules
+
+Arrays use `(x index, y index)` while Matplotlib display uses `field.T` so x is horizontal and y is vertical. The origin is the front-cylinder centre in the nondimensional Legacy coordinate convention. The frozen domain is `29 <= x/D <= 54`, `|y/D| <= 5`. Fixed partitions are geometry `29–32D`, near wake `32–40D`, middle wake `40–47D`, far wake `47–54D`.
+
+Bodies used for overlays are `(front, upper, lower)` with centres approximately `(30,0)`, `(31.3,+0.75)`, `(31.3,-0.75)` and radius `0.5D`. The mask is the exact solver-fluid common mask; body circles are reader-view geometry overlays and do not replace the solver mask.
+
+## Color and interpretation
+
+The preferred mode plot uses one robust shared signed velocity limit across all three modes, with a diverging map centered at zero. The rotation proxy has its own shared robust signed limit. A panel's color sign is meaningful only relative to the same panel's zero and the globally arbitrary mode sign; compare `coefficient × mode` for phase-specific physical sign.
+
+## Minimal NumPy redraw
+
+```python
+import numpy as np
+import matplotlib.pyplot as plt
+
+with np.load("plot_data.npz", allow_pickle=False) as z:
+    x, y = z["x_D"], z["y_D"]
+    mask = z["domain_mask"]
+    ux = z["ccd_mode_velocity"][0, 0]
+
+fig, ax = plt.subplots()
+scale = np.nanpercentile(np.abs(ux[mask]), 99.5)
+image = ax.pcolormesh(x, y, np.where(mask, ux, np.nan).T,
+                      shading="nearest", cmap="RdBu_r",
+                      vmin=-scale, vmax=scale)
+ax.set(xlim=(29,54), ylim=(-5,5), xlabel="x/D", ylabel="y/D")
+ax.set_aspect("equal")
+fig.colorbar(image, ax=ax, label="CCD mode 1 $u_x$")
+plt.show()
+```
+
+For a fresh publication, preserve the same field semantics, record the source result manifest hash, record plotting parameters, refuse accidental overwrite, and create a new sibling output. Do not copy or regenerate the dense canonical source arrays into the repository.
+
+## Provenance
+
+The publication `manifest.json` authenticates its files and records its source result manifest hash. `plot_data.json` repeats the source hash and field semantics. A third-party figure is reproducible only when its source hash, field selection, mask, coordinate orientation, color limits, and plotting parameters are retained together.
@@ -1,38 +1,47 @@
-# CCD three-part reset: final result boundary
+# Corrected Kármán validated review-fix result

-## Direct-dq facts
+The validated non-destructive publication transaction is:

-Two schema-v3, three-role, 450-control acquisitions completed with the corrected Legacy `q/U0` decoder and exact acquisition-relative timelines:
+`/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/publication-review-fix-per-cycle-20260808T1736+0800/`

- Kármán `karman_re100`: field interval 2000; retain exactly the 120 snapshots with relative lattice step `>120000` (122000–360000). The weighted vector-RMS of `e_target = q_ctl - q_target` is `0.15965716015602738`. Algebraic closure `e_target = dq_ctl - dq_tar` has maximum absolute residual `9.5367431640625e-07`.
- Illusion `illusion_1.0L`: field interval 1250; retain exactly the 144 snapshots with relative lattice step `>90000` (91250–270000). The weighted vector-RMS target error is `0.08031177071338107`. Maximum closure residual is `1.7881393432617188e-07`.
+It is bound to:

-Both results preserve full-resolution instantaneous role and difference fields, means, solver-mask intersection, exact x/D=35/40/45 profiles, mask-aware vorticity, and declared prefix/suffix convergence diagnostics. The 90/120 Kármán and 108/144 Illusion windows show decreasing mean-field deviations, but one trajectory cannot provide independent-realization uncertainty or by itself prove asymptotic convergence.
+- `/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/results/roi-mean/`
+- `/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/results/cycle-template-ccd-review-fix-per-cycle-20260808T1736+0800/`
+- `/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/phase/drl-review-fix-per-cycle-20260808T1736+0800/`
+- `/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/phase/constant_mean-review-fix-per-cycle-20260808T1736+0800/`

-These are same-time direct differences. They are not phase-conditioned results. `dq_ctl` and `dq_tar` share `-q_blk`; similarity between them is descriptive and is not causal mechanism evidence. The momentum-flux quantity is explicitly an incomplete proxy, not a complete momentum balance.
+The primary inclusive wake ROI is `34 <= x/D <= 54`, `|y/D| <= 5` on 80,200 common solver-fluid points. ROI mean target errors remain zero `0.46943522730424153`, constant mean `0.11101407390593347`, and DRL `0.08577878599374199`; full-domain metrics remain secondary.

-Authoritative immutable results:
+The primary CCD remains `N=190`, `Q=1`, `M=160400`, rank 3, with strengths `[0.017414192011156633, 0.00696948952702867, 0.001050035589643993]`. True retained-boundary per-cycle alternatives give:

- `evidence/direct-dq-karman-burn120000/`
- `evidence/direct-dq-illusion-authorized-burn90000/`
+- 8 bins: spectrum `[0.017211289396355783, 0.0067331096086703296, 0.0009696007073896092]`, rank-3 projector cosine `0.9947903592280113`
+- 12 bins: spectrum `[0.01707422700055497, 0.007134728727796222, 0.001118003774459821]`, rank-3 projector cosine `0.9977713810248983`
+- half-bin 10: spectrum `[0.017596242405215415, 0.00686633344257812, 0.0010454157443352114]`, rank-3 projector cosine `0.9982476660894443`

-A validated reload depends on the immutable live acquisition roots recorded by absolute path in each result; the result directories are not standalone portable evidence.
+The global least-favorable check is `bin_and_origin/bins8` with projector cosine `0.9947903592280113`. The global largest relative singular-value shift is `0.07660205334721537`. The alternative sample counts are N=152, N=228, and N=190; individual modes/action vectors are only the primary SVD basis representation of the stable rank-3 subspace.

-## Original CCD capability
+Fresh-process phase/ROI reload, essential CCD recomputation, and publication reload all passed. No CFD or CUDA was run. Older unsuffixed immutable artifacts remain preserved and superseded; no destructive canonical promotion was performed. Final independent review and the durable-memory update are complete.

-`original_ccd/` implements the original full-field Lyu operator `A = P U†/(N sqrt(LQ))`, weighted modes, exact lag/block handling, coefficients, lag functions, and rank-selected reconstruction/residual. It does not use preliminary POD, row standardization, or whitening. The published-scale synthetic derivation and production-reference tests pass.
+## Derived mode-interpretation extension

-The frozen real-case observable decision was executed for Kármán only. The authoritative Q=1, tau=0 result is `evidence/real-ccd-karman-q1-tau0-burn120000-v1/` (manifest SHA-256 `97d5cb300d642295bd3dffedd85d940fde06c8cc1c236fade390db78c3e3810e`). Its singular values are `0.14676751183519246`, `0.05920814023903697`, and `0.0063884442355860915`; squared cross-correlation strengths are `0.02154070253029336`, `0.003505603870565469`, and `0.00004081221975119316`. The three action means (front/upper/lower, native units) are `-0.004750394590640402`, `-0.0417400509895136`, and `0.03969304291531443`. The result has numerical rank 3, no degenerate singular blocks, modes shape `(1299184, 3)`, coefficients shape `(3, 120)`, and weighted relative residuals `0.7823805118116768`, `0.4698397615223296`, `0.45542496362972573` after complete blocks 1, 2, 3. Fresh-process live-provenance reload and essential identity recomputation passed. The singular values are cross-correlation strengths, not field energy, explained variance, canonical coefficients, causal effects, mechanism evidence, uncertainty, or a CCD-versus-POD comparison. Illusion CCD was not run and requires a separate user authorization decision.
+The validated artifact-only interpretation sibling is `/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/results/mode-interpretation-20260809T1025+0800/`, with five writing-oriented PNG/PDF figure groups at `/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/publication-mode-interpretation-20260809T1025+0800/`. It is derived from, and does not replace, the primary corrected v3 publication authority.
+
+It rebuilds the exact admitted `M=160400`, `N=190` residual through the cycle-template CCD matrix path, proves bit equality of the admitted empirical mean before centering, and applies full-fit weighted POD only as a descriptive field-energy baseline. Circular harmonics are computed per `(19,10)` cycle/phase array and retain every cycle's amplitudes/phases; they are not temporal PSD. Individual mode parity remains basis-dependent, while rank-3 reflection and CCD/POD comparisons include subspace-level diagnostics. The result manifest is `ab43cc40a60ece5a5a2d07c67fb6dd579f07195609d05d48fd9b96c37a4665ae`; the figure-publication manifest is `56208bcb5cee282797117ab59fac84a251fc36c1f311054e2afdd1f3319425aa`.


-## Karman figure subset
+## Interpretation closeout status

-The deterministic CPU artifact-derived package is `evidence/real-ccd-karman-figures-v1/`. It uses only the public verified real-CCD loader, retains the 1280 x 512 grid without downsampling, and publishes PNG/PDF figures plus quantitative Markdown/JSON interpretation. Selected reconstruction snapshots are acquisition-relative lattice steps 122000, 242000, and 360000; no full MxN reconstruction is persisted. The governing two-case figures review remains partial/in progress because Illusion was not run and awaits user authorization. The unique next entry remains that authorization decision.
+The same-object full-fit weighted rank-3 POD comparison gives principal cosines `[0.9999168, 0.9991137, 0.2639093]`. Modes 1 and 2 are strongly antisymmetric and have a descriptive streamwise shift relation of about `2.3D` with negative velocity/vorticity correlation; this is a shift-related basis observation, not a frequency, response-time, or actuator-causality result. Mode 3 is mixed-parity and near-wake/enstrophy localized, so it is not assigned a mechanism name.

-## Kármán dynamic-increment final publication
+The interpretation extension is artifact-only and full-fit: no held-out prediction, CCD superiority, independent uncertainty, causality, or mechanism claim.

-The accepted-plan final stage is complete at `data/karman-dynamic/karman-dynamic-v1-production/publication-v2/`. Six artifact-only PNG/PDF figure pairs and concise `RESULTS.json`/`RESULTS.md` report the four-role mean performance decomposition, DRL-versus-constant 10-bin phase statistics and centered phase difference, temporal negative-lag CCD spectrum/modes/left lag functions/common-support sensitivity, and an explicit phase-domain **DOWNGRADE**. Zero appears only in mean statistics; no zero phase figure or claim exists because its phase gate failed.

-Mean target errors are zero `0.2177701`, constant mean `0.1680192`, and DRL `0.1599946`. Thus zero-to-constant reduces error by `0.0497509`; constant-to-DRL contributes a further `0.00802465`, about `13.9%` of the total zero-to-DRL mean-error reduction. Constant/DRL cycle-mean 10-bin phase target errors are `0.2126461`/`0.1821272`. Temporal leading strengths are `2.77997e-3`, `1.91480e-5`, and `1.32484e-7`; common-support leading-three principal cosines remain at least `0.999998746`, while native support is sample-support-sensitive. Phase-domain CCD remains DOWNGRADE because first-harmonic rank is 2 instead of primary rank 3.
+## Pinball-inclusive convergence

-Strict code/science/claim review passed after one allowed remediation: `publication-v1` omitted the target bar despite the four-role requirement; immutable `publication-v2` fixes it. Dense fields were not deleted role-by-role: DRL is required for independent temporal recomputation, constant_mean and target for independent source-level mean/statistics recomputation, and zero has no compact field replacement after its failed phase gate. Focused Kármán tests: 26 passed; full active CCD suite: 178 passed; changed-file lints and scoped diff whitespace checks passed; old real-CCD/direct-dq loaders remained isolated from the publication process.
+The derived sibling `pinball-domain-20260809T1651+0800` extends the exact corrected residual to `29 <= x/D <= 54`, `|y/D| <= 5` on the four-role common fluid mask. It contains 99,264 fluid points and uses the same 19-by-10 samples, action observable, centering, and Q=1 operator.
+
+Raw full-fit weighted POD principal cosines are `[0.9999616609, 0.9998971005, 0.9976439804]`, and the minimum declared bin/cycle stability cosine is `0.9974651867`. The preregistered stop rule is met: the extended rank-3 residual is POD-equivalent; CCD contributes an observable/action coordinate, not a distinct field subspace. Modes 1/2 are predominantly antisymmetric and mode 3 predominantly symmetric, with solver-fluid vorticity concentrated in the geometry partition. These are descriptive diagnostics, not wall-vorticity, causality, mechanism, or response-time evidence. CCD is closed at this scope; SR and additional CFD are deferred.
+
+
+The refreshed publication sibling `publication-pinball-domain-20260809T1925+0800` is the preferred reader-facing view: mode panels use shared robust signed velocity limits and a descriptive body-local rotation proxy; plot-ready NPZ/JSON data are attached for later figure tuning.
@@ -0,0 +1,202 @@
+# CCD mathematical derivation and interpretation contract
+
+This document expands [CORRECTED_MATH_CONTRACT.md](karman_dynamic/CORRECTED_MATH_CONTRACT.md) into a writing-facing derivation. It describes exactly what is computed and what is not inferred.
+
+## 1. Coordinates, fields, and weights
+
+Let `x_D`, `y_D` be the nondimensional Cartesian grid coordinates and let `S` be the exact intersection of the DRL, constant-mean, target, and zero solver-fluid masks. For a selected domain, let `K` be the selected fluid points and `M_p=|K|`. Persisted velocity is already nondimensional; no second division by `U0` is allowed.
+
+The field vector is component-major:
+
+\[
+q = [u_x(k_1),\ldots,u_x(k_{M_p}),u_y(k_1),\ldots,u_y(k_{M_p})]^T,
+\quad M=2M_p.
+\]
+
+Coordinate quadrature is `w_x w_y`. The velocity weight matrix is diagonal,
+`W=diag(w_1,...,w_{M_p},w_1,...,w_{M_p})`; duplication gives both velocity components the same cell area.
+
+## 2. Mean estimands
+
+For role `r`, define the coordinate-weighted target error
+
+\[
+E_r=\left(\frac{\sum_{k\in K}w_k\|\bar q_r(k)-q_T(k)\|_2^2}{\sum_{k\in K}w_k}\right)^{1/2}.
+\]
+
+The primary corrections are
+
+\[
+\Delta\bar q_{0\to C}=\bar q_C-\bar q_0,
+\qquad
+\Delta\bar q_{C\to D}=\bar q_D-\bar q_C.
+\]
+
+The field identity
+
+\[
+\bar q_D-q_T=(\bar q_0-q_T)+\Delta\bar q_{0\to C}+\Delta\bar q_{C\to D}
+\]
+
+is exact up to floating-point arithmetic. Error reductions are scalar differences of `E`, not energies and not causal effects.
+
+## 3. Phase indexing and reference residual
+
+There are `C=19` independent DRL cycles and `B=10` phase bins. A sample is indexed by `(c,b)` and flattened in cycle-major/phase-minor order:
+
+\[
+j=10c+b,\qquad c=0,...,18,\quad b=0,...,9.
+\]
+
+The constant-control phase template is an ensemble mean over its 19 cycles:
+
+\[
+\hat q^C_b=\frac{1}{19}\sum_{c=0}^{18}q^C_{c,b}.
+\]
+
+Each role is centered over its own full retained field ensemble:
+
+\[
+\bar q_D=\frac{1}{190}\sum_{c,b}q^D_{c,b},
+\qquad
+\bar q_C=\frac{1}{190}\sum_{c,b}q^C_{c,b}.
+\]
+
+The field sample used by CCD is
+
+\[
+U_{c,b}=(q^D_{c,b}-\bar q_D)-(\hat q^C_b-\bar q_C).
+\]
+
+This is a phase-conditioned reference residual. DRL and constant cycles are not paired; therefore `U` is not a paired counterfactual and is not an isolated causal DRL response.
+
+## 4. Observable and centering
+
+Let `a^D_{c,b}` be the same-boundary effective action in native order `(front, upper, lower)`. Remove the exact retained DRL physical action mean `\mu_a` and then empirically center the samples:
+
+\[
+P_{c,b}=a^D_{c,b}-\mu_a,
+\qquad
+P_c=P-\frac1N P\mathbf 1\mathbf 1^T.
+\]
+
+The code uses the physical mean lineage checked against the fresh retained DRL statistic. The second centering is part of the literal CCD operator and is not optional preprocessing.
+
+## 5. Weighted Lyu CCD
+
+Stack fields and observables as columns:
+
+\[
+U=[U_0, ..., U_{N-1}]\in\mathbb R^{M\times N},
+\qquad P\in\mathbb R^{3Q\times N},
+\quad N=190, Q=1.
+\]
+
+After field and observable empirical centering, define
+
+\[
+A=\frac{P_c(W^{1/2}U_c)^T}{N\sqrt{3Q}}
+\in\mathbb R^{3Q\times M}.
+\]
+
+With `Q=1`, the normalization is `N sqrt(3)`. No whitening, standardization, POD preprojection, lag, wrap, or field-energy normalization is applied.
+
+Compute the compact SVD
+
+\[
+A=L\Sigma V^T.
+\]
+
+The left vectors `L` are observable-side CCD functions; the weighted field vectors `V` satisfy `V^T V=I` in transformed coordinates. Physical field modes are
+
+\[
+\Phi=W^{-1/2}V,
+\qquad \Phi^T W\Phi=I.
+\]
+
+The coefficient of mode `i` at sample `j` is
+
+\[
+c_{ij}=\phi_i^T W U_{c,j},
+\qquad
+\widehat U^{(r)}=\Phi_{1:r}C_{1:r}.
+\]
+
+The singular values are cross-correlation strengths. They are not field energy, explained variance, or reconstruction error fractions.
+
+## 6. Raw weighted POD comparison
+
+For the pinball-inclusive interpretation, use the same centered `U_c` and weights, but do not use `P` to select the basis. The method-of-snapshots Gram matrix is
+
+\[
+G=U_c^T W U_c\in\mathbb R^{190\times190}.
+\]
+
+Let `G r_i=\lambda_i r_i`, sorted descending. Then
+
+\[
+\sigma_i=\sqrt{\lambda_i},
+\qquad
+\psi_i=\frac{U_c r_i}{\sigma_i},
+\qquad
+U_c\approx\Psi_r(\Sigma_rR_r^T).
+\]
+
+The modes satisfy `\Psi^T W\Psi=I`. This is raw full-fit field-energy POD. It is not CCD-guided.
+
+For two rank-3 weighted-orthonormal bases `Phi` and `Psi`, principal cosines are the singular values of
+
+\[
+\Phi^T W\Psi.
+\]
+
+These are the basis-invariant subspace comparison. An aligned POD basis is obtained only by an orthogonal Procrustes rotation for visual side-by-side display; it does not alter the POD decomposition and cannot be used to claim CCD/POD identity or difference.
+
+## 7. Phase reconstruction and spatial diagnostics
+
+For each phase bin `b`, average reconstructed fields over cycles:
+
+\[
+\widehat U_b^{(r)}=\frac1{19}\sum_c\Phi_r C_{r,(c,b)}.
+\]
+
+A weighted regional residual and cosine are computed using fixed partitions:
+
+- geometry: `29–32D`;
+- near wake: `32–40D`;
+- middle wake: `40–47D`;
+- far wake: `47–54D`.
+
+These are registered coordinate partitions, not discovered wake boundaries.
+
+Reflection uses `u_x(x,y)→u_x(x,-y)` and `u_y(x,y)→-u_y(x,-y)`. Individual parity is basis-dependent; subspace reflection diagnostics are more stable.
+
+The publication's body-local rotation proxy is a visualization-only signed angular-momentum-like projection:
+
+\[
+R(x,y)=\frac{\sum_k g_k(x,y)[(x-x_k)u_y-(y-y_k)u_x]}{\sum_k g_k(x,y)},
+\]
+
+where `g_k` is a Gaussian localization around body `k`. It is not wall vorticity. The complete-stencil solver-fluid vorticity diagnostic is also auxiliary and should not be interpreted as boundary generation or actuator causality.
+
+## 8. Sign and basis caveats
+
+Every singular vector pair has the invariance
+
+\[
+(\phi_i,c_i)\equiv(-\phi_i,-c_i).
+\]
+
+Within a mode, relative signs between front and rear regions remain meaningful; the global sign is not. In near-degenerate subspaces, individual vectors may rotate while the subspace and rank-r reconstructions remain stable.
+
+## 9. What the subtraction supports
+
+Subtracting constant-control phase template removes the dominant reference backbone in an estimand-defined sense. It supports descriptive DRL-minus-reference co-variation. It does not prove that the subtraction isolates an additive physical force, a causal DRL contribution, an actuator response, or a term-wise SR mechanism. Any SR comparison must be described as term/action-coordinate association and must retain its own Legacy SR evidence contract.
+
+## 10. Code and artifact lookup
+
+- Residual/CCD: `karman_dynamic/cycle_template_ccd.py`.
+- Domain comparison/reconstruction/publication: `karman_dynamic/domain_comparison.py`.
+- Pure diagnostics: `karman_dynamic/mode_diagnostics.py`.
+- Raw POD: `karman_dynamic/pod_baseline.py` and the method-of-snapshots implementation in `domain_comparison.py`.
+- Compact arrays/config/summary/manifest: preferred external pinball-domain result indexed in [RESULTS_INDEX.md](RESULTS_INDEX.md).
@@ -1,114 +1,42 @@
-# CCD analysis: Kármán dynamic-increment campaign
+# Corrected Kármán CCD analysis

-This is the active authority. It supersedes the Illusion-next-entry checkpoint. The unique scientific question is the Legacy `karman_re100` (code Re100, physical `Re_D=50`) increment of time-varying DRL over constant control fixed to the same fresh DRL run's retained, three-channel `effective_applied_action` mean. `zero` is only the passive baseline; `target` only defines cloaking error. Neither is the dynamic CCD subtraction.
+## Five-minute entrypoint

-The frozen identity is `q_D(phi)-q_C(phi) = (mean(q_D)-mean(q_C)) + [(q_D(phi)-mean(q_D))-(q_C(phi)-mean(q_C))]`. The first term is a policy-induced statistical mean change; the second is phase-coherent unsteady change. Independent phase-conditioned trajectories are not pointwise counterfactuals, response measurements, or causal effects. Historical acquisition/direct-dq/Q=1 results below remain immutable evidence, but are superseded as the next execution direction.
+This package is the writing-ready documentation and validation layer for the corrected Legacy `karman_re100` CCD analysis. The computation is closed: this repository contains no new CFD campaign and no new scientific claim beyond the bounded results indexed here.

-`karman_dynamic/` is the campaign-specific schema and launcher. Roles are `target`, `zero`, `drl`, `constant_mean`; execution is DRL first, then constant_mean, target, zero. It wraps the unchanged active schema-v3 acquisition payload and reuses its exact clocks, q/U0 decoder, solver mask, requested/effective actions, and fresh role runtime. Campaign telemetry explicitly stores center-sensor `uy=sensors[:,3]`. Constant mean is hash-bound to the fresh DRL wrapper and retained effective actions, without hand entry or symmetry forcing.
+Read in this order:

-CFD execution requires `CONDA_DEFAULT_ENV=pycuda_3_10`, exactly one `CUDA_VISIBLE_DEVICES` token, CPU PPO inference, an exclusive campaign lease, a stable exact Optane symlink, fresh no-clobber roots, semantic reload after each child, and at least 30 s between starts (default 120 s). Failure quarantines the campaign and stops later roles. Do not run real CFD from `pinball_math`.
+1. [WRITING_HANDOFF.md](WRITING_HANDOFF.md) — manuscript-ready summary and safe wording.
+2. [RESULTS_INDEX.md](RESULTS_INDEX.md) — conclusion → number → artifact → figure lookup.
+3. [MATHEMATICAL_DERIVATION.md](MATHEMATICAL_DERIVATION.md) — continuous definitions and derivations.
+4. [CLAIMS.md](CLAIMS.md) — current, contextual, diagnostic, historical, and unsupported claims.
+5. [FIGURES_AND_DATA.md](FIGURES_AND_DATA.md) — figure inventory and third-party redraw contract.
+6. [FINAL_RESULTS.md](FINAL_RESULTS.md) — validated corrected result report.
+7. [karman_dynamic/CORRECTED_MATH_CONTRACT.md](karman_dynamic/CORRECTED_MATH_CONTRACT.md) — compact machine-facing math contract.
+8. [EXECUTION_CHECKPOINT.json](EXECUTION_CHECKPOINT.json) — machine-readable execution state.

-Current status: all four production roles were fresh-loaded. Immutable phase publications exist for every role at sibling `ROLE-phase-compact-v1` paths. DRL, constant_mean, and target pass; zero fails closed (period CV `0.0828308 > 0.05`, amplitude CV `0.101389 > 0.10`) and therefore publishes metrics only, with no compact fields. DRL–constant phase differences remain authorized because those two independent gates pass; zero phase-target metrics are unavailable. The immutable four-role result is `data/karman-dynamic/karman-dynamic-v1-production/dynamic-increment-v2`. Mean target errors are zero `0.2177701`, constant_mean `0.1680192`, and DRL `0.1599946`; corresponding target-error reductions are `0.0497509` and `0.00802465`. Dense `payload/fields.npz` remains intact for all roles and deletion is not yet authorized because zero has no sufficient compact replacement. The earlier `dynamic-increment-v1` is immutable partial-stage evidence superseded by v2.
+## Authority and path layers

+| Layer | Location | Meaning |
+|---|---|---|
+| Machine authority | `/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/` | Immutable source, phase, result, and publication parents bound by manifests and hashes. |
+| Repository reader view | `src/CCD_analysis/data/karman-dynamic/canonical/` | Tracked/visible reader-facing plots and plot data; it does not replace the external dense authority. |
+| Historical archive | `src/CCD_analysis/archive/` | Superseded code and evidence; never active imports and never current authority. |

-The authoritative primary temporal result is `data/karman-dynamic/karman-dynamic-v1-production/drl-temporal-negative-lag-ccd-v1`. It uses full-resolution mask-compressed `q_DRL(t)-mean(q_DRL)` and exact same-boundary three-channel effective-action fluctuations at 800-step cadence, with complete phase cycles as non-crossing blocks. The predeclared grid is `tau/800 = -17,...,0` (`Q=18`, about one measured shedding period), giving `N=21` complete columns and `M=1,299,184`. Leading cross-correlation strengths are `2.77997e-3`, `1.91480e-5`, and `1.32484e-7`; mode-1 left-lag energy peaks at `tau/800=-6` and has centroid `-7.969`, while modes 2/3 peak at `0/-5`. On fixed common support, dropping the oldest one/two lags changes the first three strengths by at most `1.58%`, `4.33%`, and `6.06%`, with leading-three weighted-subspace principal cosines at least `0.9999987`. Native endpoint support grows to `N=40/59` and substantially changes strengths and the third leading subspace direction, so those native-window comparisons are sample-support-sensitive and are not interpreted as timing evidence. This is closed-loop temporal co-variation only, with no causal or response-time claim.
+Preferred reader-facing artifact:

-The final deterministic artifact-only publication is `data/karman-dynamic/karman-dynamic-v1-production/publication-v2`. It contains six concise PNG/PDF figure pairs plus hash-bound JSON/Markdown, generated after fresh live-provenance reloads and essential recomputation. `publication-v1` is immutable superseded pre-review evidence; final review found that its first panel omitted the target reference, and the single remediation added the explicit fourth role in v2. No figure reports zero phase results.
+- Result: `/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/results/pinball-domain-20260809T1925+0800/`
+- Repository publication view: `data/karman-dynamic/canonical/publication-pinball-domain-20260809T1925+0800/`
+- External publication view: `/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/publication-pinball-domain-20260809T1925+0800/`

-Dense-field deletion decision: **RETAIN ALL FOUR ROLE SOURCES**. DRL dense fields remain required to independently recompute the temporal CCD and its provenance; constant_mean and target dense fields remain required to independently recompute their mean/statistical products from source; zero failed the phase gate and has no compact field artifact, so its dense source is the only source-level basis for its accepted mean/statistics. The compact products suffice for currently published downstream figures, but not for every source-level downstream recomputation/provenance contract; therefore deletion is not independently justified.
+The earlier `1651` pinball-domain sibling is retained as a superseded predecessor. Do not delete or overwrite it. The old wake-only corrected result remains the primary corrected performance authority; the pinball-inclusive result is a derived interpretation sibling.

-The exploratory phase-domain CCD is complete at `data/karman-dynamic/karman-dynamic-v1-production/drl-constant-phase-domain-ccd-v1` and is **DOWNGRADED**, not retained as a stable three-direction result. It uses only the separately centered phase-coherent difference `Δq′_phase(φ)=(q_DRL(φ)-mean(q_DRL))-(q_constant_mean(φ)-mean(q_constant_mean))` and the DRL phase-conditioned effective-action fluctuation in the literal weighted Lyu `Q=1` operator. The 10-bin strengths are `0.0420404`, `0.0157939`, and `0.00414469` (rank 3). The 8/12-bin and half-bin-origin projectors are stable (minimum leading-three cosine `0.995038`), and left functions are stable (minimum absolute cosine `0.999867` across those bin tests), but first-harmonic truncation has rank 2 rather than 3; therefore rank stability fails. Circular offsets change strengths but are phase offsets only, never time-response lags. The immutable artifact is published because its scientific contract and downgrade boundary are explicit; it supports only low-order exploratory circular co-variation.
+## Scientific bottom line

-Recommended smoke plan command (prints fresh child command; no CFD):
+The mean-control increment is dominant: constant mean accounts for about `93.42%` of the zero-to-DRL ROI target-error reduction, while the DRL dynamic increment accounts for about `6.58%`. The pinball-inclusive dynamic residual has a stable rank-3 CCD subspace, but raw full-fit weighted POD spans essentially the same field subspace (principal cosines `0.99996166`, `0.99989710`, `0.99764398`). CCD's bounded added value is observable/action-coordinate association, not a distinct field decomposition or proven mechanism.

-```bash
-PYTHONPATH="$PWD/src" conda run -n pinball_math python -m CCD_analysis.karman_dynamic orchestrate --campaign-id karman-dynamic-v1 --root src/CCD_analysis/data/karman-dynamic/karman-dynamic-v1-smoke --warmup-intervals 1 --collect-boundaries 2 --launch-delay-seconds 120 --smoke
-```
+All interpretation is descriptive phase-coherent DRL-minus-constant-template co-variation. It is not a paired counterfactual, causal response, mechanism identification, response-time estimate, independent-realization uncertainty estimate, explained-variance claim, or CCD superiority claim.

-Recommended exact DRL smoke execution command (real CUDA CFD; run only after Optane mapping/lease preflight):
+## Reproduction boundary

-```bash
-CONDA_DEFAULT_ENV=pycuda_3_10 CUDA_VISIBLE_DEVICES=0 PYTHONPATH="$PWD/src" python -m CCD_analysis.karman_dynamic role --campaign-id karman-dynamic-v1 --role drl --output src/CCD_analysis/data/karman-dynamic/karman-dynamic-v1-smoke/drl --warmup-intervals 1 --collect-boundaries 2 --launch-delay-seconds 120 --smoke
-```
-
-## Superseded historical authority (immutable evidence)
-
-This is the active authoritative surface. Historical material under `archive/2026-08-03-pre-three-part-reset/payload/` is immutable, non-authoritative, and never imported by active code.
-
-## Acquisition status
-
-`acquisition/` now defines CPU-import-safe contract/artifact schema v3 for exactly `karman_re100` and `illusion_1.0L`, each with roles `q_target`, `q_blk`, and `q_ctl`. It includes frozen geometry, signed front/upper/lower action identities, exact absolute lattice clocks, same-step telemetry, solver-derived masks, full-grid coordinates/velocity fields, controller/history identity, no-clobber publication, and fresh-process sequential orchestration.
-
-Kármán's code label Re100 uses the historical `2D` reference and is physically `Re_D=50`. Its front/rear/sensor x locations are 30/31.3/40 D. Illusion uses the strict +11D deployment geometry: front/rear/sensors/target x = 30/31.3/41/31 D, rear y = ±0.75 D, sensors y = +2/0/-2 D. This Illusion deployment differs from training geometry; replay/history smoke is mandatory before production.
-
-Velocity fields are the Legacy solver's nondimensional velocity `q/U0`, decoded exactly on solver-flagged fluid cells as `ux=(f1+f5+f8-f3-f6-f7)/u0` and `uy=(f2+f5+f6-f4-f7-f8)/u0`. This is not momentum divided by density. Nonfluid cells are exact zero and their unused populations may be ignored. The decoder schema and formula hash are mandatory artifact identities.
-
-All six artifacts under `evidence/smoke-20260804/` were produced with the incorrect density-normalized decoder. They are withdrawn, invalid under the current schema, and retained only as immutable negative evidence; they must not be overwritten. Those v2 artifacts are also structurally obsolete because they lack complete control-boundary lineage. The corrected schema-v3 smoke and the two 450-control three-role production acquisitions have completed. Schema v3 persists the exact pre-action `(150,12)` FIFO, every interval-average `(control_count,12)` boundary observation independently of field cadence, every policy source/input observation and source hash, zero-origin policy harmonic phase indices, and complete requested-action control histories. Current production roots are `evidence/production-20260804-q-over-u0-v3-karman-450-fi2000/` and `evidence/production-20260804-q-over-u0-v3-illusion-authorized-450-fi1250/`.
-
-`q_target` is a desired reference generated with different bodies. It is not a same-checkpoint counterfactual. Equal absolute time establishes same-time sampling only; `phase_reference` is evidence for later phase validation and does not itself prove same phase.
-
-### Commands
-
-CPU tests (no CFD):
-
-```bash
-PYTHONPATH="$PWD/src" conda run -n pinball_math python -m pytest src/CCD_analysis/tests/test_acquisition.py -q
-```
-
-Inspect the sequential fresh-process command plan without running CFD:
-
-```bash
-PYTHONPATH="$PWD/src" conda run -n pinball_math python -m CCD_analysis.acquisition orchestrate --case karman_re100 --output /fresh/path --control-count 450 --field-interval 2000
-```
-
-CFD entry points require `CONDA_DEFAULT_ENV=pycuda_3_10` and fail otherwise. Do not use `conda run` if it does not propagate that variable correctly. The first independent gate failed on its no-op runner/EMA lifecycle, and the second re-review failed on exact initialization and controller-reference semantics. Both sets are now remediated. Initialization is explicitly `stabilize → full current/temp/lifecycle checkpoint → zero-action normalization trajectory → exact restore → FIFO warmup`, with zero post-stabilization EMA restored before warmup. Karman's first PPO input is exactly zero. Illusion PPO consumes only the frozen training normalization (`9ec5dd…`) and two target-force harmonics (`f13566…`) from the frozen Illusion training-reference files; newly measured +11D eight-channel harmonics are phase evidence only, never controller input. The deployment/reference mismatch still requires the strict +11D compatibility replay/history smoke. The third review then failed only because solver, rollout, and policy-phase clocks were conflated. They are now explicit: solver-absolute lattice/control lineage comes from the public solver accessor; acquisition-relative lattice/control starts at zero; policy harmonic phase independently starts at zero. `lattice_steps`/`sample_ids` are solver-absolute, while `acquisition_relative_lattice_steps`, rollout-relative `control_indices`, and solver-absolute control indices are persisted separately. A final static review found that snapshots called the boundary-only clock API during an active split. Legacy now has a distinct `active_step_clock_state()` valid only after at least one completed split step; snapshots use it, while initialization and post-interval checks retain boundary-only `solver_clock_state()`. The corrected q/U0 decoder, complete lineage contract, and explicit orchestration schedule passed fresh CUDA smoke before the two production acquisitions. Kármán used 450 controls with field interval 2000; Illusion used 450 controls with field interval 1250. These runs establish successful acquisition and exact artifact lineage, not physical-phase equality or a mechanism claim.
-
-## Direct-dq status
-
-`direct_dq/` is the completed CPU-only strict same-time analysis core. It accepts three explicit completed active acquisition artifact directories, revalidates each manifest plus file/config/state-array hashes, and requires exact case/role/schema, finite float32 `(time,x,y)` fields, exact float32 coordinates, a common grid, and exact full acquisition-relative timelines before selection. Solver-absolute origins may differ by role when each lineage is internally valid. A required exclusive `--start-after-relative-step` burn-in bound and optional inclusive `--end-at-relative-step` select samples by exact integer physical-step inequalities only; no boundary lookup, nearest match, or index trimming is allowed. A supplied end must be the terminal selected sample, and the interval must be nonempty. There is no trimming, nearest-time/station substitution, phase guessing, crop, translation, or coordinate-generated mask.
-
-Only selected role columns enter the analysis. Results preserve the original common timeline, exact selection bounds, selected original indices, selected steps, and selected count. The transparent analysis domain is the intersection of the three preserved solver-derived fluid masks. Outputs include instantaneous and full-resolution time-mean `e_target`, `dq_ctl`, and `dq_tar`; all three role means; exact-station streamwise profiles; weighted vector RMS target error; signed target-relative streamwise deficit; an explicitly incomplete momentum-flux proxy; positive-deficit wake area/centroid/width; mask-aware nonuniform-coordinate vorticity; and declared nested prefix/suffix convergence diagnostics. `dq_ctl` and `dq_tar` share `-q_blk`, so agreement is not mechanism evidence. Phase-conditioned output fails closed because cross-role physical-phase equality has not been independently proven. Prefix/suffix windows are not independent-realization uncertainty.
-
-Run on completed artifacts in `pinball_math` (repeat station/window options as needed; station tokens are preserved and canonically converted to exact float32 grid values without nearest/tolerance matching; the final window must equal the full sample count):
-
-```bash
-PYTHONPATH="$PWD/src" conda run -n pinball_math python -m CCD_analysis.direct_dq \
-  --case karman_re100 --q-target /path/q_target --q-blk /path/q_blk --q-ctl /path/q_ctl \
-  --output /fresh/result --start-after-relative-step 120000 \
-  --station-x-D 35 --station-x-D 40 --station-x-D 45 \
-  --window-size 30 --window-size 60 --window-size 90 --window-size 120
-```
-
-For the completed Illusion acquisition, use burn step `90000`, stations `35,40,45`, and windows `36,72,108,144`. The canonical production grids contain each station exactly once as float32. Convergence windows are validated against the selected count.
-
-Results are immutable/no-clobber atomic directories containing `arrays.npz`, `summary.json`, `config.json`, `input_hashes.json`, and a hash manifest. Provenance-validated `load_result()` always rereads all three acquisition directories at their recorded absolute paths, revalidates their complete semantics and hashes, and requires their masks/grid/full timeline to equal the persisted result and each selected instantaneous role array to equal the exact selected live input columns. Therefore results are not portable by themselves: moving, deleting, or changing an acquisition directory makes provenance validation fail closed. `load_result_metadata_unverified()` is explicitly metadata/internal-science-only and cannot support a provenance claim; publication and CLI never use it as success validation. The completed real results are `evidence/direct-dq-karman-burn120000/` (120 retained fields) and `evidence/direct-dq-illusion-authorized-burn90000/` (144 retained fields). Their summaries report weighted vector-RMS target errors of 0.1596571602 and 0.08031177071338107, respectively. These are direct same-time estimands after the declared burn-in selections; they are not phase-conditioned, independent-realization uncertainty, or causal mechanism evidence. Tests additionally use synthetic artifacts written through the active acquisition writer.
-
-## Original CCD production status
-
-`original_ccd/ORIGINAL_CCD_MATH.md` is the active mathematical contract for the original full-field Lyu CCD. The public CPU package now implements literal full-field `U`, lag-stacked `P`, `A = P U†/(N sqrt(LQ))`, direct rectangular SVD, physical diagonal or dense complex-HPD weighting, separately declared `U`/`P` centering, exact timestamp/block lag construction, weighted modes and physical-amplitude coefficients, and selected-mode reconstruction/residuals. There is no POD pre-reduction, whitening, row standardization, nearest-time matching, or implicit centering.
-
-`build_lagged_observables()` preserves exact admitted field-column indices, and `fit()` applies that mapping to `U`. Reconstruction is the weighted projection of field snapshots onto selected CCD modes, not observable prediction. Exactly degenerate singular blocks identify subspaces/projectors rather than unique individual modes; the optional deterministic phase convention does not resolve that non-uniqueness. See `original_ccd/README.md` for API usage.
-
-Production tests compare against the private literal derivation reference for hand, random real/complex, diagonal/dense weighted, centered, multiobservable, and exact-lag cases; verify singular equations, chunk invariance, weighted orthogonality, full supported-basis reconstruction, truncation residual/projector behavior, validation failures, and a smaller Lyu equations (3.1)-(3.2) production case. The exact published-scale stochastic derivation test remains in the full suite and is not redundantly rerun in a second production test. The production API tests run no CFD. The real Kármán CCD result is reported below; Illusion CCD was not run. No observable prediction, causal result, or CCD-versus-POD claim is included.
-
-## Real-case CCD contract status
-
-`original_ccd/REAL_CASE_CCD_CONTRACT.md` freezes the first real-data estimand for both cases: centered full-resolution `dq_ctl` on the authoritative solver-mask intersection; centered q_ctl `effective_applied_action` front/upper/lower channels at each exact field time in native units; and the literal weighted original CCD with `Q=1`, `tau=0`. It fixes component-major mask flattening, coordinate quadrature, mean-action algebra, prohibited preprocessing, provenance validation, OOM-safe passes, immutable result requirements, and claim limits. The active `real_ccd/` package implements mandatory direct-dq/live-acquisition provenance validation, exact q_ctl field-time effective actions, component-major mask compression, explicit two-sided centering and coordinate weighting, a three-pass 3-by-M thin-SVD decomposition, streamed coefficients/block residuals, conservative fail-closed RAM/scratch accounting, an immutable verified result loader, and a CPU-only `preflight`/`run` CLI with explicit direct-dq root, chunk size, safe host admission budget, fresh output, and no-clobber publication.
-
-The authoritative Kármán Q=1, tau=0 decomposition is complete at `evidence/real-ccd-karman-q1-tau0-burn120000-v1/`, using the immutable `evidence/direct-dq-karman-burn120000/` input, chunk size 8, 120 samples, and 1,299,184 spatial degrees of freedom. Singular values (cross-correlation strengths) are `[0.14676751183519246, 0.05920814023903697, 0.0063884442355860915]`; their squares are `[0.02154070253029336, 0.003505603870565469, 0.00004081221975119316]`. The numerical rank is 3 with no degenerate blocks, and weighted relative residuals after complete blocks 1/2/3 are `[0.7823805118116768, 0.4698397615223296, 0.45542496362972573]`. These values are not field energy, explained variance, or canonical coefficients. The fresh-process loader passed live provenance and essential-identity recomputation. Illusion CCD was not run and is not authorized by this completion; the unique next entry is a user decision on `real-ccd-illusion` authorization.
-
-The Kármán-only figure subset is complete at `evidence/real-ccd-karman-figures-v1/`: full-resolution PNG/PDF spectrum, authoritative mean and three physical modes, zero-lag action left vectors, exact-step coefficients, complete-block residuals, and on-demand rank reconstructions at relative steps 122000/242000/360000. The two-case figures todo remains `in_progress`: Illusion figures await the same user authorization decision as Illusion CCD, which remains the unique next entry.
-
-Run CPU tests (no CFD):
-
-```bash
-PYTHONPATH="$PWD/src" conda run -n pinball_math python -m pytest src/CCD_analysis/tests -q
-```
-
-## Other active parts
-
- `original_ccd/` — production original full-field CCD API and its frozen mathematical contract.
- `tests/` — active CPU tests only.
- `evidence/` — machine-readable review evidence.
-
-`EXECUTION_CHECKPOINT.json` is the authoritative resumable execution state.
+Documentation and plotting use `pinball_math` and existing immutable artifacts only. They must not initialize CFD or CUDA. See [FIGURES_AND_DATA.md](FIGURES_AND_DATA.md) for the reader-facing redraw command and [CLAIMS.md](CLAIMS.md) for the claim boundary.
@@ -0,0 +1,44 @@
+# CCD results index
+
+This is the first lookup when a number or figure is needed. Paths are immutable identities; display copies are reader views.
+
+## Conclusion-to-evidence map
+
+### R1 — Mean control is the dominant cloaking contribution — Current
+- Estimand: inclusive wake ROI `34 <= x/D <= 54`, `|y/D| <= 5`, exact four-role common fluid mask; coordinate-weighted velocity target error.
+- Errors: zero `0.4694352273`, constant mean `0.1110140739`, DRL `0.0857787860`.
+- Reduction: zero→constant `0.3584211534`; constant→DRL `0.0252352879`; total `0.3836564413`.
+- Shares: constant mean `93.42%`; dynamic increment `6.58%`.
+- Evidence: external `.../canonical/results/roi-mean/`; reader report [FINAL_RESULTS.md](FINAL_RESULTS.md).
+- Figure: corrected publication `01_roi_mean_performance` and `02_roi_mean_target_error_fields` under the external review-fix publication.
+
+### R2 — Corrected CCD identifies a stable rank-3 dynamic residual subspace — Current, bounded
+- Residual: DRL cycle/bin field minus constant-control ensemble phase template, with separate role centering.
+- Samples: `19 cycles × 10 bins = N=190`; `Q=1`; wake-only `M=160400`.
+- Strengths: `0.0174141920`, `0.0069694895`, `0.0010500356`.
+- Stability: least favorable wake rank-3 projector cosine `0.9947903592` across persisted bin/origin alternatives.
+- Evidence: external `.../canonical/results/cycle-template-ccd-review-fix-per-cycle-20260808T1736+0800/`.
+- Math: [MATHEMATICAL_DERIVATION.md](MATHEMATICAL_DERIVATION.md).
+
+### R3 — Pinball-inclusive field and POD subspaces are practically equivalent — Current, bounded
+- Domain: `29 <= x/D <= 54`, `|y/D| <= 5`; `99,264` common fluid points; `M=198528`, `N=190`.
+- CCD/POD principal cosines: `0.9999616609`, `0.9998971005`, `0.9976439804`.
+- Stability: minimum declared bin/cycle projector cosine `0.9974651867`.
+- Evidence: external `.../canonical/results/pinball-domain-20260809T1925+0800/`.
+- Preferred figures/data: [publication-pinball-domain-20260809T1925+0800](data/karman-dynamic/canonical/publication-pinball-domain-20260809T1925+0800/).
+- Meaning: CCD adds action-conditioned interpretation to an essentially POD-equivalent residual field subspace.
+
+### R4 — Near-body mode interpretation — Diagnostic, not mechanism
+- Modes 1/2 are predominantly antisymmetric (`0.972`, `0.930`); mode 3 is predominantly symmetric (`0.833`).
+- Action coordinates use `front`, `rear-symmetric`, `rear-antisymmetric`.
+- Body-local rotation proxy and signed velocity are visualization diagnostics only; they are not wall-vorticity or causality evidence.
+- Figure: `02_ccd_mode_identity`; reconstruction: `03_phase_reconstruction`.
+
+### R5 — SR interface — Contextual cross-analysis
+- Current SR law and term ranking are documented in [SR claims](../SR_analysis/CLAIMS.md).
+- The rear constant is mean-control backbone context; SR–CCD correspondence would remain descriptive term/action association, not causal decomposition.
+- Do not merge Legacy SR and CCD estimands or claim SR term necessity from CCD.
+
+## Stop decision
+
+CCD is closed at this scope. Do not add CFD, new role arms, Re50, front-variable/rear-variable cases, or SR interventions under the CCD authority. A future SR–CCD correspondence study would require a separate plan and claim ledger.
@@ -0,0 +1,36 @@
+# CCD writing handoff
+
+## Read this first
+
+The main story is: constant mean control supplies most measured cloaking improvement; the remaining DRL-minus-constant-template dynamic residual is phase-coherent and rank-3; pinball-inclusive raw POD spans essentially the same field subspace; CCD's bounded extra information is the association between this residual and effective-action coordinates.
+
+Recommended reading order: [README.md](README.md) → [RESULTS_INDEX.md](RESULTS_INDEX.md) → [MATHEMATICAL_DERIVATION.md](MATHEMATICAL_DERIVATION.md) → [CLAIMS.md](CLAIMS.md) → [FIGURES_AND_DATA.md](FIGURES_AND_DATA.md) → [FINAL_RESULTS.md](FINAL_RESULTS.md).
+
+## Safe manuscript paragraph
+
+In the tested Legacy Kármán setting, constant-mean rotation accounts for approximately 93.42% of the zero-to-DRL reduction in the inclusive wake-ROI target error, while the remaining 6.58% is a dynamic increment. We therefore analyzed a phase-coherent residual formed by subtracting a separately centered constant-control phase template from separately centered DRL cycle/bin fields. A literal weighted rank-3 CCD of this residual against the same DRL effective-action observable was stable under the declared cycle and bin checks. Extending the domain to include the pinball, a raw full-fit weighted POD produced nearly the same rank-3 field subspace; CCD therefore serves primarily as an observable/action-coordinate interpretation of a POD-equivalent residual subspace, not as evidence of a distinct field basis or causal mechanism.
+
+## Figure choices
+
+Use the preferred publication under `data/karman-dynamic/canonical/publication-pinball-domain-20260809T1925+0800/`. Start with `01_domain_estimand`, `02_ccd_mode_identity`, `03_phase_reconstruction`, and `04_raw_pod_vs_ccd`. `02` uses shared robust signed velocity scales and a descriptive body-local rotation proxy. Always retain its caption/RESULTS context.
+
+## Common misreadings
+
+- A mode and its global sign-reversed version are identical; compare relative signs within a field or use the reconstructed product `coefficient × mode`.
+- Ten-bin harmonics are circular phase diagnostics, not temporal PSD frequencies.
+- Near/geometry/middle/far are fixed coordinate partitions, not independently discovered physical zones.
+- Raw POD is not CCD-guided; aligned POD is visual registration only.
+- DRL minus constant-template is a reference residual, not a paired counterfactual or causal DRL field.
+- The body-local rotation proxy is not wall vorticity.
+
+## Detail lookup
+
+- Definitions: [MATHEMATICAL_DERIVATION.md](MATHEMATICAL_DERIVATION.md).
+- Numbers and paths: [RESULTS_INDEX.md](RESULTS_INDEX.md).
+- Claim wording: [CLAIMS.md](CLAIMS.md).
+- Plot fields and redraw: [FIGURES_AND_DATA.md](FIGURES_AND_DATA.md).
+- Machine state: [EXECUTION_CHECKPOINT.json](EXECUTION_CHECKPOINT.json).
+
+## Blind-navigation check
+
+A reader starting at `README.md` can reach the main conclusion in `RESULTS_INDEX.md`, the continuous derivation in `MATHEMATICAL_DERIVATION.md`, the claim boundary in `CLAIMS.md`, and the preferred four-figure/data package in `FIGURES_AND_DATA.md` without inspecting `evidence/`, `canonical/`, or `archive/` directories. Machine-level details remain one hop away through `EXECUTION_CHECKPOINT.json` and `artifact_catalog.json`.
@@ -0,0 +1,46 @@
+{
+  "authority_root": "/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical",
+  "entries": [
+    {
+      "manifest_sha256": "f7f6141042678e52314a4d8f76d589108da8e844c8ce504f8f523a6f88a17e73",
+      "path": "/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/results/roi-mean",
+      "preferred": false,
+      "role": "mean performance",
+      "schema_id": "ccd-karman-roi-mean/v2",
+      "scientific_status": "current authority"
+    },
+    {
+      "manifest_sha256": "b114fbcc5f821034829ca3771ac5ba13df4b75df312b0615b662f81a9e5737d8",
+      "path": "/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/results/cycle-template-ccd-review-fix-per-cycle-20260808T1736+0800",
+      "preferred": false,
+      "role": "wake-only corrected CCD",
+      "schema_id": "ccd-karman-cycle-template/v1",
+      "scientific_status": "current authority"
+    },
+    {
+      "manifest_sha256": "7e107a895364ed3e34abf2be513a3ce9e12736766ec6383c5b831891802d5828",
+      "path": "/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/results/pinball-domain-20260809T1925+0800",
+      "preferred": true,
+      "role": "pinball-inclusive CCD/POD",
+      "schema_id": "ccd-karman-pinball-domain/v1",
+      "scientific_status": "preferred derived result"
+    },
+    {
+      "manifest_sha256": "31888f8eeb6c70dd9f0ff3449fa879f44abf381338ec261b0e45f69603691fdd",
+      "path": "/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/publication-pinball-domain-20260809T1925+0800",
+      "preferred": true,
+      "role": "reader-facing four-figure publication",
+      "schema_id": "ccd-karman-pinball-publication/v1",
+      "scientific_status": "preferred reader view"
+    },
+    {
+      "manifest_sha256": "8b952d5dd8cf9813a241da31295497e5c79db6df7dc4ed12b95a6706f99e37b6",
+      "path": "/home/frank14f/optane/DynamisLab/ccd/karman-dynamic/canonical/results/pinball-domain-20260809T1651+0800",
+      "preferred": false,
+      "role": "pinball-inclusive predecessor",
+      "schema_id": "ccd-karman-pinball-domain/v1",
+      "scientific_status": "superseded predecessor"
+    }
+  ],
+  "schema_id": "ccd-artifact-catalog/v1"
+}
@@ -0,0 +1,13 @@
+# Corrected Kármán analysis contract
+
+The primary performance domain is the inclusive wake ROI `34 <= x/D <= 54` and `|y/D| <= 5`, restricted to the exact intersection of the four solver fluid masks. Coordinate quadrature is `w_x w_y`; persisted velocity is already nondimensional and is not divided by `U0` again. For role `r`, `E_r^ROI = [sum_ROI w ||mean(q_r)-mean(q_target)||^2 / sum_ROI w]^(1/2)`.
+
+The primary mean estimands are `mean(q_constant)-mean(q_zero)` and `mean(q_DRL)-mean(q_constant)`. Full-domain errors are retained only as explicitly secondary diagnostics. Error reductions obey exact zero-to-constant plus constant-to-DRL additivity.
+
+For DRL cycle `c` and phase bin `b`, `delta q'_(c,b) = (q_DRL_(c,b)-mean(q_DRL)) - (template_constant_b-mean(q_constant))`, where the constant template is its 19-cycle ensemble phase mean. Constant and DRL cycles are never paired. The observable is the same DRL cycle/bin three-channel effective-action fluctuation about the retained DRL physical action mean. CCD empirically centers both matrices and applies `A=P_c (W^(1/2) U_c)^T / (N sqrt(3))` with `Q=1`. Component-major flattening gives `M=2*80,200=160,400`, `N=19*10=190`.
+
+Sensitivities use cycles as blocks: DRL and constant-template leave-one-cycle-out; odd/even and prefix/suffix for each source; and independently recomputed retained-boundary 8-bin, 12-bin, and 10-bin half-origin phase fields and action curves persisted by the phase artifacts. Fourier resampling of the primary 10-bin products is forbidden. No cross-role cycle pairing is allowed.
+
+Results are descriptive phase-coherent co-variation. They do not establish paired counterfactual states, causality, mechanism, response time, independent uncertainty, or superiority to POD.
+
+Only rank-3 subspace stability is supported unless mode-wise evidence is separately added. Published individual modes and left action vectors are labeled as the primary SVD basis representation, not individually stable physical objects.
@@ -3,8 +3,7 @@ from .artifacts import load_role_artifact
 from .contracts import CONTRACT,ROLES,constant_mean_provenance,contract_snapshot,effective_action_mean,verify_decomposition
 from .phase import evaluate_gate,load_phase_compact,publish_phase_compact,recover_phase
 from .dynamic_increment import load_dynamic_increment,publish_dynamic_increment
-from .temporal_ccd import TemporalConfig,TemporalTransaction,decompose_temporal,load_temporal_input,load_temporal_result
-from .phase_domain_ccd import PhaseDomainTransaction,decompose_phase_domain,load_phase_domain_result
+from .cycle_template_ccd import CycleTemplateTransaction,decompose_cycle_template,load_cycle_template_result
 from .publication import load_dynamic_publication,publish_dynamic_figures
 from .orchestration import CampaignSchedule,orchestrate,role_command
-__all__=["CONTRACT","ROLES","CampaignSchedule","constant_mean_provenance","contract_snapshot","effective_action_mean","load_role_artifact","recover_phase","evaluate_gate","publish_phase_compact","load_phase_compact","publish_dynamic_increment","load_dynamic_increment","TemporalConfig","TemporalTransaction","load_temporal_input","decompose_temporal","load_temporal_result","PhaseDomainTransaction","decompose_phase_domain","load_phase_domain_result","publish_dynamic_figures","load_dynamic_publication","orchestrate","role_command","verify_decomposition"]
+__all__=["CONTRACT","ROLES","CampaignSchedule","constant_mean_provenance","contract_snapshot","effective_action_mean","load_role_artifact","recover_phase","evaluate_gate","publish_phase_compact","load_phase_compact","publish_dynamic_increment","load_dynamic_increment","CycleTemplateTransaction","decompose_cycle_template","load_cycle_template_result","publish_dynamic_figures","load_dynamic_publication","orchestrate","role_command","verify_decomposition"]
@@ -5,27 +5,29 @@ from .contracts import ROLES
 from .orchestration import CampaignSchedule,orchestrate
 from .runtime import execute_role
 from .phase import publish_phase_compact
-from .temporal_ccd import TemporalConfig, TemporalTransaction, decompose_temporal, load_temporal_input
-from .phase_domain_ccd import PhaseDomainTransaction,decompose_phase_domain
+from .cycle_template_ccd import CycleTemplateTransaction,decompose_cycle_template
 from .publication import publish_dynamic_figures
+from .domain_comparison import publish_result as publish_domain_result, publish_figures as publish_domain_figures
 def main(argv=None):
 p=argparse.ArgumentParser(); s=p.add_subparsers(dest="command",required=True); common=argparse.ArgumentParser(add_help=False); common.add_argument("--campaign-id",required=True); common.add_argument("--warmup-intervals",type=int,default=480); common.add_argument("--collect-boundaries",type=int,default=360); common.add_argument("--launch-delay-seconds",type=float,default=120)
 role=s.add_parser("role",parents=[common]); role.add_argument("--role",choices=ROLES,required=True); role.add_argument("--output",type=Path,required=True); role.add_argument("--drl-artifact",type=Path); role.add_argument("--phase-artifact",type=Path); role.add_argument("--smoke",action="store_true")
 orch=s.add_parser("orchestrate",parents=[common]); orch.add_argument("--root",type=Path,required=True); orch.add_argument("--execute",action="store_true"); orch.add_argument("--smoke",action="store_true")
 phase=s.add_parser("phase-gate"); phase.add_argument("--role",choices=ROLES,default="drl"); phase.add_argument("--role-artifact",type=Path); phase.add_argument("--drl-artifact",type=Path); phase.add_argument("--output",type=Path,required=True)
- temporal=s.add_parser("temporal-ccd"); temporal.add_argument("--drl-artifact",type=Path,required=True); temporal.add_argument("--phase-artifact",type=Path,required=True); temporal.add_argument("--output",type=Path,required=True); temporal.add_argument("--chunk-size",type=int,default=8); temporal.add_argument("--ram-budget-bytes",type=int,required=True)
- phase_ccd=s.add_parser("phase-domain-ccd"); phase_ccd.add_argument("--drl-phase",type=Path,required=True); phase_ccd.add_argument("--constant-mean-phase",type=Path,required=True); phase_ccd.add_argument("--dynamic-increment",type=Path,required=True); phase_ccd.add_argument("--output",type=Path,required=True)
- publication=s.add_parser("publication"); publication.add_argument("--dynamic-increment",type=Path,required=True); publication.add_argument("--temporal-ccd",type=Path,required=True); publication.add_argument("--phase-domain-ccd",type=Path,required=True); publication.add_argument("--output",type=Path,required=True)
+ ccd=s.add_parser("cycle-template-ccd"); ccd.add_argument("--drl-phase",type=Path,required=True); ccd.add_argument("--constant-mean-phase",type=Path,required=True); ccd.add_argument("--roi-mean",type=Path,required=True); ccd.add_argument("--output",type=Path,required=True)
+ publication=s.add_parser("publication"); publication.add_argument("--roi-mean",type=Path,required=True); publication.add_argument("--cycle-template-ccd",type=Path,required=True); publication.add_argument("--output",type=Path,required=True)
+ domain=s.add_parser("domain-comparison"); domain.add_argument("--drl-phase",type=Path,required=True); domain.add_argument("--constant-mean-phase",type=Path,required=True); domain.add_argument("--roi-mean",type=Path,required=True); domain.add_argument("--wake-ccd",type=Path,required=True); domain.add_argument("--output",type=Path,required=True)
+ domain_fig=s.add_parser("domain-publication"); domain_fig.add_argument("--result",type=Path,required=True); domain_fig.add_argument("--output",type=Path,required=True)
 a=p.parse_args(argv)
+ if a.command=="domain-publication":
+  published=publish_domain_figures(a.result,a.output); print(json.dumps({"result":str(published.resolve())},sort_keys=True)); return 0
+ if a.command=="domain-comparison":
+  parents={"drl_phase":a.drl_phase,"constant_mean_phase":a.constant_mean_phase,"roi_mean":a.roi_mean,"wake_ccd":a.wake_ccd}
+  published=publish_domain_result(parents,a.output); print(json.dumps({"result":str(published.resolve())},sort_keys=True)); return 0
 if a.command=="publication":
-  published=publish_dynamic_figures(a.dynamic_increment,a.temporal_ccd,a.phase_domain_ccd,a.output); print(json.dumps({"result":str(published.resolve())},sort_keys=True)); return 0
- if a.command=="phase-domain-ccd":
-  result=decompose_phase_domain(a.drl_phase,a.constant_mean_phase,a.dynamic_increment)
-  with PhaseDomainTransaction(a.output) as tx: tx.write(result); published=tx.publish()
-  print(json.dumps({"result":str(published.resolve()),"summary":result.summary},sort_keys=True)); return 0
- if a.command=="temporal-ccd":
-  inp=load_temporal_input(a.drl_artifact,a.phase_artifact); result=decompose_temporal(inp,streaming_config=TemporalConfig(a.chunk_size,a.ram_budget_bytes));
-  with TemporalTransaction(a.output) as tx: tx.write(result); published=tx.publish()
+  published=publish_dynamic_figures(a.roi_mean,a.cycle_template_ccd,a.output); print(json.dumps({"result":str(published.resolve())},sort_keys=True)); return 0
+ if a.command=="cycle-template-ccd":
+  result=decompose_cycle_template(a.drl_phase,a.constant_mean_phase,a.roi_mean)
+  with CycleTemplateTransaction(a.output) as tx: tx.write(result); published=tx.publish()
  print(json.dumps({"result":str(published.resolve()),"summary":result.summary},sort_keys=True)); return 0
 if a.command=="phase-gate":
  source=a.role_artifact or a.drl_artifact
@@ -0,0 +1,142 @@
+"""Weighted Q=1 CCD of DRL cycle/bin residuals against a constant phase template."""
+from __future__ import annotations
+from dataclasses import dataclass
+import json, shutil, tempfile
+from pathlib import Path
+from typing import Any
+import numpy as np
+from CCD_analysis.acquisition.artifacts import file_sha256, rename_noreplace
+from CCD_analysis.direct_dq.schema import canonical_array_sha256
+from CCD_analysis.original_ccd import CCDConfig, decompose
+from .contracts import canonical_json
+from .dynamic_increment import load_dynamic_increment
+from .phase import load_phase_compact
+
+SCHEMA_ID="ccd-karman-cycle-template/v1"
+CHANNELS=("front","upper","lower")
+PRIMARY_BINS=10
+CLAIM_BOUNDARY="descriptive phase-coherent DRL-minus-constant-template co-variation; not paired counterfactual, causality, mechanism, response time, independent uncertainty, or CCD superiority"
+
+@dataclass(frozen=True)
+class CycleTemplateResult:
+ arrays: dict[str,np.ndarray]
+ config: dict[str,Any]
+ summary: dict[str,Any]
+ input_hashes: dict[str,Any]
+
+def _fit(u,p,w):
+ result=decompose(np.asarray(u,np.float64),np.asarray(p,np.float64),weight=np.asarray(w,np.float64),config=CCDConfig(center_snapshots=True,center_observables=True))
+ return {"cross":result.cross_correlation,"left":result.left_functions,"singular":result.singular_values,"modes":result.physical_modes,"coefficients":result.coefficients,"rank":result.identifiable_rank}
+
+def _ensemble_matrix(drl_ensemble,drl_actions,constant_ensemble,mean_drl,mean_constant,selector,weights,physical_action_mean):
+ fields=np.asarray(drl_ensemble,dtype=np.float64); actions=np.asarray(drl_actions,dtype=np.float64); template=np.asarray(constant_ensemble,dtype=np.float64)
+ if fields.ndim!=3 or fields.shape[1]!=2 or template.shape!=fields.shape or actions.shape!=(len(fields),3): raise ValueError("persisted alternative phase products malformed")
+ residual=(fields[:,:,selector]-mean_drl[None])-(template[:,:,selector]-mean_constant[None])
+ u=residual.reshape(len(fields),2*int(selector.sum())).T; physical_mean=np.asarray(physical_action_mean,dtype=np.float64)
+ if physical_mean.shape!=(3,) or not np.isfinite(physical_mean).all(): raise ValueError("exact DRL retained physical action mean must have shape (3,)")
+ p=(actions-physical_mean).T
+ return u,p,np.concatenate((weights,weights))
+
+def _projector_cosines(a,b,w):
+ r=min(a["rank"],b["rank"],3)
+ if not r: return []
+ return np.linalg.svd((a["modes"][:,:r]*np.sqrt(w)[:,None]).T@(b["modes"][:,:r]*np.sqrt(w)[:,None]),compute_uv=False).tolist()
+
+def _comparison(primary,other,w):
+ c=_projector_cosines(primary,other,w)
+ return {"rank":other["rank"],"leading_spectrum":other["singular"][:3].tolist(),"projector_principal_cosines":c,"minimum_projector_cosine":min(c) if c else None}
+
+def _matrix(drl_fields,drl_actions,constant_fields,mean_drl,mean_constant,selector,weights,physical_action_mean,*,bins=10,origin=0.0,drl_cycles=None,constant_cycles=None):
+ dc=np.arange(len(drl_fields)) if drl_cycles is None else np.asarray(drl_cycles); cc=np.arange(len(constant_fields)) if constant_cycles is None else np.asarray(constant_cycles)
+ if not len(dc) or not len(cc): raise ValueError("cycle subset is empty")
+ template=np.mean(constant_fields[cc],axis=0,dtype=np.float64)
+ actions=drl_actions[dc].astype(np.float64)
+ fields=drl_fields[dc].astype(np.float64)
+ if bins!=10 or origin: raise ValueError("primary cycle matrix does not resample; use persisted alternatives")
+ residual=(fields[:,:, :,selector]-mean_drl[None,None])-(template[None,:, :,selector]-mean_constant[None,None])
+ u=residual.transpose(0,1,2,3).reshape(len(dc)*bins,2*int(selector.sum())).T
+ physical_mean=np.asarray(physical_action_mean,dtype=np.float64)
+ if physical_mean.shape!=(3,) or not np.isfinite(physical_mean).all(): raise ValueError("exact DRL retained physical action mean must have shape (3,)")
+ p=(actions-physical_mean).reshape(len(dc)*bins,3).T
+ return u,p,np.concatenate((weights,weights)),dc,np.tile(np.arange(bins),len(dc)),physical_mean
+
+def _load(drl_phase,constant_phase,roi_mean):
+ dp,cp=load_phase_compact(drl_phase),load_phase_compact(constant_phase); dm=load_dynamic_increment(roi_mean)
+ if not dp["summary"]["gate_passed"] or dp["summary"].get("role")!="drl" or not cp["summary"]["gate_passed"] or cp["summary"].get("role")!="constant_mean": raise ValueError("passing DRL and constant_mean phase artifacts required")
+ with np.load(Path(drl_phase)/"compact.npz",allow_pickle=False) as dz,np.load(Path(constant_phase)/"compact.npz",allow_pickle=False) as cz,np.load(Path(roi_mean)/"arrays.npz",allow_pickle=False) as mz:
+  if not np.array_equal(dz["fluid_mask"],cz["fluid_mask"]): raise ValueError("phase masks differ")
+  roi=mz["roi_fluid_mask"]; selector=roi[dz["fluid_mask"]]
+  if int(selector.sum())!=len(mz["roi_quadrature_weights"]): raise ValueError("ROI mask/weights mismatch")
+  provenance=dp["summary"].get("constant_mean_provenance")
+  if not isinstance(provenance,dict) or provenance.get("source")!="fresh DRL retained effective_applied_action mean": raise ValueError("fresh DRL physical action provenance required")
+  physical_mean=np.asarray(provenance.get("constant_mean_physical_action"),dtype=np.float64)
+  dense_mean=np.asarray(dm["summary"]["statistics"]["drl"]["effective_action_mean"],dtype=np.float64)
+  if physical_mean.shape!=(3,) or not np.array_equal(physical_mean.astype(np.float32),dense_mean.astype(np.float32)): raise ValueError("phase physical mean does not match fresh retained DRL mean")
+  alternatives={}
+  for name,bins in (("bins8",8),("bins12",12),("bins10_half_shift",10)):
+   field_key=f"ensemble_{name}"; action_key=f"ensemble_{name}_effective_actions"; count_key=f"counts_{name}"
+   per_field_key=f"ensemble_{name}_per_cycle_fields"; per_action_key=f"ensemble_{name}_per_cycle_effective_actions"
+   if field_key not in dz.files or field_key not in cz.files or action_key not in dz.files or per_field_key not in dz.files or per_field_key not in cz.files or per_action_key not in dz.files or per_action_key not in cz.files or count_key not in dz.files or count_key not in cz.files: raise ValueError(f"required persisted alternative missing: {name}")
+   if dz[field_key].shape!=(bins,2,int(dz["fluid_mask"].sum())) or cz[field_key].shape!=dz[field_key].shape or dz[per_field_key].shape!=(len(dz["cycle_ids"]),bins,2,int(dz["fluid_mask"].sum())) or cz[per_field_key].shape!=dz[per_field_key].shape or dz[action_key].shape!=(bins,3) or dz[per_action_key].shape!=(len(dz["cycle_ids"]),bins,3) or cz[per_action_key].shape!=dz[per_action_key].shape or dz[count_key].shape!=(len(dz["cycle_ids"]),bins) or cz[count_key].shape!=(len(cz["cycle_ids"]),bins) or np.any(dz[count_key]<=0) or np.any(cz[count_key]<=0): raise ValueError(f"required persisted alternative malformed: {name}")
+   if not np.allclose(dz[field_key],dz[per_field_key].mean(0),rtol=1e-6,atol=1e-6) or not np.allclose(cz[field_key],cz[per_field_key].mean(0),rtol=1e-6,atol=1e-6) or not np.allclose(dz[action_key],dz[per_action_key].mean(0),rtol=1e-6,atol=1e-6): raise ValueError(f"persisted alternative ensemble mismatch: {name}")
+   alternatives[name]={"drl_fields":dz[per_field_key].copy(),"drl_actions":dz[per_action_key].copy(),"constant_fields":cz[per_field_key].copy(),"drl_counts":dz[count_key].copy(),"constant_counts":cz[count_key].copy()}
+  return dp,cp,dm,physical_mean,dz["cycle_bin_fields"].copy(),dz["cycle_bin_effective_actions"].copy(),dz["cycle_bin_counts"].copy(),dz["cycle_ids"].copy(),cz["cycle_bin_fields"].copy(),cz["cycle_bin_counts"].copy(),cz["cycle_ids"].copy(),mz["mean_drl"][:,mz["roi_selector_on_common_mask"]].copy(),mz["mean_constant_mean"][:,mz["roi_selector_on_common_mask"]].copy(),selector,mz["roi_quadrature_weights"].copy(),alternatives
+
+def decompose_cycle_template(drl_phase,constant_phase,roi_mean)->CycleTemplateResult:
+ dp,cp,dm,physical_mean,df,da,dcounts,dcycle_labels,cf,ccounts,ccycle_labels,md,mc,selector,w0,alternatives=_load(drl_phase,constant_phase,roi_mean)
+ if len(df)!=19 or len(cf)!=19: raise ValueError("exactly 19 DRL and 19 constant cycles required")
+ u,p,w,cycle_ids,bin_ids,physical_mean=_matrix(df,da,cf,md,mc,selector,w0,physical_mean); primary=_fit(u,p,w)
+ if u.shape[1]!=190 or u.shape[0]!=160400: raise ValueError(f"primary dimensions must be M=160400,N=190; got {u.shape}")
+ n=len(df); h=n//2; sensitivities={"drl_leave_one_cycle_out":[],"constant_template_leave_one_cycle_out":[],"drl_splits":{},"constant_template_splits":{},"bin_and_origin":{}}
+ for i in range(n):
+  args=_matrix(df,da,cf,md,mc,selector,w0,physical_mean,drl_cycles=np.delete(np.arange(n),i)); sensitivities["drl_leave_one_cycle_out"].append({"omitted_cycle":i,**_comparison(primary,_fit(*args[:3]),w)})
+ for i in range(len(cf)):
+  args=_matrix(df,da,cf,md,mc,selector,w0,physical_mean,constant_cycles=np.delete(np.arange(len(cf)),i)); sensitivities["constant_template_leave_one_cycle_out"].append({"omitted_cycle":i,**_comparison(primary,_fit(*args[:3]),w)})
+ split={"odd":np.arange(n)[::2],"even":np.arange(n)[1::2],"prefix":np.arange(n)[:h],"suffix":np.arange(n)[-h:]}
+ for name,ids in split.items():
+  a=_matrix(df,da,cf,md,mc,selector,w0,physical_mean,drl_cycles=ids); sensitivities["drl_splits"][name]=_comparison(primary,_fit(*a[:3]),w)
+  a=_matrix(df,da,cf,md,mc,selector,w0,physical_mean,constant_cycles=ids); sensitivities["constant_template_splits"][name]=_comparison(primary,_fit(*a[:3]),w)
+ for name,bins in (("bins8",8),("bins12",12),("bins10_half_shift",10)):
+  alt=alternatives[name]; full_points=int(alt["drl_fields"].shape[-1]); a=_ensemble_matrix(alt["drl_fields"].reshape(19*bins,2,full_points),alt["drl_actions"].reshape(19*bins,3),np.tile(alt["constant_fields"].mean(0),(19,1,1,1)).reshape(19*bins,2,full_points),md,mc,selector,w0,physical_mean); comparison=_comparison(primary,_fit(*a),w); comparison.update({"source":"persisted independent retained-boundary per-cycle re-binning","sample_count":int(a[0].shape[1]),"expected_sample_count":19*bins,"drl_boundary_count_total":int(alt["drl_counts"].sum()),"constant_boundary_count_total":int(alt["constant_counts"].sum())}); sensitivities["bin_and_origin"][name]=comparison
+ arrays={"coordinate_weights":w,"primary_cross_correlation":primary["cross"],"primary_left_functions":primary["left"],"primary_singular_values":primary["singular"],"primary_physical_modes":primary["modes"],"primary_coefficients":primary["coefficients"],"drl_cycle_ids":dcycle_labels[cycle_ids],"phase_bin_ids":bin_ids,"drl_cycle_bin_boundary_counts":dcounts.reshape(-1),"constant_template_cycle_ids":ccycle_labels,"constant_template_cycle_bin_boundary_counts":ccounts,"drl_physical_action_mean_exact":physical_mean,"admitted_field_empirical_mean":u.mean(axis=1,dtype=np.float64),"admitted_observable_after_physical_mean_empirical_mean":p.mean(axis=1,dtype=np.float64),"admitted_observable_raw_empirical_mean":p.mean(axis=1,dtype=np.float64)+physical_mean}
+ config={"schema_id":SCHEMA_ID,"N":190,"M":160400,"Q":1,"field_estimand":"DRL cycle/bin separately centered minus constant ensemble phase template separately centered","observable_estimand":"same DRL cycle/bin effective-action fluctuation about exact fresh retained DRL physical mean, followed by CCD empirical centering","operator":"A=P_c(W^(1/2)U_c)^T/(190*sqrt(3))","component_order":"ux ROI points, then uy ROI points","cycle_pairing":False,"center_snapshots":True,"center_observables":True,"standardization":False,"whitening":False,"pod":False,"physical_action_mean_source":"phase summary constant_mean_provenance, exactly checked against fresh retained DRL dense statistic","alternative_phase_products":"persisted independent retained-boundary re-binnings; no Fourier resampling","column_bootstrap_api":False,"claim_boundary":CLAIM_BOUNDARY}
+ summary={"schema_id":SCHEMA_ID,"N":190,"M":160400,"Q":1,"drl_cycle_count":n,"constant_template_cycle_count":len(cf),"numerical_rank":primary["rank"],"primary_singular_values":primary["singular"].tolist(),"sensitivities":sensitivities,"spectrum_label":"cross-correlation strength; not field energy or explained variance","claim_boundary":CLAIM_BOUNDARY}
+ parents={"drl_phase":{"path":str(Path(drl_phase).resolve()),"manifest_sha256":file_sha256(Path(drl_phase)/"manifest.json")},"constant_mean_phase":{"path":str(Path(constant_phase).resolve()),"manifest_sha256":file_sha256(Path(constant_phase)/"manifest.json")},"roi_mean":{"path":str(Path(roi_mean).resolve()),"manifest_sha256":file_sha256(Path(roi_mean)/"manifest.json")}}
+ hashes={"parents":parents,"canonical_arrays":{k:canonical_array_sha256(v) for k,v in arrays.items()}}
+ validate_cycle_template(arrays,config,summary,hashes); return CycleTemplateResult(arrays,config,summary,hashes)
+
+def validate_cycle_template(arrays,config,summary,hashes):
+ if config.get("schema_id")!=SCHEMA_ID or summary.get("schema_id")!=SCHEMA_ID or config.get("N")!=190 or config.get("M")!=160400 or config.get("Q")!=1 or config.get("cycle_pairing") is not False: raise ValueError("cycle-template scientific contract contradicted")
+ if config.get("alternative_phase_products")!="persisted independent retained-boundary re-binnings; no Fourier resampling" or config.get("column_bootstrap_api") is not False or config.get("operator")!="A=P_c(W^(1/2)U_c)^T/(190*sqrt(3))" or any(config.get(k) is not False for k in ("standardization","whitening","pod")): raise ValueError("cycle-template operator contract contradicted")
+ d={k:np.asarray(v) for k,v in arrays.items()}
+ if d["primary_cross_correlation"].shape!=(3,160400) or d["primary_left_functions"].shape!=(3,3) or d["primary_coefficients"].shape[1]!=190: raise ValueError("cycle-template dimensions invalid")
+ if set(hashes)=={"parents","canonical_arrays"} and all(hashes["canonical_arrays"].get(k)==canonical_array_sha256(v) for k,v in d.items()): canonical_json(config); canonical_json(summary); canonical_json(hashes); return
+ raise ValueError("cycle-template hashes invalid")
+
+class CycleTemplateTransaction:
+ def __init__(self,destination): self.destination=Path(destination); self.stage=None
+ def __enter__(self):
+  if self.destination.exists(): raise FileExistsError(self.destination)
+  self.destination.parent.mkdir(parents=True,exist_ok=True); self.stage=Path(tempfile.mkdtemp(prefix=f".{self.destination.name}.partial-",dir=self.destination.parent)); return self
+ def write(self,result):
+  validate_cycle_template(result.arrays,result.config,result.summary,result.input_hashes); np.savez_compressed(self.stage/"arrays.npz",**result.arrays)
+  for name,value in (("config.json",result.config),("summary.json",result.summary),("input_hashes.json",result.input_hashes)): (self.stage/name).write_bytes(canonical_json(value))
+  files={q.name:file_sha256(q) for q in self.stage.iterdir()}; (self.stage/"manifest.json").write_bytes(canonical_json({"schema_id":SCHEMA_ID,"complete":True,"files":files}))
+ def publish(self):
+  load_cycle_template_result(self.stage,recompute=False); rename_noreplace(self.stage,self.destination); self.stage=None; load_cycle_template_result(self.destination,recompute=True); return self.destination
+ def __exit__(self,*args):
+  if self.stage is not None: shutil.rmtree(self.stage,ignore_errors=True)
+
+def load_cycle_template_result(path,recompute=True):
+ root=Path(path); manifest=json.loads((root/"manifest.json").read_text())
+ if manifest.get("schema_id")!=SCHEMA_ID or not manifest.get("complete") or set(manifest.get("files",{}))!={"arrays.npz","config.json","summary.json","input_hashes.json"}: raise ValueError("cycle-template manifest invalid")
+ for name,digest in manifest["files"].items():
+  if file_sha256(root/name)!=digest: raise ValueError("cycle-template file hash mismatch")
+ with np.load(root/"arrays.npz",allow_pickle=False) as z: arrays={k:z[k].copy() for k in z.files}
+ config=json.loads((root/"config.json").read_text()); summary=json.loads((root/"summary.json").read_text()); hashes=json.loads((root/"input_hashes.json").read_text()); validate_cycle_template(arrays,config,summary,hashes)
+ for parent in hashes["parents"].values():
+  if file_sha256(Path(parent["path"])/"manifest.json")!=parent["manifest_sha256"]: raise ValueError("cycle-template live parent identity changed")
+ if recompute:
+  p=hashes["parents"]; fresh=decompose_cycle_template(p["drl_phase"]["path"],p["constant_mean_phase"]["path"],p["roi_mean"]["path"])
+  for key in arrays: np.testing.assert_allclose(arrays[key],fresh.arrays[key],rtol=2e-11,atol=2e-12)
+ return {"arrays":arrays,"config":config,"summary":summary,"input_hashes":hashes,"manifest":manifest,"provenance_validation":"VERIFIED live parents and essential recomputation" if recompute else "VERIFIED_HASHES_AND_LIVE_PARENTS"}
@@ -0,0 +1,567 @@
+"""Pinball-inclusive, artifact-only CCD/POD comparison and publication."""
+from __future__ import annotations
+
+import json
+import shutil
+import tempfile
+from pathlib import Path
+
+import matplotlib
+matplotlib.use("Agg")
+import matplotlib.pyplot as plt
+import numpy as np
+
+from CCD_analysis.acquisition.artifacts import file_sha256, rename_noreplace
+from CCD_analysis.direct_dq.schema import canonical_array_sha256
+from .contracts import canonical_json
+from .cycle_template_ccd import _comparison, _ensemble_matrix, _fit, _matrix, load_cycle_template_result
+from .dynamic_increment import load_dynamic_increment
+from .mode_diagnostics import (
+    action_coordinates, bounded_streamwise_shift_correlation,
+    component_major_velocity, masked_vorticity, phase_harmonics,
+    reshape_cycle_phase, subspace_symmetry_diagnostics,
+    symmetry_diagnostics, x_localization,
+)
+from .phase import load_phase_compact
+from .pod_baseline import fit_weighted_pod, weighted_principal_cosines
+
+SCHEMA_ID = "ccd-karman-pinball-domain/v1"
+PUBLICATION_SCHEMA_ID = "ccd-karman-pinball-publication/v1"
+DOMAIN = {"x_min_D": 29.0, "x_max_D": 54.0, "abs_y_max_D": 5.0}
+REGIONS = (
+    ("geometry", 29.0, 32.0), ("near_wake", 32.0, 40.0),
+    ("middle_wake", 40.0, 47.0), ("far_wake", 47.0, 54.000001),
+)
+BODIES = (
+    ("front", 30.0, 0.0, 0.5),
+    ("upper", 31.3, 0.75, 0.5),
+    ("lower", 31.3, -0.75, 0.5),
+)
+CLAIMS = (
+    "artifact-only descriptive phase-coherent DRL-minus-constant-template "
+    "co-variation; no counterfactual, causality, mechanism, response time, "
+    "independent uncertainty, or CCD superiority"
+)
+REPRESENTATIVE_BINS = np.array([0, 2, 5, 7], dtype=int)
+
+
+def _jsonable(value):
+    if isinstance(value, np.ndarray):
+        return value.tolist()
+    if isinstance(value, np.generic):
+        return value.item()
+    if isinstance(value, dict):
+        return {key: _jsonable(item) for key, item in value.items()}
+    if isinstance(value, (list, tuple)):
+        return [_jsonable(item) for item in value]
+    return value
+
+
+def domain_mask(x, y, common_mask, bounds=DOMAIN):
+    """Build the frozen inclusive pinball-domain mask on the common fluid mask."""
+    x, y = np.asarray(x, float), np.asarray(y, float)
+    common = np.asarray(common_mask, bool)
+    if common.shape != (x.size, y.size) or np.any(np.diff(x) <= 0) or np.any(np.diff(y) <= 0):
+        raise ValueError("coordinates and common mask are inconsistent")
+    mask = common & (x[:, None] >= bounds["x_min_D"]) & (x[:, None] <= bounds["x_max_D"]) & (np.abs(y[None, :]) <= bounds["abs_y_max_D"])
+    if not mask.any():
+        raise ValueError("frozen domain has no common fluid points")
+    for name, cx, cy, radius in BODIES:
+        vicinity = (x[:, None] - cx) ** 2 + (y[None, :] - cy) ** 2 <= (radius + 0.15) ** 2
+        if not np.any(mask & vicinity):
+            raise ValueError(f"domain does not contain common fluid around {name} body")
+    return mask
+
+
+def _load_primary(drl_phase, constant_phase, roi_mean):
+    """Load only primary 10-bin products and full common-grid means/weights."""
+    drl_root, constant_root, mean_root = map(Path, (drl_phase, constant_phase, roi_mean))
+    dp, cp = load_phase_compact(drl_root), load_phase_compact(constant_root)
+    mean = load_dynamic_increment(mean_root)
+    if not dp["summary"]["gate_passed"] or dp["summary"].get("role") != "drl":
+        raise ValueError("passing DRL phase artifact required")
+    if not cp["summary"]["gate_passed"] or cp["summary"].get("role") != "constant_mean":
+        raise ValueError("passing constant_mean phase artifact required")
+    with np.load(drl_root / "compact.npz", allow_pickle=False) as dz, np.load(constant_root / "compact.npz", allow_pickle=False) as cz, np.load(mean_root / "arrays.npz", allow_pickle=False) as mz:
+        phase_mask = dz["fluid_mask"].copy()
+        if not np.array_equal(phase_mask, cz["fluid_mask"]) or not np.array_equal(phase_mask, mz["four_role_fluid_mask"]):
+            raise ValueError("phase and four-role common masks differ")
+        if not np.array_equal(dz["x_D"], mz["x_D"]) or not np.array_equal(dz["y_D"], mz["y_D"]):
+            raise ValueError("phase and mean coordinates differ")
+        x, y = mz["x_D"].copy(), mz["y_D"].copy()
+        mask = domain_mask(x, y, phase_mask)
+        selector = mask[phase_mask]
+        full_weights = mz["quadrature_weights"].copy()
+        if full_weights.shape != (int(phase_mask.sum()),) or np.any(full_weights <= 0):
+            raise ValueError("full common-mask quadrature weights invalid")
+        physical = np.asarray(dp["summary"]["constant_mean_provenance"]["constant_mean_physical_action"], float)
+        dense_mean = np.asarray(mean["summary"]["statistics"]["drl"]["effective_action_mean"], float)
+        if physical.shape != (3,) or not np.array_equal(physical.astype(np.float32), dense_mean.astype(np.float32)):
+            raise ValueError("physical action mean lineage mismatch")
+        loaded = {
+            "x": x, "y": y, "common_mask": phase_mask, "mask": mask,
+            "selector": selector, "point_weights": full_weights[selector],
+            "mean_drl": mz["mean_drl"][:, selector].copy(),
+            "mean_constant": mz["mean_constant_mean"][:, selector].copy(),
+            "drl_fields": dz["cycle_bin_fields"].copy(),
+            "drl_actions": dz["cycle_bin_effective_actions"].copy(),
+            "drl_counts": dz["cycle_bin_counts"].copy(),
+            "drl_cycle_ids": dz["cycle_ids"].copy(),
+            "constant_fields": cz["cycle_bin_fields"].copy(),
+            "constant_counts": cz["cycle_bin_counts"].copy(),
+            "constant_cycle_ids": cz["cycle_ids"].copy(),
+            "physical_mean": physical,
+        }
+    if loaded["drl_fields"].shape[:3] != (19, 10, 2) or loaded["constant_fields"].shape[:3] != (19, 10, 2):
+        raise ValueError("exact 19-cycle by 10-bin primary fields required")
+    return loaded, mean
+
+
+def _selected_alternative(root, name, selector):
+    bins = {"bins8": 8, "bins12": 12, "bins10_half_shift": 10}[name]
+    with np.load(Path(root) / "compact.npz", allow_pickle=False) as source:
+        fields = source[f"ensemble_{name}_per_cycle_fields"][:, :, :, selector].copy()
+        actions = source[f"ensemble_{name}_per_cycle_effective_actions"].copy()
+        counts = source[f"counts_{name}"].copy()
+    if fields.shape[:3] != (19, bins, 2) or actions.shape != (19, bins, 3) or counts.shape != (19, bins):
+        raise ValueError(f"persisted {name} products malformed")
+    return fields, actions, counts
+
+
+def _grid_weights(mask, point_weights):
+    result = np.zeros(mask.shape, float)
+    result[mask] = point_weights
+    return result
+
+
+def _region_selector(mask, x, low, high):
+    grid_x = np.broadcast_to(np.asarray(x)[:, None], mask.shape)
+    return mask & (grid_x >= low) & (grid_x < high)
+
+
+def _weighted_metrics(truth, reconstruction, weights, region_points):
+    """Return per-sample residual and cosine for one component-major region."""
+    points = weights.size // 2
+    point_select = np.asarray(region_points, bool)
+    dof_select = np.r_[point_select, point_select]
+    w = weights[dof_select]
+    a, b = truth[dof_select], reconstruction[dof_select]
+    numerator = np.sum(w[:, None] * (a - b) ** 2, axis=0)
+    denominator = np.sum(w[:, None] * a ** 2, axis=0)
+    residual = np.sqrt(np.divide(numerator, denominator, out=np.zeros_like(numerator), where=denominator > 0))
+    dot = np.sum(w[:, None] * a * b, axis=0)
+    norm = np.sqrt(np.sum(w[:, None] * a * a, axis=0) * np.sum(w[:, None] * b * b, axis=0))
+    cosine = np.divide(dot, norm, out=np.zeros_like(dot), where=norm > 0)
+    return residual, np.clip(cosine, -1, 1)
+
+
+def _reconstruction_diagnostics(centered, weights, mask, x, ccd, pod):
+    points = int(mask.sum())
+    region_point_selectors = {"full": np.ones(points, bool)}
+    for name, low, high in REGIONS:
+        region_point_selectors[name] = _region_selector(mask, x, low, high)[mask]
+    diagnostics = {}
+    representative = {"raw": centered[:, REPRESENTATIVE_BINS].copy()}
+    for method, modes, coefficients in (
+        ("ccd", ccd["modes"][:, :3], ccd["coefficients"][:3]),
+        ("pod", pod["modes"][:, :3], pod["coefficients"][:3]),
+    ):
+        diagnostics[method] = {}
+        for rank in (1, 2, 3):
+            reconstruction = modes[:, :rank] @ coefficients[:rank]
+            diagnostics[method][f"rank_{rank}"] = {}
+            for region, selector in region_point_selectors.items():
+                residual, cosine = _weighted_metrics(centered, reconstruction, weights, selector)
+                phase_residual = residual.reshape(19, 10).mean(axis=0)
+                phase_cosine = cosine.reshape(19, 10).mean(axis=0)
+                diagnostics[method][f"rank_{rank}"][region] = {
+                    "residual_mean": float(residual.mean()),
+                    "residual_cycle_std": float(residual.reshape(19, 10).mean(1).std()),
+                    "phase_residual_mean": phase_residual,
+                    "cosine_mean": float(cosine.mean()),
+                    "phase_cosine_mean": phase_cosine,
+                }
+            if rank in (2, 3):
+                phase_mean = reconstruction.reshape(2 * points, 19, 10).mean(axis=1)
+                representative[f"{method}_rank_{rank}"] = phase_mean[:, REPRESENTATIVE_BINS].copy()
+    raw_phase_mean = centered.reshape(2 * points, 19, 10).mean(axis=1)
+    representative["raw"] = raw_phase_mean[:, REPRESENTATIVE_BINS].copy()
+    representative["ccd_mode_1"] = (ccd["modes"][:, :1] @ ccd["coefficients"][:1]).reshape(2 * points, 19, 10).mean(1)[:, REPRESENTATIVE_BINS]
+    representative["ccd_mode_2"] = (ccd["modes"][:, 1:2] @ ccd["coefficients"][1:2]).reshape(2 * points, 19, 10).mean(1)[:, REPRESENTATIVE_BINS]
+    return diagnostics, representative
+
+
+def _mode_diagnostics(modes, mask, x, y, point_weights):
+    grid_weights = _grid_weights(mask, point_weights)
+    metrics, velocities, vorticities, valid_masks = [], [], [], []
+    for index in range(3):
+        velocity = component_major_velocity(modes[:, index], mask)
+        omega, valid = masked_vorticity(*velocity, x, y, mask)
+        velocities.append(velocity); vorticities.append(omega); valid_masks.append(valid)
+        metrics.append({
+            "mode": index + 1,
+            "symmetry": symmetry_diagnostics(velocity, grid_weights, mask, y),
+            "velocity_localization": x_localization(np.sum(velocity ** 2, axis=0), x, grid_weights, mask, REGIONS),
+            "vorticity_localization": x_localization(omega ** 2, x, grid_weights, valid, REGIONS),
+        })
+    shifts = {}
+    for first, second in ((0, 1), (0, 2), (1, 2)):
+        key = f"mode_{first + 1}_{second + 1}"
+        shifts[f"{key}_velocity"] = bounded_streamwise_shift_correlation(velocities[first], velocities[second], x, grid_weights, mask, 0.05)
+        shifts[f"{key}_vorticity"] = bounded_streamwise_shift_correlation(vorticities[first], vorticities[second], x, grid_weights, valid_masks[first] & valid_masks[second], 0.05)
+    return {
+        "metrics": metrics,
+        "subspace_symmetry": subspace_symmetry_diagnostics(modes[:, :3], grid_weights, mask, y),
+        "shifts": {key: {name: value[name] for name in ("best_shift_x_D", "correlation", "paired_point_count", "index_offset")} for key, value in shifts.items()},
+        "velocities": np.asarray(velocities), "vorticities": np.asarray(vorticities),
+        "vorticity_valid": np.asarray(valid_masks),
+    }
+
+
+def compute(parents):
+    """Compute the frozen pinball-domain result and all declared sensitivities."""
+    roots = {name: Path(path) for name, path in parents.items()}
+    required = {"drl_phase", "constant_mean_phase", "roi_mean", "wake_ccd"}
+    if set(roots) != required:
+        raise ValueError("all four declared parents are required")
+    loaded, mean_result = _load_primary(roots["drl_phase"], roots["constant_mean_phase"], roots["roi_mean"])
+    u, p, weights, cycle_ids, bin_ids, physical_mean = _matrix(
+        loaded["drl_fields"], loaded["drl_actions"], loaded["constant_fields"],
+        loaded["mean_drl"], loaded["mean_constant"], loaded["selector"],
+        loaded["point_weights"], loaded["physical_mean"],
+    )
+    empirical_mean = u.mean(axis=1, dtype=np.float64)
+    centered = u - empirical_mean[:, None]
+    ccd = _fit(u, p, weights)
+    # Exact method-of-snapshots POD: diagonalize the 190x190 weighted Gram matrix
+    # instead of asking LAPACK for an SVD of the 160400x190 field matrix.
+    gram = centered.T @ (weights[:, None] * centered)
+    gram = 0.5 * (gram + gram.T)
+    eigenvalues, right = np.linalg.eigh(gram)
+    order = np.argsort(eigenvalues)[::-1]
+    eigenvalues, right = np.maximum(eigenvalues[order], 0.0), right[:, order]
+    singular = np.sqrt(eigenvalues)
+    pod_modes = centered @ right[:, :3] / np.maximum(singular[:3], np.finfo(float).tiny)
+    pod_coefficients = singular[:3, None] * right[:, :3].T
+    pod_energy = singular ** 2
+    pod_total = float(pod_energy.sum())
+    pod_retained = np.r_[0.0, np.cumsum(pod_energy[:3])]
+    pod = {"modes": pod_modes, "coefficients": pod_coefficients, "singular_values": singular, "weighted_residual_fraction": np.maximum(1 - pod_retained / pod_total, 0), "label": "raw full-fit weighted POD field-energy spectrum; descriptive only"}
+    principal = weighted_principal_cosines(ccd["modes"], pod["modes"], weights)
+
+    wake = load_cycle_template_result(roots["wake_ccd"], recompute=False)
+    wake_mask = (loaded["common_mask"] & (loaded["x"][:, None] >= 34) & (loaded["x"][:, None] <= 54) & (np.abs(loaded["y"][None, :]) <= 5))
+    wake_selector = wake_mask[loaded["common_mask"]]
+    wake_u, wake_p, wake_weights, *_ = _matrix(
+        loaded["drl_fields"], loaded["drl_actions"], loaded["constant_fields"],
+        np.load(roots["roi_mean"] / "arrays.npz")["mean_drl"][:, wake_selector],
+        np.load(roots["roi_mean"] / "arrays.npz")["mean_constant_mean"][:, wake_selector],
+        wake_selector, np.load(roots["roi_mean"] / "arrays.npz")["quadrature_weights"][wake_selector],
+        loaded["physical_mean"],
+    )
+    wake_fit = _fit(wake_u, wake_p, wake_weights)
+    np.testing.assert_allclose(wake_fit["singular"][:3], wake["arrays"]["primary_singular_values"][:3], rtol=2e-11, atol=2e-12)
+    wake_cosines = np.linalg.svd((wake_fit["modes"][:, :3] * np.sqrt(wake_weights)[:, None]).T @ (wake["arrays"]["primary_physical_modes"][:, :3] * np.sqrt(wake_weights)[:, None]), compute_uv=False)
+    if min(wake_cosines) < 1 - 1e-10 or not np.array_equal(wake_u.mean(axis=1, dtype=np.float64), wake["arrays"]["admitted_field_empirical_mean"]):
+        raise ValueError("generic selector failed canonical wake-domain regression")
+
+    def fast_subset_fit(field_matrix, observable_matrix, ids, corrections=None):
+        ids = np.asarray(ids, dtype=int)
+        uu = field_matrix[:, ids]
+        pp = observable_matrix[:, ids]
+        if corrections is not None:
+            uu = uu + corrections[:, ids]
+        uc = uu - uu.mean(axis=1, keepdims=True)
+        pc = pp - pp.mean(axis=1, keepdims=True)
+        cross = np.zeros((pc.shape[0], uc.shape[0]), dtype=np.float64)
+        # One 3-vector outer product per admitted sample avoids repeated MxN products.
+        for column in range(uc.shape[1]):
+            cross += np.outer(pc[:, column], weights * uc[:, column])
+        cross /= uc.shape[1] * np.sqrt(3)
+        left, singular, vh = np.linalg.svd(cross, full_matrices=False)
+        modes = vh.T / np.sqrt(weights)[:, None]
+        return {"rank": int(singular.size), "singular": singular, "modes": modes, "left": left}
+
+    def fast_comparison(fit):
+        return _comparison(ccd, fit, weights)
+
+    n = 19
+    sensitivities = {"drl_leave_one_cycle_out": [], "constant_template_leave_one_cycle_out": [], "drl_splits": {}, "constant_template_splits": {}, "bin_and_origin": {}}
+    sample_ids = np.arange(190)
+    for index in range(n):
+        keep = np.delete(np.arange(n), index)
+        ids = np.concatenate([10 * keep + bin for bin in range(10)])
+        # _matrix uses cycle-major, phase-minor ordering.
+        ids = np.sort(ids)
+        sensitivities["drl_leave_one_cycle_out"].append({"omitted_cycle": index, **fast_comparison(fast_subset_fit(u, p, ids))})
+        template_full = loaded["constant_fields"].mean(axis=0)
+        template_loo = loaded["constant_fields"][np.delete(np.arange(n), index)].mean(axis=0)
+        delta = -(template_loo - template_full)[:, :, loaded["selector"]].transpose(1, 0, 2).reshape(2 * int(loaded["selector"].sum()), 10)
+        correction = np.tile(delta, (1, n))
+        sensitivities["constant_template_leave_one_cycle_out"].append({"omitted_cycle": index, **fast_comparison(fast_subset_fit(u, p, sample_ids, correction))})
+    splits = {"odd": np.arange(n)[::2], "even": np.arange(n)[1::2], "prefix": np.arange(n)[: n // 2], "suffix": np.arange(n)[-(n // 2):]}
+    for name, cycles in splits.items():
+        ids = np.sort(np.concatenate([10 * cycles + bin for bin in range(10)]))
+        sensitivities["drl_splits"][name] = fast_comparison(fast_subset_fit(u, p, ids))
+        template_full = loaded["constant_fields"].mean(axis=0)
+        template_split = loaded["constant_fields"][cycles].mean(axis=0)
+        delta = -(template_split - template_full)[:, :, loaded["selector"]].transpose(1, 0, 2).reshape(2 * int(loaded["selector"].sum()), 10)
+        correction = np.tile(delta, (1, n))
+        sensitivities["constant_template_splits"][name] = fast_comparison(fast_subset_fit(u, p, ids, correction))
+    for name, bins in (("bins8", 8), ("bins12", 12), ("bins10_half_shift", 10)):
+        df, da, dc = _selected_alternative(roots["drl_phase"], name, loaded["selector"])
+        cf, _ca, cc = _selected_alternative(roots["constant_mean_phase"], name, loaded["selector"])
+        template = np.tile(cf.mean(0), (19, 1, 1, 1)).reshape(19 * bins, 2, -1)
+        args = _ensemble_matrix(df.reshape(19 * bins, 2, -1), da.reshape(19 * bins, 3), template, loaded["mean_drl"], loaded["mean_constant"], np.ones(int(loaded["selector"].sum()), bool), loaded["point_weights"], physical_mean)
+        comparison = _comparison(ccd, _fit(*args), weights)
+        comparison.update(sample_count=19 * bins, drl_boundary_count_total=int(dc.sum()), constant_boundary_count_total=int(cc.sum()))
+        sensitivities["bin_and_origin"][name] = comparison
+
+    def all_comparisons(value):
+        if isinstance(value, list):
+            return [item for entry in value for item in all_comparisons(entry)]
+        if isinstance(value, dict):
+            if "minimum_projector_cosine" in value:
+                return [value]
+            return [item for entry in value.values() for item in all_comparisons(entry)]
+        return []
+    minimum_sensitivity = min(item["minimum_projector_cosine"] for item in all_comparisons(sensitivities))
+
+    ccd_diag = _mode_diagnostics(ccd["modes"][:, :3], loaded["mask"], loaded["x"], loaded["y"], loaded["point_weights"])
+    pod_diag = _mode_diagnostics(pod["modes"][:, :3], loaded["mask"], loaded["x"], loaded["y"], loaded["point_weights"])
+    reconstruction, representative = _reconstruction_diagnostics(centered, weights, loaded["mask"], loaded["x"], ccd, pod)
+    actions_native = loaded["drl_actions"].reshape(190, 3).T
+    left_physical = action_coordinates(ccd["left"][:, :3])
+    harmonics = phase_harmonics(reshape_cycle_phase(ccd["coefficients"][:3]), 5)
+    region_counts = {name: int(_region_selector(loaded["mask"], loaded["x"], low, high).sum()) for name, low, high in REGIONS}
+    minimum_sensitivity = min(item["minimum_projector_cosine"] for category in sensitivities.values() for item in (category.values() if isinstance(category, dict) else category))
+
+    arrays = {
+        "x_D": loaded["x"], "y_D": loaded["y"], "domain_mask": loaded["mask"],
+        "point_weights": loaded["point_weights"], "ccd_modes": ccd["modes"][:, :3],
+        "ccd_coefficients": ccd["coefficients"][:3], "ccd_left_functions": ccd["left"][:, :3],
+        "left_action_physical_coordinates": left_physical,
+        "ccd_singular_values": ccd["singular"][:3], "pod_modes_raw": pod["modes"][:, :3],
+        "pod_coefficients": pod["coefficients"][:3], "pod_field_energy_singular_values": pod["singular_values"],
+        "pod_weighted_residual_fraction": pod["weighted_residual_fraction"],
+        "ccd_pod_principal_cosines": principal, "empirical_field_mean": empirical_mean,
+        "representative_phase_bins": REPRESENTATIVE_BINS,
+        "ccd_mode_velocity": ccd_diag["velocities"], "ccd_mode_vorticity": ccd_diag["vorticities"],
+        "ccd_mode_rotation_proxy": np.asarray([_rotation_proxy(value, loaded["x"], loaded["y"], loaded["mask"]) for value in ccd_diag["velocities"]]),
+        "ccd_vorticity_valid": ccd_diag["vorticity_valid"],
+        "pod_mode_velocity": pod_diag["velocities"], "pod_mode_vorticity": pod_diag["vorticities"],
+        "pod_vorticity_valid": pod_diag["vorticity_valid"],
+        "actions_native": actions_native,
+    }
+    arrays.update({f"representative_{key}": value for key, value in representative.items()})
+    for key, value in harmonics.items():
+        arrays[f"ccd_harmonic_{key}"] = value
+    summary = {
+        "schema_id": SCHEMA_ID, "M": int(u.shape[0]), "N": 190, "Q": 1, "rank": ccd["rank"],
+        "domain_fluid_point_count": int(loaded["mask"].sum()), "region_fluid_point_counts": region_counts,
+        "primary_singular_values": ccd["singular"][:3], "ccd_pod_principal_cosines": principal,
+        "pod_weighted_residual_fraction_rank_0_to_3": pod["weighted_residual_fraction"],
+        "minimum_declared_sensitivity_projector_cosine": minimum_sensitivity,
+        "sensitivities": sensitivities, "ccd_mode_diagnostics": {"metrics": ccd_diag["metrics"], "subspace_symmetry": ccd_diag["subspace_symmetry"], "shifts": ccd_diag["shifts"]},
+        "pod_mode_diagnostics": {"metrics": pod_diag["metrics"], "subspace_symmetry": pod_diag["subspace_symmetry"], "shifts": pod_diag["shifts"]},
+        "reconstruction": reconstruction,
+        "wake_domain_regression": {"M": int(wake_u.shape[0]), "singular_values": wake_fit["singular"][:3], "subspace_cosines": wake_cosines, "admitted_mean_bit_equal": True},
+        "raw_U_canonical_sha256": canonical_array_sha256(u), "centered_U_canonical_sha256": canonical_array_sha256(centered),
+        "contribution_context": {"constant_mean_share": 0.9342242558837475, "dynamic_share": 0.06577574411625256},
+        "decision_rule": "if all CCD/POD principal cosines >=0.95, close as POD-equivalent; otherwise only stable distinct directions are observable-selected candidates",
+        "claim_boundary": CLAIMS,
+    }
+    config = {
+        "schema_id": SCHEMA_ID, "domain": DOMAIN, "regions": REGIONS, "bodies_front_upper_lower": BODIES,
+        "sensor_x_D": 40.0, "component_order": "ux selected points, then uy selected points",
+        "field_estimand": "DRL cycle/bin separately centered minus constant ensemble phase template separately centered",
+        "operator": "A=P_c(W^(1/2)U_c)^T/(190*sqrt(3))", "phase_shape": [19, 10],
+        "representative_phase_bins": REPRESENTATIVE_BINS, "pod": "raw full-fit weighted rank-3 field-energy POD",
+        "standardization": False, "whitening": False, "cycle_pairing": False,
+        "execution": "pinball_math CPU-only; no CFD/CUDA/SR", "claim_boundary": CLAIMS,
+    }
+    hashes = {name: {"path": str(root.resolve()), "manifest_sha256": file_sha256(root / "manifest.json")} for name, root in roots.items()}
+    return arrays, _jsonable(config), _jsonable(summary), hashes
+
+
+def _validate_files(root, schema):
+    manifest = json.loads((root / "manifest.json").read_text())
+    if manifest.get("schema_id") != schema or not manifest.get("complete"):
+        raise ValueError("domain artifact manifest invalid")
+    for name, digest in manifest.get("files", {}).items():
+        if file_sha256(root / name) != digest:
+            raise ValueError("domain artifact file hash mismatch")
+    return manifest
+
+
+def publish_result(parents, destination):
+    destination = Path(destination)
+    if destination.exists():
+        raise FileExistsError(destination)
+    destination.parent.mkdir(parents=True, exist_ok=True)
+    stage = Path(tempfile.mkdtemp(prefix=f".{destination.name}.partial-", dir=destination.parent))
+    try:
+        arrays, config, summary, hashes = compute(parents)
+        np.savez_compressed(stage / "arrays.npz", **arrays)
+        for name, value in (("config.json", config), ("summary.json", summary), ("input_hashes.json", hashes)):
+            (stage / name).write_bytes(canonical_json(value))
+        (stage / "RESULTS.md").write_text(_results_text(summary))
+        files = {path.name: file_sha256(path) for path in stage.iterdir()}
+        (stage / "manifest.json").write_bytes(canonical_json({"schema_id": SCHEMA_ID, "complete": True, "files": files}))
+        load_result(stage, recompute=False)
+        rename_noreplace(stage, destination); stage = None
+        return destination
+    finally:
+        if stage is not None:
+            shutil.rmtree(stage, ignore_errors=True)
+
+
+def load_result(path, recompute=False):
+    root = Path(path); manifest = _validate_files(root, SCHEMA_ID)
+    expected = {"arrays.npz", "config.json", "summary.json", "input_hashes.json", "RESULTS.md"}
+    if set(manifest["files"]) != expected:
+        raise ValueError("domain result inventory invalid")
+    hashes = json.loads((root / "input_hashes.json").read_text())
+    summary = json.loads((root / "summary.json").read_text())
+    for value in hashes.values():
+        if file_sha256(Path(value["path"]) / "manifest.json") != value["manifest_sha256"]:
+            raise ValueError("domain live parent identity changed")
+    if recompute:
+        fresh = compute({name: value["path"] for name, value in hashes.items()})[2]
+        for key in ("raw_U_canonical_sha256", "centered_U_canonical_sha256"):
+            if fresh[key] != summary[key]:
+                raise ValueError(f"fresh {key} changed")
+        np.testing.assert_allclose(fresh["ccd_pod_principal_cosines"], summary["ccd_pod_principal_cosines"], rtol=1e-10, atol=1e-12)
+    return {"path": root, "summary": summary, "manifest": manifest}
+
+
+def _results_text(summary):
+    principal = summary["ccd_pod_principal_cosines"]
+    decision = "POD-equivalent rank-3 residual subspace" if min(principal) >= 0.95 else "only non-overlapping stable directions remain observable-selected candidates"
+    return (
+        "# Pinball-inclusive corrected Karman CCD/POD\n\n"
+        f"The frozen domain is `29 <= x/D <= 54, |y/D| <= 5` on {summary['domain_fluid_point_count']} common fluid points. "
+        f"The exact residual has `M={summary['M']}, N=190, Q=1`.\n\n"
+        f"CCD/POD principal cosines are `{principal}`; the preregistered decision is: **{decision}**. "
+        f"The least declared CCD sensitivity projector cosine is `{summary['minimum_declared_sensitivity_projector_cosine']}`.\n\n"
+        "This is an artifact-only descriptive extension. It does not establish actuator causality, mechanism, response time, independent uncertainty, or CCD superiority.\n"
+    )
+
+
+def _field(values, mask):
+    result = np.full(mask.shape, np.nan); result[mask] = values; return result
+
+
+
+def _rotation_proxy(velocity, x, y, mask):
+    """Return a signed body-local angular-velocity proxy for visualization only."""
+    ux, uy = velocity
+    proxy = np.zeros(mask.shape, float)
+    total = np.zeros(mask.shape, float)
+    for _name, cx, cy, radius in BODIES:
+        distance2 = (x[:, None] - cx) ** 2 + (y[None, :] - cy) ** 2
+        weight = np.exp(-distance2 / (2.0 * (1.25 * radius) ** 2))
+        proxy += weight * ((x[:, None] - cx) * uy - (y[None, :] - cy) * ux)
+        total += weight
+    return np.divide(proxy, total, out=np.zeros_like(proxy), where=total > 1e-15) * mask
+
+def _panel(ax, x, y, field, title, limit=None, cmap="RdBu_r"):
+    opts = {} if limit is None else {"vmin": -limit, "vmax": limit}
+    image = ax.pcolormesh(x, y, field.T, shading="nearest", cmap=cmap, rasterized=True, **opts)
+    for _name, cx, cy, radius in BODIES:
+        ax.add_patch(plt.Circle((cx, cy), radius, facecolor="none", edgecolor="black", lw=.7))
+    ax.axvline(40, color="black", ls=":", lw=.6)
+    ax.set(xlim=(29, 54), ylim=(-5, 5), xlabel="x/D", ylabel="y/D", title=title); ax.set_aspect("equal")
+    return image
+
+
+def _save(fig, stage, stem):
+    for ext in ("png", "pdf"):
+        fig.savefig(stage / f"{stem}.{ext}", dpi=250 if ext == "png" else None, bbox_inches="tight", metadata={"Creator": "CCD pinball-domain publication"})
+    plt.close(fig)
+
+
+def publish_figures(result_root, output):
+    result_root, output = Path(result_root), Path(output)
+    if output.exists():
+        raise FileExistsError(output)
+    load_result(result_root); output.parent.mkdir(parents=True, exist_ok=True)
+    stage = Path(tempfile.mkdtemp(prefix=f".{output.name}.partial-", dir=output.parent))
+    source = np.load(result_root / "arrays.npz", allow_pickle=False)
+    summary = json.loads((result_root / "summary.json").read_text())
+    try:
+        x, y, mask = source["x_D"], source["y_D"], source["domain_mask"]
+        fig, ax = plt.subplots(figsize=(12, 4), layout="constrained")
+        grid = np.where(mask, 1.0, np.nan); _panel(ax, x, y, grid, "Frozen pinball-inclusive common-fluid domain", None, "Greys")
+        for name, low, high in REGIONS:
+            ax.axvline(low, color="black", lw=.5); ax.text((low + high) / 2, 4.4, name.replace("_", " "), ha="center", fontsize=8)
+        fig.suptitle("Corrected residual: DRL cycle/bin fluctuation minus constant ensemble phase template")
+        _save(fig, stage, "01_domain_estimand")
+
+        fig, axes = plt.subplots(3, 3, figsize=(15, 7.5), sharex=True, sharey=True, layout="constrained")
+        left = source["left_action_physical_coordinates"]
+        velocity_values = source["ccd_mode_velocity"][:, :, mask]
+        velocity_limit = float(np.percentile(np.abs(velocity_values), 99.5))
+        rotation_values = source["ccd_mode_rotation_proxy"][:, mask]
+        rotation_limit = float(np.percentile(np.abs(rotation_values), 99.5))
+        velocity_limit = max(velocity_limit, 1e-12); rotation_limit = max(rotation_limit, 1e-12)
+        for mode in range(3):
+            velocity = source["ccd_mode_velocity"][mode]
+            rotation = source["ccd_mode_rotation_proxy"][mode]
+            _panel(axes[mode, 0], x, y, np.where(mask, velocity[0], np.nan), f"CCD M{mode+1} $u_x$", velocity_limit)
+            _panel(axes[mode, 1], x, y, np.where(mask, velocity[1], np.nan), f"CCD M{mode+1} $u_y$", velocity_limit)
+            action = left[:, mode]
+            _panel(axes[mode, 2], x, y, np.where(mask, rotation, np.nan), f"body-local rotation proxy; F/Rs/Ra={action.round(2)}", rotation_limit)
+        fig.suptitle("Pinball-inclusive CCD modes; shared robust signed color scales\n"
+                      "rotation proxy is descriptive, not wall vorticity or causality")
+        _save(fig, stage, "02_ccd_mode_identity")
+
+        bins = source["representative_phase_bins"]
+        fig, axes = plt.subplots(4, 4, figsize=(15, 8), sharex=True, sharey=True, layout="constrained")
+        fields = (("raw", source["representative_raw"]), ("CCD r2", source["representative_ccd_rank_2"]), ("CCD r3", source["representative_ccd_rank_3"]), ("POD r3", source["representative_pod_rank_3"]))
+        points = int(mask.sum())
+        global_limit = max(np.max(np.abs(values)) for _, values in fields)
+        for row, (label, values) in enumerate(fields):
+            for col, phase_bin in enumerate(bins):
+                speed_signed = values[:points, col]
+                _panel(axes[row, col], x, y, _field(speed_signed, mask), f"{label}; bin {phase_bin}", global_limit)
+        fig.suptitle("Four preregistered phase-bin mean residual reconstructions (ux; common scale)")
+        _save(fig, stage, "03_phase_reconstruction")
+
+        fig, axes = plt.subplots(3, 3, figsize=(15, 7.5), sharex=True, sharey=True, layout="constrained")
+        for mode in range(3):
+            ccd_v, pod_v = source["ccd_mode_velocity"][mode], source["pod_mode_velocity"][mode]
+            ccd_lim = np.max(np.abs(ccd_v[:, mask])); pod_lim = np.max(np.abs(pod_v[:, mask]))
+            _panel(axes[mode, 0], x, y, np.where(mask, ccd_v[0], np.nan), f"CCD M{mode+1} ux", ccd_lim)
+            _panel(axes[mode, 1], x, y, np.where(mask, pod_v[0], np.nan), f"raw POD E{mode+1} ux", pod_lim)
+            ccd_phase = summary["reconstruction"]["ccd"][f"rank_{mode+1}"]["full"]["phase_residual_mean"]
+            pod_phase = summary["reconstruction"]["pod"][f"rank_{mode+1}"]["full"]["phase_residual_mean"]
+            axes[mode, 2].plot(range(10), ccd_phase, "o-", label="CCD"); axes[mode, 2].plot(range(10), pod_phase, "s-", label="POD")
+            axes[mode, 2].set(title=f"rank {mode+1} residual", xlabel="phase bin", ylabel="relative residual"); axes[mode, 2].legend(); axes[mode, 2].grid(alpha=.2)
+        fig.suptitle(f"Raw energy-ranked POD versus CCD; principal cosines={np.round(summary['ccd_pod_principal_cosines'],3)}")
+        _save(fig, stage, "04_raw_pod_vs_ccd")
+
+        plot_keys = ("x_D", "y_D", "domain_mask", "point_weights", "ccd_mode_velocity",
+                     "ccd_mode_rotation_proxy", "ccd_modes", "ccd_left_functions",
+                     "left_action_physical_coordinates", "representative_phase_bins")
+        np.savez_compressed(stage / "plot_data.npz", **{key: source[key] for key in plot_keys})
+        (stage / "plot_data.json").write_bytes(canonical_json({
+            "schema_id": "ccd-karman-pinball-plot-data/v1",
+            "source_result_manifest_sha256": file_sha256(result_root / "manifest.json"),
+            "field_semantics": {
+                "ccd_mode_velocity": "(mode, component ux/uy, x, y), shared signed velocity scale",
+                "ccd_mode_rotation_proxy": "(mode, x, y), Gaussian body-local signed angular-momentum proxy; visualization only",
+                "ccd_modes": "component-major vectors, ux selected points then uy selected points",
+            },
+            "geometry": BODIES, "domain": DOMAIN, "regions": REGIONS,
+            "plot_files": ["02_ccd_mode_identity.png", "02_ccd_mode_identity.pdf"],
+        }))
+        (stage / "RESULTS.md").write_text(_results_text(summary) + "\nFour PNG/PDF groups provide domain, CCD mode, phase reconstruction, and raw-POD comparison views. `plot_data.npz` and `plot_data.json` are attached for reproducible later figure adjustment.\n")
+        files = {path.name: file_sha256(path) for path in stage.iterdir()}
+        (stage / "manifest.json").write_bytes(canonical_json({"schema_id": PUBLICATION_SCHEMA_ID, "complete": True, "source_manifest_sha256": file_sha256(result_root / "manifest.json"), "files": files}))
+        rename_noreplace(stage, output); stage = None
+        return output
+    finally:
+        source.close()
+        if stage is not None:
+            shutil.rmtree(stage, ignore_errors=True)
@@ -9,7 +9,9 @@ from .artifacts import load_role_artifact
 from .contracts import canonical_json, verify_decomposition
 from .phase import load_phase_compact

-SCHEMA_ID="ccd-karman-dynamic-increment/v1"
+SCHEMA_ID="ccd-karman-roi-mean/v2"
+ROI_X_D=(34.0,54.0)
+ROI_ABS_Y_D_MAX=5.0
 ROLES=("drl","constant_mean","target","zero")

 def _wrms(field, weights):
@@ -43,15 +45,23 @@ def publish_dynamic_increment(role_paths,phase_paths,output):
  if not np.array_equal(a["x_D"],x) or not np.array_equal(a["y_D"],y): raise ValueError(f"{r} grid mismatch")
 common=np.logical_and.reduce([loaded[r]["legacy_arrays"]["fluid_mask"] for r in ROLES])
 if common.sum()<4: raise ValueError("four-mask intersection too small")
- weights=(coordinate_weights(x)[:,None]*coordinate_weights(y)[None,:])[common]
+ full_weights=(coordinate_weights(x)[:,None]*coordinate_weights(y)[None,:])[common]
+ roi_grid=(x[:,None]>=ROI_X_D[0])&(x[:,None]<=ROI_X_D[1])&(np.abs(y[None,:])<=ROI_ABS_Y_D_MAX)
+ roi_mask=common&roi_grid
+ if roi_mask.sum()<4: raise ValueError("wake ROI intersection too small")
+ common_flat_roi=roi_grid[common]
+ roi_weights=full_weights[common_flat_roi]
 means={}; stats={}
- for r in ROLES: means[r],stats[r]=_dense_statistics(Path(role_paths[r]),common,weights)
+ for r in ROLES: means[r],stats[r]=_dense_statistics(Path(role_paths[r]),common,full_weights)
 for r in ("zero","constant_mean","drl"):
-  stats[r]["mean_target_error_weighted_vector_rms"]=_wrms(means[r]-means["target"],weights)
- stats["target"]["mean_target_error_weighted_vector_rms"]=0.0
+  delta=means[r]-means["target"]
+  stats[r]["mean_target_error_roi_weighted_vector_rms"]=_wrms(delta[:,common_flat_roi],roi_weights)
+  stats[r]["mean_target_error_full_domain_weighted_vector_rms_secondary"]=_wrms(delta,full_weights)
+ stats["target"]["mean_target_error_roi_weighted_vector_rms"]=0.0
+ stats["target"]["mean_target_error_full_domain_weighted_vector_rms_secondary"]=0.0
 zero_constant=means["constant_mean"]-means["zero"]; constant_drl=means["drl"]-means["constant_mean"]
 phase={r:_phase_mean(phase_paths[r],common) for r in ROLES}; phase_available=phase["drl"] is not None and phase["constant_mean"] is not None
- arrays={"x_D":x,"y_D":y,"four_role_fluid_mask":common,"quadrature_weights":weights,"mean_drl":means["drl"],"mean_constant_mean":means["constant_mean"],"mean_target":means["target"],"mean_zero":means["zero"],"mean_increment_zero_to_constant":zero_constant,"mean_increment_constant_to_drl":constant_drl}
+ arrays={"x_D":x,"y_D":y,"four_role_fluid_mask":common,"quadrature_weights":full_weights,"roi_fluid_mask":roi_mask,"roi_selector_on_common_mask":common_flat_roi,"roi_quadrature_weights":roi_weights,"mean_drl":means["drl"],"mean_constant_mean":means["constant_mean"],"mean_target":means["target"],"mean_zero":means["zero"],"mean_increment_zero_to_constant":zero_constant,"mean_increment_constant_to_drl":constant_drl}
 phase_metrics={}
 if phase_available:
  total=phase["drl"]-phase["constant_mean"]; centered=(phase["drl"]-means["drl"])-(phase["constant_mean"]-means["constant_mean"])
@@ -61,14 +71,14 @@ def publish_dynamic_increment(role_paths,phase_paths,output):
  if phase["target"] is not None:
   for r in ("zero","constant_mean","drl"):
    if phase[r] is not None:
-     vals=[_wrms(phase[r][b]-phase["target"][b],weights) for b in range(10)]; phase_metrics[r]={"target_error_by_bin_weighted_vector_rms":vals,"target_error_cycle_mean_weighted_vector_rms":float(np.mean(vals))}
+     vals=[_wrms((phase[r][b]-phase["target"][b])[:,common_flat_roi],roi_weights) for b in range(10)]; phase_metrics[r]={"target_error_by_bin_weighted_vector_rms":vals,"target_error_cycle_mean_weighted_vector_rms":float(np.mean(vals))}
  closure=total-((means["drl"]-means["constant_mean"])[None]+centered); closure_max=float(np.max(np.abs(closure)))
  if closure_max>2e-6 or not verify_decomposition(phase["drl"],phase["constant_mean"]): raise ValueError("phase decomposition closure failed")
 else: closure_max=None
- benefits={"zero_to_constant_overall_mean_control_benefit_target_error_reduction":stats["zero"]["mean_target_error_weighted_vector_rms"]-stats["constant_mean"]["mean_target_error_weighted_vector_rms"],"constant_to_drl_dynamic_increment_target_error_reduction":stats["constant_mean"]["mean_target_error_weighted_vector_rms"]-stats["drl"]["mean_target_error_weighted_vector_rms"],"zero_to_drl_total_target_error_reduction":stats["zero"]["mean_target_error_weighted_vector_rms"]-stats["drl"]["mean_target_error_weighted_vector_rms"]}
+ benefits={"zero_to_constant_overall_mean_control_benefit_target_error_reduction":stats["zero"]["mean_target_error_roi_weighted_vector_rms"]-stats["constant_mean"]["mean_target_error_roi_weighted_vector_rms"],"constant_to_drl_dynamic_increment_target_error_reduction":stats["constant_mean"]["mean_target_error_roi_weighted_vector_rms"]-stats["drl"]["mean_target_error_roi_weighted_vector_rms"],"zero_to_drl_total_target_error_reduction":stats["zero"]["mean_target_error_roi_weighted_vector_rms"]-stats["drl"]["mean_target_error_roi_weighted_vector_rms"]}
 benefit_closure=benefits["zero_to_constant_overall_mean_control_benefit_target_error_reduction"]+benefits["constant_to_drl_dynamic_increment_target_error_reduction"]-benefits["zero_to_drl_total_target_error_reduction"]
 parents={r:{"role_path":str(Path(role_paths[r]).resolve()),"role_campaign_manifest_sha256":file_sha256(Path(role_paths[r])/"campaign_manifest.json"),"phase_path":str(Path(phase_paths[r]).resolve()),"phase_manifest_sha256":file_sha256(Path(phase_paths[r])/"manifest.json"),"phase_gate_passed":phase[r] is not None} for r in ROLES}
- summary={"schema_id":SCHEMA_ID,"complete":True,"drl_constant_phase_differences_available":phase_available,"phase_gate_passed_by_role":{r:phase[r] is not None for r in ROLES},"phase_blockers":[r for r in ROLES if phase[r] is None],"mask_definition":"exact intersection of all four solver-derived fluid masks","analysis_fluid_point_count":int(common.sum()),"quadrature_rule":"coordinate_weights(x_D)*coordinate_weights(y_D) on common mask; not area-normalized","statistics":stats,"phase_target_error_metrics":phase_metrics,"benefits":benefits,"closure":{"benefit_additivity_absolute_residual":float(abs(benefit_closure)),"mean_increment_max_absolute_residual":float(np.max(np.abs((means['drl']-means['zero'])-(zero_constant+constant_drl)))),"phase_decomposition_max_absolute_residual":closure_max},"claims":"independent phase-conditioned trajectory means; not pointwise counterfactual, response, or causal effect","dense_fields_deleted":False,"parents":parents}
+ summary={"schema_id":SCHEMA_ID,"complete":True,"drl_constant_phase_differences_available":phase_available,"phase_gate_passed_by_role":{r:phase[r] is not None for r in ROLES},"phase_blockers":[r for r in ROLES if phase[r] is None],"mask_definition":"exact intersection of all four solver-derived fluid masks","analysis_fluid_point_count":int(common.sum()),"primary_domain":"wake ROI inclusive 34<=x/D<=54 and |y/D|<=5","roi_fluid_point_count":int(roi_mask.sum()),"full_domain_metrics_role":"secondary","quadrature_rule":"coordinate_weights(x_D)*coordinate_weights(y_D) on common mask; not area-normalized","statistics":stats,"phase_target_error_metrics":phase_metrics,"benefits":benefits,"closure":{"benefit_additivity_absolute_residual":float(abs(benefit_closure)),"mean_increment_max_absolute_residual":float(np.max(np.abs((means['drl']-means['zero'])-(zero_constant+constant_drl)))),"phase_decomposition_max_absolute_residual":closure_max},"claims":"independent phase-conditioned trajectory means; not pointwise counterfactual, response, or causal effect","dense_fields_deleted":False,"parents":parents}
 output.parent.mkdir(parents=True,exist_ok=True); stage=Path(tempfile.mkdtemp(prefix=f".{output.name}.partial-",dir=output.parent))
 try:
  np.savez_compressed(stage/"arrays.npz",**arrays); (stage/"summary.json").write_bytes(canonical_json(summary)); files={p.name:file_sha256(p) for p in stage.iterdir() if p.is_file()}; (stage/"manifest.json").write_bytes(canonical_json({"schema_id":SCHEMA_ID,"complete":True,"files":files})); rename_noreplace(stage,output); return load_dynamic_increment(output)
@@ -80,10 +90,10 @@ def load_dynamic_increment(path):
 for n,h in manifest["files"].items():
  if file_sha256(path/n)!=h: raise ValueError("dynamic increment file hash mismatch")
 for r,p in summary["parents"].items():
-  load_role_artifact(p["role_path"],expected_role=r); load_phase_compact(p["phase_path"])
+  load_role_artifact(p["role_path"],expected_role=r); load_phase_compact(p["phase_path"],allow_legacy=True)
  if file_sha256(Path(p["role_path"])/"campaign_manifest.json")!=p["role_campaign_manifest_sha256"] or file_sha256(Path(p["phase_path"])/"manifest.json")!=p["phase_manifest_sha256"]: raise ValueError("dynamic increment live provenance mismatch")
 with np.load(path/"arrays.npz",allow_pickle=False) as z:
-  required={"x_D","y_D","four_role_fluid_mask","quadrature_weights","mean_drl","mean_constant_mean","mean_target","mean_zero","mean_increment_zero_to_constant","mean_increment_constant_to_drl"}
+  required={"x_D","y_D","four_role_fluid_mask","quadrature_weights","roi_fluid_mask","roi_selector_on_common_mask","roi_quadrature_weights","mean_drl","mean_constant_mean","mean_target","mean_zero","mean_increment_zero_to_constant","mean_increment_constant_to_drl"}
  if not required.issubset(z.files): raise ValueError("dynamic increment arrays incomplete")
  if np.max(np.abs((z["mean_drl"]-z["mean_zero"])-(z["mean_increment_zero_to_constant"]+z["mean_increment_constant_to_drl"])))>2e-6: raise ValueError("mean increment closure mismatch")
 return {"path":path,"summary":summary}
@@ -0,0 +1,320 @@
+"""No-clobber publisher for compact CCD mode and full-fit POD interpretation."""
+from __future__ import annotations
+
+import json
+import shutil
+import tempfile
+from pathlib import Path
+
+import matplotlib
+matplotlib.use("Agg")
+import matplotlib.pyplot as plt
+import numpy as np
+
+from CCD_analysis.acquisition.artifacts import file_sha256, rename_noreplace
+from CCD_analysis.direct_dq.schema import canonical_array_sha256
+from .contracts import canonical_json
+from .cycle_template_ccd import load_cycle_template_result
+from .dynamic_increment import load_dynamic_increment
+from .mode_diagnostics import (
+    action_coordinates, bounded_streamwise_shift_correlation,
+    component_major_velocity, masked_vorticity, phase_harmonics,
+    reshape_cycle_phase, subspace_symmetry_diagnostics,
+    symmetry_diagnostics, x_localization,
+)
+from .phase import load_phase_compact
+from .pod_baseline import (
+    align_basis_to_reference, fit_weighted_pod,
+    reconstruct_residual_from_parents, weighted_principal_cosines,
+)
+
+SCHEMA_ID = "ccd-karman-mode-interpretation/v2"
+PUBLICATION_SCHEMA_ID = "ccd-karman-mode-publication/v2"
+EXPECTED = {
+    "drl_phase": "b6c8a8b716a384cb3225298523aad775cff977124f8c1d54a79a52c5438260c6",
+    "constant_mean_phase": "1d608df3a658e334b79e20a4f690c4978418d30e229e6cf0f244d6484bcf84b4",
+    "roi_mean": "f7f6141042678e52314a4d8f76d589108da8e844c8ce504f8f523a6f88a17e73",
+    "cycle_template_ccd": "b114fbcc5f821034829ca3771ac5ba13df4b75df312b0615b662f81a9e5737d8",
+    "publication": "a6f0c55515b0e4a3b8980174988bf6bbea3617d60f8648aab25a7ffa0a2a24cd",
+}
+CLAIMS = "full-fit artifact-only descriptive diagnostics; no held-out test, prediction, causality, mechanism, independent uncertainty, or CCD superiority"
+REGIONS = (("near", 34.0, 40.0), ("middle", 40.0, 47.0), ("far", 47.0, 54.000001))
+
+
+def _jsonable(value):
+    if isinstance(value, np.ndarray):
+        return value.tolist()
+    if isinstance(value, np.generic):
+        return value.item()
+    if isinstance(value, dict):
+        return {key: _jsonable(item) for key, item in value.items()}
+    return value
+
+
+def _field(values, mask):
+    field = np.full(mask.shape, np.nan)
+    field[mask] = values
+    return field
+
+
+def _panel(ax, x, y, field, title, limit=None, cmap="RdBu_r"):
+    options = {} if limit is None else {"vmin": -limit, "vmax": limit}
+    image = ax.pcolormesh(x, y, field.T, shading="nearest", cmap=cmap, rasterized=True, **options)
+    ax.set(xlim=(34, 54), ylim=(-5, 5), xlabel="x/D", ylabel="y/D", title=title)
+    ax.set_aspect("equal")
+    return image
+
+
+def compute(parents):
+    """Compute the compact interpretation while binding all live parents."""
+    roots = {name: Path(path) for name, path in parents.items()}
+    if set(roots) != set(EXPECTED):
+        raise ValueError("all five declared parents are required")
+    for name, digest in EXPECTED.items():
+        if file_sha256(roots[name] / "manifest.json") != digest:
+            raise ValueError(f"{name} parent manifest identity mismatch")
+
+    load_phase_compact(roots["drl_phase"])
+    load_phase_compact(roots["constant_mean_phase"])
+    mean_result = load_dynamic_increment(roots["roi_mean"])
+    ccd = load_cycle_template_result(roots["cycle_template_ccd"], recompute=False)
+    admitted_mean = ccd["arrays"]["admitted_field_empirical_mean"]
+    residual = reconstruct_residual_from_parents(
+        roots["drl_phase"], roots["constant_mean_phase"], roots["roi_mean"], admitted_mean
+    )
+    centered, weights = residual["centered"], residual["weights"]
+    pod = fit_weighted_pod(centered, weights, rank=3, center=False)
+    ccd_modes = ccd["arrays"]["primary_physical_modes"][:, :3]
+    ccd_coefficients = ccd["arrays"]["primary_coefficients"][:3]
+    projected_coefficients = ccd_modes.T @ (weights[:, None] * centered)
+    np.testing.assert_allclose(projected_coefficients, ccd_coefficients, rtol=2e-12, atol=2e-12)
+    principal_cosines = weighted_principal_cosines(ccd_modes, pod["modes"], weights)
+    pod_aligned, pod_rotation = align_basis_to_reference(ccd_modes, pod["modes"], weights)
+    matched_cosines = np.diag(ccd_modes.T @ (weights[:, None] * pod_aligned))
+
+    with np.load(roots["roi_mean"] / "arrays.npz", allow_pickle=False) as source:
+        x, y, mask = source["x_D"].copy(), source["y_D"].copy(), source["roi_fluid_mask"].copy()
+        grid_weights = np.zeros(mask.shape)
+        grid_weights[mask] = source["roi_quadrature_weights"]
+    with np.load(roots["drl_phase"] / "compact.npz", allow_pickle=False) as source:
+        actions_native = source["cycle_bin_effective_actions"].reshape(190, 3).T
+    actions_transformed = action_coordinates(actions_native)
+    left_action_coordinates = action_coordinates(ccd["arrays"]["primary_left_functions"][:, :3])
+
+    velocities, vorticities, vorticity_valid, mode_metrics = [], [], [], []
+    for index in range(3):
+        velocity = component_major_velocity(ccd_modes[:, index], mask)
+        vorticity, valid = masked_vorticity(*velocity, x, y, mask)
+        velocities.append(velocity)
+        vorticities.append(vorticity)
+        vorticity_valid.append(valid)
+        mode_metrics.append({
+            "mode": index + 1,
+            "symmetry": symmetry_diagnostics(velocity, grid_weights, mask, y),
+            "velocity_localization": x_localization(np.sum(velocity**2, axis=0), x, grid_weights, mask, REGIONS),
+            "vorticity_localization": x_localization(vorticity**2, x, grid_weights, valid, REGIONS),
+        })
+    velocities, vorticities, vorticity_valid = map(np.asarray, (velocities, vorticities, vorticity_valid))
+    subspace_symmetry = subspace_symmetry_diagnostics(ccd_modes, grid_weights, mask, y)
+
+    shifts = {}
+    for first, second in ((0, 1), (0, 2), (1, 2)):
+        key = f"mode_{first + 1}_{second + 1}"
+        shifts[f"{key}_velocity"] = bounded_streamwise_shift_correlation(
+            velocities[first], velocities[second], x, grid_weights, mask, 3.0
+        )
+        shifts[f"{key}_vorticity"] = bounded_streamwise_shift_correlation(
+            vorticities[first], vorticities[second], x, grid_weights,
+            vorticity_valid[first] & vorticity_valid[second], 3.0,
+        )
+
+    coefficient_harmonics = phase_harmonics(reshape_cycle_phase(ccd_coefficients), 5)
+    action_harmonics_native = phase_harmonics(reshape_cycle_phase(actions_native), 5)
+    action_harmonics_transformed = phase_harmonics(reshape_cycle_phase(actions_transformed), 5)
+    arrays = {
+        "ccd_modes": ccd_modes, "ccd_coefficients": ccd_coefficients,
+        "mode_velocity": velocities, "mode_vorticity": vorticities,
+        "vorticity_valid_mask": vorticity_valid,
+        "ccd_subspace_reflection_principal_cosines": subspace_symmetry["reflection_principal_cosines"],
+        "actions_native_front_upper_lower": actions_native,
+        "actions_front_rear_symmetric_antisymmetric": actions_transformed,
+        "left_action_physical_coordinates": left_action_coordinates,
+        "pod_modes": pod["modes"], "pod_modes_aligned_to_ccd": pod_aligned,
+        "pod_alignment_rotation": pod_rotation, "pod_coefficients": pod["coefficients"],
+        "pod_field_energy_singular_values": pod["singular_values"],
+        "pod_weighted_residual_fraction": pod["weighted_residual_fraction"],
+        "ccd_pod_principal_cosines": principal_cosines,
+        "ccd_pod_aligned_mode_cosines": matched_cosines,
+        "x_D": x, "y_D": y, "roi_mask": mask, "roi_weights_grid": grid_weights,
+    }
+    for prefix, harmonic in (("ccd_coefficient", coefficient_harmonics), ("action_native", action_harmonics_native), ("action_transformed", action_harmonics_transformed)):
+        for name, values in harmonic.items():
+            arrays[f"{prefix}_harmonic_{name}"] = values
+
+    benefits = mean_result["summary"]["benefits"]
+    total = benefits["zero_to_drl_total_target_error_reduction"]
+    summary = {
+        "schema_id": SCHEMA_ID, "M": 160400, "N": 190, "rank": 3,
+        "raw_U_canonical_sha256": canonical_array_sha256(residual["raw"]),
+        "centered_U_canonical_sha256": canonical_array_sha256(centered),
+        "admitted_mean_bit_equal": True,
+        "ccd_coefficient_projection_max_absolute_residual": float(np.max(np.abs(projected_coefficients - ccd_coefficients))),
+        "mode_metrics": mode_metrics, "subspace_symmetry": _jsonable(subspace_symmetry),
+        "shift_metrics": {name: {key: _jsonable(value[key]) for key in ("best_shift_x_D", "correlation", "paired_point_count", "index_offset")} for name, value in shifts.items()},
+        "ccd_pod_principal_cosines": principal_cosines.tolist(),
+        "ccd_pod_principal_angles_degrees": np.degrees(np.arccos(principal_cosines)).tolist(),
+        "ccd_pod_aligned_mode_cosines": matched_cosines.tolist(),
+        "pod_weighted_residual_fraction_rank_0_to_3": pod["weighted_residual_fraction"].tolist(),
+        "contribution": {
+            "constant_mean_share_of_zero_to_drl_improvement": benefits["zero_to_constant_overall_mean_control_benefit_target_error_reduction"] / total,
+            "drl_dynamic_share_of_zero_to_drl_improvement": benefits["constant_to_drl_dynamic_increment_target_error_reduction"] / total,
+        },
+        "phase_harmonic_semantics": "10-bin circular harmonics computed separately for each of 19 cycles; coherent complex mean, mean cycle amplitude, phase locking, and every cycle amplitude/phase are persisted; not temporal PSD",
+        "mode_semantics": "individual modes are basis-dependent primary SVD vectors; reflection-subspace invariance is basis-invariant",
+        "pod_semantics": pod["label"], "claim_boundary": CLAIMS,
+    }
+    config = {
+        "schema_id": SCHEMA_ID, "component_order": "ux ROI points, then uy ROI points",
+        "residual": "same raw U rebuilt through cycle_template_ccd._load/_matrix, bit-matched admitted empirical mean, then exactly CCD-centered",
+        "quadrature": "ROI coordinate weights duplicated for ux and uy",
+        "vorticity": "axis 0=x and axis 1=y; omega_z=d(uy)/dx-d(ux)/dy; nonuniform three-point stencils requiring complete fluid support",
+        "reflection": "explicit y->-y; vector parity ux even, uy odd; individual and subspace diagnostics",
+        "phase_shape": [19, 10],
+        "action_coordinates": ["front", "rear_symmetric=(upper+lower)/sqrt(2)", "rear_antisymmetric=(upper-lower)/sqrt(2)"],
+        "shift_sign": "positive shift correlates a(x) with b(x+shift); geometric-mean pair weights",
+        "pod": "full-fit weighted rank-3 field-energy baseline, descriptive only",
+        "claim_boundary": CLAIMS,
+    }
+    hashes = {name: {"path": str(roots[name].resolve()), "manifest_sha256": EXPECTED[name]} for name in EXPECTED}
+    return arrays, config, summary, hashes
+
+
+def publish_result(parents, destination):
+    destination = Path(destination)
+    if destination.exists():
+        raise FileExistsError(destination)
+    destination.parent.mkdir(parents=True, exist_ok=True)
+    stage = Path(tempfile.mkdtemp(prefix=f".{destination.name}.partial-", dir=destination.parent))
+    try:
+        arrays, config, summary, hashes = compute(parents)
+        np.savez_compressed(stage / "arrays.npz", **arrays)
+        for name, value in (("config.json", config), ("summary.json", summary), ("input_hashes.json", hashes)):
+            (stage / name).write_bytes(canonical_json(_jsonable(value)))
+        (stage / "RESULTS.md").write_text(
+            "# Kármán mode interpretation extension\n\n"
+            f"This derived artifact analyzes the exact corrected residual U (M=160400, N=190). {CLAIMS}.\n\n"
+            f"CCD/POD rank-3 principal cosines are {summary['ccd_pod_principal_cosines']}; full-fit POD residual fractions rank 0..3 are {summary['pod_weighted_residual_fraction_rank_0_to_3']}. Harmonics are circular 10-bin phase diagnostics, never a flattened 190-point FFT or temporal PSD.\n"
+        )
+        files = {path.name: file_sha256(path) for path in stage.iterdir()}
+        (stage / "manifest.json").write_bytes(canonical_json({"schema_id": SCHEMA_ID, "complete": True, "files": files}))
+        rename_noreplace(stage, destination)
+        stage = None
+        load_result(destination, recompute=True)
+        return destination
+    finally:
+        if stage is not None:
+            shutil.rmtree(stage, ignore_errors=True)
+
+
+def load_result(path, recompute=False):
+    root = Path(path)
+    manifest = json.loads((root / "manifest.json").read_text())
+    expected = {"arrays.npz", "config.json", "summary.json", "input_hashes.json", "RESULTS.md"}
+    if manifest.get("schema_id") != SCHEMA_ID or set(manifest.get("files", {})) != expected:
+        raise ValueError("interpretation manifest invalid")
+    for name, digest in manifest["files"].items():
+        if file_sha256(root / name) != digest:
+            raise ValueError("interpretation file hash mismatch")
+    hashes = json.loads((root / "input_hashes.json").read_text())
+    summary = json.loads((root / "summary.json").read_text())
+    for name, value in hashes.items():
+        if file_sha256(Path(value["path"]) / "manifest.json") != value["manifest_sha256"]:
+            raise ValueError(f"live parent changed: {name}")
+    if recompute:
+        fresh = compute({name: value["path"] for name, value in hashes.items()})[2]
+        for key in ("raw_U_canonical_sha256", "centered_U_canonical_sha256"):
+            if summary[key] != fresh[key]:
+                raise ValueError(f"fresh {key} changed")
+        np.testing.assert_allclose(summary["ccd_pod_principal_cosines"], fresh["ccd_pod_principal_cosines"], rtol=1e-11, atol=1e-12)
+        np.testing.assert_allclose(summary["pod_weighted_residual_fraction_rank_0_to_3"], fresh["pod_weighted_residual_fraction_rank_0_to_3"], rtol=1e-11, atol=1e-12)
+    return {"path": root, "summary": summary, "manifest": manifest, "provenance_validation": "VERIFIED live parents and fresh key recomputation" if recompute else "VERIFIED"}
+
+
+def _save(fig, stage, stem):
+    for extension in ("png", "pdf"):
+        fig.savefig(stage / f"{stem}.{extension}", dpi=250 if extension == "png" else None, bbox_inches="tight", metadata={"Creator": "CCD mode interpretation"})
+    plt.close(fig)
+
+
+def publish_figures(result_root, old_publication, output):
+    """Publish four actual PNG/PDF diagnostic groups as a no-clobber sibling."""
+    result, root, old, output = load_result(result_root), Path(result_root), Path(old_publication), Path(output)
+    if output.exists():
+        raise FileExistsError(output)
+    output.parent.mkdir(parents=True, exist_ok=True)
+    stage = Path(tempfile.mkdtemp(prefix=f".{output.name}.partial-", dir=output.parent))
+    arrays = np.load(root / "arrays.npz", allow_pickle=False)
+    summary = result["summary"]
+    try:
+        shutil.copy2(old / "01_roi_mean_performance.png", stage / "01_contribution_performance.png")
+        shutil.copy2(old / "01_roi_mean_performance.pdf", stage / "01_contribution_performance.pdf")
+        shutil.copy2(old / "02_roi_mean_target_error_fields.png", stage / "02_target_error_fields.png")
+        shutil.copy2(old / "02_roi_mean_target_error_fields.pdf", stage / "02_target_error_fields.pdf")
+        x, y, mask = arrays["x_D"], arrays["y_D"], arrays["roi_mask"]
+
+        fig, axes = plt.subplots(2, 2, figsize=(11, 7), layout="constrained")
+        amplitudes = arrays["ccd_coefficient_harmonic_coherent_amplitude"][:, 1:]
+        for mode in range(3): axes[0, 0].plot(range(1, 6), amplitudes[mode], "o-", label=f"mode {mode+1}")
+        axes[0, 0].set(title="CCD coefficient coherent phase harmonics", xlabel="10-bin harmonic", ylabel="amplitude"); axes[0, 0].legend()
+        transformed = arrays["action_transformed_harmonic_coherent_amplitude"][:, 1:]
+        for channel, label in enumerate(("front", "rear symmetric", "rear antisymmetric")): axes[0, 1].plot(range(1, 6), transformed[channel], "o-", label=label)
+        axes[0, 1].set(title="Physical action-coordinate harmonics", xlabel="10-bin harmonic", ylabel="amplitude"); axes[0, 1].legend()
+        symmetric = [item["symmetry"]["symmetric_fraction"] for item in summary["mode_metrics"]]
+        antisymmetric = [item["symmetry"]["antisymmetric_fraction"] for item in summary["mode_metrics"]]
+        axes[1, 0].bar(np.arange(3)+1, symmetric, label="symmetric"); axes[1, 0].bar(np.arange(3)+1, antisymmetric, bottom=symmetric, label="antisymmetric")
+        axes[1, 0].set(title=f"Vector reflection parity; subspace invariance={summary['subspace_symmetry']['invariance_fraction']:.3f}", xlabel="CCD basis mode", ylabel="paired weighted fraction", ylim=(0, 1)); axes[1, 0].legend()
+        region_names = ("near", "middle", "far"); bottom = np.zeros(3)
+        for region in region_names:
+            values = np.array([item["vorticity_localization"]["region_fractions"][region] for item in summary["mode_metrics"]]); axes[1, 1].bar(np.arange(3)+1, values, bottom=bottom, label=region); bottom += values
+        axes[1, 1].set(title="Vorticity localization", xlabel="CCD basis mode", ylabel="weighted fraction", ylim=(0, 1)); axes[1, 1].legend()
+        for axis in axes.flat: axis.grid(alpha=.2)
+        fig.suptitle("Mode physical diagnostics (basis-dependent vectors; subspace reflection also reported)")
+        _save(fig, stage, "03_mode_identity_diagnostics")
+
+        fig, axes = plt.subplots(3, 2, figsize=(11, 8), sharex=True, sharey=True, layout="constrained")
+        for mode in range(3):
+            velocity, omega = arrays["mode_velocity"][mode], arrays["mode_vorticity"][mode]
+            speed = np.sqrt(np.sum(velocity**2, axis=0)); speed_limit = np.nanmax(speed[mask]); omega_limit = np.nanmax(np.abs(omega[arrays["vorticity_valid_mask"][mode]]))
+            _panel(axes[mode, 0], x, y, np.where(mask, speed, np.nan), f"CCD mode {mode+1}: speed", None, "viridis").set_clim(0, speed_limit)
+            _panel(axes[mode, 1], x, y, np.where(arrays["vorticity_valid_mask"][mode], omega, np.nan), f"vorticity; centroid {summary['mode_metrics'][mode]['vorticity_localization']['centroid_x_D']:.2f}D", omega_limit)
+        action_labels = []
+        for mode in range(3):
+            coords = arrays["left_action_physical_coordinates"][:, mode]
+            harmonic = arrays["ccd_coefficient_harmonic_coherent_amplitude"][mode, 1]
+            action_labels.append(f"M{mode + 1}: action=({coords[0]:.2f},{coords[1]:.2f},{coords[2]:.2f}), k=1 amp={harmonic:.3f}")
+        fig.suptitle("Corrected residual CCD spatial modes and masked vorticity\n" + "  |  ".join(action_labels), fontsize=10)
+        _save(fig, stage, "04_ccd_spatial_vorticity")
+
+        fig, axes = plt.subplots(3, 3, figsize=(14, 8), sharex=True, sharey=True, layout="constrained")
+        for mode in range(3):
+            ccd_velocity = component_major_velocity(arrays["ccd_modes"][:, mode], mask)
+            pod_velocity = component_major_velocity(arrays["pod_modes_aligned_to_ccd"][:, mode], mask)
+            difference = ccd_velocity - pod_velocity
+            limit = max(np.max(np.abs(ccd_velocity[:, mask])), np.max(np.abs(pod_velocity[:, mask])))
+            _panel(axes[mode, 0], x, y, np.where(mask, ccd_velocity[0], np.nan), f"CCD {mode+1} ux", limit)
+            _panel(axes[mode, 1], x, y, np.where(mask, pod_velocity[0], np.nan), f"aligned POD {mode+1} ux; cos={summary['ccd_pod_aligned_mode_cosines'][mode]:.3f}", limit)
+            difference_limit = max(np.max(np.abs(difference[0, mask])), 1e-15)
+            _panel(axes[mode, 2], x, y, np.where(mask, difference[0], np.nan), f"CCD - aligned POD {mode+1}", difference_limit)
+        fig.suptitle(f"Rank-3 spatial comparison; principal cosines {np.round(summary['ccd_pod_principal_cosines'],3)}; full-fit POD descriptive only")
+        _save(fig, stage, "05_ccd_vs_full_fit_pod_spatial")
+
+        (stage / "RESULTS.md").write_text((root / "RESULTS.md").read_text() + "\nFive PNG/PDF groups show contribution, target-error fields, mode identity diagnostics, CCD spatial/vorticity, and aligned CCD-versus-full-fit-POD spatial evidence.\n")
+        files = {path.name: file_sha256(path) for path in stage.iterdir()}
+        (stage / "manifest.json").write_bytes(canonical_json({"schema_id": PUBLICATION_SCHEMA_ID, "complete": True, "source_manifest_sha256": file_sha256(root / "manifest.json"), "parent_publication_manifest_sha256": EXPECTED["publication"], "files": files}))
+        rename_noreplace(stage, output)
+        stage = None
+        return output
+    finally:
+        arrays.close()
+        if stage is not None:
+            shutil.rmtree(stage, ignore_errors=True)
@@ -0,0 +1,251 @@
+"""Pure spatial and phase diagnostics for component-major Kármán modes."""
+from __future__ import annotations
+
+import numpy as np
+
+
+def component_major_velocity(mode, mask):
+    """Expand [ux(mask), uy(mask)] into a zero-filled (2, nx, ny) grid."""
+    mask = np.asarray(mask, dtype=bool)
+    mode = np.asarray(mode, dtype=float)
+    points = int(mask.sum())
+    if mode.shape != (2 * points,) or not np.isfinite(mode).all():
+        raise ValueError("a finite component-major mode is required")
+    velocity = np.zeros((2, *mask.shape))
+    velocity[0, mask], velocity[1, mask] = mode[:points], mode[points:]
+    return velocity
+
+
+def _derivative(field, coordinate, mask, axis):
+    """Apply explicit nonuniform three-point Lagrange derivative stencils."""
+    field = np.asarray(field, dtype=float)
+    coordinate = np.asarray(coordinate, dtype=float)
+    mask = np.asarray(mask, dtype=bool)
+    if field.shape != mask.shape or axis not in (0, 1):
+        raise ValueError("invalid field, mask, or derivative axis")
+    if coordinate.ndim != 1 or coordinate.size != field.shape[axis]:
+        raise ValueError("coordinate does not match derivative axis")
+    if coordinate.size < 3 or np.any(np.diff(coordinate) <= 0):
+        raise ValueError("at least three strictly increasing coordinates are required")
+
+    values = np.moveaxis(field, axis, 0)
+    fluid = np.moveaxis(mask, axis, 0)
+    output = np.zeros_like(values)
+    valid = np.zeros_like(fluid)
+    for index in range(coordinate.size):
+        if index == 0:
+            stencil = np.array([0, 1, 2])
+        elif index == coordinate.size - 1:
+            stencil = np.arange(coordinate.size - 3, coordinate.size)
+        else:
+            stencil = np.arange(index - 1, index + 2)
+        nodes = coordinate[stencil]
+        at = coordinate[index]
+        coefficients = []
+        for local in range(3):
+            other = [item for item in range(3) if item != local]
+            numerator = 2 * at - nodes[other[0]] - nodes[other[1]]
+            denominator = (nodes[local] - nodes[other[0]]) * (nodes[local] - nodes[other[1]])
+            coefficients.append(numerator / denominator)
+        stencil_values = values[stencil]
+        admitted = np.all(fluid[stencil], axis=0) & np.all(np.isfinite(stencil_values), axis=0)
+        evaluated = np.tensordot(coefficients, stencil_values, axes=(0, 0))
+        output[index] = np.where(admitted, evaluated, 0.0)
+        valid[index] = admitted
+    return np.moveaxis(output, 0, axis), np.moveaxis(valid, 0, axis)
+
+
+def masked_vorticity(ux, uy, x, y, mask):
+    """Return omega_z=d(uy)/dx-d(ux)/dy on an axis-0=x, axis-1=y grid."""
+    ux, uy, mask = np.asarray(ux), np.asarray(uy), np.asarray(mask, dtype=bool)
+    if ux.shape != mask.shape or uy.shape != mask.shape:
+        raise ValueError("velocity and mask shapes must match")
+    duy_dx, x_valid = _derivative(uy, x, mask, axis=0)
+    dux_dy, y_valid = _derivative(ux, y, mask, axis=1)
+    valid = mask & x_valid & y_valid
+    return np.where(valid, duy_dx - dux_dy, 0.0), valid
+
+
+def reflection_indices(y, atol=1e-10):
+    """Map each y index to its unique -y reflection index."""
+    y = np.asarray(y, dtype=float)
+    reflected = []
+    for value in y:
+        matches = np.flatnonzero(np.isclose(y, -value, rtol=0, atol=atol))
+        if matches.size != 1:
+            raise ValueError("the grid does not provide a unique y reflection")
+        reflected.append(matches[0])
+    return np.asarray(reflected)
+
+
+def reflect_velocity(velocity, y):
+    """Reflect a vector field: ux is even and uy is odd under y -> -y."""
+    velocity = np.asarray(velocity, dtype=float)
+    indices = reflection_indices(y)
+    if velocity.shape[0] != 2 or velocity.shape[2] != indices.size:
+        raise ValueError("velocity must have shape (2, nx, len(y))")
+    return np.stack((velocity[0][:, indices], -velocity[1][:, indices]))
+
+
+def symmetry_diagnostics(velocity, weights, mask, y):
+    """Compute fixed-parity fractions for one basis-dependent mode."""
+    velocity = np.asarray(velocity, dtype=float)
+    weights, mask = np.asarray(weights, dtype=float), np.asarray(mask, dtype=bool)
+    indices = reflection_indices(y)
+    paired = mask & mask[:, indices]
+    pair_weights = 0.5 * (weights + weights[:, indices])
+    reflected = reflect_velocity(velocity, y)
+    symmetric = 0.5 * (velocity + reflected)
+    antisymmetric = 0.5 * (velocity - reflected)
+
+    def energy(field):
+        return float(np.sum(pair_weights[paired] * np.sum(field[:, paired] ** 2, axis=0)))
+
+    total = energy(velocity)
+    if total <= 0 or not np.isfinite(total):
+        raise ValueError("mode has no finite paired energy")
+    return {
+        "symmetric_fraction": energy(symmetric) / total,
+        "antisymmetric_fraction": energy(antisymmetric) / total,
+        "paired_weight_fraction": float(weights[paired].sum() / weights[mask].sum()),
+    }
+
+
+def subspace_symmetry_diagnostics(modes, weights, mask, y):
+    """Return basis-invariant reflection overlap and parity energy fractions."""
+    modes = np.asarray(modes, dtype=float)
+    weights, mask = np.asarray(weights, dtype=float), np.asarray(mask, dtype=bool)
+    component_weights = np.tile(weights[mask], 2)
+    gram = modes.T @ (component_weights[:, None] * modes)
+    if not np.allclose(gram, np.eye(modes.shape[1]), rtol=1e-8, atol=1e-10):
+        raise ValueError("subspace modes must be W-orthonormal")
+    reflected = []
+    symmetric_energy = antisymmetric_energy = 0.0
+    for mode in modes.T:
+        grid = reflect_velocity(component_major_velocity(mode, mask), y)
+        column = np.r_[grid[0, mask], grid[1, mask]]
+        reflected.append(column)
+        symmetric_energy += 0.25 * np.sum(component_weights * (mode + column) ** 2)
+        antisymmetric_energy += 0.25 * np.sum(component_weights * (mode - column) ** 2)
+    reflected = np.column_stack(reflected)
+    overlap = modes.T @ (component_weights[:, None] * reflected)
+    cosines = np.clip(np.linalg.svd(overlap, compute_uv=False), 0, 1)
+    rank = modes.shape[1]
+    return {
+        "reflection_principal_cosines": cosines,
+        "invariance_fraction": float(np.sum(cosines**2) / rank),
+        "symmetric_fraction": float(symmetric_energy / rank),
+        "antisymmetric_fraction": float(antisymmetric_energy / rank),
+    }
+
+
+def x_localization(density, x, weights, mask, regions=()):
+    """Summarize a nonnegative weighted density along x."""
+    density, weights = np.asarray(density, dtype=float), np.asarray(weights, dtype=float)
+    mask, x = np.asarray(mask, dtype=bool), np.asarray(x, dtype=float)
+    mass = np.where(mask, np.maximum(density, 0) * weights, 0)
+    total = float(mass.sum())
+    if total <= 0 or not np.isfinite(total):
+        raise ValueError("positive finite localization mass is required")
+    x_grid = np.broadcast_to(x[:, None], mask.shape)
+    centroid = float(np.sum(mass * x_grid) / total)
+    spread = float(np.sqrt(np.sum(mass * (x_grid - centroid) ** 2) / total))
+    fractions = {name: float(mass[(x_grid >= low) & (x_grid < high)].sum() / total) for name, low, high in regions}
+    return {"centroid_x_D": centroid, "spread_x_D": spread, "region_fractions": fractions}
+
+
+def reshape_cycle_phase(values, cycles=19, bins=10):
+    """Expose the cycle-major sample axis as (..., 19, 10)."""
+    values = np.asarray(values, dtype=float)
+    if values.shape[-1] != cycles * bins:
+        raise ValueError("cycle-major phase-minor sample axis required")
+    return values.reshape(*values.shape[:-1], cycles, bins)
+
+
+def phase_harmonics(values, max_harmonic=5):
+    """Compute 10-bin harmonics per cycle and an intentional coherent mean."""
+    values = np.asarray(values, dtype=float)
+    if values.ndim < 2 or values.shape[-2:] != (19, 10):
+        raise ValueError("phase harmonics require explicit (..., 19, 10) values")
+    if not np.isfinite(values).all() or not 0 <= max_harmonic <= 5:
+        raise ValueError("invalid phase values or harmonic order")
+    complex_cycle = np.fft.rfft(values, axis=-1)[..., : max_harmonic + 1] / 10
+    cycle_amplitude = np.abs(complex_cycle)
+    cycle_amplitude[..., 1:] *= 2
+    coherent = complex_cycle.mean(axis=-2)
+    coherent_amplitude = np.abs(coherent)
+    coherent_amplitude[..., 1:] *= 2
+    mean_amplitude = cycle_amplitude.mean(axis=-2)
+    locking = np.divide(coherent_amplitude, mean_amplitude, out=np.zeros_like(mean_amplitude), where=mean_amplitude > 0)
+    return {
+        "coherent_complex": coherent,
+        "coherent_amplitude": coherent_amplitude,
+        "coherent_phase_radians": np.angle(coherent),
+        "mean_cycle_amplitude": mean_amplitude,
+        "phase_locking": np.clip(locking, 0, 1),
+        "per_cycle_complex": complex_cycle,
+        "per_cycle_amplitude": cycle_amplitude,
+        "per_cycle_phase_radians": np.angle(complex_cycle),
+    }
+
+
+def action_coordinates(values, inverse=False):
+    """Map front/upper/lower to front/rear-symmetric/rear-antisymmetric."""
+    values = np.asarray(values, dtype=float)
+    if values.shape[0] != 3 or not np.isfinite(values).all():
+        raise ValueError("finite front/upper/lower channels are required")
+    scale = 1 / np.sqrt(2)
+    transform = np.array([[1, 0, 0], [0, scale, scale], [0, scale, -scale]])
+    result = (transform.T if inverse else transform) @ values.reshape(3, -1)
+    return result.reshape(values.shape)
+
+
+def bounded_streamwise_shift_correlation(a, b, x, weights, mask, max_shift_D=3.0):
+    """Correlate a(x) with b(x+shift); positive shift places b downstream."""
+    a, b = np.asarray(a, dtype=float), np.asarray(b, dtype=float)
+    x, weights, mask = np.asarray(x), np.asarray(weights), np.asarray(mask, dtype=bool)
+    if a.shape != b.shape or a.shape[-2:] != mask.shape or weights.shape != mask.shape:
+        raise ValueError("field, weight, and mask shapes are inconsistent")
+    if x.size != mask.shape[0] or np.any(np.diff(x) <= 0) or max_shift_D < 0:
+        raise ValueError("invalid x coordinates or shift bound")
+    component_axes = tuple(range(a.ndim - 2))
+    records = []
+    for offset in range(-x.size + 1, x.size):
+        if offset >= 0:
+            ia = np.arange(x.size - offset)
+            ib = ia + offset
+        else:
+            ib = np.arange(x.size + offset)
+            ia = ib - offset
+        displacement = x[ib] - x[ia]
+        shift = float(np.median(displacement))
+        if abs(shift) > max_shift_D + 1e-12:
+            continue
+        # Matching is by one exact integer x-index offset. On grids whose
+        # stored coordinates carry tiny spacing roundoff, report the median
+        # physical displacement and retain its spread instead of inventing
+        # interpolation or rejecting an otherwise exact index pairing.
+        shift_spread = float(np.ptp(displacement))
+        va, vb = a[..., ia, :], b[..., ib, :]
+        finite = np.all(np.isfinite(va), axis=component_axes) & np.all(np.isfinite(vb), axis=component_axes)
+        paired = mask[ia] & mask[ib] & finite
+        pair_weights = np.where(paired, np.sqrt(weights[ia] * weights[ib]), 0)
+        va = np.where(np.isfinite(va), va, 0.0)
+        vb = np.where(np.isfinite(vb), vb, 0.0)
+        dot = np.sum(va * vb, axis=component_axes) if component_axes else va * vb
+        na = np.sum(va**2, axis=component_axes) if component_axes else va**2
+        nb = np.sum(vb**2, axis=component_axes) if component_axes else vb**2
+        denominator = np.sqrt(np.sum(pair_weights * na) * np.sum(pair_weights * nb))
+        correlation = np.sum(pair_weights * dot) / denominator if denominator > 0 else np.nan
+        if np.isfinite(correlation):
+            records.append((shift, float(np.clip(correlation, -1, 1)), int(paired.sum()), offset, shift_spread))
+    if not records:
+        raise ValueError("no finite nonzero shifted correlation is available")
+    best = max(records, key=lambda record: abs(record[1]))
+    return {
+        "best_shift_x_D": best[0], "correlation": best[1],
+        "paired_point_count": best[2], "index_offset": best[3], "shift_spread_x_D": best[4],
+        "shifts_x_D": np.array([record[0] for record in records]),
+        "correlations": np.array([record[1] for record in records]),
+        "paired_point_counts": np.array([record[2] for record in records]),
+    }
@@ -6,7 +6,8 @@ import numpy as np
 from CCD_analysis.acquisition.artifacts import file_sha256, rename_noreplace
 from .artifacts import load_role_artifact
 from .contracts import canonical_json, constant_mean_provenance
-PHASE_SCHEMA_ID="ccd-karman-dynamic-phase-compact/v1"
+PHASE_SCHEMA_ID="ccd-karman-dynamic-phase-compact/v2"
+LEGACY_PHASE_SCHEMA_ID="ccd-karman-dynamic-phase-compact/v1"
 LIMITS={"minimum_complete_cycles":10,"maximum_period_cv":0.05,"maximum_amplitude_cv":0.10,"minimum_cycle_amplitude_fraction":0.25,"maximum_prefix_suffix_period_shift":0.05,"maximum_prefix_suffix_amplitude_shift":0.10,"maximum_double_crossing_fraction":0.0,"maximum_field_sensitivity_relative_rms":0.15}
 def _crossings(s,rising):
 s=np.asarray(s,np.float64); q=(s[:-1]<0)&(s[1:]>=0) if rising else (s[:-1]>0)&(s[1:]<=0); i=np.flatnonzero(q); return i-s[i]/(s[i+1]-s[i])
@@ -46,8 +47,8 @@ def publish_phase_compact(role_path,output,*,role="drl"):
  if not blockers:
   ux,uy=a["ux"][start:stop],a["uy"][start:stop]; actions=d["telemetry"]["effective_applied_action"][start:stop]; mask=a["fluid_mask"]; primary,pacts,counts=_bins(ux,uy,actions,mask,r["phase"],r["cycle_id"],10,0.); ensembles={}; all_counts={}
   for bins,origin,name in ((8,0.,"bins8"),(12,0.,"bins12"),(10,np.pi/10,"bins10_half_shift")):
-    q,_,c=_bins(ux,uy,actions,mask,r["phase"],r["cycle_id"],bins,origin); ensembles[name]=np.mean(q,axis=0,dtype=np.float64).astype(np.float32); all_counts[name]=c
-   n=len(primary); h=n//2; base=np.mean(primary,axis=0,dtype=np.float64); splits={"odd_cycles":np.mean(primary[::2],axis=0,dtype=np.float64),"even_cycles":np.mean(primary[1::2],axis=0,dtype=np.float64),"prefix_cycles":np.mean(primary[:h],axis=0,dtype=np.float64),"suffix_cycles":np.mean(primary[-h:],axis=0,dtype=np.float64)}; sensitivity={k:_curve_distance(base,v,.05 if k=="bins10_half_shift" else 0.) for k,v in {**ensembles,**splits}.items()}; m["field_sensitivity_relative_rms"]=sensitivity
+    q,qa,c=_bins(ux,uy,actions,mask,r["phase"],r["cycle_id"],bins,origin); ensembles[name]=np.mean(q,axis=0,dtype=np.float64).astype(np.float32); ensembles[f"{name}_per_cycle_fields"]=q; ensembles[f"{name}_effective_actions"]=np.mean(qa,axis=0,dtype=np.float64).astype(np.float32); ensembles[f"{name}_per_cycle_effective_actions"]=qa; all_counts[name]=c
+   n=len(primary); h=n//2; base=np.mean(primary,axis=0,dtype=np.float64); splits={"odd_cycles":np.mean(primary[::2],axis=0,dtype=np.float64),"even_cycles":np.mean(primary[1::2],axis=0,dtype=np.float64),"prefix_cycles":np.mean(primary[:h],axis=0,dtype=np.float64),"suffix_cycles":np.mean(primary[-h:],axis=0,dtype=np.float64)}; sensitivity={k:_curve_distance(base,v,.05 if k=="bins10_half_shift" else 0.) for k,v in {**{k:v for k,v in ensembles.items() if not (k.endswith("_effective_actions") or k.endswith("_per_cycle_fields") or k.endswith("_per_cycle_effective_actions"))},**splits}.items()}; m["field_sensitivity_relative_rms"]=sensitivity
   unstable={k:v for k,v in sensitivity.items() if v>LIMITS["maximum_field_sensitivity_relative_rms"]}
   if unstable: blockers.append(f"field sensitivity relative RMS exceeds 0.15: {unstable}")
   if not blockers:
@@ -56,17 +57,21 @@ def publish_phase_compact(role_path,output,*,role="drl"):
     provenance=constant_mean_provenance(drl_manifest_sha256=parent,effective=d["telemetry"]["effective_applied_action"],retained_start=start); (stage/"constant_mean_provenance.json").write_bytes(canonical_json(provenance))
  summary={"schema_id":PHASE_SCHEMA_ID,"complete":True,"gate_passed":not blockers,"blockers":blockers,"quality_limits":LIMITS,"metrics":m,"phase_definition":"linear phase between independent rising zero crossings of retained center sensor uy","primary_bins":10,"sensitivity":{"bin_counts":[8,10,12],"origin":"zero and half-bin for 10","cycle_splits":["odd/even","prefix/suffix"]},"role":role,"source":{"absolute_path":str(role_path.resolve()),"campaign_manifest_sha256":parent,"payload_manifest_sha256":file_sha256(role_path/"payload/manifest.json"),"retained_slice":[start,stop]},"dense_fields_deleted":False,"constant_mean_provenance":provenance}; (stage/"summary.json").write_bytes(canonical_json(summary)); files={p.name:file_sha256(p) for p in stage.iterdir() if p.is_file()}; (stage/"manifest.json").write_bytes(canonical_json({"schema_id":PHASE_SCHEMA_ID,"complete":True,"files":files})); rename_noreplace(stage,output); return load_phase_compact(output)
 except Exception: shutil.rmtree(stage,ignore_errors=True); raise
-def load_phase_compact(path):
+def load_phase_compact(path,*,allow_legacy=False):
 path=Path(path); manifest=json.loads((path/"manifest.json").read_text()); summary=json.loads((path/"summary.json").read_text())
- if manifest.get("schema_id")!=PHASE_SCHEMA_ID or not manifest.get("complete") or summary.get("schema_id")!=PHASE_SCHEMA_ID: raise ValueError("phase compact schema/incomplete")
+ schema=manifest.get("schema_id")
+ if schema!=summary.get("schema_id") or schema not in ({PHASE_SCHEMA_ID,LEGACY_PHASE_SCHEMA_ID} if allow_legacy else {PHASE_SCHEMA_ID}) or not manifest.get("complete"): raise ValueError("phase compact schema/incomplete")
 for n,h in manifest["files"].items():
  if file_sha256(path/n)!=h: raise ValueError("phase compact file hash mismatch")
 role=summary.get("role","drl"); source=Path(summary["source"]["absolute_path"]); load_role_artifact(source,expected_role=role)
 if file_sha256(source/"campaign_manifest.json")!=summary["source"]["campaign_manifest_sha256"] or file_sha256(source/"payload/manifest.json")!=summary["source"]["payload_manifest_sha256"]: raise ValueError("phase compact live source provenance mismatch")
 if summary["gate_passed"]:
  with np.load(path/"compact.npz",allow_pickle=False) as z:
-   required={"fluid_mask","x_D","y_D","cycle_bin_fields","cycle_bin_effective_actions","cycle_bin_counts","cycle_ids","rising_crossings","periods","cycle_amplitudes","ensemble_bins8","ensemble_bins12","ensemble_bins10_half_shift","counts_bins8","counts_bins12","counts_bins10_half_shift"}
+   base={"fluid_mask","x_D","y_D","cycle_bin_fields","cycle_bin_effective_actions","cycle_bin_counts","cycle_ids","rising_crossings","periods","cycle_amplitudes","ensemble_bins8","ensemble_bins12","ensemble_bins10_half_shift","counts_bins8","counts_bins12","counts_bins10_half_shift"}; action_alternatives={"ensemble_bins8_effective_actions","ensemble_bins12_effective_actions","ensemble_bins10_half_shift_effective_actions","ensemble_bins8_per_cycle_fields","ensemble_bins12_per_cycle_fields","ensemble_bins10_half_shift_per_cycle_fields","ensemble_bins8_per_cycle_effective_actions","ensemble_bins12_per_cycle_effective_actions","ensemble_bins10_half_shift_per_cycle_effective_actions"}; required=base|action_alternatives if schema==PHASE_SCHEMA_ID else base
   if set(z.files)!=required or z["cycle_bin_fields"].dtype!=np.float32 or z["cycle_bin_fields"].shape[1:3]!=(10,2) or np.any(z["cycle_bin_counts"]<=0): raise ValueError("phase compact arrays invalid")
+   for name,bins in (("bins8",8),("bins12",12),("bins10_half_shift",10)):
+    action_invalid=schema==PHASE_SCHEMA_ID and (z[f"ensemble_{name}_effective_actions"].shape!=(bins,3) or z[f"ensemble_{name}_per_cycle_fields"].shape!=(len(z["cycle_ids"]),bins,2,int(z["fluid_mask"].sum())) or z[f"ensemble_{name}_per_cycle_effective_actions"].shape!=(len(z["cycle_ids"]),bins,3) or not np.isfinite(z[f"ensemble_{name}_effective_actions"]).all() or not np.isfinite(z[f"ensemble_{name}_per_cycle_fields"]).all() or not np.isfinite(z[f"ensemble_{name}_per_cycle_effective_actions"]).all())
+    if z[f"ensemble_{name}"].shape!=(bins,2,int(z["fluid_mask"].sum())) or not np.isfinite(z[f"ensemble_{name}"]).all() or action_invalid or z[f"counts_{name}"].shape!=(len(z["cycle_ids"]),bins) or np.any(z[f"counts_{name}"]<=0): raise ValueError("phase compact alternative arrays invalid")
 elif "constant_mean_provenance.json" in manifest["files"] or (path/"compact.npz").exists(): raise ValueError("failed gate published forbidden products")
 if role!="drl" and (summary.get("constant_mean_provenance") is not None or "constant_mean_provenance.json" in manifest["files"]): raise ValueError("non-DRL phase artifact contains DRL mean provenance")
 return {"path":path,"summary":summary}
@@ -1,119 +0,0 @@
-"""Exploratory circular phase-domain CCD for DRL minus constant-mean dynamics."""
-from __future__ import annotations
-from dataclasses import dataclass
-from pathlib import Path
-from typing import Any
-import json, shutil, tempfile
-import numpy as np
-from CCD_analysis.acquisition.artifacts import file_sha256, rename_noreplace
-from CCD_analysis.direct_dq.schema import canonical_array_sha256
-from .contracts import canonical_json
-from .dynamic_increment import load_dynamic_increment
-from .phase import load_phase_compact
-
-SCHEMA_ID="ccd-karman-phase-domain-exploratory/v1"
-CHANNELS=("front","upper","lower")
-VARIANTS=(("bins8",8,0.0),("primary",10,0.0),("bins12",12,0.0),("bins10_half_shift",10,0.5))
-HARMONIC_ORDERS=(1,2,3)
-CLAIM_BOUNDARY="exploratory circular phase co-variation only; phase offsets are not time-response lags, causality, mechanism, uncertainty, or observable prediction"
-
-@dataclass(frozen=True)
-class PhaseDomainResult:
- arrays:dict[str,np.ndarray]; config:dict[str,Any]; summary:dict[str,Any]; input_hashes:dict[str,Any]
-
-def _periodic_resample(curve,n,origin_fraction=0.0,order=None):
- x=np.asarray(curve,np.float64); old_n=len(x); old=(np.arange(old_n)+.5)/old_n; new=(np.arange(n)+.5+origin_fraction)/n
- if order is None: order=min((old_n-1)//2,4)
- out=np.broadcast_to(x.mean(0),(n,)+x.shape[1:]).copy(); flat=x.reshape(old_n,-1); outf=out.reshape(n,-1)
- for k in range(1,order+1):
-  c=(2/old_n)*np.sum(flat*np.cos(2*np.pi*k*old)[:,None],axis=0); s=(2/old_n)*np.sum(flat*np.sin(2*np.pi*k*old)[:,None],axis=0)
-  outf+=np.cos(2*np.pi*k*new)[:,None]*c+np.sin(2*np.pi*k*new)[:,None]*s
- return out
-
-def _decompose(u,p,w):
- u=np.asarray(u,np.float64); p=np.asarray(p,np.float64); w=np.asarray(w,np.float64); n=u.shape[1]
- uc=u-u.mean(1,keepdims=True); pc=p-p.mean(1,keepdims=True); roots=np.sqrt(w)
- a=pc@(uc*roots[:,None]).T/(n*np.sqrt(3.0)); left,s,vh=np.linalg.svd(a,full_matrices=False); weighted=vh.T
- for k in range(weighted.shape[1]):
-  pivot=int(np.argmax(np.abs(weighted[:,k])))
-  if weighted[pivot,k]<0: weighted[:,k]*=-1; left[:,k]*=-1
- rank=int(np.sum(s>1e-10*s[0])) if len(s) and s[0]>0 else 0
- return {"cross":a,"left":left,"singular":s,"modes":weighted/roots[:,None],"rank":rank}
-
-def _comparison(primary,other,w):
- r=min(primary["rank"],other["rank"],3); out={"rank":other["rank"],"leading_spectrum":other["singular"][:3].tolist()}
- if r:
-  root=np.sqrt(w)[:,None]; vp=primary["modes"][:,:r]*root; vo=other["modes"][:,:r]*root
-  out["projector_principal_cosines"]=np.linalg.svd(vp.T@vo,compute_uv=False).tolist(); out["minimum_projector_cosine"]=float(min(out["projector_principal_cosines"])); out["left_function_absolute_cosines"]=[float(abs(primary["left"][:,k]@other["left"][:,k])) for k in range(r)]
- else: out.update(projector_principal_cosines=[],minimum_projector_cosine=None,left_function_absolute_cosines=[])
- return out
-
-def _load_inputs(drl_phase,constant_phase,dynamic):
- dp,cp=load_phase_compact(drl_phase),load_phase_compact(constant_phase); di=load_dynamic_increment(dynamic)
- if not dp["summary"]["gate_passed"] or dp["summary"].get("role","drl")!="drl" or not cp["summary"]["gate_passed"] or cp["summary"].get("role","drl")!="constant_mean": raise ValueError("passing DRL and constant_mean phase artifacts required")
- with np.load(Path(drl_phase)/"compact.npz",allow_pickle=False) as dz, np.load(Path(constant_phase)/"compact.npz",allow_pickle=False) as cz, np.load(Path(dynamic)/"arrays.npz",allow_pickle=False) as iz:
-  if not np.array_equal(dz["fluid_mask"],cz["fluid_mask"]): raise ValueError("phase masks differ")
-  common=iz["four_role_fluid_mask"].copy(); selector=common[dz["fluid_mask"]]
-  if int(selector.sum())!=len(iz["quadrature_weights"]): raise ValueError("dynamic weights/common mask mismatch")
-  means={r:iz[f"mean_{r}"].copy() for r in ("drl","constant_mean")}; weights=np.concatenate((iz["quadrature_weights"],iz["quadrature_weights"])).astype(np.float64); fields={}
-  for name,_,_ in VARIANTS:
-   kd="cycle_bin_fields" if name=="primary" else f"ensemble_{name}"; dd=np.mean(dz[kd],axis=0,dtype=np.float64) if name=="primary" else dz[kd].astype(np.float64); cc=np.mean(cz[kd],axis=0,dtype=np.float64) if name=="primary" else cz[kd].astype(np.float64)
-   fields[name]=(dd[:,:,selector]-means["drl"][None])-(cc[:,:,selector]-means["constant_mean"][None])
-  actions=np.mean(dz["cycle_bin_effective_actions"],axis=0,dtype=np.float64)
- return dp,cp,di,fields,actions,weights
-
-def decompose_phase_domain(drl_phase,constant_phase,dynamic_increment)->PhaseDomainResult:
- dp,cp,di,fields,actions,w=_load_inputs(drl_phase,constant_phase,dynamic_increment); results={}; arrays={"coordinate_weights":w,"primary_phase_action_curve":actions}
- for name,n,origin in VARIANTS:
-  u=fields[name].transpose(1,2,0).reshape(len(w),n); p=_periodic_resample(actions,n,origin_fraction=origin).T; r=_decompose(u,p,w); results[name]=r
-  for key in ("cross","left","singular","modes"): arrays[f"{name}_{key}"]=r[key]
- primary=results["primary"]; comparisons={name:_comparison(primary,r,w) for name,r in results.items() if name!="primary"}; harmonic={}
- for order in HARMONIC_ORDERS:
-  u0=fields["primary"].transpose(1,2,0).reshape(len(w),10).T; r=_decompose(_periodic_resample(u0,10,order=order).T,_periodic_resample(actions,10,order=order).T,w); harmonic[str(order)]=_comparison(primary,r,w)
-  for key in ("cross","left","singular","modes"): arrays[f"harmonic_{order}_{key}"]=r[key]
- offsets=[]
- for off in range(10):
-  r=_decompose(fields["primary"].transpose(1,2,0).reshape(len(w),10),np.roll((actions-actions.mean(0)).T,off,axis=1),w); signed=off if off<=5 else off-10
-  offsets.append({"phase_offset_bins":signed,"phase_offset_radians":float(2*np.pi*signed/10),"leading_spectrum":r["singular"][:3].tolist(),"rank":r["rank"]})
- checks=list(comparisons.values())+list(harmonic.values()); stable=[x["minimum_projector_cosine"] for x in checks if x.get("minimum_projector_cosine") is not None]; left=[min(x["left_function_absolute_cosines"]) for x in checks if x["left_function_absolute_cosines"]]; rank_stable=all(x["rank"]==primary["rank"] for x in checks)
- decision="PASS_EXPLORATORY" if primary["rank"]>0 and rank_stable and min(stable,default=0)>=0.9 and min(left,default=0)>=0.8 else "DOWNGRADE"
- config={"schema_id":SCHEMA_ID,"field_estimand":"separately_centered_phase_coherent_difference_delta_q_prime_phase(phi)=(q_DRL(phi)-mean_q_DRL)-(q_constant_mean(phi)-mean_q_constant_mean)","observable_estimand":"DRL phase-conditioned three-channel effective_applied_action fluctuation","operator":"A=P(W^(1/2)U)^T/(N*sqrt(3Q)); Q=1","primary_bins":10,"sensitivity_bins":[8,12],"half_bin_origin":True,"harmonic_orders":list(HARMONIC_ORDERS),"action_grid_sensitivity":"periodic Fourier interpolation of immutable primary 10-bin DRL action curve, maximum order 4; harmonic tests truncate both U and P","phase_offset_semantics":"circular phase offsets only; not time-response lags","center_snapshots":True,"center_observables":True,"standardization":False,"whitening":False,"pod":False,"claim_boundary":CLAIM_BOUNDARY}
- summary={"schema_id":SCHEMA_ID,"decision":decision,"N":10,"Q":1,"M":len(w),"numerical_rank":primary["rank"],"primary_singular_values":primary["singular"].tolist(),"primary_squared_singular_values":(primary["singular"]**2).tolist(),"primary_left_functions":primary["left"].tolist(),"sensitivity":{"bin_and_origin":comparisons,"harmonic_order":harmonic},"circular_phase_offsets":offsets,"stability_gate":{"rank_stable":rank_stable,"minimum_projector_principal_cosine":min(stable,default=None),"minimum_matched_left_function_absolute_cosine":min(left,default=None),"thresholds":{"projector_cosine":0.9,"left_function_absolute_cosine":0.8}},"spectrum_label":"cross-correlation strength; not field energy or explained variance","claim_boundary":CLAIM_BOUNDARY,"scientific_contract_explicit":True}
- parents={"drl_phase":{"path":str(Path(drl_phase).resolve()),"manifest_sha256":file_sha256(Path(drl_phase)/"manifest.json")},"constant_mean_phase":{"path":str(Path(constant_phase).resolve()),"manifest_sha256":file_sha256(Path(constant_phase)/"manifest.json")},"dynamic_increment":{"path":str(Path(dynamic_increment).resolve()),"manifest_sha256":file_sha256(Path(dynamic_increment)/"manifest.json")}}; hashes={"parents":parents,"canonical_arrays":{k:canonical_array_sha256(v) for k,v in arrays.items()}}
- validate_phase_domain(arrays,config,summary,hashes); return PhaseDomainResult(arrays,config,summary,hashes)
-
-def validate_phase_domain(arrays,config,summary,hashes):
- if config.get("schema_id")!=SCHEMA_ID or summary.get("schema_id")!=SCHEMA_ID or config.get("field_estimand")!="separately_centered_phase_coherent_difference_delta_q_prime_phase(phi)=(q_DRL(phi)-mean_q_DRL)-(q_constant_mean(phi)-mean_q_constant_mean)" or config.get("operator")!="A=P(W^(1/2)U)^T/(N*sqrt(3Q)); Q=1": raise ValueError("phase-domain scientific contract contradicted")
- if config.get("phase_offset_semantics")!="circular phase offsets only; not time-response lags" or any(config.get(k) is not False for k in ("standardization","whitening","pod")): raise ValueError("phase-domain preprocessing/offset contract contradicted")
- d={k:np.asarray(v) for k,v in arrays.items()}; m=len(d["coordinate_weights"])
- if summary.get("N")!=10 or summary.get("Q")!=1 or summary.get("M")!=m or d["primary_cross"].shape!=(3,m) or d["primary_left"].shape!=(3,3): raise ValueError("phase-domain dimensions invalid")
- if set(hashes)!={"parents","canonical_arrays"} or set(hashes["canonical_arrays"])!=set(d) or any(hashes["canonical_arrays"][k]!=canonical_array_sha256(v) for k,v in d.items()): raise ValueError("phase-domain hashes invalid")
- canonical_json(config); canonical_json(summary); canonical_json(hashes)
-
-class PhaseDomainTransaction:
- def __init__(self,destination): self.destination=Path(destination); self.stage=None
- def __enter__(self):
-  if self.destination.exists(): raise FileExistsError(self.destination)
-  self.destination.parent.mkdir(parents=True,exist_ok=True); self.stage=Path(tempfile.mkdtemp(prefix=f".{self.destination.name}.partial-",dir=self.destination.parent)); return self
- def write(self,result):
-  validate_phase_domain(result.arrays,result.config,result.summary,result.input_hashes); np.savez_compressed(self.stage/"arrays.npz",**result.arrays)
-  for n,v in (("config.json",result.config),("summary.json",result.summary),("input_hashes.json",result.input_hashes)): (self.stage/n).write_bytes(canonical_json(v))
-  files={p.name:file_sha256(p) for p in self.stage.iterdir()}; (self.stage/"manifest.json").write_bytes(canonical_json({"schema_id":SCHEMA_ID,"complete":True,"files":files}))
- def publish(self):
-  load_phase_domain_result(self.stage,recompute=False); rename_noreplace(self.stage,self.destination); self.stage=None; load_phase_domain_result(self.destination,recompute=True); return self.destination
- def __exit__(self,*args):
-  if self.stage is not None: shutil.rmtree(self.stage,ignore_errors=True)
-
-def load_phase_domain_result(path,recompute=True):
- root=Path(path); manifest=json.loads((root/"manifest.json").read_text())
- if manifest.get("schema_id")!=SCHEMA_ID or not manifest.get("complete") or set(manifest.get("files",{}))!={"arrays.npz","config.json","summary.json","input_hashes.json"}: raise ValueError("phase-domain manifest invalid")
- for n,h in manifest["files"].items():
-  if file_sha256(root/n)!=h: raise ValueError("phase-domain file hash mismatch")
- with np.load(root/"arrays.npz",allow_pickle=False) as z: arrays={k:z[k].copy() for k in z.files}
- config=json.loads((root/"config.json").read_text()); summary=json.loads((root/"summary.json").read_text()); hashes=json.loads((root/"input_hashes.json").read_text()); validate_phase_domain(arrays,config,summary,hashes)
- for parent in hashes["parents"].values():
-  if file_sha256(Path(parent["path"])/"manifest.json")!=parent["manifest_sha256"]: raise ValueError("phase-domain live parent identity changed")
- if recompute:
-  p=hashes["parents"]; fresh=decompose_phase_domain(p["drl_phase"]["path"],p["constant_mean_phase"]["path"],p["dynamic_increment"]["path"])
-  for k in arrays: np.testing.assert_allclose(arrays[k],fresh.arrays[k],rtol=2e-11,atol=2e-12)
- return {"arrays":arrays,"config":config,"summary":summary,"input_hashes":hashes,"manifest":manifest,"provenance_validation":"VERIFIED: live compact and dynamic-increment parents reread; essential decomposition recomputed" if recompute else "VERIFIED_HASHES_AND_LIVE_PARENTS"}
@@ -0,0 +1,72 @@
+"""Full-fit weighted POD baseline for the exact corrected CCD residual."""
+from __future__ import annotations
+
+import numpy as np
+
+from .cycle_template_ccd import _load, _matrix
+
+
+def reconstruct_residual_from_parents(drl_phase, constant_phase, roi_mean, admitted_empirical_mean=None):
+    """Rebuild raw U with the exact CCD loader/matrix path, then center it."""
+    loaded = _load(drl_phase, constant_phase, roi_mean)
+    physical_mean = loaded[3]
+    drl_fields, drl_actions, constant_fields = loaded[4], loaded[5], loaded[8]
+    mean_drl, mean_constant, selector, roi_weights = loaded[11:15]
+    raw, _actions, weights, _cycles, _bins, _mean = _matrix(
+        drl_fields, drl_actions, constant_fields, mean_drl, mean_constant,
+        selector, roi_weights, physical_mean,
+    )
+    empirical_mean = raw.mean(axis=1, dtype=np.float64)
+    if admitted_empirical_mean is not None:
+        admitted = np.asarray(admitted_empirical_mean, dtype=np.float64)
+        if admitted.shape != empirical_mean.shape or not np.array_equal(empirical_mean, admitted):
+            raise ValueError("reconstructed U does not bit-match the admitted CCD mean")
+        empirical_mean = admitted
+    if raw.shape != (160400, 190) or weights.shape != (160400,):
+        raise ValueError(f"corrected residual dimensions changed: {raw.shape}, {weights.shape}")
+    return {"raw": raw, "centered": raw - empirical_mean[:, None], "weights": weights, "empirical_mean": empirical_mean}
+
+
+def fit_weighted_pod(snapshots, weights, rank=3, center=True):
+    """Fit a descriptive full-sample POD in weighted field coordinates."""
+    snapshots, weights = np.asarray(snapshots, dtype=float), np.asarray(weights, dtype=float)
+    if snapshots.ndim != 2 or weights.shape != (snapshots.shape[0],):
+        raise ValueError("POD snapshots/weights have inconsistent dimensions")
+    if not np.isfinite(snapshots).all() or np.any(~np.isfinite(weights)) or np.any(weights <= 0):
+        raise ValueError("finite snapshots and positive finite weights are required")
+    if not 0 <= rank <= min(snapshots.shape):
+        raise ValueError("invalid POD rank")
+    mean = snapshots.mean(axis=1, keepdims=True) if center else np.zeros((snapshots.shape[0], 1))
+    centered = snapshots - mean
+    weighted_modes, singular_values, right_vectors = np.linalg.svd(np.sqrt(weights)[:, None] * centered, full_matrices=False)
+    modes = weighted_modes[:, :rank] / np.sqrt(weights)[:, None]
+    coefficients = singular_values[:rank, None] * right_vectors[:rank]
+    energy, total = singular_values**2, float(np.sum(singular_values**2))
+    retained = np.r_[0, np.cumsum(energy[:rank])]
+    residual = np.zeros(rank + 1) if total == 0 else np.maximum(1 - retained / total, 0)
+    if not np.allclose(modes.T @ (weights[:, None] * modes), np.eye(rank), rtol=1e-10, atol=1e-12):
+        raise ValueError("computed POD modes are not W-orthonormal")
+    return {"mean": mean[:, 0], "modes": modes, "coefficients": coefficients, "singular_values": singular_values, "weighted_residual_fraction": residual, "label": "full-fit weighted POD field-energy spectrum; descriptive only"}
+
+
+def weighted_principal_cosines(first, second, weights, rank=3):
+    """Return principal cosines after verifying both bases are W-orthonormal."""
+    first, second, weights = np.asarray(first), np.asarray(second), np.asarray(weights)
+    selected = min(rank, first.shape[1], second.shape[1])
+    identity = np.eye(selected)
+    for name, basis in (("first", first[:, :selected]), ("second", second[:, :selected])):
+        gram = basis.T @ (weights[:, None] * basis)
+        if not np.allclose(gram, identity, rtol=1e-8, atol=1e-10):
+            raise ValueError(f"{name} basis is not W-orthonormal")
+    overlap = first[:, :selected].T @ (weights[:, None] * second[:, :selected])
+    return np.clip(np.linalg.svd(overlap, compute_uv=False), 0, 1)
+
+
+def align_basis_to_reference(reference, candidate, weights):
+    """Orthogonally align a candidate basis to a reference for visual comparison."""
+    reference, candidate, weights = np.asarray(reference), np.asarray(candidate), np.asarray(weights)
+    if reference.shape != candidate.shape or weights.shape != (reference.shape[0],):
+        raise ValueError("basis-alignment dimensions are inconsistent")
+    left, _singular, right_t = np.linalg.svd(reference.T @ (weights[:, None] * candidate), full_matrices=False)
+    rotation = right_t.T @ left.T
+    return candidate @ rotation, rotation
@@ -1,6 +1,6 @@
-"""Deterministic artifact-only publication for the Karman dynamic campaign."""
+"""Concise artifact-only publication for ROI means and cycle-template CCD."""
 from __future__ import annotations
-import json, os, shutil, uuid
+import json, shutil, tempfile
 from pathlib import Path
 import matplotlib
 matplotlib.use("Agg")
@@ -9,64 +9,89 @@ import numpy as np
 from CCD_analysis.acquisition.artifacts import file_sha256, rename_noreplace
 from .contracts import canonical_json
 from .dynamic_increment import load_dynamic_increment
-from .temporal_ccd import load_temporal_result, CHANNELS
-from .phase_domain_ccd import load_phase_domain_result
-SCHEMA_ID="ccd-karman-dynamic-publication/v1"
-STEMS=("01_four_role_mean_performance","02_drl_constant_phase_difference","03_temporal_ccd_spectrum_sensitivity","04_temporal_ccd_leading_modes","05_temporal_ccd_left_lag_functions","06_phase_domain_ccd_downgrade")
-def _save(fig,root,stem):
- names=[]
- for ext in ("png","pdf"):
-  q=root/f"{stem}.{ext}"; fig.savefig(q,dpi=300 if ext=="png" else None,bbox_inches="tight",metadata={"Creator":"CCD_analysis.karman_dynamic.publication"}); names.append(q.name)
- plt.close(fig); return names
-def _wrms(field,w):
- a=np.asarray(field,float); return float(np.sqrt(np.sum(w*np.sum(a*a,axis=0))/np.sum(w)))
-def _components(v,mask):
- n=int(mask.sum()); out=[]
- for a in (v[:n],v[n:]):
-  f=np.full(mask.shape,np.nan); f[mask]=a; out.append(f)
- return out
-def _panel(ax,f,x,y,mask,lim,title):
- ax.pcolormesh(x,y,f.T,shading="nearest",cmap="RdBu_r",vmin=-lim,vmax=lim,rasterized=True); solid=np.ma.masked_where(mask,np.ones(mask.shape)); ax.pcolormesh(x,y,solid.T,shading="nearest",cmap="Greys",vmin=0,vmax=1); ax.set(title=title,xlabel="x/D",ylabel="y/D"); ax.set_aspect("equal")
-def _source(root): return {"path":str(root.resolve()),"manifest_sha256":file_sha256(root/"manifest.json")}
-def publish_dynamic_figures(dynamic_root,temporal_root,phase_domain_root,output):
- dynamic_root,temporal_root,phase_domain_root=map(Path,(dynamic_root,temporal_root,phase_domain_root)); destination=Path(output)
+from .cycle_template_ccd import CHANNELS, load_cycle_template_result
+SCHEMA_ID="ccd-karman-corrected-publication/v3"
+FIGURE_STEMS=("01_roi_mean_performance","02_roi_mean_target_error_fields","03_roi_mean_corrections","04_cycle_template_spectrum_stability","05_cycle_template_leading_modes","06_cycle_template_left_action_vectors")
+FIGURE_FILES=tuple(f"{stem}.{ext}" for stem in FIGURE_STEMS for ext in ("png","pdf"))
+
+def _save(fig,stage,stem):
+ for ext in ("png","pdf"): fig.savefig(stage/f"{stem}.{ext}",dpi=300 if ext=="png" else None,bbox_inches="tight",metadata={"Creator":"CCD_analysis.karman_dynamic.publication"})
+ plt.close(fig)
+
+def _roi_field(values,roi_mask):
+ out=np.full(roi_mask.shape,np.nan); out[roi_mask]=values; return out
+
+def _panel(ax,x,y,field,title,lim=None,cmap="RdBu_r"):
+ kw={} if lim is None else {"vmin":-lim,"vmax":lim}; ax.pcolormesh(x,y,field.T,shading="nearest",cmap=cmap,rasterized=True,**kw); ax.set(xlim=(34,54),ylim=(-5,5),xlabel="x/D",ylabel="y/D",title=title); ax.set_aspect("equal")
+
+def _stability(summary):
+ sensitivities=summary["sensitivities"]; categories={
+  "drl_leave_one_cycle_out":sensitivities["drl_leave_one_cycle_out"],
+  "constant_template_leave_one_cycle_out":sensitivities["constant_template_leave_one_cycle_out"],
+  "drl_splits":[{"check":name,**value} for name,value in sensitivities["drl_splits"].items()],
+  "constant_template_splits":[{"check":name,**value} for name,value in sensitivities["constant_template_splits"].items()],
+  "bin_and_origin":[{"check":name,**value} for name,value in sensitivities["bin_and_origin"].items()],
+ }
+ primary=np.asarray(summary["primary_singular_values"][:3],dtype=float); disclosed={}; global_items=[]
+ for category,items in categories.items():
+  records=[]
+  for item in items:
+   spectrum=np.asarray(item["leading_spectrum"],dtype=float); relative=(spectrum/primary[:len(spectrum)]).tolist(); record={"label":str(item.get("check",item.get("omitted_cycle"))),"rank":item["rank"],"minimum_projector_cosine":item.get("minimum_projector_cosine"),"leading_spectrum":spectrum.tolist(),"singular_value_ratio_to_primary":relative,"maximum_absolute_relative_singular_shift":float(np.max(np.abs(np.asarray(relative)-1))) if relative else None}
+   records.append(record); global_items.append((category,record))
+  disclosed[category]=records
+ valid=[(category,record) for category,record in global_items if record["minimum_projector_cosine"] is not None]
+ least=min(valid,key=lambda item:item[1]["minimum_projector_cosine"]) if valid else None; shifts=[(category,record) for category,record in global_items if record["maximum_absolute_relative_singular_shift"] is not None]; worst=max(shifts,key=lambda item:item[1]["maximum_absolute_relative_singular_shift"]) if shifts else None
+ summary_out={"scope":"all declared DRL/template cycle-block, split, and persisted bin/origin checks","checks":disclosed,"global_least_favorable_projector_check":None if least is None else {"category":least[0],**least[1]},"global_largest_singular_shift_check":None if worst is None else {"category":worst[0],**worst[1]}}
+ plot={category:[record["minimum_projector_cosine"] for record in records if record["minimum_projector_cosine"] is not None] for category,records in disclosed.items()}
+ return plot,summary_out
+
+def publish_dynamic_figures(roi_mean_root,cycle_template_root,output):
+ roi_mean_root,cycle_template_root,destination=map(Path,(roi_mean_root,cycle_template_root,output)); destination.parent.mkdir(parents=True,exist_ok=True)
 if destination.exists(): raise FileExistsError(destination)
- dynamic=load_dynamic_increment(dynamic_root); temporal=load_temporal_result(temporal_root,recompute=True); phase=load_phase_domain_result(phase_domain_root,recompute=True)
- if phase["summary"]["decision"]!="DOWNGRADE": raise ValueError("frozen phase-domain decision must be DOWNGRADE")
- partial=destination.with_name(f".{destination.name}.partial.{os.getpid()}.{uuid.uuid4().hex}"); partial.mkdir(parents=True)
+ mean=load_dynamic_increment(roi_mean_root); ccd=load_cycle_template_result(cycle_template_root,recompute=True); stage=Path(tempfile.mkdtemp(prefix=f".{destination.name}.partial-",dir=destination.parent))
 try:
-  with np.load(dynamic_root/"arrays.npz",allow_pickle=False) as z: inc={k:z[k].copy() for k in z.files}
-  ds=dynamic["summary"]; ta,ts=temporal["arrays"],temporal["summary"]; ps=phase["summary"]; files=[]
-  roles=("target","zero","constant_mean","drl"); errors=[ds["statistics"][r]["mean_target_error_weighted_vector_rms"] for r in roles]; b=ds["benefits"]
-  fig,ax=plt.subplots(figsize=(6.4,3.8),layout="constrained"); bars=ax.bar(("Target","Zero","Constant mean","DRL"),errors,color=("#59A14F","#777777","#4C78A8","#E45756")); ax.bar_label(bars,fmt="%.4f",padding=3); ax.annotate(f"-{b['zero_to_constant_overall_mean_control_benefit_target_error_reduction']:.4f}",xy=(2,errors[2]),xytext=(1,errors[1]+.012),arrowprops={"arrowstyle":"->"},ha="center"); ax.annotate(f"-{b['constant_to_drl_dynamic_increment_target_error_reduction']:.4f}",xy=(3,errors[3]),xytext=(2,errors[2]+.012),arrowprops={"arrowstyle":"->"},ha="center"); ax.set(ylabel="Mean-field target error (weighted vector RMS)",ylim=(0,max(errors)*1.25),title="Mean performance decomposition"); ax.grid(axis="y",alpha=.25); files+=_save(fig,partial,STEMS[0])
-  bins=np.arange(1,11); cp=ds["phase_target_error_metrics"]["constant_mean"]["target_error_by_bin_weighted_vector_rms"]; dp=ds["phase_target_error_metrics"]["drl"]["target_error_by_bin_weighted_vector_rms"]; centered=inc["phase_difference_centered_drl_minus_constant"]; w=inc["quadrature_weights"]; cr=[_wrms(centered[i],w) for i in range(10)]
-  fig,axs=plt.subplots(1,2,figsize=(10,3.7),layout="constrained"); axs[0].plot(bins,cp,"o-",label="Constant mean"); axs[0].plot(bins,dp,"o-",label="DRL"); axs[0].set(ylabel="Phase target error (weighted vector RMS)",title="Independent 10-bin phase means"); axs[0].legend(); axs[1].plot(bins,cr,"o-",color="#7A5195"); axs[1].set(ylabel="Centered DRL-constant difference (weighted vector RMS)",title="Phase-coherent increment"); [a.set(xlabel="Phase bin",xticks=bins) for a in axs]; [a.grid(alpha=.25) for a in axs]; fig.suptitle("DRL versus constant mean; zero omitted because its phase gate failed"); files+=_save(fig,partial,STEMS[1])
-  sigma=ta["primary_singular_values"]; sens=ts["sensitivity"]; fig,axs=plt.subplots(1,2,figsize=(10,3.7),layout="constrained"); axs[0].semilogy(np.arange(1,len(sigma)+1),sigma,"o-",ms=3); axs[0].set(xlabel="Mode",ylabel="Cross-correlation strength",title="Negative-lag CCD spectrum"); spectra=[sigma[:3]]+[np.asarray(v["common_support_vs_primary"]["leading_spectrum"]) for v in sens]; labels=["-17...0"]+[f"{v['lags'][0]}...0" for v in sens]
-  for vals,label in zip(spectra,labels): axs[1].plot((1,2,3),vals/sigma[:3],"o-",label=label)
-  axs[1].axhline(1,color="black",lw=.7); axs[1].set(xlabel="Leading mode",ylabel="Strength / primary strength",xticks=(1,2,3),title="Common-support sensitivity (N=21)"); axs[1].legend(title="Lag window"); [a.grid(alpha=.25) for a in axs]; fig.suptitle("Closed-loop temporal co-variation; not causality or response time"); files+=_save(fig,partial,STEMS[2])
-  mask,x,y=ta["fluid_mask"],ta["x_D"],ta["y_D"]; comps=[_components(ta["primary_physical_modes"][:,j],mask) for j in range(3)]; fig,axs=plt.subplots(3,2,figsize=(12,7),sharex=True,sharey=True,layout="constrained")
-  for j,c in enumerate(comps):
-   lim=float(np.percentile(np.abs(np.concatenate([q[np.isfinite(q)] for q in c])),99)) or 1
-   for ax,q,label in zip(axs[j],c,("ux","uy")): _panel(ax,q,x,y,mask,lim,f"Mode {j+1} {label}")
-  fig.suptitle("Leading temporal CCD physical modes (full grid; per-mode symmetric scale)"); files+=_save(fig,partial,STEMS[3])
-  left=ta["primary_left_functions"].reshape(3,len(ta["primary_lags"]),-1); fig,axs=plt.subplots(1,3,figsize=(11,3.4),sharey=True,layout="constrained")
-  for j,ax in enumerate(axs):
-   for ch,name in enumerate(CHANNELS): ax.plot(ta["primary_lags"],left[ch,:,j],"o-",ms=3,label=name)
-   peak=ts["primary_left_lag_metrics"]["modes"][j]["peak_lag_boundaries"]; ax.axvline(peak,color="black",ls="--",lw=.8); ax.set(title=f"Mode {j+1}; energy peak {peak}",xlabel="tau / 800 lattice steps"); ax.grid(alpha=.25)
-  axs[0].set_ylabel("Left lag-function component"); axs[-1].legend(); fig.suptitle("tau < 0 means action precedes field; lag structure is descriptive"); files+=_save(fig,partial,STEMS[4])
-  primary=np.asarray(ps["primary_singular_values"][:3]); variants={"10 bins":primary,"8 bins":ps["sensitivity"]["bin_and_origin"]["bins8"]["leading_spectrum"],"12 bins":ps["sensitivity"]["bin_and_origin"]["bins12"]["leading_spectrum"],"half-bin":ps["sensitivity"]["bin_and_origin"]["bins10_half_shift"]["leading_spectrum"],"harmonic 1":ps["sensitivity"]["harmonic_order"]["1"]["leading_spectrum"]}; fig,axs=plt.subplots(1,2,figsize=(10,3.7),layout="constrained"); axs[0].bar((1,2,3),primary,color="#F58518"); axs[0].set(xlabel="Mode",ylabel="Cross-correlation strength",title="Primary 10-bin spectrum")
-  for name,vals in variants.items(): axs[1].plot((1,2,3),np.asarray(vals)/primary,"o-",label=name)
-  axs[1].set(xlabel="Mode",ylabel="Strength / primary strength",xticks=(1,2,3),title="Resolution and harmonic sensitivity"); axs[1].legend(fontsize=8); [a.grid(alpha=.25) for a in axs]; fig.suptitle("PHASE-DOMAIN CCD - DOWNGRADE: first-harmonic rank is 2, not 3"); files+=_save(fig,partial,STEMS[5])
-  report={"schema_id":SCHEMA_ID,"artifact_only":True,"sources":{"dynamic_increment":_source(dynamic_root),"temporal_ccd":_source(temporal_root),"phase_domain_ccd":_source(phase_domain_root)},"source_reload":{"dynamic_increment":"VERIFIED live four-role and phase parents","temporal_ccd":temporal["provenance_validation"],"phase_domain_ccd":phase["provenance_validation"]},"figure_files":files,"mean_target_errors":dict(zip(roles,errors)),"mean_target_error_reductions":b,"phase_target_error_cycle_means":{r:ds["phase_target_error_metrics"][r]["target_error_cycle_mean_weighted_vector_rms"] for r in ("constant_mean","drl")},"centered_phase_difference_weighted_vector_rms_by_bin":cr,"zero_phase_output":"PROHIBITED: zero failed the phase gate; no zero phase result is plotted or claimed","temporal_leading_singular_values":sigma[:3].tolist(),"temporal_common_support_sensitivity":sens,"phase_domain_decision":ps["decision"],"phase_domain_leading_singular_values":primary.tolist(),"phase_domain_stability_gate":ps["stability_gate"],"claim_boundary":"Independent trajectory statistics and closed-loop co-variation only; no pointwise counterfactual, causal, mechanism, response-time, uncertainty, CCD>POD, or explained-variance claim."}; (partial/"RESULTS.json").write_bytes(canonical_json(report)); fraction=100*b["constant_to_drl_dynamic_increment_target_error_reduction"]/b["zero_to_drl_total_target_error_reduction"]
-  lines=["# Karman dynamic-increment publication results","","All outputs were generated from fresh verified reloads of immutable artifacts; no CFD or CUDA was used.","","## Mean performance","",f"Mean target error decreases from zero `{errors[1]:.7f}` to constant mean `{errors[2]:.7f}` (reduction `{b['zero_to_constant_overall_mean_control_benefit_target_error_reduction']:.7f}`), then to DRL `{errors[3]:.7f}` (additional reduction `{b['constant_to_drl_dynamic_increment_target_error_reduction']:.7f}`). The latter is about {fraction:.1f}% of the total zero-to-DRL reduction.","","## Phase and CCD boundary","",f"The cycle-mean 10-bin target error is `{report['phase_target_error_cycle_means']['constant_mean']:.7f}` for constant mean and `{report['phase_target_error_cycle_means']['drl']:.7f}` for DRL. Zero is absent because it failed the phase gate.","",f"Temporal leading strengths are `{report['temporal_leading_singular_values']}`. Common-support comparisons preserve the leading-three subspace; native-support changes are not timing evidence.","",f"Phase-domain CCD is **{ps['decision']}**. Primary rank is 3, but first-harmonic rank is `{ps['sensitivity']['harmonic_order']['1']['rank']}`; only exploratory circular co-variation is supported.","",f"Claim boundary: {report['claim_boundary']}",""]; (partial/"RESULTS.md").write_text("\n".join(lines)); hashes={q.name:file_sha256(q) for q in partial.iterdir() if q.is_file()}; (partial/"manifest.json").write_bytes(canonical_json({"schema_id":SCHEMA_ID,"complete":True,"files":hashes})); rename_noreplace(partial,destination); return destination
- except Exception: shutil.rmtree(partial,ignore_errors=True); raise
+  s,cs=mean["summary"],ccd["summary"]; roles=("target","zero","constant_mean","drl"); errors=[s["statistics"][r]["mean_target_error_roi_weighted_vector_rms"] for r in roles]
+  with np.load(roi_mean_root/"arrays.npz",allow_pickle=False) as z: a={k:z[k].copy() for k in z.files}
+  roi,x,y=a["roi_fluid_mask"],a["x_D"],a["y_D"]; selector=a["roi_selector_on_common_mask"]
+  fig,ax=plt.subplots(figsize=(6.5,3.8),layout="constrained"); bars=ax.bar(("Target","Zero","Constant mean","DRL"),errors,color=("#59A14F","#777777","#4C78A8","#E45756")); ax.bar_label(bars,fmt="%.5f"); ax.set(ylabel="ROI weighted vector RMS",title="Mean target error: 34 <= x/D <= 54, |y/D| <= 5"); ax.grid(axis="y",alpha=.25); _save(fig,stage,FIGURE_STEMS[0])
+  target=a["mean_target"][:,selector]; role_errors={r:a[f"mean_{r}"][:,selector]-target for r in ("zero","constant_mean","drl")}; limits=[max(float(np.max(np.abs(v[c]))) for v in role_errors.values()) or 1 for c in range(2)]
+  fig,axs=plt.subplots(3,2,figsize=(10,8),sharex=True,sharey=True,layout="constrained")
+  for row,r in enumerate(("zero","constant_mean","drl")):
+   for c,name in enumerate(("ux","uy")): _panel(axs[row,c],x,y,_roi_field(role_errors[r][c],roi),f"{r} - target: {name}",limits[c])
+  fig.suptitle("ROI mean target-error fields; shared scale within each component"); _save(fig,stage,FIGURE_STEMS[1])
+  incs=(a["mean_increment_zero_to_constant"][:,selector],a["mean_increment_constant_to_drl"][:,selector]); component_limits=[max(float(np.max(np.abs(v[c]))) for v in incs) or 1 for c in range(2)]; mag_limit=max(float(np.max(np.sqrt(np.sum(v.astype(float)**2,axis=0)))) for v in incs) or 1
+  fig,axs=plt.subplots(2,3,figsize=(13,6),sharex=True,sharey=True,layout="constrained")
+  for row,(inc,label) in enumerate(zip(incs,("constant - zero","DRL - constant"))):
+   for c,name in enumerate(("ux","uy")): _panel(axs[row,c],x,y,_roi_field(inc[c],roi),f"{label}: {name}",component_limits[c])
+   _panel(axs[row,2],x,y,_roi_field(np.sqrt(np.sum(inc.astype(float)**2,axis=0)),roi),f"{label}: magnitude",None,"viridis"); axs[row,2].collections[0].set_clim(0,mag_limit)
+  fig.suptitle("ROI mean corrections; shared scales across estimands"); _save(fig,stage,FIGURE_STEMS[2])
+  sigma=ccd["arrays"]["primary_singular_values"]; stability_values,stability_summary=_stability(cs); fig,axs=plt.subplots(1,2,figsize=(10,3.8),layout="constrained"); axs[0].bar(range(1,len(sigma)+1),sigma); axs[0].set(xlabel="Mode",ylabel="Cross-correlation strength",title="N=190, Q=1 spectrum")
+  for name,values in stability_values.items(): axs[1].plot(range(len(values)),values,"o",ms=3,label=name)
+  axs[1].set(xlabel="Check within category",ylabel="Rank-3 projector minimum cosine",title="All declared stability checks"); axs[1].legend(); [q.grid(alpha=.25) for q in axs]; _save(fig,stage,FIGURE_STEMS[3])
+  modes=ccd["arrays"]["primary_physical_modes"]; n=int(roi.sum()); mode_count=min(3,modes.shape[1]); fig,axs=plt.subplots(mode_count,2,figsize=(10,2.7*mode_count),sharex=True,sharey=True,layout="constrained")
+  if mode_count==1: axs=np.asarray([axs])
+  for k in range(mode_count):
+   limit=max(float(np.max(np.abs(modes[:n,k]))),float(np.max(np.abs(modes[n:,k])))) or 1
+   for c,name in enumerate(("ux","uy")): _panel(axs[k,c],x,y,_roi_field(modes[c*n:(c+1)*n,k],roi),f"Mode {k+1}: {name}",limit)
+  fig.suptitle("Primary-basis representation of stable rank-3 ROI subspace; individual-vector stability not claimed"); _save(fig,stage,FIGURE_STEMS[4])
+  left=ccd["arrays"]["primary_left_functions"]; fig,ax=plt.subplots(figsize=(7,4),layout="constrained"); width=.22; positions=np.arange(len(CHANNELS))
+  for k in range(min(3,left.shape[1])): ax.bar(positions+(k-1)*width,left[:,k],width,label=f"Mode {k+1}")
+  ax.set(xticks=positions,xticklabels=CHANNELS,ylabel="Unit left-vector component",title="Primary-basis Q=1 action vectors; individual-vector stability not claimed"); ax.axhline(0,color="black",lw=.7); ax.legend(); ax.grid(axis="y",alpha=.25); _save(fig,stage,FIGURE_STEMS[5])
+  full={r:s["statistics"][r]["mean_target_error_full_domain_weighted_vector_rms_secondary"] for r in roles}; sources={"roi_mean":{"path":str(roi_mean_root.resolve()),"manifest_sha256":file_sha256(roi_mean_root/"manifest.json")},"cycle_template_ccd":{"path":str(cycle_template_root.resolve()),"manifest_sha256":file_sha256(cycle_template_root/"manifest.json")}}
+  report={"schema_id":SCHEMA_ID,"sources":sources,"figure_files":list(FIGURE_FILES),"primary_domain":s["primary_domain"],"roi_fluid_point_count":s["roi_fluid_point_count"],"roi_grid_shape":[int(np.count_nonzero((x>=34)&(x<=54))),int(np.count_nonzero(np.abs(y)<=5))],"roi_mean_target_errors":dict(zip(roles,errors)),"full_domain_mean_target_errors_secondary":full,"mean_target_error_reductions":s["benefits"],"closure":s["closure"],"cycle_template":{"N":cs["N"],"M":cs["M"],"Q":cs["Q"],"numerical_rank":cs["numerical_rank"],"singular_values":sigma.tolist(),"left_action_vectors":left.tolist(),"exact_physical_action_mean":ccd["arrays"]["drl_physical_action_mean_exact"].tolist(),"admitted_observable_raw_empirical_mean":ccd["arrays"]["admitted_observable_raw_empirical_mean"].tolist(),"admitted_observable_after_physical_mean_empirical_mean":ccd["arrays"]["admitted_observable_after_physical_mean_empirical_mean"].tolist(),"stability":stability_summary,"basis_interpretation":"individual modes and action vectors are the primary SVD basis representation of the rank-3 subspace; only subspace stability is supported"},"full_domain_metrics":"secondary only; retained numerically","claim_boundary":cs["claim_boundary"]}; (stage/"RESULTS.json").write_bytes(canonical_json(report))
+  b=s["benefits"]; (stage/"RESULTS.md").write_text(f"# Corrected Karman results\n\nThe primary domain is the fixed downstream wake ROI `34 <= x/D <= 54`, `|y/D| <= 5`, containing `{s['roi_fluid_point_count']}` common solver-fluid cells.\n\nMean target errors are zero `{errors[1]:.8g}`, constant mean `{errors[2]:.8g}`, and DRL `{errors[3]:.8g}`. Constant control improves the zero baseline by `{b['zero_to_constant_overall_mean_control_benefit_target_error_reduction']:.8g}`; DRL adds `{b['constant_to_drl_dynamic_increment_target_error_reduction']:.8g}`. Full-domain values are retained only as secondary diagnostics.\n\nThe CCD uses 19 DRL cycles by 10 phase bins (`N=190`), `Q=1`, and `M=160400`. Constant-control cycles form one ensemble phase template and are never paired with DRL cycles. The globally least-favorable rank-3 projector check across all declared LOO, split, and persisted bin/origin variants is `{stability_summary['global_least_favorable_projector_check']}`; the largest relative singular-value shift is `{stability_summary['global_largest_singular_shift_check']}`.\n\nThe displayed ROI modes and action vectors are the primary SVD basis representation of the stable rank-3 subspace; individual-vector stability is not claimed. The figures report descriptive co-variation only. {cs['claim_boundary']}\n")
+  files={q.name:file_sha256(q) for q in stage.iterdir()}; (stage/"manifest.json").write_bytes(canonical_json({"schema_id":SCHEMA_ID,"complete":True,"files":files})); rename_noreplace(stage,destination); fresh=load_dynamic_publication(destination); return fresh["path"]
+ except Exception: shutil.rmtree(stage,ignore_errors=True); raise
+
 def load_dynamic_publication(path):
 root=Path(path); manifest=json.loads((root/"manifest.json").read_text())
 if manifest.get("schema_id")!=SCHEMA_ID or not manifest.get("complete"): raise ValueError("publication manifest invalid")
- for name,h in manifest.get("files",{}).items():
-  if file_sha256(root/name)!=h: raise ValueError("publication file hash mismatch")
+ if set(manifest.get("files",{}))!=set(FIGURE_FILES)|{"RESULTS.json","RESULTS.md"}: raise ValueError("publication inventory invalid")
+ for name,digest in manifest["files"].items():
+  if file_sha256(root/name)!=digest: raise ValueError("publication file hash mismatch")
 report=json.loads((root/"RESULTS.json").read_text())
- if report.get("schema_id")!=SCHEMA_ID or report.get("phase_domain_decision")!="DOWNGRADE" or not report.get("zero_phase_output","").startswith("PROHIBITED"): raise ValueError("publication claim contract invalid")
+ if report.get("schema_id")!=SCHEMA_ID or report.get("figure_files")!=list(FIGURE_FILES) or set(report.get("sources",{}))!={"roi_mean","cycle_template_ccd"} or not report.get("full_domain_metrics","").startswith("secondary"): raise ValueError("publication claim contract invalid")
 for source in report["sources"].values():
  if file_sha256(Path(source["path"])/"manifest.json")!=source["manifest_sha256"]: raise ValueError("publication source identity changed")
- return {"manifest":manifest,"report":report,"provenance_validation":"VERIFIED publication hashes and immutable source manifest identities"}
+ load_dynamic_increment(report["sources"]["roi_mean"]["path"]); load_cycle_template_result(report["sources"]["cycle_template_ccd"]["path"],recompute=True)
+ return {"path":root,"manifest":manifest,"report":report,"provenance_validation":"VERIFIED publication, live parents, and cycle-template recomputation"}
@@ -1,161 +0,0 @@
-"""Immutable block-local temporal negative-lag CCD for the Karman DRL role."""
-from __future__ import annotations
-from dataclasses import dataclass
-from pathlib import Path
-from typing import Any, Iterator
-import json, os, shutil, tempfile
-import numpy as np
-from CCD_analysis.acquisition.artifacts import file_sha256, rename_noreplace
-from CCD_analysis.acquisition.contracts import ACTION_IDENTITIES, canonical_json
-from CCD_analysis.direct_dq.analysis import coordinate_weights
-from CCD_analysis.direct_dq.schema import canonical_array_sha256
-from .artifacts import load_role_artifact
-from .phase import load_phase_compact, recover_phase
-
-SCHEMA_ID="ccd-karman-temporal-lagged/v1"
-ROW_ORDER="channel-major_delay-minor"
-CHANNELS=("front","upper","lower")
-PRIMARY_LAGS=tuple(range(-17,1))
-NEIGHBOR_LAG_WINDOWS=(tuple(range(-16,1)),tuple(range(-15,1)))
-CLAIM_BOUNDARY="closed-loop temporal co-variation only; no causal, response-time, mechanism, uncertainty, CCD>POD, or observable-prediction claim"
-
-@dataclass(frozen=True)
-class TemporalInput:
- role_root:Path; phase_root:Path; role_manifest_sha256:str; phase_manifest_sha256:str
- x_D:np.ndarray; y_D:np.ndarray; mask:np.ndarray; fields:np.ndarray; actions:np.ndarray
- relative_steps:np.ndarray; absolute_steps:np.ndarray; cycle_ids:np.ndarray
-
-@dataclass(frozen=True)
-class TemporalConfig:
- chunk_size:int; ram_budget_bytes:int; safety_margin:float=1.25
- def __post_init__(self):
-  if type(self.chunk_size) is not int or self.chunk_size<=0 or type(self.ram_budget_bytes) is not int or self.ram_budget_bytes<=0 or not np.isfinite(self.safety_margin) or self.safety_margin<1: raise ValueError("positive chunk/RAM and safety_margin>=1 required")
-
-@dataclass(frozen=True)
-class TemporalResult:
- arrays:dict[str,np.ndarray]; config:dict[str,Any]; summary:dict[str,Any]; input_hashes:dict[str,Any]
-
-def _memory(inp:TemporalInput,cfg:TemporalConfig)->dict[str,Any]:
- n=len(inp.actions); m=2*int(inp.mask.sum()); q=max(len(PRIMARY_LAGS),*(len(x) for x in NEIGHBOR_LAG_WINDOWS)); c=min(n,cfg.chunk_size)
- terms={"loaded_role_fields_float32":int(inp.fields.nbytes),"loaded_actions_clocks_masks":int(inp.actions.nbytes+inp.relative_steps.nbytes+inp.absolute_steps.nbytes+inp.cycle_ids.nbytes+inp.mask.nbytes),"field_mean_weights_float64":3*m*8,"largest_cross_modes_float64":3*(3*q)*m*8,"chunk_float64":m*c*8,"small_factors":(3*q*n+3*q*3*q)*8}
- raw=sum(terms.values()); peak=int(np.ceil(raw*cfg.safety_margin)); out={"terms_bytes":terms,"raw_peak_ram_bytes":raw,"estimated_peak_ram_bytes":peak,"ram_budget_bytes":cfg.ram_budget_bytes,"safety_margin":cfg.safety_margin,"scratch_bytes":0,"decision":"PASS" if peak<=cfg.ram_budget_bytes else "FAIL","formula":"ceil(safety_margin*declared terms); no MxM or full float64 MxN"}
- if peak>cfg.ram_budget_bytes: raise MemoryError(f"temporal CCD memory admission failed: {out}")
- return out
-
-def load_temporal_input(role_root:str|Path,phase_root:str|Path)->TemporalInput:
- role_root,phase_root=Path(role_root).resolve(),Path(phase_root).resolve(); d=load_role_artifact(role_root,expected_role="drl"); phase=load_phase_compact(phase_root)
- if not phase["summary"]["gate_passed"]: raise ValueError("passing DRL phase artifact required")
- if phase["summary"]["source"]["campaign_manifest_sha256"]!=file_sha256(role_root/"campaign_manifest.json"): raise ValueError("phase/role parent mismatch")
- start,stop=d["metadata"]["retained_slice"]; a=d["legacy_arrays"]; t=d["telemetry"]
- if d["metadata"]["contract"]["contract"]["control_interval"]!=800: raise ValueError("exact 800-lattice-step cadence required")
- rel=t["acquisition_relative_lattice_steps"][start:stop].copy(); absolute=a["lattice_steps"][start:stop].copy()
- if rel.dtype!=np.int64 or absolute.dtype!=np.int64 or not np.all(np.diff(rel)==800) or not np.all(np.diff(absolute)==800): raise ValueError("exact 800-step clocks required")
- r=recover_phase(t["center_sensor_uy"][start:stop]); ids=r["cycle_id"].copy()
- with np.load(phase_root/"compact.npz",allow_pickle=False) as z:
-  if not np.array_equal(z["cycle_ids"],np.unique(ids[ids>=0])) or not np.array_equal(z["fluid_mask"],a["fluid_mask"]): raise ValueError("phase blocks/mask differ from live role")
- fields=np.stack((a["ux"][start:stop],a["uy"][start:stop]),axis=1)
- actions=t["effective_applied_action"][start:stop,-3:].copy()
- if fields.dtype!=np.float32 or actions.dtype!=np.float32 or fields.shape[0]!=len(ids) or actions.shape!=(len(ids),3): raise ValueError("full-resolution DRL fields/actions invalid")
- return TemporalInput(role_root,phase_root,file_sha256(role_root/"campaign_manifest.json"),file_sha256(phase_root/"manifest.json"),a["x_D"].copy(),a["y_D"].copy(),a["fluid_mask"].copy(),fields,actions,rel,absolute,ids)
-
-def _admit(ids:np.ndarray,lags:tuple[int,...],support:np.ndarray|None=None)->tuple[np.ndarray,np.ndarray]:
- if not lags or any(type(x) is not int or x>0 for x in lags) or tuple(sorted(set(lags)))!=lags: raise ValueError("lags must be unique increasing nonpositive integers")
- local={}; positions={}
- for i,b in enumerate(ids):
-  if b>=0: positions[i]=len(local.setdefault(int(b),[])); local[int(b)].append(i)
- field=[]; obs=[]
- for i in range(len(ids)):
-  if i not in positions: continue
-  seq=local[int(ids[i])]; pos=positions[i]; targets=[pos+lag for lag in lags]
-  if min(targets)<0: continue
-  if support is not None and i not in support: continue
-  field.append(i); obs.append([seq[j] for j in targets])
- if not field: raise ValueError("no complete block-local lag columns")
- return np.asarray(field,np.int64),np.asarray(obs,np.int64)
-
-def _chunks(inp:TemporalInput,indices:np.ndarray,chunk:int)->Iterator[tuple[slice,np.ndarray]]:
- mask=inp.mask
- for s in range(0,len(indices),chunk):
-  e=min(s+chunk,len(indices)); raw=inp.fields[indices[s:e]]; u=np.concatenate((raw[:,0][:,mask],raw[:,1][:,mask]),axis=1).T.astype(np.float64); yield slice(s,e),u
-
-def _decompose(inp:TemporalInput,lags:tuple[int,...],cfg:TemporalConfig,field_mean:np.ndarray,weights:np.ndarray,support:np.ndarray|None=None)->dict[str,np.ndarray]:
- fi,oi=_admit(inp.cycle_ids,lags,support); p=inp.actions[oi].transpose(2,1,0).reshape(3*len(lags),len(fi)).astype(np.float64); pm=p.mean(1); roots=np.sqrt(weights); cross=np.zeros((len(pm),len(weights)))
- for sl,u in _chunks(inp,fi,cfg.chunk_size): cross+=(p[:,sl]-pm[:,None])@((u-field_mean[:,None])*roots[:,None]).T
- cross/=len(fi)*np.sqrt(3*len(lags)); left,s,vh=np.linalg.svd(cross,full_matrices=False); weighted=vh.T
- for k in range(weighted.shape[1]):
-  pivot=int(np.argmax(np.abs(weighted[:,k])))
-  if weighted[pivot,k]<0: weighted[:,k]*=-1; left[:,k]*=-1
- modes=weighted/roots[:,None]; coeff=np.empty((len(s),len(fi))); total=0.
- for sl,u in _chunks(inp,fi,cfg.chunk_size): x=(u-field_mean[:,None])*roots[:,None]; coeff[:,sl]=weighted.T@x; total+=float(np.sum(x*x))
- residual=np.sqrt(np.maximum(total-np.cumsum(np.sum(coeff*coeff,axis=1)),0)/max(total,np.finfo(float).tiny))
- return {"lags":np.asarray(lags,np.int64),"field_indices":fi,"observable_indices":oi,"observable_mean":pm,"cross_correlation":cross,"left_functions":left,"singular_values":s,"physical_modes":modes,"coefficients":coeff,"weighted_relative_residuals":residual}
-
-def _left_metrics(result:dict[str,np.ndarray])->dict[str,Any]:
- lags=result["lags"]; left=result["left_functions"].reshape(3,len(lags),-1); out=[]
- for k in range(left.shape[2]):
-  lag_energy=np.sum(left[:,:,k]**2,axis=0); lag_energy/=lag_energy.sum(); channel=np.sum(left[:,:,k]**2,axis=1)
-  out.append({"mode":k+1,"peak_lag_boundaries":int(lags[int(np.argmax(lag_energy))]),"peak_lag_lattice_steps":800*int(lags[int(np.argmax(lag_energy))]),"lag_energy_centroid_boundaries":float(np.sum(lags*lag_energy)),"channel_squared_norms":{CHANNELS[j]:float(channel[j]) for j in range(3)},"coefficient_rms":float(np.sqrt(np.mean(result["coefficients"][k]**2)))})
- return {"modes":out,"lag_energy_by_mode":np.sum(left*left,axis=0).T.tolist()}
-
-def _compare(a:dict[str,np.ndarray],b:dict[str,np.ndarray],weights:np.ndarray)->dict[str,Any]:
- r=min(3,len(a["singular_values"]),len(b["singular_values"])); roots=np.sqrt(weights)[:,None]; va=a["physical_modes"][:,:r]*roots; vb=b["physical_modes"][:,:r]*roots
- overlap=np.linalg.svd(va.T@vb,compute_uv=False)
- return {"leading_spectrum":b["singular_values"][:r].tolist(),"leading_spectrum_relative_change":((b["singular_values"][:r]-a["singular_values"][:r])/np.maximum(a["singular_values"][:r],np.finfo(float).tiny)).tolist(),"leading_weighted_subspace_principal_cosines":overlap.tolist()}
-
-def decompose_temporal(inp:TemporalInput,*,streaming_config:TemporalConfig)->TemporalResult:
- mem=_memory(inp,streaming_config); mask=inp.mask; point=(coordinate_weights(inp.x_D)[:,None]*coordinate_weights(inp.y_D)[None,:])[mask]; weights=np.concatenate((point,point)); all_idx=np.arange(len(inp.actions),dtype=np.int64)
- fsum=np.zeros(len(weights)); count=0
- for _,u in _chunks(inp,all_idx,streaming_config.chunk_size): fsum+=u.sum(1); count+=u.shape[1]
- field_mean=fsum/count
- primary=_decompose(inp,PRIMARY_LAGS,streaming_config,field_mean,weights); common=primary["field_indices"]
- neighbors=[]; common_neighbors=[]
- for lags in NEIGHBOR_LAG_WINDOWS:
-  neighbors.append(_decompose(inp,lags,streaming_config,field_mean,weights)); common_neighbors.append(_decompose(inp,lags,streaming_config,field_mean,weights,support=common))
- arrays={"x_D":inp.x_D,"y_D":inp.y_D,"fluid_mask":mask,"coordinate_weights":weights,"full_run_field_mean":field_mean,"full_run_action_mean":inp.actions.astype(np.float64).mean(0),"effective_actions":inp.actions,"relative_steps":inp.relative_steps,"absolute_steps":inp.absolute_steps,"cycle_ids":inp.cycle_ids}
- for prefix,r in [("primary",primary)]+[(f"neighbor_{i}",v) for i,v in enumerate(neighbors)]+[(f"common_neighbor_{i}",v) for i,v in enumerate(common_neighbors)]:
-  for k,v in r.items(): arrays[f"{prefix}_{k}"]=v
- sensitivity=[]
- for i,(n,c) in enumerate(zip(neighbors,common_neighbors)):
-  sensitivity.append({"lags":n["lags"].tolist(),"native_N":int(len(n["field_indices"])),"common_support_N":int(len(c["field_indices"])),"native_vs_primary":_compare(primary,n,weights),"common_support_vs_primary":_compare(primary,c,weights)})
- config={"schema_id":"ccd-karman-temporal-lagged-config/v1","case_id":"karman_re100","role":"drl","primary_lags_boundaries":list(PRIMARY_LAGS),"neighbor_lag_windows_boundaries":[list(x) for x in NEIGHBOR_LAG_WINDOWS],"lag_sign":"tau<0 means action precedes field","cadence_lattice_steps":800,"row_order":ROW_ORDER,"observable_channels":list(CHANNELS),"action_identities":list(ACTION_IDENTITIES),"field_estimand":"full-resolution mask-compressed q_DRL(t)-mean_over_all_360_retained_q_DRL","observable_estimand":"exact same-boundary three-channel effective_applied_action fluctuation; per-lag-row admitted-support mean","operator":"A=P(W^(1/2)U)^T/(N*sqrt(3Q))","center_snapshots":True,"center_observables":True,"standardization":False,"whitening":False,"interpolation":False,"nearest":False,"wrap":False,"block_definition":"independent complete rising-zero-crossing phase cycles","chunk_size":streaming_config.chunk_size,"memory":mem,"claim_boundary":CLAIM_BOUNDARY}
- summary={"schema_id":"ccd-karman-temporal-lagged-summary/v1","N":int(len(primary["field_indices"])),"Q":len(PRIMARY_LAGS),"M":len(weights),"cycle_count":int(len(np.unique(inp.cycle_ids[inp.cycle_ids>=0]))),"numerical_rank":int(np.sum(primary["singular_values"]>1e-10*primary["singular_values"][0])),"spectrum_label":"cross-correlation strength; not field energy or explained variance","primary_singular_values":primary["singular_values"].tolist(),"primary_squared_singular_values":(primary["singular_values"]**2).tolist(),"primary_left_lag_metrics":_left_metrics(primary),"sensitivity":sensitivity,"claim_boundary":CLAIM_BOUNDARY,"provenance_status":"VERIFIED_LIVE_ROLE_AND_PHASE_REQUIRED_ON_LOAD"}
- hashes={"parents":{"role":{"path":str(inp.role_root),"campaign_manifest_sha256":inp.role_manifest_sha256},"phase":{"path":str(inp.phase_root),"manifest_sha256":inp.phase_manifest_sha256}},"canonical_arrays":{k:canonical_array_sha256(v) for k,v in arrays.items()}}
- validate_temporal_result(arrays,config,summary,hashes); return TemporalResult(arrays,config,summary,hashes)
-
-def validate_temporal_result(arrays,config,summary,hashes):
- if config.get("schema_id")!="ccd-karman-temporal-lagged-config/v1" or summary.get("schema_id")!="ccd-karman-temporal-lagged-summary/v1" or config.get("primary_lags_boundaries")!=list(PRIMARY_LAGS) or config.get("row_order")!=ROW_ORDER or config.get("claim_boundary")!=CLAIM_BOUNDARY or summary.get("claim_boundary")!=CLAIM_BOUNDARY: raise ValueError("temporal CCD frozen schema contradicted")
- if config.get("field_estimand")!="full-resolution mask-compressed q_DRL(t)-mean_over_all_360_retained_q_DRL" or config.get("observable_channels")!=list(CHANNELS) or config.get("cadence_lattice_steps")!=800 or any(config.get(k) is not False for k in ("standardization","whitening","interpolation","nearest","wrap")): raise ValueError("temporal CCD estimand contradicted")
- d={k:np.asarray(v) for k,v in arrays.items()}; n=len(d["effective_actions"]); m=2*int(d["fluid_mask"].sum())
- if summary.get("N")!=len(d["primary_field_indices"]) or summary.get("Q")!=18 or summary.get("M")!=m or d["primary_cross_correlation"].shape!=(54,m) or d["primary_left_functions"].shape[0]!=54 or d["full_run_field_mean"].shape!=(m,): raise ValueError("temporal CCD dimensions invalid")
- if d["effective_actions"].shape!=(n,3) or d["cycle_ids"].shape!=(n,) or not np.all(np.diff(d["relative_steps"])==800): raise ValueError("temporal CCD clocks/actions invalid")
- if set(hashes)!={"parents","canonical_arrays"} or set(hashes["canonical_arrays"])!=set(d) or any(hashes["canonical_arrays"][k]!=canonical_array_sha256(v) for k,v in d.items()): raise ValueError("temporal CCD hashes invalid")
- canonical_json(config); canonical_json(summary); canonical_json(hashes); return d
-
-class TemporalTransaction:
- def __init__(self,destination): self.destination=Path(destination); self.stage=None
- def __enter__(self):
-  if self.destination.exists(): raise FileExistsError(self.destination)
-  self.destination.parent.mkdir(parents=True,exist_ok=True); self.stage=Path(tempfile.mkdtemp(prefix=f".{self.destination.name}.partial-",dir=self.destination.parent)); return self
- def write(self,result):
-  validate_temporal_result(result.arrays,result.config,result.summary,result.input_hashes); np.savez_compressed(self.stage/"arrays.npz",**result.arrays)
-  for n,v in (("config.json",result.config),("summary.json",result.summary),("input_hashes.json",result.input_hashes)): (self.stage/n).write_bytes(canonical_json(v))
-  files={p.name:file_sha256(p) for p in self.stage.iterdir()}; (self.stage/"manifest.json").write_bytes(canonical_json({"schema_id":SCHEMA_ID,"complete":True,"files":files}))
- def publish(self):
-  load_temporal_result(self.stage,recompute=False); rename_noreplace(self.stage,self.destination); self.stage=None; load_temporal_result(self.destination,recompute=True); return self.destination
- def __exit__(self,*args):
-  if self.stage is not None: shutil.rmtree(self.stage,ignore_errors=True)
-
-def load_temporal_result(path,recompute=True):
- root=Path(path); manifest=json.loads((root/"manifest.json").read_text())
- if manifest.get("schema_id")!=SCHEMA_ID or not manifest.get("complete") or set(manifest.get("files",{}))!={"arrays.npz","config.json","summary.json","input_hashes.json"}: raise ValueError("temporal CCD manifest invalid")
- for n,h in manifest["files"].items():
-  if file_sha256(root/n)!=h: raise ValueError("temporal CCD file hash mismatch")
- with np.load(root/"arrays.npz",allow_pickle=False) as z: arrays={k:z[k].copy() for k in z.files}
- config=json.loads((root/"config.json").read_text()); summary=json.loads((root/"summary.json").read_text()); hashes=json.loads((root/"input_hashes.json").read_text()); validate_temporal_result(arrays,config,summary,hashes)
- inp=load_temporal_input(hashes["parents"]["role"]["path"],hashes["parents"]["phase"]["path"])
- if inp.role_manifest_sha256!=hashes["parents"]["role"]["campaign_manifest_sha256"] or inp.phase_manifest_sha256!=hashes["parents"]["phase"]["manifest_sha256"]: raise ValueError("temporal CCD live parent identity changed")
- if not np.array_equal(inp.actions,arrays["effective_actions"]) or not np.array_equal(inp.relative_steps,arrays["relative_steps"]) or not np.array_equal(inp.cycle_ids,arrays["cycle_ids"]): raise ValueError("temporal CCD live telemetry/blocks changed")
- if recompute:
-  cfg=TemporalConfig(config["chunk_size"],config["memory"]["ram_budget_bytes"],config["memory"]["safety_margin"]); fresh=decompose_temporal(inp,streaming_config=cfg)
-  for k in arrays: np.testing.assert_allclose(arrays[k],fresh.arrays[k],rtol=2e-11,atol=2e-12)
- return {"arrays":arrays,"config":config,"summary":summary,"input_hashes":hashes,"manifest":manifest,"provenance_validation":"VERIFIED: live DRL role and passing phase parents reread; essential decomposition recomputed" if recompute else "VERIFIED_HASHES_AND_LIVE_PARENTS"}
@@ -0,0 +1,36 @@
+#!/usr/bin/env python3
+"""Redraw one CCD mode field from an authenticated plot_data.npz view."""
+from __future__ import annotations
+import argparse
+from pathlib import Path
+import matplotlib
+matplotlib.use("Agg")
+import matplotlib.pyplot as plt
+import numpy as np
+
+def main(argv=None):
+    parser=argparse.ArgumentParser()
+    parser.add_argument("--plot-data", type=Path, required=True)
+    parser.add_argument("--mode", type=int, choices=(1,2,3), default=1)
+    parser.add_argument("--component", choices=("ux","uy","rotation"), default="ux")
+    parser.add_argument("--output", type=Path, required=True)
+    parser.add_argument("--limit", type=float)
+    args=parser.parse_args(argv)
+    if args.output.exists(): raise FileExistsError(args.output)
+    with np.load(args.plot_data, allow_pickle=False) as data:
+        x,y,mask=data["x_D"],data["y_D"],data["domain_mask"]
+        if args.component == "rotation": values=data["ccd_mode_rotation_proxy"][args.mode-1]
+        else: values=data["ccd_mode_velocity"][args.mode-1, 0 if args.component == "ux" else 1]
+        values=np.where(mask, values, np.nan)
+    limit=args.limit or float(np.nanpercentile(np.abs(values[mask]),99.5))
+    args.output.parent.mkdir(parents=True, exist_ok=True)
+    fig,ax=plt.subplots(figsize=(9,4.5),layout="constrained")
+    image=ax.pcolormesh(x,y,values.T,shading="nearest",cmap="RdBu_r",vmin=-limit,vmax=limit)
+    ax.set(xlim=(29,54),ylim=(-5,5),xlabel="x/D",ylabel="y/D",title=f"CCD mode {args.mode} {args.component}")
+    ax.set_aspect("equal"); fig.colorbar(image,ax=ax)
+    for ext in ("png","pdf"):
+        target=args.output.with_suffix("."+ext)
+        fig.savefig(target,dpi=250 if ext=="png" else None,bbox_inches="tight")
+    plt.close(fig)
+    return 0
+if __name__ == "__main__": raise SystemExit(main())
@@ -0,0 +1,106 @@
+import json
+from pathlib import Path
+import numpy as np
+import pytest
+from CCD_analysis.direct_dq.schema import canonical_array_sha256
+from CCD_analysis.karman_dynamic.cycle_template_ccd import (
+ CycleTemplateResult,CycleTemplateTransaction,_ensemble_matrix,_fit,_load,_matrix,
+ load_cycle_template_result,validate_cycle_template,SCHEMA_ID,
+)
+from CCD_analysis.original_ccd import CCDConfig,decompose
+
+
+def _synthetic(points=5,seed=4):
+ rng=np.random.default_rng(seed); cycles,bins=19,10
+ df=rng.normal(size=(cycles,bins,2,points)); cf=rng.normal(size=(cycles,bins,2,points)); da=rng.normal(size=(cycles,bins,3))
+ return df,cf,da,np.ones(points,bool),np.linspace(.5,1.5,points),np.zeros((2,points)),np.array([.11,-.07,.03])
+
+
+def test_shapes_component_order_literal_operator_exact_physical_mean_and_no_pairing():
+ df,cf,da,selector,w0,means,physical_mean=_synthetic()
+ u,p,w,cycle_ids,bin_ids,stored_mean=_matrix(df,da,cf,means,means,selector,w0,physical_mean)
+ assert u.shape==(10,190) and p.shape==(3,190) and w.shape==(10,)
+ np.testing.assert_array_equal(u[:5],(df-cf.mean(0)[None]).reshape(190,2,5).transpose(1,2,0).reshape(10,190)[:5])
+ np.testing.assert_array_equal(stored_mean,physical_mean)
+ np.testing.assert_allclose(p.mean(1)+physical_mean,da.reshape(190,3).mean(0),rtol=1e-15,atol=1e-15)
+ np.testing.assert_array_equal(np.unique(cycle_ids),np.arange(19)); np.testing.assert_array_equal(bin_ids,np.tile(np.arange(10),19))
+ got=_fit(u,p,w); ref=decompose(u,p,weight=w,config=CCDConfig(center_snapshots=True,center_observables=True))
+ expected=(p-p.mean(1,keepdims=True))@((u-u.mean(1,keepdims=True))*np.sqrt(w)[:,None]).T/(190*np.sqrt(3))
+ np.testing.assert_allclose(got['cross'],expected,rtol=1e-14,atol=1e-14); np.testing.assert_allclose(got['cross'],ref.cross_correlation,rtol=1e-14,atol=1e-14)
+ u_reordered,*_=_matrix(df,da,cf[::-1],means,means,selector,w0,physical_mean)
+ np.testing.assert_allclose(u,u_reordered,rtol=1e-14,atol=1e-14)
+
+
+def test_cycle_block_and_constant_template_leave_one_out_are_distinct():
+ df,cf,da,selector,w,means,physical_mean=_synthetic(seed=2)
+ primary=_matrix(df,da,cf,means,means,selector,w,physical_mean)
+ kept=np.delete(np.arange(19),9); drl_loo=_matrix(df,da,cf,means,means,selector,w,physical_mean,drl_cycles=kept)
+ template_loo=_matrix(df,da,cf,means,means,selector,w,physical_mean,constant_cycles=np.arange(1,19))
+ assert primary[0].shape[1]==190 and drl_loo[0].shape[1]==180 and template_loo[0].shape[1]==190
+ assert not np.allclose(primary[0],template_loo[0]); np.testing.assert_allclose(drl_loo[0],np.concatenate((primary[0][:,:90],primary[0][:,100:]),axis=1)); np.testing.assert_array_equal(drl_loo[3],kept)
+
+
+def test_persisted_alternatives_are_used_without_resampling():
+ df,cf,da,selector,w,means,physical_mean=_synthetic(seed=7); drl_alt=df.mean(0)[:8].copy(); const_alt=cf.mean(0)[:8].copy(); action_alt=da.mean(0)[:8].copy(); drl_alt[0]+=123.; action_alt[0]+=17.
+ u,p,weights=_ensemble_matrix(drl_alt,action_alt,const_alt,means,means,selector,w,physical_mean)
+ expected=((drl_alt-const_alt).reshape(8,10).T); np.testing.assert_array_equal(u,expected); np.testing.assert_array_equal(p,(action_alt-physical_mean).T); assert weights.shape==(10,)
+ with pytest.raises(ValueError,match='malformed'): _ensemble_matrix(drl_alt,action_alt[:7],const_alt,means,means,selector,w,physical_mean)
+
+
+
+def test_loader_returns_persisted_alternatives_and_fails_closed_when_missing(tmp_path,monkeypatch):
+ import CCD_analysis.karman_dynamic.cycle_template_ccd as mod
+ drl,constant,roi=[tmp_path/name for name in ('drl','constant','roi')]; [q.mkdir() for q in (drl,constant,roi)]
+ mask=np.ones((2,2),bool); cycles=19; sentinel={}
+ def phase_values(role):
+  values={'fluid_mask':mask,'cycle_bin_fields':np.zeros((cycles,10,2,4),np.float32),'cycle_bin_effective_actions':np.zeros((cycles,10,3),np.float32),'cycle_bin_counts':np.ones((cycles,10),np.int64),'cycle_ids':np.arange(cycles)}
+  for name,bins in (('bins8',8),('bins12',12),('bins10_half_shift',10)):
+   marker=(11 if role=='drl' else 23)+bins; values[f'ensemble_{name}']=np.full((bins,2,4),marker,np.float32); values[f'ensemble_{name}_per_cycle_fields']=np.full((cycles,bins,2,4),marker,np.float32); values[f'counts_{name}']=np.ones((cycles,bins),np.int64)
+   values[f'ensemble_{name}_effective_actions']=np.full((bins,3),marker+100,np.float32); values[f'ensemble_{name}_per_cycle_effective_actions']=np.full((cycles,bins,3),marker+100,np.float32); sentinel[name]=(marker,marker+100) if role=='drl' else sentinel.get(name)
+  return values
+ np.savez(drl/'compact.npz',**phase_values('drl')); np.savez(constant/'compact.npz',**phase_values('constant'))
+ np.savez(roi/'arrays.npz',roi_fluid_mask=mask,roi_quadrature_weights=np.ones(4),roi_selector_on_common_mask=np.ones(4,bool),mean_drl=np.zeros((2,4)),mean_constant_mean=np.zeros((2,4)))
+ provenance={'source':'fresh DRL retained effective_applied_action mean','constant_mean_physical_action':[.1,.2,.3]}; summaries={str(drl):{'gate_passed':True,'role':'drl','constant_mean_provenance':provenance},str(constant):{'gate_passed':True,'role':'constant_mean'}}
+ monkeypatch.setattr(mod,'load_phase_compact',lambda path:{'summary':summaries[str(Path(path))]}); monkeypatch.setattr(mod,'load_dynamic_increment',lambda path:{'summary':{'statistics':{'drl':{'effective_action_mean':[.1,.2,.3]}}}})
+ loaded=_load(drl,constant,roi); alternatives=loaded[-1]
+ for name,(field_marker,action_marker) in sentinel.items(): assert np.all(alternatives[name]['drl_fields']==field_marker) and np.all(alternatives[name]['drl_actions']==action_marker)
+ with np.load(drl/'compact.npz') as z: broken={k:z[k] for k in z.files if k!='ensemble_bins8_effective_actions'}
+ np.savez(drl/'compact.npz',**broken)
+ with pytest.raises(ValueError,match='required persisted alternative missing'): _load(drl,constant,roi)
+
+def _valid_large_result(tmp_path):
+ m,n=160400,190; arrays={
+  'coordinate_weights':np.ones(m), 'primary_cross_correlation':np.zeros((3,m)),
+  'primary_left_functions':np.eye(3), 'primary_singular_values':np.ones(3),
+  'primary_physical_modes':np.zeros((m,3)), 'primary_coefficients':np.zeros((3,n)),
+  'drl_cycle_ids':np.repeat(np.arange(19),10), 'phase_bin_ids':np.tile(np.arange(10),19),
+  'drl_cycle_bin_boundary_counts':np.ones(n,np.int64), 'constant_template_cycle_ids':np.arange(19),
+  'constant_template_cycle_bin_boundary_counts':np.ones((19,10),np.int64),
+  'drl_physical_action_mean_exact':np.zeros(3), 'admitted_field_empirical_mean':np.zeros(m),
+  'admitted_observable_after_physical_mean_empirical_mean':np.zeros(3), 'admitted_observable_raw_empirical_mean':np.zeros(3),
+ }
+ parents={}
+ for name in ('drl_phase','constant_mean_phase','roi_mean'):
+  root=tmp_path/name; root.mkdir(); (root/'manifest.json').write_text(name); parents[name]={'path':str(root),'manifest_sha256':name}
+ config={'schema_id':SCHEMA_ID,'N':190,'M':m,'Q':1,'cycle_pairing':False,'operator':'A=P_c(W^(1/2)U_c)^T/(190*sqrt(3))','standardization':False,'whitening':False,'pod':False,'alternative_phase_products':'persisted independent retained-boundary re-binnings; no Fourier resampling','column_bootstrap_api':False}
+ summary={'schema_id':SCHEMA_ID,'N':190,'M':m,'Q':1}; hashes={'parents':parents,'canonical_arrays':{k:canonical_array_sha256(v) for k,v in arrays.items()}}
+ return CycleTemplateResult(arrays,config,summary,hashes)
+
+
+def test_no_column_bootstrap_contract_and_hash_bound_reload(tmp_path,monkeypatch):
+ result=_valid_large_result(tmp_path); validate_cycle_template(result.arrays,result.config,result.summary,result.input_hashes)
+ assert result.config['column_bootstrap_api'] is False and 'bootstrap' not in {name for name in dir(result) if not name.startswith('_')}
+ import CCD_analysis.karman_dynamic.cycle_template_ccd as mod
+ real_hash=mod.file_sha256
+ def hash_or_label(path):
+  path=Path(path)
+  if path.name=='manifest.json' and path.parent.name in ('drl_phase','constant_mean_phase','roi_mean'): return path.parent.name
+  return real_hash(path)
+ monkeypatch.setattr(mod,'file_sha256',hash_or_label); monkeypatch.setattr(mod,'decompose_cycle_template',lambda *args:result)
+ out=tmp_path/'nested'/'result'
+ with CycleTemplateTransaction(out) as tx: tx.write(result); assert tx.publish()==out
+ loaded=load_cycle_template_result(out,recompute=True); assert loaded['provenance_validation'].startswith('VERIFIED live')
+ with pytest.raises(FileExistsError):
+  with CycleTemplateTransaction(out): pass
+ (out/'summary.json').write_text('{}')
+ with pytest.raises(ValueError,match='hash'): load_cycle_template_result(out,recompute=False)
@@ -0,0 +1,77 @@
+import json
+from pathlib import Path
+
+import numpy as np
+import pytest
+
+from CCD_analysis.karman_dynamic.domain_comparison import (
+    DOMAIN, _weighted_metrics, domain_mask, load_result,
+)
+from CCD_analysis.karman_dynamic.pod_baseline import fit_weighted_pod
+
+
+def test_domain_selector_is_inclusive_common_mask_and_contains_body_vicinities():
+    x = np.linspace(28, 55, 271)
+    y = np.linspace(-6, 6, 121)
+    common = np.ones((x.size, y.size), bool)
+    common[np.argmin(abs(x - 30)), np.argmin(abs(y))] = False
+    selected = domain_mask(x, y, common)
+    assert selected.shape == common.shape
+    assert not selected[x < DOMAIN["x_min_D"]].any()
+    assert not selected[x > DOMAIN["x_max_D"]].any()
+    assert not selected[:, abs(y) > DOMAIN["abs_y_max_D"]].any()
+    assert selected[np.isclose(x, 29)].any() and selected[np.isclose(x, 54)].any()
+    assert not selected[np.argmin(abs(x - 30)), np.argmin(abs(y))]
+
+
+def test_domain_selector_rejects_malformed_and_missing_geometry():
+    x, y = np.arange(4.0), np.arange(-2.0, 3.0)
+    with pytest.raises(ValueError, match="inconsistent"):
+        domain_mask(x, y, np.ones((4, 4), bool))
+    with pytest.raises(ValueError, match="no common fluid"):
+        domain_mask(np.linspace(29, 54, 30), np.linspace(-5, 5, 20), np.zeros((30, 20), bool))
+
+
+def test_weighted_reconstruction_identity_and_region_metrics():
+    truth = np.array([[1., 2.], [2., 1.], [3., 4.], [4., 3.]])
+    weights = np.array([1., 2., 1., 2.])
+    residual, cosine = _weighted_metrics(truth, truth, weights, np.array([True, True]))
+    np.testing.assert_array_equal(residual, 0)
+    np.testing.assert_allclose(cosine, 1)
+    zero_residual, zero_cosine = _weighted_metrics(truth, np.zeros_like(truth), weights, np.array([True, False]))
+    np.testing.assert_allclose(zero_residual, 1)
+    np.testing.assert_array_equal(zero_cosine, 0)
+
+
+def test_raw_pod_is_weight_orthonormal_and_energy_ranked():
+    rng = np.random.default_rng(7)
+    snapshots = rng.normal(size=(12, 8))
+    weights = np.linspace(.5, 1.5, 12)
+    pod = fit_weighted_pod(snapshots, weights, rank=3)
+    np.testing.assert_allclose(pod["modes"].T @ (weights[:, None] * pod["modes"]), np.eye(3), atol=1e-12)
+    assert np.all(np.diff(pod["singular_values"]) <= 0)
+    assert np.all(np.diff(pod["weighted_residual_fraction"]) <= 0)
+
+
+def test_result_loader_is_hash_bound_and_no_clobber_inventory(tmp_path, monkeypatch):
+    root = tmp_path / "result"; root.mkdir()
+    parent = tmp_path / "parent"; parent.mkdir(); (parent / "manifest.json").write_text("parent")
+    for name, content in {
+        "arrays.npz": b"arrays", "config.json": b"{}", "summary.json": b"{}",
+        "input_hashes.json": json.dumps({"p": {"path": str(parent), "manifest_sha256": "parent-hash"}}).encode(),
+        "RESULTS.md": b"result",
+    }.items():
+        (root / name).write_bytes(content)
+    import CCD_analysis.karman_dynamic.domain_comparison as mod
+    real_hash = mod.file_sha256
+    def fake_hash(path):
+        path = Path(path)
+        if path == parent / "manifest.json": return "parent-hash"
+        return real_hash(path)
+    monkeypatch.setattr(mod, "file_sha256", fake_hash)
+    files = {p.name: fake_hash(p) for p in root.iterdir()}
+    (root / "manifest.json").write_bytes(mod.canonical_json({"schema_id": mod.SCHEMA_ID, "complete": True, "files": files}))
+    assert load_result(root)["path"] == root
+    (root / "RESULTS.md").write_text("changed")
+    with pytest.raises(ValueError, match="hash"):
+        load_result(root)
@@ -6,7 +6,7 @@ from CCD_analysis.karman_dynamic.dynamic_increment import publish_dynamic_increm

 def test_synthetic_four_role_increment_closure_and_provenance(tmp_path,monkeypatch):
 from CCD_analysis.karman_dynamic import dynamic_increment as mod
- roles=('drl','constant_mean','target','zero'); nx,ny,n=4,3,12; mask=np.ones((nx,ny),bool); x=np.arange(nx,dtype=np.float32); y=np.arange(ny,dtype=np.float32)
+ roles=('drl','constant_mean','target','zero'); nx,ny,n=4,3,12; mask=np.ones((nx,ny),bool); x=np.asarray([34,40,54,60],np.float32); y=np.asarray([-5,0,5],np.float32)
 offsets={'target':1.,'zero':4.,'constant_mean':3.,'drl':2.}
 role_data={}; role_paths={}; phase_paths={}
 for role in roles:
@@ -19,9 +19,9 @@ def test_synthetic_four_role_increment_closure_and_provenance(tmp_path,monkeypat
 def phase(path):
  role=Path(path).name.replace('-phase',''); fields=np.full((10,2,mask.sum()),offsets[role],np.float32)
  return {'summary':{'gate_passed':True},'fields':fields}
- monkeypatch.setattr(mod,'load_phase_compact',lambda path:phase(path))
+ monkeypatch.setattr(mod,'load_phase_compact',lambda path,**kwargs:phase(path))
 monkeypatch.setattr(mod,'_phase_mean',lambda path,common:phase(path)['fields'])
 out=tmp_path/'result'; result=publish_dynamic_increment(role_paths,phase_paths,out)
- s=result['summary']; assert s['drl_constant_phase_differences_available']; assert s['benefits']['zero_to_constant_overall_mean_control_benefit_target_error_reduction']>0; assert s['benefits']['constant_to_drl_dynamic_increment_target_error_reduction']>0
+ s=result['summary']; assert s['drl_constant_phase_differences_available']; assert s['primary_domain'].startswith('wake ROI'); assert s['roi_fluid_point_count']==9; assert s['full_domain_metrics_role']=='secondary'; assert s['benefits']['zero_to_constant_overall_mean_control_benefit_target_error_reduction']>0; assert s['benefits']['constant_to_drl_dynamic_increment_target_error_reduction']>0
 assert s['closure']['benefit_additivity_absolute_residual']<1e-12 and not s['dense_fields_deleted']
 with pytest.raises(FileExistsError): publish_dynamic_increment(role_paths,phase_paths,out)
@@ -0,0 +1,54 @@
+import numpy as np
+import pytest
+from CCD_analysis.karman_dynamic.mode_diagnostics import (action_coordinates,bounded_streamwise_shift_correlation,component_major_velocity,masked_vorticity,phase_harmonics,reshape_cycle_phase,subspace_symmetry_diagnostics,symmetry_diagnostics,x_localization)
+from CCD_analysis.karman_dynamic.pod_baseline import align_basis_to_reference,fit_weighted_pod,weighted_principal_cosines
+
+
+def test_component_major_vorticity_sign_and_xy_orientation():
+ x=np.array([-1.,-.4,.2,1.,2.1]); y=np.array([-2.,-.8,0.,.9,2.]); xx,yy=np.meshgrid(x,y,indexing="ij"); mask=np.ones(xx.shape,bool)
+ velocity=component_major_velocity(np.r_[(-yy)[mask],xx[mask]],mask); omega,valid=masked_vorticity(*velocity,x,y,mask)
+ assert np.all(valid); np.testing.assert_allclose(omega,2.,rtol=1e-13,atol=1e-13)
+ with pytest.raises(ValueError,match="coordinate"): masked_vorticity(*velocity,y[:-1],x,mask)
+
+
+def test_masked_vorticity_requires_complete_three_point_stencils():
+ x=np.array([0.,.2,.7,1.5,2.6]); y=np.array([-1.,-.2,.5,1.4,2.5]); xx,yy=np.meshgrid(x,y,indexing="ij"); mask=np.ones(xx.shape,bool); mask[2,2]=False
+ omega,valid=masked_vorticity(np.zeros_like(xx),xx**2,x,y,mask)
+ assert np.isfinite(omega).all() and not valid[1:4,2].any() and not valid[2,1:4].any(); np.testing.assert_allclose(omega[valid],2*xx[valid],rtol=1e-13,atol=1e-13)
+
+
+def test_vector_parity_and_subspace_invariance():
+ x=np.linspace(0,4,9); y=np.linspace(-2,2,7); xx,yy=np.meshgrid(x,y,indexing="ij"); mask=np.ones(xx.shape,bool); weights=np.ones_like(xx)
+ velocity=np.stack((np.exp(-(xx-2)**2),np.zeros_like(xx))); diagnostics=symmetry_diagnostics(velocity,weights,mask,y)
+ assert diagnostics["symmetric_fraction"]==pytest.approx(1) and diagnostics["antisymmetric_fraction"]==pytest.approx(0)
+ columns=np.column_stack((np.r_[np.exp(-(xx-1.5)**2)[mask],(yy*np.exp(-(xx-2.5)**2))[mask]],np.r_[(yy*np.exp(-(xx-1.5)**2))[mask],np.exp(-(xx-2.5)**2)[mask]]))
+ component_weights=np.tile(weights[mask],2); q,_=np.linalg.qr(np.sqrt(component_weights)[:,None]*columns); modes=q/np.sqrt(component_weights)[:,None]
+ subspace=subspace_symmetry_diagnostics(modes,weights,mask,y); np.testing.assert_allclose(subspace["reflection_principal_cosines"],1); assert subspace["invariance_fraction"]==pytest.approx(1)
+ rotation=np.array([[1.,1.],[-1.,1.]])/np.sqrt(2); rotated=subspace_symmetry_diagnostics(modes@rotation,weights,mask,y)
+ assert rotated["symmetric_fraction"]==pytest.approx(subspace["symmetric_fraction"]); assert x_localization(np.exp(-8*(xx-2)**2),x,weights,mask)["centroid_x_D"]==pytest.approx(2)
+
+
+def test_harmonics_require_19_by_10_and_retain_each_cycle():
+ phase=2*np.pi*np.arange(10)/10; drift=np.arange(19)[:,None]*.02; signal=3*np.cos(phase[None,:]-.4+drift)+2*np.cos(2*phase[None,:]+.2)
+ harmonics=phase_harmonics(reshape_cycle_phase(signal.reshape(-1)),2)
+ assert harmonics["per_cycle_amplitude"].shape==(19,3) and harmonics["mean_cycle_amplitude"][1]==pytest.approx(3) and harmonics["coherent_amplitude"][1]<3 and harmonics["coherent_amplitude"][2]==pytest.approx(2)
+ with pytest.raises(ValueError,match="explicit"): phase_harmonics(np.zeros(190))
+ with pytest.raises(ValueError): reshape_cycle_phase(np.zeros(189))
+
+
+def test_action_transform_literal_definitions_and_inverse():
+ native=np.array([[2.,4.],[3.,5.],[1.,-1.]]); transformed=action_coordinates(native)
+ np.testing.assert_allclose(transformed[0],native[0]); np.testing.assert_allclose(transformed[1],(native[1]+native[2])/np.sqrt(2)); np.testing.assert_allclose(transformed[2],(native[1]-native[2])/np.sqrt(2)); np.testing.assert_allclose(action_coordinates(transformed,inverse=True),native)
+
+
+def test_shift_sign_weights_and_finite_guard():
+ x=np.arange(8.); a=np.zeros((2,8,3)); b=np.zeros_like(a); a[:,2]=np.array([[1.],[2.]]); b[:,4]=a[:,2]; weights=np.linspace(1,2,24).reshape(8,3); mask=np.ones((8,3),bool)
+ result=bounded_streamwise_shift_correlation(a,b,x,weights,mask,3); assert result["best_shift_x_D"]==pytest.approx(2) and result["correlation"]==pytest.approx(1) and result["index_offset"]==2
+ b[:,4,0]=np.nan; b[:,4,1:]=a[:,2,1:]; assert np.isfinite(bounded_streamwise_shift_correlation(a,b,x,weights,mask,3)["correlation"])
+
+
+def test_weighted_pod_residual_cosines_and_alignment():
+ rng=np.random.default_rng(4); snapshots=rng.normal(size=(12,3))@rng.normal(size=(3,20)); weights=np.linspace(.5,2,12); pod=fit_weighted_pod(snapshots,weights,3)
+ np.testing.assert_allclose(pod["modes"].T@(weights[:,None]*pod["modes"]),np.eye(3),atol=1e-12); assert np.all(np.diff(pod["weighted_residual_fraction"])<=0) and pod["weighted_residual_fraction"][-1]<1e-25; np.testing.assert_allclose(weighted_principal_cosines(pod["modes"],pod["modes"],weights),1)
+ rotation,_=np.linalg.qr(rng.normal(size=(3,3))); aligned,recovered=align_basis_to_reference(pod["modes"],pod["modes"]@rotation,weights); np.testing.assert_allclose(aligned,pod["modes"],atol=1e-12); np.testing.assert_allclose(recovered.T@recovered,np.eye(3),atol=1e-12)
+ with pytest.raises(ValueError,match="orthonormal"): weighted_principal_cosines(2*pod["modes"],pod["modes"],weights)
@@ -18,7 +18,7 @@ def test_gate_fails_closed_for_too_few_and_low_amplitude_cycle():
 def test_compact_role_artifact_cycle_balances_and_freezes_mean(tmp_path,monkeypatch):
 root,role=fake_role(tmp_path,monkeypatch,signal()); out=tmp_path/'derived'/'phase-v1'; result=publish_phase_compact(root,out); assert result['summary']['gate_passed'] and not result['summary']['dense_fields_deleted']; assert result['summary']['constant_mean_provenance']['symmetrized'] is False
 with np.load(out/'compact.npz') as z:
-  assert z['cycle_bin_fields'].shape[1:3]==(10,2); assert np.all(z['cycle_bin_counts']>0); assert z['ensemble_bins8'].shape[0]==8 and z['ensemble_bins12'].shape[0]==12
+  assert z['cycle_bin_fields'].shape[1:3]==(10,2); assert np.all(z['cycle_bin_counts']>0); assert z['ensemble_bins8'].shape[0]==8 and z['ensemble_bins12'].shape[0]==12; assert z['ensemble_bins8_effective_actions'].shape==(8,3) and z['ensemble_bins12_effective_actions'].shape==(12,3) and z['ensemble_bins10_half_shift_effective_actions'].shape==(10,3); assert z['ensemble_bins8_per_cycle_fields'].shape[0]==15 and z['ensemble_bins12_per_cycle_effective_actions'].shape==(15,12,3)
 with pytest.raises(FileExistsError): publish_phase_compact(root,out)
 def test_failed_gate_publishes_metrics_but_no_mean_or_fields(tmp_path,monkeypatch):
 root,_=fake_role(tmp_path,monkeypatch,signal(cycles=7)); out=tmp_path/'derived'/'failed'; r=publish_phase_compact(root,out); assert not r['summary']['gate_passed'] and r['summary']['constant_mean_provenance'] is None; assert not (out/'compact.npz').exists()
@@ -1,40 +0,0 @@
-import json
-from pathlib import Path
-import numpy as np
-import pytest
-from CCD_analysis.karman_dynamic.phase_domain_ccd import _decompose,_periodic_resample,decompose_phase_domain,PhaseDomainTransaction,load_phase_domain_result
-from CCD_analysis.original_ccd import CCDConfig,decompose
-
-def test_literal_weighted_reference_and_harmonic_resample():
- rng=np.random.default_rng(2); u=rng.normal(size=(9,10)); p=rng.normal(size=(3,10)); w=np.linspace(.5,1.5,9); got=_decompose(u,p,w); ref=decompose(u,p,weight=w,config=CCDConfig(center_snapshots=True,center_observables=True))
- np.testing.assert_allclose(got['cross'],ref.cross_correlation,rtol=1e-14,atol=1e-14); np.testing.assert_allclose(got['singular'],ref.singular_values,rtol=1e-14,atol=1e-14)
- theta=(np.arange(10)+.5)/10; curve=np.stack((np.cos(2*np.pi*theta),np.sin(4*np.pi*theta)),axis=1); np.testing.assert_allclose(_periodic_resample(curve,10,order=2),curve,atol=1e-14)
-
-def _synthetic(monkeypatch,tmp_path):
- from CCD_analysis.karman_dynamic import phase_domain_ccd as mod
- roots={k:tmp_path/k for k in ('drl','constant','dynamic')}
- for root in roots.values(): root.mkdir(); (root/'manifest.json').write_text('{}')
- mask=np.ones((2,2),bool); theta=(np.arange(10)+.5)*2*np.pi/10; spatial=np.arange(8).reshape(2,4)+1; delta=np.cos(theta)[:,None,None]*spatial[None]
- drl=np.zeros((2,10,2,4)); const=np.zeros_like(drl); drl[:]=delta[None]; actions=np.stack((np.cos(theta),np.sin(theta),np.cos(2*theta)),axis=1)
- def write_phase(root,role,fields):
-  vals={'fluid_mask':mask,'cycle_bin_fields':fields.astype(np.float32),'cycle_bin_effective_actions':np.broadcast_to(actions,(2,10,3)).astype(np.float32)}
-  for n in (8,12): vals[f'ensemble_bins{n}']=_periodic_resample(fields.mean(0),n).astype(np.float32)
-  vals['ensemble_bins10_half_shift']=_periodic_resample(fields.mean(0),10,.5).astype(np.float32); np.savez(root/'compact.npz',**vals)
- write_phase(roots['drl'],'drl',drl); write_phase(roots['constant'],'constant_mean',const)
- means=np.zeros((2,4),np.float32); np.savez(roots['dynamic']/'arrays.npz',four_role_fluid_mask=mask,quadrature_weights=np.ones(4),mean_drl=means,mean_constant_mean=means)
- monkeypatch.setattr(mod,'load_phase_compact',lambda path:{'summary':{'gate_passed':True,'role':'drl' if Path(path)==roots['drl'] else 'constant_mean'}}); monkeypatch.setattr(mod,'load_dynamic_increment',lambda path:{'summary':{}}); monkeypatch.setattr(mod,'file_sha256',lambda path:'a'*64)
- return roots
-
-def test_synthetic_estimand_sensitivity_and_reference(monkeypatch,tmp_path):
- roots=_synthetic(monkeypatch,tmp_path); r=decompose_phase_domain(roots['drl'],roots['constant'],roots['dynamic']); assert r.config['field_estimand'].startswith('separately_centered'); assert r.config['phase_offset_semantics'].startswith('circular phase offsets'); assert r.summary['numerical_rank']>=1; assert len(r.summary['circular_phase_offsets'])==10
- ref=decompose(r.arrays['primary_modes'].T*0 if False else np.zeros((1,1)),np.zeros((1,1)),config=CCDConfig()) if False else None
- np.testing.assert_allclose(r.arrays['primary_singular'],np.linalg.svd(r.arrays['primary_cross'],compute_uv=False))
-
-def test_immutable_reload_recomputes_and_tamper_rejected(monkeypatch,tmp_path):
- roots=_synthetic(monkeypatch,tmp_path); result=decompose_phase_domain(roots['drl'],roots['constant'],roots['dynamic']); out=tmp_path/'result'
- with PhaseDomainTransaction(out) as tx: tx.write(result); tx.publish()
- loaded=load_phase_domain_result(out,recompute=True); assert loaded['summary']['scientific_contract_explicit']
- with pytest.raises(FileExistsError):
-  with PhaseDomainTransaction(out): pass
- (out/'summary.json').write_text('{}')
- with pytest.raises(ValueError,match='hash|contract'): load_phase_domain_result(out)
@@ -1,10 +1,44 @@
 import json
 from pathlib import Path
 import pytest
-from CCD_analysis.karman_dynamic.publication import SCHEMA_ID,load_dynamic_publication
-def test_publication_loader_hash_source_and_claim_contract(tmp_path,monkeypatch):
- root=tmp_path/"publication"; root.mkdir(); source=tmp_path/"source"; source.mkdir(); (source/"manifest.json").write_text("source")
- report={"schema_id":SCHEMA_ID,"phase_domain_decision":"DOWNGRADE","zero_phase_output":"PROHIBITED: failed gate","sources":{"source":{"path":str(source),"manifest_sha256":"source-hash"}}}; (root/"RESULTS.json").write_text(json.dumps(report)); (root/"manifest.json").write_text(json.dumps({"schema_id":SCHEMA_ID,"complete":True,"files":{"RESULTS.json":"result-hash"}}))
- monkeypatch.setattr("CCD_analysis.karman_dynamic.publication.file_sha256",lambda path:"result-hash" if Path(path).name=="RESULTS.json" else "source-hash"); assert load_dynamic_publication(root)["report"]["phase_domain_decision"]=="DOWNGRADE"
- report["zero_phase_output"]="available"; (root/"RESULTS.json").write_text(json.dumps(report))
- with pytest.raises(ValueError,match="claim"): load_dynamic_publication(root)
+from CCD_analysis.karman_dynamic.publication import FIGURE_FILES,SCHEMA_ID,_stability,load_dynamic_publication
+
+
+def _write_publication(root,sources):
+ root.mkdir(); report={'schema_id':SCHEMA_ID,'figure_files':list(FIGURE_FILES),'full_domain_metrics':'secondary only; retained numerically','sources':sources}
+ files={name:f'hash-{name}' for name in (*FIGURE_FILES,'RESULTS.md')}
+ for name in files: (root/name).write_text(name)
+ (root/'RESULTS.json').write_text(json.dumps(report)); files['RESULTS.json']='hash-RESULTS.json'
+ (root/'manifest.json').write_text(json.dumps({'schema_id':SCHEMA_ID,'complete':True,'files':files}))
+
+
+def test_publication_loader_requires_complete_inventory_and_live_corrected_sources(tmp_path,monkeypatch):
+ sources={}
+ for name in ('roi_mean','cycle_template_ccd'):
+  q=tmp_path/name; q.mkdir(); (q/'manifest.json').write_text(name); sources[name]={'path':str(q),'manifest_sha256':f'source-{name}'}
+ root=tmp_path/'publication'; _write_publication(root,sources)
+ def fake_hash(path):
+  path=Path(path)
+  if path.name=='manifest.json' and path.parent.name in sources: return f'source-{path.parent.name}'
+  return f'hash-{path.name}'
+ monkeypatch.setattr('CCD_analysis.karman_dynamic.publication.file_sha256',fake_hash)
+ monkeypatch.setattr('CCD_analysis.karman_dynamic.publication.load_dynamic_increment',lambda path:{'summary':{}})
+ monkeypatch.setattr('CCD_analysis.karman_dynamic.publication.load_cycle_template_result',lambda path,recompute=True:{'summary':{},'recompute':recompute})
+ loaded=load_dynamic_publication(root); assert set(loaded['report']['sources'])=={'roi_mean','cycle_template_ccd'}; assert loaded['report']['figure_files']==list(FIGURE_FILES); assert len(FIGURE_FILES)==12
+ manifest=json.loads((root/'manifest.json').read_text()); removed=next(iter(FIGURE_FILES)); manifest['files'].pop(removed); (root/'manifest.json').write_text(json.dumps(manifest))
+ with pytest.raises(ValueError,match='inventory'): load_dynamic_publication(root)
+
+
+def test_publication_loader_rejects_obsolete_source_contract(tmp_path,monkeypatch):
+ source=tmp_path/'source'; source.mkdir(); (source/'manifest.json').write_text('source')
+ root=tmp_path/'publication'; _write_publication(root,{'roi_mean':{'path':str(source),'manifest_sha256':'source-hash'}})
+ monkeypatch.setattr('CCD_analysis.karman_dynamic.publication.file_sha256',lambda path:f'hash-{Path(path).name}')
+ with pytest.raises(ValueError,match='claim'): load_dynamic_publication(root)
+
+
+def test_stability_discloses_all_categories_global_minimum_and_singular_shift():
+ def item(cosine,spectrum): return {'rank':3,'minimum_projector_cosine':cosine,'leading_spectrum':spectrum}
+ summary={'primary_singular_values':[4.,2.,1.],'sensitivities':{'drl_leave_one_cycle_out':[{'omitted_cycle':0,**item(.99,[4.,2.,1.])}],'constant_template_leave_one_cycle_out':[{'omitted_cycle':1,**item(.98,[4.1,2.,1.])}],'drl_splits':{'odd':item(.97,[3.8,2.,1.])},'constant_template_splits':{'suffix':item(.96,[4.,1.8,1.])},'bin_and_origin':{'bins8':item(.91,[3.,2.,1.])}}}
+ plot,report=_stability(summary); assert set(plot)=={'drl_leave_one_cycle_out','constant_template_leave_one_cycle_out','drl_splits','constant_template_splits','bin_and_origin'}
+ assert report['global_least_favorable_projector_check']['category']=='bin_and_origin'; assert report['global_least_favorable_projector_check']['minimum_projector_cosine']==.91
+ assert report['global_largest_singular_shift_check']['category']=='bin_and_origin'; assert report['scope'].startswith('all declared')
@@ -1,39 +0,0 @@
-import numpy as np
-import pytest
-from CCD_analysis.karman_dynamic.temporal_ccd import _admit, _compare, PRIMARY_LAGS, ROW_ORDER
-from CCD_analysis.original_ccd import CCDConfig, decompose
-
-def test_block_local_negative_lags_endpoint_sign_and_no_crossing():
- ids=np.array([-1,0,0,0,1,1,1,-1]); fi,oi=_admit(ids,(-2,-1,0))
- np.testing.assert_array_equal(fi,[3,6]); np.testing.assert_array_equal(oi,[[1,2,3],[4,5,6]])
- assert np.all(oi<=fi[:,None]) and ROW_ORDER=='channel-major_delay-minor'
- with pytest.raises(ValueError): _admit(ids,(0,1))
-
-def test_common_support_and_row_order_dense_original_reference():
- rng=np.random.default_rng(4); ids=np.repeat(np.arange(3),6); lags=(-2,-1,0); fi,oi=_admit(ids,lags); u=rng.normal(size=(7,len(ids))); actions=rng.normal(size=(len(ids),3)); p=actions[oi].transpose(2,1,0).reshape(9,len(fi)); w=np.linspace(.5,1.5,7)
- dense=decompose(u[:,fi],p.reshape(3,3,len(fi)),weight=w,config=CCDConfig(center_snapshots=True,center_observables=True))
- uc=u[:,fi]-u[:,fi].mean(1,keepdims=True); pc=p-p.mean(1,keepdims=True); expected=pc@(uc*np.sqrt(w)[:,None]).T/(len(fi)*np.sqrt(9))
- np.testing.assert_allclose(dense.cross_correlation,expected,rtol=1e-14,atol=1e-14)
- fi2,_=_admit(ids,(-1,0),support=fi); np.testing.assert_array_equal(fi2,fi)
-
-def test_weighted_projector_comparison_is_identity_for_same_modes():
- rng=np.random.default_rng(1); w=np.linspace(.5,1.5,6); q,_=np.linalg.qr(rng.normal(size=(6,3))*np.sqrt(w)[:,None]); modes=q/np.sqrt(w)[:,None]; result={'singular_values':np.array([3.,2.,1.]),'physical_modes':modes}
- got=_compare(result,result,w); np.testing.assert_allclose(got['leading_weighted_subspace_principal_cosines'],np.ones(3),rtol=1e-14,atol=1e-14)
-
-
-def _synthetic_input(tmp_path):
- from CCD_analysis.karman_dynamic.temporal_ccd import TemporalInput
- rng=np.random.default_rng(8); nx,ny=4,3; n=57; actions=rng.normal(size=(n,3)).astype(np.float32); fields=rng.normal(size=(n,2,nx,ny)).astype(np.float32)
- ids=np.repeat(np.arange(3),19).astype(np.int64); rel=np.arange(n,dtype=np.int64)*800
- return TemporalInput(tmp_path/'role',tmp_path/'phase','r','p',np.arange(nx,dtype=np.float32),np.arange(ny,dtype=np.float32),np.ones((nx,ny),bool),fields,actions,rel,rel+1000,ids)
-
-def test_streaming_chunk_invariance_memory_and_literal_dense(tmp_path):
- from CCD_analysis.karman_dynamic.temporal_ccd import TemporalConfig,decompose_temporal
- inp=_synthetic_input(tmp_path); r1=decompose_temporal(inp,streaming_config=TemporalConfig(1,10_000_000)); r2=decompose_temporal(inp,streaming_config=TemporalConfig(7,10_000_000))
- for key in ('primary_cross_correlation','primary_singular_values'):
-  np.testing.assert_allclose(r1.arrays[key],r2.arrays[key],rtol=1e-13,atol=1e-13)
- identifiable=int(np.sum(r1.arrays['primary_singular_values']>1e-10*r1.arrays['primary_singular_values'][0])); w=r1.arrays['coordinate_weights']; v1=r1.arrays['primary_physical_modes'][:,:identifiable]*np.sqrt(w)[:,None]; v2=r2.arrays['primary_physical_modes'][:,:identifiable]*np.sqrt(w)[:,None]
- np.testing.assert_allclose(v1@v1.T,v2@v2.T,rtol=1e-12,atol=1e-12); np.testing.assert_allclose(r1.arrays['primary_coefficients'][:identifiable],r2.arrays['primary_coefficients'][:identifiable],rtol=1e-11,atol=1e-11)
- fi,oi=_admit(inp.cycle_ids,PRIMARY_LAGS); u=np.concatenate((inp.fields[:,0].reshape(len(inp.fields),-1),inp.fields[:,1].reshape(len(inp.fields),-1)),axis=1).T.astype(float); p=inp.actions[oi].transpose(2,1,0).reshape(54,len(fi)); w=r1.arrays['coordinate_weights']; dense=decompose(u[:,fi]-u.mean(1,keepdims=True),p,weight=w,config=CCDConfig(center_snapshots=False,center_observables=True))
- np.testing.assert_allclose(r1.arrays['primary_cross_correlation'],dense.cross_correlation,rtol=1e-13,atol=1e-13)
- with pytest.raises(MemoryError): decompose_temporal(inp,streaming_config=TemporalConfig(2,100))