-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathspm_read_netcdf.m
132 lines (118 loc) · 4.25 KB
/
spm_read_netcdf.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
function cdf = spm_read_netcdf(fname)
% Read the header information from a NetCDF file into a data structure.
% FORMAT cdf = spm_read_netcdf(fname)
% fname - name of NetCDF file
% cdf - data structure
%
% See: http://www.unidata.ucar.edu/packages/netcdf/
% _________________________________________________________________________
% Copyright (C) 1999-2011 Wellcome Trust Centre for Neuroimaging
% John Ashburner
% $Id: spm_read_netcdf.m 4182 2011-02-01 12:29:09Z guillaume $
dsiz = [1 1 2 4 4 8];
fp=fopen(fname,'r','ieee-be');
if fp==-1
cdf = [];
return;
end
% Return null if not a CDF file.
%--------------------------------------------------------------------------
mgc = fread(fp,4,'uchar')';
if ~all(['CDF' 1] == mgc)
cdf = [];
fclose(fp);
if all(mgc==[137,72,68,70])
fprintf(['"%s" appears to be based around HDF.\n',...
'This is a newer version of MINC that SPM can not yet read.\n'],...
fname);
end
return;
end
% I've no idea what this is for
numrecs = fread(fp,1,'uint32');
cdf = struct('numrecs', numrecs,...
'dim_array', [], ...
'gatt_array', [], ...
'var_array', []);
dt = fread(fp,1,'uint32');
if dt == 10
% Dimensions
nelem = fread(fp,1,'uint32');
for j=1:nelem
str = readname(fp);
dim_length = fread(fp,1,'uint32');
cdf.dim_array(j).name = str;
cdf.dim_array(j).dim_length = dim_length;
end
dt = fread(fp,1,'uint32');
end
while ~dt, dt = fread(fp,1,'uint32'); end
if dt == 12
% Attributes
nelem = fread(fp,1,'uint32');
for j=1:nelem
str = readname(fp);
nc_type = fread(fp,1,'uint32');
nnelem = fread(fp,1,'uint32');
val = fread(fp,nnelem,dtypestr(nc_type));
if nc_type == 2, val = deblank([val' ' ']); end
padding= fread(fp,...
ceil(nnelem*dsiz(nc_type)/4)*4-nnelem*dsiz(nc_type),'uchar');
cdf.gatt_array(j).name = str;
cdf.gatt_array(j).nc_type = nc_type;
cdf.gatt_array(j).val = val;
end
dt = fread(fp,1,'uint32');
end
while ~dt, dt = fread(fp,1,'uint32'); end
if dt == 11
% Variables
nelem = fread(fp,1,'uint32');
for j=1:nelem
str = readname(fp);
nnelem = fread(fp,1,'uint32');
val = fread(fp,nnelem,'uint32');
cdf.var_array(j).name = str;
cdf.var_array(j).dimid = val+1;
cdf.var_array(j).nc_type = 0;
cdf.var_array(j).vsize = 0;
cdf.var_array(j).begin = 0;
dt0 = fread(fp,1,'uint32');
if dt0 == 12
nelem0 = fread(fp,1,'uint32');
for jj=1:nelem0
str = readname(fp);
nc_type= fread(fp,1,'uint32');
nnelem = fread(fp,1,'uint32');
val = fread(fp,nnelem,dtypestr(nc_type));
if nc_type == 2, val = deblank([val' ' ']); end
padding= fread(fp,...
ceil(nnelem*dsiz(nc_type)/4)*4-nnelem*dsiz(nc_type),'uchar');
cdf.var_array(j).vatt_array(jj).name = str;
cdf.var_array(j).vatt_array(jj).nc_type = nc_type;
cdf.var_array(j).vatt_array(jj).val = val;
end
dt0 = fread(fp,1,'uint32');
end
cdf.var_array(j).nc_type = dt0;
cdf.var_array(j).vsize = fread(fp,1,'uint32');
cdf.var_array(j).begin = fread(fp,1,'uint32');
end
dt = fread(fp,1,'uint32');
end
fclose(fp);
%==========================================================================
% function str = dtypestr(i)
%==========================================================================
function str = dtypestr(i)
% Returns a string appropriate for reading or writing the CDF data-type.
types = char('uint8','uint8','int16','int32','float32','float64');
str = deblank(types(i,:));
%==========================================================================
% function name = readname(fp)
%==========================================================================
function name = readname(fp)
% Extracts a name from a CDF file pointed to at the right location by fp.
stlen = fread(fp,1,'uint32');
name = deblank([fread(fp,stlen,'uchar')' ' ']);
padding = fread(fp,ceil(stlen/4)*4-stlen,'uchar');