1 |
function [madt]=aviso_get_madt(loc,grid_aviso) |
2 |
% function [madt]=aviso_get_madt(loc,grid) |
3 |
% extract madt value (mapped absolute dynamic topography) at 1 location |
4 |
% linear interpolation in space and time at each location. |
5 |
% |
6 |
% loc is a structure with three array of size (N,1) : |
7 |
% loc.lat : latitude [deg_north] |
8 |
% loc.lon : longitude [deg_east] |
9 |
% loc.time : date in Matlab julian date format |
10 |
% |
11 |
% grid is the structure returned by aviso_get_grid |
12 |
% |
13 |
% madt is the output array [unit: m], of size (N,1). |
14 |
% |
15 |
|
16 |
N=length(loc.lat); |
17 |
madt=NaN*zeros(N,1); |
18 |
|
19 |
loc.lon(loc.lon<0)=loc.lon(loc.lon<0)+360; |
20 |
|
21 |
Itime=interp1(grid_aviso.time,(1:grid_aviso.Ntime)'-1,loc.time); |
22 |
Ilon=interp1(grid_aviso.NbLongitudes,(1:grid_aviso.Nlon)'-1,loc.lon); |
23 |
Ilat=interp1(grid_aviso.NbLatitudes,(1:grid_aviso.Nlat)'-1,loc.lat); |
24 |
|
25 |
[X,Y,T]=meshgrid([0 1],[0 1],[0 1]); |
26 |
|
27 |
for kk=1:N, |
28 |
|
29 |
if isnan(Itime(kk)*Ilon(kk)*Ilat(kk)), continue, end |
30 |
it=floor(Itime(kk)); |
31 |
il=floor(Ilon(kk)); |
32 |
iL=floor(Ilat(kk)); |
33 |
URL=sprintf('http://opendap.aviso.oceanobs.com/thredds/dodsC/dataset-duacs-dt-ref-global-merged-madt-h?Grid_0001[%d:%d][%d:%d][%d:%d]',... |
34 |
it,it+1,il,il+1,iL,iL+1); |
35 |
a=loaddap(URL);data=a.Grid_0001; % data.Grid_0001 : lon x lat x time |
36 |
if any(data.Grid_0001(:)>1e10), continue, end |
37 |
madt(kk)=interp3(X,Y,T,data.Grid_0001/100,Ilon(kk)-il,Ilat(kk)-iL,Itime(kk)-it); |
38 |
|
39 |
end |