1 |
function [data,xax,yax,zax,months,time,dim] = ... |
2 |
DiagLoadGradsData(Grads,dad,fln,DiagDebug) |
3 |
|
4 |
|
5 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
6 |
% Parse table file % |
7 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
8 |
ShiftData = 0; |
9 |
format = 'NONSEQUENTIAL'; |
10 |
|
11 |
tablfile = [dad,'/',Grads,'.tabl']; |
12 |
file = textread(tablfile,'%s','delimiter','\n','whitespace',''); |
13 |
|
14 |
for iline = 1:length(file) |
15 |
rem = file{iline}; tokens = {}; itoken = 0; |
16 |
while isstr(rem) |
17 |
[token,rem] = strtok(rem); itoken = itoken + 1; |
18 |
if ~isempty(token), tokens{itoken} = token; end |
19 |
end |
20 |
|
21 |
if isequal(tokens{1},'XDEF') |
22 |
if ~isequal(tokens{3},'LINEAR') |
23 |
error('Unable to deal with nonlinear grads x-axis data!'); |
24 |
end |
25 |
num = str2num(tokens{2}); |
26 |
ini = str2num(tokens{4}); |
27 |
inc = str2num(tokens{5}); |
28 |
xax = [ini:inc:ini+(num-1)*inc]; |
29 |
if min(xax) >= 0 && max(xax) > 180 |
30 |
ShiftData = 1; |
31 |
end |
32 |
nx = length(xax); |
33 |
end |
34 |
|
35 |
if isequal(tokens{1},'YDEF') |
36 |
if ~isequal(tokens{3},'LINEAR') |
37 |
error('Unable to deal with nonlinear grads y-axis data!'); |
38 |
end |
39 |
num = str2num(tokens{2}); |
40 |
ini = str2num(tokens{4}); |
41 |
inc = str2num(tokens{5}); |
42 |
yax = [ini:inc:ini+(num-1)*inc]; ny = length(yax); |
43 |
end |
44 |
|
45 |
if isequal(tokens{1},'ZDEF') |
46 |
if isequal(tokens{2},'1') |
47 |
zax = [0]; zax_found = 1; dim = 2; |
48 |
else |
49 |
if ~isequal(tokens{3},'LINEAR') |
50 |
error('Unable to deal with nonlinear grads z-axis data!'); |
51 |
end |
52 |
num = str2num(tokens{2}); |
53 |
ini = str2num(tokens{4}); |
54 |
inc = str2num(tokens{5}); |
55 |
zax = [ini:inc:ini+(num-1)*inc]; dim = 3; |
56 |
end |
57 |
nz = length(zax); |
58 |
end |
59 |
|
60 |
if isequal(tokens{1},'TDEF'), ini=tokens{4}; |
61 |
if ~isequal(tokens{3},'LINEAR') |
62 |
error('Currently unable to deal with nonlinear grads z-axis data!'); |
63 |
end |
64 |
if ~isequal(tokens{5},'1mo') |
65 |
error('Currently unable to deal with non monthly mean grads data!'); |
66 |
end |
67 |
if length(ini) == 13 |
68 |
monchar=ini(9:11); |
69 |
end |
70 |
if isequal(monchar,'JAN'), inimonth = 1; |
71 |
elseif isequal(monchar,'FEB'), inimonth = 2; |
72 |
elseif isequal(monchar,'MAR'), inimonth = 3; |
73 |
elseif isequal(monchar,'APR'), inimonth = 4; |
74 |
elseif isequal(monchar,'MAY'), inimonth = 5; |
75 |
elseif isequal(monchar,'JUN'), inimonth = 6; |
76 |
elseif isequal(monchar,'JUL'), inimonth = 7; |
77 |
elseif isequal(monchar,'AUG'), inimonth = 8; |
78 |
elseif isequal(monchar,'SEP'), inimonth = 9; |
79 |
elseif isequal(monchar,'OCT'), inimonth = 10; |
80 |
elseif isequal(monchar,'NOV'), inimonth = 11; |
81 |
elseif isequal(monchar,'DEC'), inimonth = 12; end |
82 |
num = str2num(tokens{2}); |
83 |
months=[inimonth:num+inimonth-1]; |
84 |
time=months/12; nt = length(months); |
85 |
|
86 |
end |
87 |
|
88 |
if isequal(tokens{1},'VARS') |
89 |
nv = str2num(tokens{2}); |
90 |
VARSline = iline; |
91 |
end |
92 |
if isequal(tokens{1},fln) |
93 |
FIELDline = iline; |
94 |
ivar = iline - VARSline; |
95 |
end |
96 |
|
97 |
if isequal(tokens{1},'UNDEF') |
98 |
undef = str2num(tokens{2}); |
99 |
end |
100 |
|
101 |
if isequal(tokens{1},'FORMAT') |
102 |
if isequal(tokens{2},'SEQUENTIAL') |
103 |
format = 'SEQUENTIAL'; |
104 |
else |
105 |
disp(['Unrecognized grads FORMAT: ',tokens{2}]); |
106 |
end |
107 |
end |
108 |
|
109 |
if isequal(tokens{1},'DSET') |
110 |
datafile = [dad,tokens{2}(2:end)]; |
111 |
end |
112 |
|
113 |
end |
114 |
|
115 |
|
116 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
117 |
% Verification and debugging % |
118 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
119 |
|
120 |
if DiagDebug, disp([' Debug -- grads data format: ',format]); end |
121 |
|
122 |
% Verify that things are there. |
123 |
parameters = {'xax' ,'yax' ,'zax' ,'time','nv' ,'undef','ivar'}; |
124 |
GRADSnames = {'XDEF','YDEF','ZDEF','TDEF','VARS','UNDEF',fln }; |
125 |
for ii = 1:length(parameters) |
126 |
try |
127 |
eval([parameters{ii},';']); |
128 |
catch |
129 |
error(['GRADS information not found: ',GRADSnames{ii}]); |
130 |
end |
131 |
end |
132 |
|
133 |
|
134 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
135 |
% Read data file % |
136 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
137 |
|
138 |
tol = 1e-5; |
139 |
fid=fopen(datafile,'r','b'); |
140 |
data = fread(fid,'real*4'); |
141 |
fclose(fid); |
142 |
if isequal(format,'SEQUENTIAL') |
143 |
index=true([nx*ny*nz*nv*nt+2*nv*nt,1]); |
144 |
index([1:nx*ny*nz+2:end])=false; |
145 |
index([2:nx*ny*nz+2:end])=false; |
146 |
data = data(index); |
147 |
end |
148 |
data = reshape(data,[nx,ny,nz,nv,nt]); |
149 |
data = squeeze(data(:,:,:,ivar,:)); |
150 |
data( abs((data-undef)/undef) < tol ) = NaN; |
151 |
if ShiftData |
152 |
indexWestHemi = xax>=180; |
153 |
indexEastHemi = xax<180; |
154 |
data = cat(1,data(indexWestHemi,:,:),data(indexEastHemi,:,:)); |
155 |
xax = cat(2,xax(indexWestHemi)-360,xax(indexEastHemi)); |
156 |
end |