Skip to content

Commit 7e522a2

Browse files
committed
Merge pull request #269 from matthew-brett/parrec-arrayproxy-api
MRG: test and extend PARREC array proxy I forgot to add the array proxy test for PARRECArrayProxy. To make the tests pass, add an interface to return slices from the proxied array. Optimize this a little using fileslice and some fancy broadcasting tricks.
2 parents 981db47 + 842fb0d commit 7e522a2

File tree

7 files changed

+242
-1
lines changed

7 files changed

+242
-1
lines changed

nibabel/fileslice.py

Lines changed: 23 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -747,3 +747,26 @@ def fileslice(fileobj, sliceobj, shape, dtype, offset=0, order='C',
747747
bytes = read_segments(fileobj, segments, n_bytes)
748748
sliced = np.ndarray(sliced_shape, dtype, buffer=bytes, order=order)
749749
return sliced[post_slicers]
750+
751+
752+
def strided_scalar(shape, scalar=0.):
753+
""" Return array shape `shape` where all entries point to value `scalar`
754+
755+
Parameters
756+
----------
757+
shape : sequence
758+
Shape of output array.
759+
scalar : scalar
760+
Scalar value with which to fill array.
761+
762+
Returns
763+
-------
764+
strided_arr : array
765+
Array of shape `shape` for which all values == `scalar`, built by
766+
setting all strides of `strided_arr` to 0, so the scalar is broadcast
767+
out to the full array `shape`.
768+
"""
769+
shape = tuple(shape)
770+
scalar = np.array(scalar)
771+
strides = [0] * len(shape)
772+
return np.lib.stride_tricks.as_strided(scalar, shape, strides)

nibabel/parrec.py

Lines changed: 17 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -102,6 +102,7 @@
102102
from .volumeutils import Recoder, array_from_file, BinOpener
103103
from .affines import from_matvec, dot_reduce, apply_affine
104104
from .nifti1 import unit_codes
105+
from .fileslice import fileslice, strided_scalar
105106

106107
# PSL to RAS affine
107108
PSL_TO_RAS = np.array([[0, 0, -1, 0], # L -> R
@@ -561,6 +562,22 @@ def __array__(self):
561562
scalings=self._slice_scaling,
562563
mmap=self._mmap)
563564

565+
def __getitem__(self, slicer):
566+
indices = self._slice_indices
567+
if indices[0] != 0 or np.any(np.diff(indices) != 1):
568+
# We can't load direct from REC file, use inefficient slicing
569+
return np.asanyarray(self)[slicer]
570+
# Slices all sequential from zero, can use fileslice
571+
# This gives more efficient volume by volume loading, for example
572+
with BinOpener(self.file_like) as fileobj:
573+
raw_data = fileslice(fileobj, slicer, self._shape, self._dtype, 0, 'F')
574+
# Broadcast scaling to shape of original data
575+
slopes, inters = self._slice_scaling
576+
fake_data = strided_scalar(self._shape)
577+
_, slopes, inters = np.broadcast_arrays(fake_data, slopes, inters)
578+
# Slice scaling to give output shape
579+
return raw_data * slopes[slicer] + inters[slicer]
580+
564581

565582
class PARRECHeader(Header):
566583
"""PAR/REC header"""
Lines changed: 129 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,129 @@
1+
# === DATA DESCRIPTION FILE ======================================================
2+
#
3+
# CAUTION - Investigational device.
4+
# Limited by Federal Law to investigational use.
5+
#
6+
# Dataset name: E:\\Export\phantom_EPI_asc_CLEAR_2_1
7+
#
8+
# CLINICAL TRYOUT Research image export tool V4.2
9+
#
10+
# === GENERAL INFORMATION ========================================================
11+
#
12+
. Patient name : phantom
13+
. Examination name : Konvertertest
14+
. Protocol name : EPI_asc CLEAR
15+
. Examination date/time : 2014.02.14 / 09:00:57
16+
. Series Type : Image MRSERIES
17+
. Acquisition nr : 2
18+
. Reconstruction nr : 1
19+
. Scan Duration [sec] : 14
20+
. Max. number of cardiac phases : 1
21+
. Max. number of echoes : 1
22+
. Max. number of slices/locations : 9
23+
. Max. number of dynamics : 3
24+
. Max. number of mixes : 1
25+
. Patient position : Head First Supine
26+
. Preparation direction : Anterior-Posterior
27+
. Technique : FEEPI
28+
. Scan resolution (x, y) : 64 39
29+
. Scan mode : MS
30+
. Repetition time [ms] : 2000.000
31+
. FOV (ap,fh,rl) [mm] : 240.000 70.000 240.000
32+
. Water Fat shift [pixels] : 11.050
33+
. Angulation midslice(ap,fh,rl)[degr]: -13.265 0.000 0.000
34+
. Off Centre midslice(ap,fh,rl) [mm] : 2.508 30.339 -16.032
35+
. Flow compensation <0=no 1=yes> ? : 0
36+
. Presaturation <0=no 1=yes> ? : 0
37+
. Phase encoding velocity [cm/sec] : 0.000000 0.000000 0.000000
38+
. MTC <0=no 1=yes> ? : 0
39+
. SPIR <0=no 1=yes> ? : 1
40+
. EPI factor <0,1=no EPI> : 39
41+
. Dynamic scan <0=no 1=yes> ? : 1
42+
. Diffusion <0=no 1=yes> ? : 0
43+
. Diffusion echo time [ms] : 0.0000
44+
. Max. number of diffusion values : 1
45+
. Max. number of gradient orients : 1
46+
. Number of label types <0=no ASL> : 0
47+
#
48+
# === PIXEL VALUES =============================================================
49+
# PV = pixel value in REC file, FP = floating point value, DV = displayed value on console
50+
# RS = rescale slope, RI = rescale intercept, SS = scale slope
51+
# DV = PV * RS + RI FP = DV / (RS * SS)
52+
#
53+
# === IMAGE INFORMATION DEFINITION =============================================
54+
# The rest of this file contains ONE line per image, this line contains the following information:
55+
#
56+
# slice number (integer)
57+
# echo number (integer)
58+
# dynamic scan number (integer)
59+
# cardiac phase number (integer)
60+
# image_type_mr (integer)
61+
# scanning sequence (integer)
62+
# index in REC file (in images) (integer)
63+
# image pixel size (in bits) (integer)
64+
# scan percentage (integer)
65+
# recon resolution (x y) (2*integer)
66+
# rescale intercept (float)
67+
# rescale slope (float)
68+
# scale slope (float)
69+
# window center (integer)
70+
# window width (integer)
71+
# image angulation (ap,fh,rl in degrees ) (3*float)
72+
# image offcentre (ap,fh,rl in mm ) (3*float)
73+
# slice thickness (in mm ) (float)
74+
# slice gap (in mm ) (float)
75+
# image_display_orientation (integer)
76+
# slice orientation ( TRA/SAG/COR ) (integer)
77+
# fmri_status_indication (integer)
78+
# image_type_ed_es (end diast/end syst) (integer)
79+
# pixel spacing (x,y) (in mm) (2*float)
80+
# echo_time (float)
81+
# dyn_scan_begin_time (float)
82+
# trigger_time (float)
83+
# diffusion_b_factor (float)
84+
# number of averages (integer)
85+
# image_flip_angle (in degrees) (float)
86+
# cardiac frequency (bpm) (integer)
87+
# minimum RR-interval (in ms) (integer)
88+
# maximum RR-interval (in ms) (integer)
89+
# TURBO factor <0=no turbo> (integer)
90+
# Inversion delay (in ms) (float)
91+
# diffusion b value number (imagekey!) (integer)
92+
# gradient orientation number (imagekey!) (integer)
93+
# contrast type (string)
94+
# diffusion anisotropy type (string)
95+
# diffusion (ap, fh, rl) (3*float)
96+
# label type (ASL) (imagekey!) (integer)
97+
#
98+
# === IMAGE INFORMATION ==========================================================
99+
# sl ec dyn ph ty idx pix scan% rec size (re)scale window angulation offcentre thick gap info spacing echo dtime ttime diff avg flip freq RR-int turbo delay b grad cont anis diffusion L.ty
100+
101+
1 1 1 1 0 2 0 16 62 64 64 -0.69352 0.65184 4.78462e-003 1070 1860 -13.26 -0.00 -0.00 2.51 -0.81 -8.69 6.000 2.000 0 1 0 2 3.750 3.750 30.00 0.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
102+
2 1 1 1 0 2 1 16 62 64 64 -0.28395 0.18528 5.01287e-003 1122 1951 -13.26 -0.00 -0.00 2.51 6.98 -10.53 6.000 2.000 0 1 0 2 3.750 3.750 30.00 0.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
103+
3 1 1 1 0 2 2 16 62 64 64 0.60964 0.48545 3.20504e-003 1137 1977 -13.26 -0.00 -0.00 2.51 14.77 -12.36 6.000 2.000 0 1 0 2 3.750 3.750 30.00 0.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
104+
4 1 1 1 0 2 3 16 62 64 64 0.23940 0.58869 5.70364e-003 1217 2116 -13.26 -0.00 -0.00 2.51 22.55 -14.20 6.000 2.000 0 1 0 2 3.750 3.750 30.00 0.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
105+
5 1 1 1 0 2 4 16 62 64 64 -0.56774 0.54425 4.19357e-003 1216 2113 -13.26 -0.00 -0.00 2.51 30.34 -16.03 6.000 2.000 0 1 0 2 3.750 3.750 30.00 0.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
106+
6 1 1 1 0 2 5 16 62 64 64 0.84933 0.61451 4.51471e-003 1141 1983 -13.26 -0.00 -0.00 2.51 38.13 -17.87 6.000 2.000 0 1 0 2 3.750 3.750 30.00 0.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
107+
7 1 1 1 0 2 6 16 62 64 64 -0.23580 2.07194 4.81163e-003 1119 1945 -13.26 -0.00 -0.00 2.51 45.91 -19.70 6.000 2.000 0 1 0 2 3.750 3.750 30.00 0.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
108+
8 1 1 1 0 2 7 16 62 64 64 -0.45610 1.75542 3.01393e-003 1097 1907 -13.26 -0.00 -0.00 2.51 53.70 -21.54 6.000 2.000 0 1 0 2 3.750 3.750 30.00 0.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
109+
9 1 1 1 0 2 8 16 62 64 64 -0.33958 0.41352 3.31281e-003 1146 1991 -13.26 -0.00 -0.00 2.51 61.49 -23.37 6.000 2.000 0 1 0 2 3.750 3.750 30.00 0.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
110+
1 1 2 1 0 2 9 16 62 64 64 -0.50054 0.02978 3.02487e-003 1071 1863 -13.26 -0.00 -0.00 2.51 -0.81 -8.69 6.000 2.000 0 1 0 2 3.750 3.750 30.00 2.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
111+
2 1 2 1 0 2 10 16 62 64 64 1.09595 2.68028 4.32900e-003 1123 1953 -13.26 -0.00 -0.00 2.51 6.98 -10.53 6.000 2.000 0 1 0 2 3.750 3.750 30.00 2.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
112+
3 1 2 1 0 2 11 16 62 64 64 0.66260 -0.01198 5.25907e-003 1135 1973 -13.26 -0.00 -0.00 2.51 14.77 -12.36 6.000 2.000 0 1 0 2 3.750 3.750 30.00 2.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
113+
4 1 2 1 0 2 12 16 62 64 64 0.49692 0.62467 3.65960e-003 1209 2101 -13.26 -0.00 -0.00 2.51 22.55 -14.20 6.000 2.000 0 1 0 2 3.750 3.750 30.00 2.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
114+
5 1 2 1 0 2 13 16 62 64 64 -0.81913 0.05637 2.96656e-003 1215 2113 -13.26 -0.00 -0.00 2.51 30.34 -16.03 6.000 2.000 0 1 0 2 3.750 3.750 30.00 2.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
115+
6 1 2 1 0 2 14 16 62 64 64 -0.46862 0.88504 5.71452e-003 1145 1990 -13.26 -0.00 -0.00 2.51 38.13 -17.87 6.000 2.000 0 1 0 2 3.750 3.750 30.00 2.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
116+
7 1 2 1 0 2 15 16 62 64 64 0.24446 3.78608 4.86407e-003 1119 1945 -13.26 -0.00 -0.00 2.51 45.91 -19.70 6.000 2.000 0 1 0 2 3.750 3.750 30.00 2.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
117+
8 1 2 1 0 2 16 16 62 64 64 0.07351 -1.06477 2.57365e-003 1093 1899 -13.26 -0.00 -0.00 2.51 53.70 -21.54 6.000 2.000 0 1 0 2 3.750 3.750 30.00 2.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
118+
9 1 2 1 0 2 17 16 62 64 64 -1.52431 2.45537 5.97195e-003 1150 1999 -13.26 -0.00 -0.00 2.51 61.49 -23.37 6.000 2.000 0 1 0 2 3.750 3.750 30.00 2.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
119+
1 1 3 1 0 2 18 16 62 64 64 0.74747 0.98210 4.98308e-003 1070 1860 -13.26 -0.00 -0.00 2.51 -0.81 -8.69 6.000 2.000 0 1 0 2 3.750 3.750 30.00 4.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
120+
2 1 3 1 0 2 19 16 62 64 64 -0.15021 0.52770 3.56176e-003 1125 1955 -13.26 -0.00 -0.00 2.51 6.98 -10.53 6.000 2.000 0 1 0 2 3.750 3.750 30.00 4.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
121+
3 1 3 1 0 2 20 16 62 64 64 1.68669 0.46589 4.71663e-003 1135 1973 -13.26 -0.00 -0.00 2.51 14.77 -12.36 6.000 2.000 0 1 0 2 3.750 3.750 30.00 4.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
122+
4 1 3 1 0 2 21 16 62 64 64 0.51405 0.88194 4.60270e-003 1211 2105 -13.26 -0.00 -0.00 2.51 22.55 -14.20 6.000 2.000 0 1 0 2 3.750 3.750 30.00 4.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
123+
5 1 3 1 0 2 22 16 62 64 64 -1.47991 2.03474 4.01615e-003 1218 2118 -13.26 -0.00 -0.00 2.51 30.34 -16.03 6.000 2.000 0 1 0 2 3.750 3.750 30.00 4.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
124+
6 1 3 1 0 2 23 16 62 64 64 -1.66364 3.21846 4.54139e-003 1143 1987 -13.26 -0.00 -0.00 2.51 38.13 -17.87 6.000 2.000 0 1 0 2 3.750 3.750 30.00 4.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
125+
7 1 3 1 0 2 24 16 62 64 64 0.60690 -0.41266 2.81800e-003 1120 1947 -13.26 -0.00 -0.00 2.51 45.91 -19.70 6.000 2.000 0 1 0 2 3.750 3.750 30.00 4.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
126+
8 1 3 1 0 2 25 16 62 64 64 1.01725 0.91765 3.19469e-003 1093 1901 -13.26 -0.00 -0.00 2.51 53.70 -21.54 6.000 2.000 0 1 0 2 3.750 3.750 30.00 4.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
127+
9 1 3 1 0 2 26 16 62 64 64 0.79204 1.65725 3.95777e-003 1151 2001 -13.26 -0.00 -0.00 2.51 61.49 -23.37 6.000 2.000 0 1 0 2 3.750 3.750 30.00 4.00 0.00 0.00 0 90.00 0 0 0 39 0.0 1 1 8 0 0.000 0.000 0.000 1
128+
129+
# === END OF DATA DESCRIPTION FILE ===============================================
216 KB
Binary file not shown.

nibabel/tests/test_fileslice.py

Lines changed: 19 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -14,7 +14,8 @@
1414
predict_shape, read_segments, _positive_slice,
1515
threshold_heuristic, optimize_slicer, slice2len,
1616
fill_slicer, optimize_read_slicers, slicers2segments,
17-
calc_slicedefs, _simple_fileslice, slice2outax)
17+
calc_slicedefs, _simple_fileslice, slice2outax,
18+
strided_scalar)
1819

1920
from nose.tools import assert_true, assert_false, assert_equal, assert_raises
2021

@@ -627,6 +628,23 @@ def test_predict_shape():
627628
assert_equal(predict_shape((1, slice(None), None), (2, 3)), (3, 1))
628629

629630

631+
def test_strided_scalar():
632+
# Utility to make numpy array of given shape from scalar using striding
633+
for shape, scalar in product(
634+
((2,), (2, 3,), (2, 3, 4)),
635+
(1, 2, np.int16(3))):
636+
expected = np.zeros(shape, dtype=np.array(scalar).dtype) + scalar
637+
observed = strided_scalar(shape, scalar)
638+
assert_array_equal(observed, expected)
639+
assert_equal(observed.shape, shape)
640+
assert_equal(observed.dtype, expected.dtype)
641+
assert_array_equal(observed.strides, 0)
642+
observed[..., 0] = 99
643+
assert_array_equal(observed, expected * 0 + 99)
644+
# Default scalar value is 0
645+
assert_array_equal(strided_scalar((2, 3, 4)), np.zeros((2, 3, 4)))
646+
647+
630648
def _check_bytes(bytes, arr):
631649
barr = np.ndarray(arr.shape, arr.dtype, buffer=bytes)
632650
assert_array_equal(barr, arr)

nibabel/tests/test_parrec.py

Lines changed: 29 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -15,6 +15,7 @@
1515
vol_is_full, PARRECImage, PARRECArrayProxy)
1616
from ..openers import Opener
1717
from ..fileholders import FileHolder
18+
from ..volumeutils import array_from_file
1819

1920
from numpy.testing import (assert_almost_equal,
2021
assert_array_equal)
@@ -40,6 +41,9 @@
4041
V4_PAR = pjoin(DATA_PATH, 'phantom_fake_v4.PAR')
4142
# Fake V4.1
4243
V41_PAR = pjoin(DATA_PATH, 'phantom_fake_v4_1.PAR')
44+
# Fake varying scaling
45+
VARY_PAR = pjoin(DATA_PATH, 'phantom_varscale.PAR')
46+
VARY_REC = pjoin(DATA_PATH, 'phantom_varscale.REC')
4347
# Affine as we determined it mid-2014
4448
AN_OLD_AFFINE = np.array(
4549
[[-3.64994708, 0., 1.83564171, 123.66276611],
@@ -554,3 +558,28 @@ def test_bitpix():
554558
for pix_size in (24, 32):
555559
hdr_defs['image pixel size'] = pix_size
556560
assert_raises(PARRECError, PARRECHeader, HDR_INFO, hdr_defs)
561+
562+
563+
def test_varying_scaling():
564+
# Check the algorithm works as expected for varying scaling
565+
img = PARRECImage.load(VARY_REC)
566+
rec_shape = (64, 64, 27)
567+
with open(VARY_REC, 'rb') as fobj:
568+
arr = array_from_file(rec_shape, '<i2', fobj)
569+
img_defs = img.header.image_defs
570+
slopes = img_defs['rescale slope']
571+
inters = img_defs['rescale intercept']
572+
sc_slopes = img_defs['scale slope']
573+
# Check dv scaling
574+
scaled_arr = arr.astype(np.float64)
575+
for i in range(arr.shape[2]):
576+
scaled_arr[:, :, i] *= slopes[i]
577+
scaled_arr[:, :, i] += inters[i]
578+
assert_almost_equal(np.reshape(scaled_arr, img.shape, order='F'),
579+
img.get_data(), 9)
580+
# Check fp scaling
581+
for i in range(arr.shape[2]):
582+
scaled_arr[:, :, i] /= (slopes[i] * sc_slopes[i])
583+
dv_img = PARRECImage.load(VARY_REC, scaling='fp')
584+
assert_almost_equal(np.reshape(scaled_arr, img.shape, order='F'),
585+
dv_img.get_data(), 9)

0 commit comments

Comments
 (0)