Skip to content

Commit a95b552

Browse files
MathieuMorlighemMathieuMorlighem
committed
CHG: cleaned up function so that it does not rely on try/catch (and sped up the code here and there)
1 parent f720e70 commit a95b552

File tree

1 file changed

+46
-44
lines changed

1 file changed

+46
-44
lines changed

‎src/m/contrib/badgeley/interpFromMEaSUREsGeotiffs.m‎

Lines changed: 46 additions & 44 deletions
Original file line numberDiff line numberDiff line change
@@ -52,6 +52,7 @@
5252
else
5353
error(['The velocity data for ', product, ' is not available, either download it first or double check the ID number.']);
5454
end
55+
5556
% get the time info from file names
5657
if ismember(p,[725,727,731,766,478])
5758
templist = dir([foldername,'*_vv_v05.0.tif']);
@@ -75,7 +76,7 @@
7576
dataVers(ii) = join(dataVers_temp(1:2),'.');
7677
elseif p == 646
7778
dataPrefix(ii) = join(tempConv(1:3), '_');
78-
startdatenum = datenum(tempConv{3});
79+
startdatenum = datenum(tempConv{3},'yyyy-mm');
7980
startdatetime = datetime(tempConv{3},'InputFormat','yyyy-MM');
8081
enddatenum = datenum(year(startdatetime),month(startdatetime)+1,1) - 1;
8182
dataTstart(ii) = date2decyear(startdatenum);
@@ -94,53 +95,54 @@
9495
dataTend(ii) = decyear(datetime('31-10-2015','InputFormat','dd-MM-yyyy'));
9596
end
9697
end
97-
disp([' Found ', num2str(Ndata), ' records in ', foldername]);
98-
disp([' from ', datestr(decyear2date(min(dataTstart)),'yyyy-mm-dd'), ' to ', datestr(decyear2date(max(dataTend)),'yyyy-mm-dd') ]);
98+
fprintf([' == Found ', num2str(Ndata), ' records in ', foldername(numel(src_dir)+1:end)]);
99+
fprintf([' (from ', datestr(decyear2date(min(dataTstart)),'yyyy-mm-dd'), ' to ', datestr(decyear2date(max(dataTend)),'yyyy-mm-dd'), ')\n']);
99100

100101
%Tstart = Tstarts(p_ind);
101102
%Tend = Tends(p_ind);
102103

103104
dataInd = find((dataTend>=Tstarts(p_ind)) & (dataTstart<=Tends(p_ind)));
104-
disp([' For the selected period: ', datestr(decyear2date((Tstarts(p_ind))),'yyyy-mm-dd'), ' to ', datestr(decyear2date((Tends(p_ind))),'yyyy-mm-dd'), ', there are ', num2str(length(dataInd)), ' records' ]);
105+
disp([' For the selected period: ', datestr(decyear2date((Tstarts(p_ind))),'yyyy-mm-dd'), ' to ', datestr(decyear2date((Tends(p_ind))),'yyyy-mm-dd'), ', there are ', num2str(length(dataInd)), ' records' ]);
105106

106107
dataToLoad = dataPrefix(dataInd);
107108
TstartToload = dataTstart(dataInd);
108109
TendToload = dataTend(dataInd);
109110
VersToload = dataVers(dataInd);
110111

111112
jj = 1;
113+
skip = 0;
112114
for ii = 1:length(dataToLoad)
113-
try
114-
if contains('vx',output) || (p==646 && TendToload(ii) <= 2016) || (p==670)
115-
vx(jj,:) = interpFromGeotiff([foldername, dataToLoad{ii}, '_vx_', VersToload{ii}, '.tif'], X, Y, -99999);
116-
end
117-
if contains('vy',output) || (p==646 && TendToload(ii) <= 2016) || (p==670)
118-
vy(jj,:) = interpFromGeotiff([foldername, dataToLoad{ii}, '_vy_', VersToload{ii}, '.tif'], X, Y, -99999);
119-
end
120-
if ((p==646) && (TendToload(ii) <= 2016)) || (p==670)
121-
vx_temp = vx(jj,:);
122-
vy_temp = vy(jj,:);
123-
vx_temp(vx_temp==-99999)=nan;
124-
vy_temp(vy_temp==-99999)=nan;
125-
vel(jj,:) = (vx(jj,:).^2 + vy(jj,:).^2).^(1/2);
126-
else
127-
vel(jj,:) = interpFromGeotiff([foldername, dataToLoad{ii}, '_vv_', VersToload{ii}, '.tif'], X, Y, -99999);
128-
end
129-
if contains('ex',output)
130-
ex(jj,:) = interpFromGeotiff([foldername, dataToLoad{ii}, '_ex_', VersToload{ii}, '.tif'], X, Y, -99999);
115+
if contains('vx',output) || (p==646 && TendToload(ii) <= 2016) || (p==670)
116+
temp = interpFromGeotiff([foldername, dataToLoad{ii}, '_vx_', VersToload{ii}, '.tif'], X, Y, -99999);
117+
if all(isnan(temp), 'all')
118+
skip=skip+1;
119+
continue;
131120
end
132-
if contains('ey',output)
133-
ey(jj,:) = interpFromGeotiff([foldername, dataToLoad{ii}, '_ey_', VersToload{ii}, '.tif'], X, Y, -99999);
134-
end
135-
Tstart(jj) = TstartToload(ii);
136-
Tend(jj) = TendToload(ii);
137-
jj = jj + 1;
138-
catch
139-
% disp('error - continuing onto next observation');
140-
continue
141-
end
121+
vx(jj,:) = temp;
122+
end
123+
if contains('vy',output) || (p==646 && TendToload(ii) <= 2016) || (p==670)
124+
vy(jj,:) = interpFromGeotiff([foldername, dataToLoad{ii}, '_vy_', VersToload{ii}, '.tif'], X, Y, -99999);
125+
end
126+
if ((p==646) && (TendToload(ii) <= 2016)) || (p==670)
127+
vx_temp = vx(jj,:);
128+
vy_temp = vy(jj,:);
129+
vx_temp(vx_temp==-99999)=NaN;
130+
vy_temp(vy_temp==-99999)=NaN;
131+
vel(jj,:) = (vx(jj,:).^2 + vy(jj,:).^2).^(1/2);
132+
else
133+
vel(jj,:) = interpFromGeotiff([foldername, dataToLoad{ii}, '_vv_', VersToload{ii}, '.tif'], X, Y, -99999);
134+
end
135+
if contains('ex',output)
136+
ex(jj,:) = interpFromGeotiff([foldername, dataToLoad{ii}, '_ex_', VersToload{ii}, '.tif'], X, Y, -99999);
137+
end
138+
if contains('ey',output)
139+
ey(jj,:) = interpFromGeotiff([foldername, dataToLoad{ii}, '_ey_', VersToload{ii}, '.tif'], X, Y, -99999);
140+
end
141+
Tstart(jj) = TstartToload(ii);
142+
Tend(jj) = TendToload(ii);
143+
jj = jj + 1;
142144
end
143-
disp([' For the given domain and product ' num2str(p) ', there are ', num2str(jj-1), ' records' ]);
145+
disp([' For the given domain and product (' num2str(p) '), there are ', num2str(jj-1), ' records (skipped ' num2str(skip) ')']);
144146

145147
if jj > 1
146148
% find unique and sort by the mean time
@@ -149,12 +151,12 @@
149151
%Tend = Tend(ia);
150152

151153
% average
152-
disp([' Before averaging: ', num2str(length(ic)), ' sets']);
154+
disp([' Before averaging: ', num2str(length(ic)), ' sets']);
153155
kk = 0;
154156
for ii = 1:length(ic)
155157
ids = find(ic == ii);
156158
vels = vel(ids,:);
157-
vels(vels<0) = nan;
159+
vels(vels<0) = NaN;
158160
vel_mean = mean(double(vels), 1, 'omitnan');
159161
if all(isnan(vel_mean))
160162
continue
@@ -163,30 +165,30 @@
163165
dataout(d_ind+kk-1).vel = vel_mean;
164166
if contains('vx',output)
165167
vxs = vx(ids,:);
166-
vxs(vxs==-99999) = nan;
168+
vxs(vxs==-99999) = NaN;
167169
dataout(d_ind+kk-1).vx = mean(double(vxs), 1, 'omitnan');
168170
end
169171
if contains('vy',output)
170172
vys = vy(ids,:);
171-
vys(vys==-99999) = nan;
173+
vys(vys==-99999) = NaN;
172174
dataout(d_ind+kk-1).vy = mean(double(vys), 1, 'omitnan');
173175
end
174176
if contains('ex',output)
175177
exs = ex(ids,:);
176-
exs(exs==-99999) = nan;
178+
exs(exs==-99999) = NaN;
177179
dataout(d_ind+kk-1).ex = mean(double(exs), 1, 'omitnan');
178180
end
179181
if contains('ey',output)
180182
eys = ey(ids,:);
181-
eys(eys==-99999) = nan;
183+
eys(eys==-99999) = NaN;
182184
dataout(d_ind+kk-1).ey = mean(double(eys), 1, 'omitnan');
183185
end
184186
dataout(d_ind+kk-1).Tstart = mean(Tstart(ids));
185187
dataout(d_ind+kk-1).Tend = mean(Tend(ids));
186188
dataout(d_ind+kk-1).Tmean = mean(Tstart(ids) + Tend(ids))./2;
187189
dataout(d_ind+kk-1).product = p;
188190
end
189-
disp([' After averaging: ', num2str(kk), ' sets']);
191+
disp([' After averaging: ', num2str(kk), ' sets']);
190192
else
191193
kk = 0;
192194
end
@@ -214,8 +216,8 @@
214216
if contains('ey',output)
215217
vel_obs(ii).ey = reshape(dataout(pos(ii)).ey,[],1);
216218
end
217-
vel_obs(ii).Tstart = dataout(pos(ii)).Tstart;
218-
vel_obs(ii).Tend = dataout(pos(ii)).Tend;
219-
vel_obs(ii).Tmean = dataout(pos(ii)).Tmean;
220-
vel_obs(ii).product = dataout(pos(ii)).product;
219+
vel_obs(ii).Tstart = dataout(pos(ii)).Tstart;
220+
vel_obs(ii).Tend = dataout(pos(ii)).Tend;
221+
vel_obs(ii).Tmean = dataout(pos(ii)).Tmean;
222+
vel_obs(ii).product = dataout(pos(ii)).product;
221223
end

0 commit comments

Comments
 (0)