-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathtest-read_fs_mgh.R
177 lines (139 loc) · 7.67 KB
/
test-read_fs_mgh.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
test_that("The dimensions of our demo 3D MGZ file are read correctly", {
brain_image = system.file("extdata", "brain.mgz", package = "freesurferformats", mustWork = TRUE);
vd = read.fs.mgh(brain_image);
expect_equal(class(vd), "array");
expect_equal(length(dim(vd)), 4); # It has 4 dimensions
expect_equal(dim(vd), c(256, 256, 256, 1));
})
test_that("The data values in the demo 3D MGZ file are read as in the reference implementation", {
# Tests that the data read is identical to the values returned by the matlab function $FREESURFER_HOME/matlab/MRIread.m
# You can also check with the `mri_info` command line tool, just be sure to shift the
# vertex indices (zero-based for mri_info, one-based in R). See comment below for an example.
brain_image = system.file("extdata", "brain.mgz", package = "freesurferformats", mustWork = TRUE);
vd = read.fs.mgh(brain_image);
expect_equal(class(vd), "array");
expect_equal(length(dim(vd)), 4); # It has 4 dimensions
expect_equal(dim(vd), c(256, 256, 256, 1));
expect_equal(vd[100, 100, 100, 1], 77); # try on command line: mri_info --voxel 99 99 99 inst/extdata/brain.mgz
expect_equal(vd[110, 110, 110, 1], 71);
expect_equal(vd[1, 1, 1, 1], 0);
expect_equal(sum(vd), 121035479);
})
test_that("The demo 1D MGZ file is read correctly", {
morph_file = system.file("extdata", "lh.curv.fwhm10.fsaverage.mgz", package = "freesurferformats", mustWork = TRUE);
vd = read.fs.mgh(morph_file);
expect_equal(class(vd), "array");
expect_equal(length(dim(vd)), 4); # It has 4 dimensions
num_vers_fsaverage = 163842
expect_equal(dim(vd), c(num_vers_fsaverage, 1, 1, 1));
vd2 = read.fs.volume(morph_file);
expect_equal(class(vd2), "array");
expect_equal(length(dim(vd2)), 4); # It has 4 dimensions
expect_equal(dim(vd2), c(num_vers_fsaverage, 1, 1, 1));
})
test_that("The header can be read", {
brain_image = system.file("extdata", "brain.mgz", package = "freesurferformats", mustWork = TRUE);
ret = read.fs.mgh(brain_image, with_header=TRUE);
expect_true(is.fs.volume(ret));
ret2 = read.fs.volume(brain_image, with_header=TRUE);
expect_true(is.fs.volume(ret2));
header = ret$header;
expect_true(is.mghheader(header));
expect_equal(header$dtype, 0); # MRI_UCHAR
expect_equal(header$dof, 0);
expect_equal(header$ras_good_flag, 1);
expect_equal(length(header$internal$delta), 3);
expect_equal(header$internal$delta, c(1, 1, 1), tolerance=1e-2);
expect_equal(length(header$internal$Mdc), 9);
expect_equal(header$internal$Mdc, matrix(c(-1,0,0,0,0,-1,0,1,0), nrow=3), tolerance=1e-2);
expect_equal(length(header$internal$Pxyz_c), 3);
expect_equal(header$internal$Pxyz_c, c(-0.5, 29.4, -48.9), tolerance=1e-2);
expect_equal(header$voldim, c(256, 256, 256, 1));
expect_equal(length(header$has_mr_params), 1);
expect_equal(length(header$mr_params), 5);
expect_equal(header$mr_params, c(2300.000000, 0.157080, 2.010000, 900.000000, 256.000), tolerance=1e-2);
expect_equal(length(header$internal$D), 9); # 3x3
expect_equal(length(header$internal$Pcrs_c), 3); # 3x1
expect_equal(length(header$internal$Pxyz_0), 3); # 3x1
expect_equal(length(header$vox2ras_matrix), 16);
expect_equal(header$vox2ras_matrix, matrix(c(-1,0,0,0, 0,0,-1,0, 0,1,0,0, 127.5,-98.6273,79.0953,1.000), nrow=4, byrow = FALSE), tolerance=1e-2);
header2 = ret2$header;
expect_true(is.mghheader(header2));
expect_equal(header2$dtype, 0); # MRI_UCHAR
expect_equal(header2$dof, 0);
expect_equal(header2$ras_good_flag, 1);
expect_equal(length(header2$internal$delta), 3);
vd = ret$data;
expect_equal(class(vd), "array");
expect_equal(length(dim(vd)), 4); # It has 4 dimensions
expect_equal(typeof(vd), "integer");
expect_equal(dim(vd), c(256, 256, 256, 1));
expect_equal(vd[80,80,80,1], 110); # value of voxel 80,80,80 must be 110 (known from ref. implementation)
expect_equal(vd[100,100,100,1], 77);
expect_equal(vd[100,100,80,1], 105);
vd2 = ret2$data;
expect_equal(class(vd2), "array");
expect_equal(length(dim(vd2)), 4); # It has 4 dimensions
expect_equal(typeof(vd2), "integer");
expect_equal(dim(vd2), c(256, 256, 256, 1));
expect_equal(vd2[80,80,80,1], 110); # value of voxel 80,80,80 must be 110 (known from ref. implementation)
expect_equal(vd2[100,100,100,1], 77);
expect_equal(vd2[100,100,80,1], 105);
})
test_that("Data is flattened when requested", {
brain_image = system.file("extdata", "brain.mgz", package = "freesurferformats", mustWork = TRUE);
vd = read.fs.mgh(brain_image, flatten=TRUE);
expect_equal(class(vd), "integer");
expect_equal(length((vd)), 16777216);
})
test_that("The gzip status is guessed as expected from a filename", {
expect_equal(guess.filename.is.gzipped("noway"), FALSE);
expect_equal(guess.filename.is.gzipped("/there/is/noway"), FALSE);
expect_equal(guess.filename.is.gzipped("file.tar.gz"), TRUE);
expect_equal(guess.filename.is.gzipped("relative/path/to/file.tar.gz"), TRUE);
expect_equal(guess.filename.is.gzipped("file.tar"), FALSE);
expect_equal(guess.filename.is.gzipped("brain.mgz"), TRUE);
expect_equal(guess.filename.is.gzipped("brain.mgh"), FALSE);
})
test_that("A real MGH can be read, rewritten, read again, and the data and header are preserved as expected.", {
brain_image = system.file("extdata", "brain.mgz", package = "freesurferformats", mustWork = TRUE);
ret_orig = read.fs.mgh(brain_image, with_header=TRUE);
expect_true(is.fs.volume(ret_orig));
new_copy = tempfile(fileext="mgz");
write.fs.mgh(new_copy, ret_orig$data, vox2ras_matrix = ret_orig$header$vox2ras_matrix, mr_params=ret_orig$header$mr_params)
ret = read.fs.mgh(new_copy, with_header=TRUE);
header = ret$header;
expect_true(is.mghheader(header));
expect_equal(header$dtype, 1); # IMPORTANT: The data type will have changed from MRI_UCHAR to MRI_INTEGER. This is fine with us for now.
expect_equal(header$dof, 0);
expect_equal(header$ras_good_flag, 1);
expect_equal(length(header$internal$delta), 3);
expect_equal(header$internal$delta, c(1, 1, 1), tolerance=1e-2);
expect_equal(length(header$internal$Mdc), 9);
expect_equal(header$internal$Mdc, matrix(c(-1,0,0,0,0,-1,0,1,0), nrow=3), tolerance=1e-2);
expect_equal(length(header$internal$D), 9); # 3x3
expect_equal(length(header$internal$Pcrs_c), 3); # 3x1
expect_equal(length(header$internal$Pxyz_0), 3); # 3x1
expect_equal(header$internal$Pxyz_c, c(-0.5, 29.4, -48.9), tolerance=1e-2);
expect_equal(header$voldim, c(256, 256, 256, 1));
expect_equal(length(header$has_mr_params), 1);
expect_equal(length(header$mr_params), 5);
expect_equal(header$mr_params, c(2300.000000, 0.157080, 2.010000, 900.000000, 256.000), tolerance=1e-2);
expect_equal(header$vox2ras_matrix, matrix(c(-1,0,0,0, 0,0,-1,0, 0,1,0,0, 127.5,-98.6273,79.0953,1.000), nrow=4, byrow = FALSE), tolerance=1e-2);
vd = ret$data;
expect_equal(class(vd), "array");
expect_equal(length(dim(vd)), 4); # It has 4 dimensions
expect_equal(dim(vd), c(256, 256, 256, 1));
expect_equal(vd[80,80,80,1], 110); # value of voxel 80,80,80 must be 110 (known from ref. implementation)
expect_equal(vd[100,100,100,1], 77);
expect_equal(vd[100,100,80,1], 105);
})
test_that("Invalid volume format leads to errors", {
brain_image = system.file("extdata", "brain.mgz", package = "freesurferformats", mustWork = TRUE);
expect_error(read.fs.volume(brain_image, format = 'invalid format')); # invalid format
})
test_that("Trying to read a non-existent file that does not match any of the file extension checks with read.fs.volume throws an error (issue #25)", {
non_existent_file = tempfile(, fileext = "");
testthat::expect_false(file.exists(non_existent_file)); # ensure it does not exist.
testthat::expect_error(read.fs.volume(non_existent_file));
})