-
Notifications
You must be signed in to change notification settings - Fork 0
/
rp_spm_vol.m
173 lines (156 loc) · 5.44 KB
/
rp_spm_vol.m
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
function V = rp_spm_vol(P)
% Get header information etc for images.
% FORMAT V = spm_vol(P)
% P - a matrix of filenames.
% V - a vector of structures containing image volume information.
% The elements of the structures are:
% V.fname - the filename of the image.
% V.dim - the x, y and z dimensions of the volume
% V.dt - A 1x2 array. First element is datatype (see spm_type).
% The second is 1 or 0 depending on the endian-ness.
% V.mat - a 4x4 affine transformation matrix mapping from
% voxel coordinates to real world coordinates.
% V.pinfo - plane info for each plane of the volume.
% V.pinfo(1,:) - scale for each plane
% V.pinfo(2,:) - offset for each plane
% The true voxel intensities of the jth image are given
% by: val*V.pinfo(1,j) + V.pinfo(2,j)
% V.pinfo(3,:) - offset into image (in bytes).
% If the size of pinfo is 3x1, then the volume is assumed
% to be contiguous and each plane has the same scalefactor
% and offset.
%____________________________________________________________________________
%
% The fields listed above are essential for the mex routines, but other
% fields can also be incorporated into the structure.
%
% The images are not memory mapped at this step, but are mapped when
% the mex routines using the volume information are called.
%
% Note that spm_vol can also be applied to the filename(s) of 4-dim
% volumes. In that case, the elements of V will point to a series of 3-dim
% images.
%
% This is a replacement for the spm_map_vol and spm_unmap_vol stuff of
% MatLab4 SPMs (SPM94-97), which is now obsolete.
%_______________________________________________________________________
% Copyright (C) 2005 Wellcome Department of Imaging Neuroscience
% John Ashburner
% $Id: spm_vol.m 982 2007-10-26 14:01:54Z john $
if nargin==0,
V = struct('fname', {},...
'dim', {},...
'dt', {},...
'pinfo', {},...
'mat', {},...
'n', {},...
'descrip', {},...
'private',{});
return;
end;
% If is already a vol structure then just return;
if isstruct(P), V = P; return; end;
V = subfunc2(P);
return;
%_______________________________________________________________________
%_______________________________________________________________________
function V = subfunc2(P)
if iscell(P),
V = cell(size(P));
for j=1:numel(P),
if iscell(P{j}),
V{j} = subfunc2(P{j});
else
V{j} = subfunc1(P{j});
end;
end;
else
V = subfunc1(P);
end;
return;
%_______________________________________________________________________
%_______________________________________________________________________
function V = subfunc1(P)
if isempty(P),
V = [];
return;
end;
counter = 0;
for i=1:size(P,1),
v = subfunc(P(i,:));
[V(counter+1:counter+size(v, 2),1).fname] = deal('');
[V(counter+1:counter+size(v, 2),1).mat] = deal([0 0 0 0]);
[V(counter+1:counter+size(v, 2),1).mat] = deal(eye(4));
[V(counter+1:counter+size(v, 2),1).mat] = deal([1 0 0]');
if isempty(v),
hread_error_message(P(i,:));
error(['Can''t get volume information for ''' P(i,:) '''']);
end
f = fieldnames(v);
for j=1:size(f,1)
eval(['[V(counter+1:counter+size(v,2),1).' f{j} '] = deal(v.' f{j} ');']);
end
counter = counter + size(v,2);
end
return;
%_______________________________________________________________________
%_______________________________________________________________________
function V = subfunc(p)
[pth,nam,ext,n1] = rp_spm_fileparts(deblank(p));
p = fullfile(pth,[nam ext]);
n = str2num(n1);
if ~exist(p,'file'),
existance_error_message(p);
error('File "%s" does not exist.', p);
end
switch ext,
case {'.nii','.NII'},
% Do nothing
case {'.img','.IMG'},
if ~exist(fullfile(pth,[nam '.hdr']),'file') && ~exist(fullfile(pth,[nam '.HDR']),'file'),
existance_error_message(fullfile(pth,[nam '.hdr'])),
error('File "%s" does not exist.', fullfile(pth,[nam '.hdr']));
end
case {'.hdr','.HDR'},
ext = '.img';
p = fullfile(pth,[nam ext]);
if ~exist(p,'file'),
existance_error_message(p),
error('File "%s" does not exist.', p);
end
otherwise,
error('File "%s" is not of a recognised type.', p);
end
if isempty(n),
V = rp_spm_vol_nifti(p);
else
V = rp_spm_vol_nifti(p,n);
end;
if isempty(n) && length(V.private.dat.dim) > 3
V0(1) = V;
for i = 2:V.private.dat.dim(4)
V0(i) = rp_spm_vol_nifti(p, i);
end
V = V0;
end
if ~isempty(V), return; end;
return;
%_______________________________________________________________________
%_______________________________________________________________________
function hread_error_message(q)
str = {...
'Error reading information on:',...
[' ',rp_spm_str_manip(q,'k40d')],...
' ',...
'Please check that it is in the correct format.'};
rp_spm('alert*',str,mfilename,sqrt(-1));
return;
%_______________________________________________________________________
function existance_error_message(q)
str = {...
'Unable to find file:',...
[' ',rp_spm_str_manip(q,'k40d')],...
' ',...
'Please check that it exists.'};
rp_spm('alert*',str,mfilename,sqrt(-1));
return;