diff --git a/src/highdicom/spatial.py b/src/highdicom/spatial.py index 043f6eb0..97d01d37 100644 --- a/src/highdicom/spatial.py +++ b/src/highdicom/spatial.py @@ -3859,19 +3859,24 @@ def get_volume_positions( (origin_distances - origin_distances.min()) / spacing ) - is_regular = np.allclose( - origin_distance_multiples, - origin_distance_multiples.round(), - rtol=rtol, - atol=atol, + # Deviations are measured in mm, against the position each slice would + # occupy if it sat exactly on its integer multiple of the spacing. + # Comparing the multiples themselves would make the tolerances mean + # "fractions of a slice" here but "mm" in the branch below, so the same + # 'atol' would accept a larger physical error the wider the spacing. + expected_distances = ( + origin_distance_multiples.round() * spacing + + origin_distances.min() ) + deviations = np.abs(origin_distances - expected_distances) + + # 'rtol' is relative to the slice spacing, matching its meaning in the + # regularly-spaced branch below, where each spacing is compared with + # the nominal spacing. + tolerance = atol + rtol * abs(spacing) + is_regular = bool((deviations <= tolerance).all()) if not is_regular: - max_deviation = float( - np.abs( - origin_distance_multiples - - origin_distance_multiples.round() - ).max() - ) + max_deviation = float(deviations.max()) logger.info("Frame positions are not regularly spaced.") logger.debug( f"Maximum spacing deviation from regular: {max_deviation:.2e} " diff --git a/tests/test_spatial.py b/tests/test_spatial.py index 8b46b411..72c95696 100644 --- a/tests/test_spatial.py +++ b/tests/test_spatial.py @@ -1047,6 +1047,46 @@ def test_get_volume_positions_missing(): assert volume_positions is None +def test_get_volume_positions_missing_tolerance_units(): + # 'atol' is documented in mm and must mean the same thing whatever the + # slice spacing is, including when missing slices are allowed. Previously + # the deviation was measured in multiples of the spacing, so the same + # physical error was accepted at wide spacings and rejected at narrow ones. + orientation = [1.0, 0.0, 0.0, 0.0, 1.0, 0.0] + error = 0.5 # mm, out of place by more than atol below + + for spacing in [1.0, 5.0, 10.0, 20.0]: + positions = [ + [0.0, 0.0, 0.0], + [0.0, 0.0, spacing], + [0.0, 0.0, 2 * spacing + error], + ] + out_spacing, volume_positions = get_volume_positions( + positions, + orientation, + allow_missing_positions=True, + atol=0.1, + ) + assert out_spacing is None + assert volume_positions is None + + # A deviation smaller than 'atol' remains acceptable at every spacing + for spacing in [1.0, 5.0, 10.0, 20.0]: + positions = [ + [0.0, 0.0, 0.0], + [0.0, 0.0, spacing], + [0.0, 0.0, 2 * spacing + 0.05], + ] + out_spacing, volume_positions = get_volume_positions( + positions, + orientation, + allow_missing_positions=True, + atol=0.1, + ) + assert out_spacing is not None + assert volume_positions is not None + + def test_get_volume_positions_missing_duplicates(): # Test ability to determine spacing and volume positions with missing # slices and duplicate positions