Derive the acquisition-side geometry: acquisition_affine(s) for fid/rawdata

[gdevenyi]

Fixes #166. Investigated against FILE_FORMAT.md §5.4 / §5.6 / §12 (and §5.2, §7.2, §4.2) and against the full local corpus — the spec's coordinate-system section was revised from the same evidence, and this PR follows it; the maintainer's request for documented method-level fallbacks is in.

API

raw = Dataset("…/5/fid", load=LOAD_STAGES["properties"])
raw.acquisition_affines()      # [4x4] one per slice (NSLICES), slice-id order
raw.acquisition_affine(3)      # slice 3

Same Visu/DICOM patient frame as Dataset.affine of the 2dseq, so the two are directly comparable: ACQ_grad_matrix (read/phase/slice direction cosines, C-order per spec 2.3), ACQ_read_offset / ACQ_phase1_offset / ACQ_slice_offset (slice centre), ACQ_fov (extent), ACQ_slice_sepn/ACQ_slice_thick (lone-slice step). Voxel (i, j[, k]) of the Fourier-transformed k-space matrix (PVM_EncMatrix points per axis) is at affine @ (i, j, k, 1); index N/2 is the field-of-view centre and the in-plane image axes run against the gradient directions, as ParaVision's own reconstruction has it. A 3-D volume gets one affine whose third column is the partition step. The slice column of a 2-D slice is the step to its neighbour in the same package (physical progression); for a lone slice it is the slice normal — which Visu instead completes right-handedly, so that sign can differ from VisuCoreOrientation[2] while the positions agree.

acqp first, method as fallback (per your comment): PVM_SPackArrGradOrient (repeated per package by PVM_SPackArrNSlices), PVM_SPackArr{Read,Phase1}Offset, PVM_SPackArrSliceOffset spread by PVM_SPackArrSliceDistance, PVM_Fov, PVM_SliceThick — each with a RuntimeWarning naming the substitution. Caveat, stated in the code: in the corpus PVM_SPackArrGradOrient equals ACQ_grad_matrix on every PV6/PV7 scan, but differs by a read-axis sign on 28/105 PV5.1 and 57/57 PV360 scans, so the orientation fallback is trustworthy only where acqp would have agreed. In 1,543 corpus experiments no acqp geometry parameter other than ACQ_slice_sepn (PV360 MSME/MGE/RARE) was ever missing, so the fallback is there for the custom-sequence case you described, not because the corpus needed it.

Spectroscopic acquisitions (ACQ_dim 1, or any ACQ_dim_desc ≠ Spatial) raise UnsupportedDatasetType. Reports are unchanged (not added to COMPUTED_REPORT_PROPERTIES — every fid reference would change; say so if you want it there).

Per slice, not per object. The issue proposed acquisition_affine(object=k). The NI objects map onto slices by a method-dependent nesting that only the reconstruction's VisuFGOrderDesc records: MGE is slice-major (obj = slice·NE + echo, ACQ_obj_order = [63 64 … 71, 45 …]), T1_EPI slice-minor (obj = movie·NSLICES + slice, [0 6 12 …]), and DtiStandard's sequential [0 1 2 … 13] is slice-minor too — indistinguishable from slice-major by acqp alone. So the API is per slice; for the common NI == NSLICES object k is slice k after brukerapi's ACQ_obj_order reordering.

What the corpus says — now written into FILE_FORMAT.md §5.4 / §5.6 / §12

Fitted and validated on 1,507 acquisitions with a reconstruction — PV5.1, 6.0, 6.0.1, 7.0.0, 360.3.4–3.7; 1,341 Head_Supine, 221 Head_Prone; the full local corpus — by searching, per dataset, the 48 signed-permutation maps P for which VisuCoreOrientation = S · ACQ_grad_matrix · Pᵀ with S a signed permutation, then resolving signs with VisuCorePosition.

1. Which frame the matrix is in depends on the version (§5.4, §12). On PV5.1/6/7 ACQ_grad_matrix is written in the ParaVision subject frame (L→R, P→A, F→H): Head_Supine and Head_Prone admit exactly the same map in every generation (PV6: 1150 and 97 scans; PV5: 86/19; PV7: 69/29) — P = VISU_DICOM_PV_MATRIX = diag(−1,−1,1), Visu only converts to DICOM and the position cancels. On PV360 the matrix is in the magnet frame ("transforms the logical gradient orientation directly to the physical gradient vector") and the declared position's magnet→subject map comes first, P = VISU_DICOM_PV_MATRIX · M_pos with M_pos the manual's table read as subject = M_pos · magnet (ParaVision's GTB_ObjPosMatrix, §5.6): Head_Prone → diag(−1,1,−1), which fits every PV360 Head_Prone scan; Head_Supine → diag(1,−1,−1), which fits the one PV360 Head_Supine scan (a compressed-sensing 3-D FLASH; it was the single outlier of the first version of this PR, and its VisuCorePosition is reproduced to 1e-5 mm once its own position is used). So a reader applies ACQ_patient_pos once on PV360 and never on PV5.1–7. Head_Left/Right and Foot_* do not occur in 1,562 acqp files; their rows of the table are implemented from the manual but untested, and a PV360 scan with an unknown position raises rather than guessing. The frame is that of the declared position, as Visu's is. ACQ_GradientMatrix == ACQ_grad_matrix in all 57 PV360 scans, so one code path.
2. ACQ_grad_matrix is indexed by acquisition position, the offset arrays by slice id (§5.4). The three-package scouts (FLASHAngio/tripilot, 15 slices, ACQ_obj_order = [0 2 4 6 8 10 12 14 1 3 5 7 9 11 13]) prove it: read per slice id, the sagittal/coronal slices land 5.66 mm off (|(4,4)|); re-indexed with argsort(ACQ_obj_order) all 15 are exact (0.000 mm). Once NI > NSLICES, ACQ_obj_order counts objects; PVM_ObjOrderList (method) is the slice order — equal to ACQ_obj_order in all 421 NI == NSLICES scans, and in the 57 NI > NSLICES 2-D scans the matrices are all identical anyway.
3. Image-axis convention (§12). S = diag(−1,−1,1) for every untransposed reconstruction and [[0,−1,0],[−1,0,0],[0,0,±1]] for RECO_transposition = 1: the read and phase image axes run against the mapped gradient directions, and the first pixel is at exactly −FOV/2 (index N/2 at the centre, N-independent) — so VisuCorePosition = P·c + (FOV_r/2)·P r̂ + (FOV_p/2)·P p̂, c the offsets along ACQ_grad_matrix. PV360's ACQ_AtsCenterDistance (356/369 mm in the corpus) does not enter: the acqp offsets reproduce VisuCorePosition without it.
4. 3-D partition grid (§7.2, §12). PV5.1 places partition N/2 at the centre like the in-plane axes; PV6/7/360 centre the grid between partitions — every PV6+ 3-D volume was off by exactly half a partition under the DFT convention (0.0625 mm at 12 mm/96, 0.625 mm at 40 mm/32, …), PV5.1 by none. Modelled on pv_version.

Result with the implemented method (Dataset.acquisition_affines() vs the 2dseq, same checks as the new test): 2-D — PV5.1 87/87, PV6 1203/1203, PV7 75/75, PV360 47/47 slice positions exact at 1e-3 mm, in-plane axes equal up to RECO_transposition; 3-D — 91/91 within half a partition (PV5 18, PV6 40, PV7 23, PV360 10). Positions compare after accounting for ParaVision's cropping of an anti-aliased field of view (PVM_AntiAlias: the k-space image is the acquired ACQ_fov, the 2dseq the cropped VisuCoreExtent; the centre is shared).

Tests

test_geometry.py: synthetic acquisitions built from real scans' numbers (write_acquisition helper) — PV5.1 0.2H2/13 (interleaved EPI: every slice on its VisuCorePosition), PV6 lego phantom/3 (axial/sagittal/coronal acquired 0,2,1 — fails if the matrix is read per slice id), PV5.1 0.2H2/1 package offsets, PV5.1 0.2H2/10 vs PV6 3-D grid, PV360 T2_TurboRARE Head_Prone (tilted, offset: (10.1949, 10.9375, −5.4053) ✓; declared Head_Supine it comes out rotated by π about the bore, an unknown position raises), the PV360 Head_Supine CS volume (VisuCorePosition at the far end of the reversed slab ✓), the method fallback with its warning, and the spectroscopic refusal. Plus test_acquisition_affine_agrees_with_the_reconstruction, parametrized over every fid/rawdata.job0 of the CI corpus that has a pdata/1/2dseq (129 checked, 22 spectroscopic skips) — the issue's acceptance test, run in CI.

Suite: 2358 passed, 26 skipped (master: 2251 / 4).

🤖 Generated with Claude Code

https://claude.ai/code/session_01NHhfFXMNZKMuPm3ppdzFXe