/[MITgcm]/MITgcm_contrib/gael/profilesMatlabProcessing/profiles_misc/aviso_get_madt.m
ViewVC logotype

Annotation of /MITgcm_contrib/gael/profilesMatlabProcessing/profiles_misc/aviso_get_madt.m

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.1 - (hide annotations) (download)
Thu Apr 14 23:37:11 2011 UTC (14 years, 3 months ago) by roquet
Branch: MAIN
New capability: extract madt from aviso opendap server at the position of profiles
and add it as a field prof_madt_aviso in the MITprof netcdf file.

1 roquet 1.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     Itime=interp1(grid_aviso.time,(1:grid_aviso.Ntime)'-1,loc.time);
20     Ilon=interp1(grid_aviso.NbLongitudes,(1:grid_aviso.Nlon)'-1,loc.lon);
21     Ilat=interp1(grid_aviso.NbLatitudes,(1:grid_aviso.Nlat)'-1,loc.lat);
22    
23     [X,Y,T]=meshgrid([0 1],[0 1],[0 1]);
24    
25     for kk=1:N,
26    
27     if isnan(Itime(kk)*Ilon(kk)*Ilat(kk)), continue, end
28     it=floor(Itime(kk));
29     il=floor(Ilon(kk));
30     iL=floor(Ilat(kk));
31     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]',...
32     it,it+1,il,il+1,iL,iL+1);
33     a=loaddap(URL);data=a.Grid_0001; % data.Grid_0001 : lon x lat x time
34     if any(data.Grid_0001>1e10), continue, end
35     madt(kk)=interp3(X,Y,T,data.Grid_0001/100,Ilon(kk)-il,Ilat(kk)-iL,Itime(kk)-it);
36    
37     end

  ViewVC Help
Powered by ViewVC 1.1.22