Skip to content
Snippets Groups Projects
Commit a9c3e176 authored by Paul McCarthy's avatar Paul McCarthy :mountain_bicyclist:
Browse files

ENH: Rewrote fslcreatehd test in python, expanded tests to handle 4th dim of (0,

1), and to handle modification of existing image
parent 5e8220fa
No related branches found
No related tags found
No related merge requests found
#!/bin/bash #!/usr/bin/env python
set -e import os
import sys
import subprocess as sp
import numpy as np
import nibabel as nib
function neq() {
n1=$1
n2=$2
result=`echo "$n1 == $n2" | bc`
if [[ "$result" == "0" ]]; then # 2=char, 4=short, 8=int, 16=float, 64=double
echo "$n1 != $n2" DTYPE_MAPPING = {
return 0; 2 : np.uint8,
fi 4 : np.int16,
8 : np.int32,
return 1; 16 : np.float32,
} 64 : np.float64
function fslcreatehd_test() {
dim0=$1; shift
dim1=$1; shift
dim2=$1; shift
dim3=$1; shift
dim4=$1; shift
pixdim1=$1; shift
pixdim2=$1; shift
pixdim3=$1; shift
pixdim4=$1; shift
originx=$1; shift
originy=$1; shift
originz=$1; shift
dtype=$1; shift
fname=$1; shift
args="$dim0 $dim1 $dim2 $dim3 $pixdim1 $pixdim2 $pixdim3 $pixdim4 $originx $originy $originz $dtype $fname"
originx=`echo " $originx * $pixdim1" | bc`
originy=`echo "-$originy * $pixdim2" | bc`
originz=`echo "-$originz * $pixdim3" | bc`
got_dim0=` fslhd $fname | egrep "^dim0" | tr -s '\t' | cut -f 2`
got_dim1=` fslhd $fname | egrep "^dim1" | tr -s '\t' | cut -f 2`
got_dim2=` fslhd $fname | egrep "^dim2" | tr -s '\t' | cut -f 2`
got_dim3=` fslhd $fname | egrep "^dim3" | tr -s '\t' | cut -f 2`
got_dim4=` fslhd $fname | egrep "^dim4" | tr -s '\t' | cut -f 2`
got_pixdim1=`fslhd $fname | egrep "^pixdim1" | tr -s '\t' | cut -f 2`
got_pixdim2=`fslhd $fname | egrep "^pixdim2" | tr -s '\t' | cut -f 2`
got_pixdim3=`fslhd $fname | egrep "^pixdim3" | tr -s '\t' | cut -f 2`
got_pixdim4=`fslhd $fname | egrep "^pixdim4" | tr -s '\t' | cut -f 2`
got_dtype=` fslhd $fname | egrep "^datatype" | tr -s '\t' | cut -f 2`
got_originx=`fslhd $fname | egrep "^sto_xyz:1" | cut -d ' ' -f 4`
got_originy=`fslhd $fname | egrep "^sto_xyz:2" | cut -d ' ' -f 4`
got_originz=`fslhd $fname | egrep "^sto_xyz:3" | cut -d ' ' -f 4`
if neq $dim0 $got_dim0; then echo $args; return 1; fi
if neq $dim1 $got_dim1; then echo $args; return 1; fi
if neq $dim2 $got_dim2; then echo $args; return 1; fi
if neq $dim3 $got_dim3; then echo $args; return 1; fi
if neq $dim4 $got_dim4; then echo $args; return 1; fi
if neq $pixdim1 $got_pixdim1; then echo $args; return 1; fi
if neq $pixdim2 $got_pixdim2; then echo $args; return 1; fi
if neq $pixdim3 $got_pixdim3; then echo $args; return 1; fi
if neq $pixdim4 $got_pixdim4; then echo $args; return 1; fi
if neq $dtype $got_dtype; then echo $args; return 1; fi
if neq $originx $got_originx; then echo $args; return 1; fi
if neq $originy $got_originy; then echo $args; return 1; fi
if neq $originz $got_originz; then echo $args; return 1; fi
return 0
} }
tests=" def validate_result(imgfile, expshape, exppixdim, exporigin, expdtype):
5 5 5 1 2 2 2 1 0 0 0 2
5 5 5 1 2 2 2 1 0 0 0 4 img = nib.load(imgfile)
5 5 5 1 2 2 2 1 0 0 0 8 hdr = img.header
5 5 5 1 2 2 2 1 0 0 0 16 expshape = list(expshape)
5 5 5 1 2 2 2 1 1 2 3 64 exppixdim = list(exppixdim)
5 5 5 1 0.5 1.5 1.25 1 0 0 0 2 exporigin = list(exporigin)
5 5 5 1 0.5 1.5 1.25 1.5 0 0 0 2
5 5 5 1 0.5 1.5 1.25 0.5 0 0 0 2 if (len(expshape) == 3) or \
expshape[3] in (0, 1):
30 30 30 5 5 10 3 5 10 20 10 2 expndims = 3
" else:
expndims = 4
echo "$tests" | while IFS= read -r test; do
if [[ $test =~ ^\ *$ ]]; then expaffine = np.diag(exppixdim[:3] + [1])
continue expaffine[0, 0] *= -1
fi expaffine[:3, 3] = exporigin
expaffine[0, 3] *= exppixdim[0]
rm `imglob -extension fslcreatehd_test_img` &> /dev/null || true expaffine[1, 3] *= -exppixdim[1]
fslcreatehd $test fslcreatehd_test_img expaffine[2, 3] *= -exppixdim[2]
fslcreatehd_test 4 $test fslcreatehd_test_img
assert len(img.shape) == expndims
done assert list(img.shape) == expshape[ :expndims]
assert list(hdr.get_zooms()) == exppixdim[:expndims]
assert img.get_data_dtype() == expdtype
assert np.all(np.isclose(img.affine, expaffine))
def create_image(shape, pixdim, origin, dtype):
pixdim = list(pixdim)
affine = np.diag(pixdim[:3] + [1])
affine[:3, 3] = origin
data = np.random.randint(1, 100, shape).astype(dtype)
hdr = nib.Nifti1Header()
hdr.set_data_dtype(dtype)
hdr.set_data_shape(shape)
hdr.set_zooms(pixdim[:len(shape)])
return nib.Nifti1Image(data, affine, hdr)
def test_new_file():
tests = [
# 4th dim of size 0 or 1 should
# result in a 3D image (this
# is coded in validate_result)
' 5 5 5 0 2 2 2 0 0 0 0 2',
' 5 5 5 0 2 2 2 1 0 0 0 2',
' 5 5 5 1 2 2 2 1 0 0 0 2',
' 5 5 5 1 2 2 2 1 0 0 0 4',
' 5 5 5 1 2 2 2 1 0 0 0 8',
' 5 5 5 1 2 2 2 1 0 0 0 16',
' 5 5 5 1 2 2 2 1 1 2 3 64',
' 5 5 5 1 0.5 1.5 1.25 1 0 0 0 2',
' 5 5 5 1 0.5 1.5 1.25 1.5 0 0 0 2',
' 5 5 5 1 0.5 1.5 1.25 0.5 0 0 0 2',
'30 30 30 5 5 10 3 5 10 20 10 2',
]
for test in tests:
args = list(test.split())
imgfile = 'image.nii.gz'
sp.run(['fslcreatehd'] + args + [imgfile])
args = [float(a) for a in args]
expshape = args[:4]
exppixdim = args[4:8]
exporigin = args[8:11]
expdtype = DTYPE_MAPPING[args[11]]
try:
validate_result(
imgfile, expshape, exppixdim, exporigin, expdtype)
finally:
os.remove(imgfile)
def test_existing_file():
# each test is a tuple containing:
# - dtype_of_existing_image
# - shape_of_existing_image
# - exp_shape (set to None if image should not be overwritten)
# - fslcreatehd_args
tests = [
# 4th dim of 0 or 1 should result in a 3D image
(2, (5, 5, 5), None, '5 5 5 0 2 2 2 1 0 0 0 2'),
(2, (5, 5, 5), None, '5 5 5 0 2 2 2 1 1 2 3 2'),
(2, (5, 5, 5), None, '5 5 5 1 2 2 2 2 0 0 0 2'),
# dtype arg should be ignored for existing images
(2, (5, 5, 5), None, '5 5 5 0 2 2 2 1 0 0 0 4'),
(4, (5, 5, 5), None, '5 5 5 0 2 2 2 1 0 0 0 2'),
# data should be overwritten if any
# of the first 3 dims are different,
# or if a 4th dim is specified
(2, (5, 5, 5), (5, 3, 5), '5 3 5 0 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5), (5, 3, 5), '5 3 5 0 2 2 2 2 0 0 0 4'),
(2, (5, 5, 5), (5, 3, 5), '5 3 5 1 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5), (5, 3, 5), '5 3 5 1 2 2 2 2 1 2 3 2'),
(2, (5, 5, 5), (5, 3, 5, 5), '5 3 5 5 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5), (5, 3, 5, 5), '5 3 5 5 2 2 2 2 1 2 3 2'),
# 4D - same rules apply
(2, (5, 5, 5, 5), None, '5 5 5 5 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5, 5), None, '5 5 5 5 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5, 5), None, '5 5 5 5 2 2 2 2 0 0 0 4'),
(4, (5, 5, 5, 5), None, '5 5 5 5 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5, 5), None, '5 5 5 5 2 2 2 2 1 2 3 2'),
# data overwritten if nelements
# change
(2, (5, 5, 5, 5), (5, 5, 5), '5 5 5 0 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5, 5), (5, 5, 5), '5 5 5 1 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5, 5), (5, 3, 5), '5 3 5 0 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5, 5), (5, 3, 5), '5 3 5 1 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5, 5), (5, 3, 5), '5 3 5 1 2 2 2 2 1 2 3 2'),
(2, (5, 5, 5, 5), (5, 5, 5, 2), '5 5 5 2 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5, 5), (5, 3, 5, 5), '5 3 5 5 2 2 2 2 0 0 0 2'),
(2, (5, 5, 5, 5), (5, 3, 5, 5), '5 3 5 5 2 2 2 2 1 2 3 2'),
]
for (dtype, shape, expshape, args) in tests:
imgfile = 'image.nii.gz'
dtype = DTYPE_MAPPING[dtype]
args = list(args.split())
img = create_image(shape, (1, 1, 1, 1), (0, 0, 0), dtype)
img.to_filename(imgfile)
sp.run(['fslcreatehd'] + args + [imgfile])
checkdata = expshape is None
if expshape is None:
expshape = shape
args = [float(a) for a in args]
exppixdim = args[4:8]
exporigin = args[8:11]
try:
validate_result(imgfile, expshape, exppixdim, exporigin, dtype)
if checkdata:
data = np.asanyarray(img.dataobj)
datacopy = np.asanyarray(nib.load(imgfile).dataobj)
assert np.all(data == datacopy)
finally:
os.remove(imgfile)
def main():
test_new_file()
test_existing_file()
if __name__ == '__main__':
sys.exit(main())
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment