[BUG] cpp backend returns sensor `p` transposed when n_sensor < Nt

[waltsims]

Summary

For the C++ backend (backend="cpp"), kspaceFirstOrder() returns sensor pressure with axes swapped relative to the docstring contract and the Python backend.

Expected

Per the unified API docstring (kwave/kspaceFirstOrder.py:218), result['p'] is shaped (n_sensor, Nt) regardless of backend.

Actual

On the cpp backend, when n_sensor < Nt the returned p is (Nt, n_sensor). The Python backend returns the documented (n_sensor, Nt).

Discovered while running a cpp-vs-python numerical comparison: same kgrid, medium, source, sensor produced sensor arrays of shape (200, 10) from cpp and (10, 200) from python. Both backends agree numerically once one is transposed — see #759 (sub-thread) for the diff results — but the orientation mismatch will silently corrupt any downstream code that broadcasts, indexes, or plots the cpp result assuming the documented (n_sensor, Nt) layout.

Reproduction

from kwave.kspaceFirstOrder import kspaceFirstOrder
# ... build kgrid (Nx=128, Ny=128), homogeneous medium, smoothed-disk p0, 10-pt line sensor ...
# Nt = 200 (kgrid.makeTime → setTime)

res_cpp = kspaceFirstOrder(kgrid, medium, source, sensor, backend="cpp", device="cpu", quiet=True)
res_py  = kspaceFirstOrder(kgrid, medium, source, sensor, backend="python", device="cpu", quiet=True)

print(res_cpp['p'].shape)  # (200, 10)  ← wrong; should be (10, 200)
print(res_py['p'].shape)   # (10, 200)  ✓ matches docstring

Full comparison script saved at /tmp/cpp-vs-py-diff.py during the test that found this; happy to attach if useful.

Root cause

kwave/solvers/cpp_simulation.py:90-113 (_fix_output_order) attempts to detect F-order sensor data and reorder rows from F-indexed to C-indexed. The trigger condition keys off val.shape[0] == n_sensor:

if is_sensor and val.ndim == 2 and val.shape[0] == n_sensor:
    # reorder rows ...

This is always false when n_sensor < Nt because the cpp binary writes the array as (Nt, n_sensor) (HDF5 row-major from the binary's column-major F-order writes), so shape[0] is Nt, not n_sensor. The reorder is then skipped and the array is returned as (Nt, n_sensor). The reorder logic also doesn't transpose — only permutes rows — so even if the trigger fired correctly it wouldn't fix the axis order.

Fix sketch

The condition should detect cpp's actual layout (likely val.shape[1] == n_sensor when the array is (Nt, n_sensor)) and transpose, then apply the row-reorder. Roughly:

if is_sensor and val.ndim == 2:
    if val.shape[1] == n_sensor:    # cpp's (Nt, n_sensor) layout
        val = val.T                  # → (n_sensor, Nt)
    if val.shape[0] == n_sensor:
        # existing F→C row reorder
        ...
    result[key] = val

Needs confirmation that the F-order row-reorder is still needed after the transpose; the two corrections may compose differently than the original.

Scope

Latent bug (predates v0.6.3 — not a regression). Doesn't block the v0.6.3rc1 → v0.6.3 promotion (#759); should be tracked as a separate fix-forward.

Impacts: any user of backend="cpp" whose sensor has fewer points than timesteps — i.e. essentially every real simulation. Downstream code that indexes p[i, t] for sensor i at time t would see p[t, i] instead, with no error message.

Test plan for the fix

Add a 2D smoke test that runs both backends on a small grid with n_sensor=10, Nt=200, and asserts result['p'].shape == (n_sensor, Nt) for both. The cpp-vs-py numerical diff (after the fix) should agree to ~1e-6 L∞-relative (float32 noise floor — see #759 sub-thread).