ecomvel

 ECOMVEL extracts a velocity time series from a 4D netCDF ECOM output file.

   USAGE:  [u,jd,depth]=ecomvel(cdf,i,j,lev,[tind]);

 cdf = cdf file e.g. 'ecomsi.cdf'
 (i,j) = station location in grid coordinates
 lev = sigma level from which to extract time series (0 is surface)
 w = complex vector of velocity with dimension time
 jd = time vector  (digital Julian day)
 depth = depth of time series (meters)

0001 function [w,jd,depth]=ecomvel(cdf,i,j,lev,iind)
0002 % ECOMVEL extracts a velocity time series from a 4D netCDF ECOM output file.
0003 %
0004 %   USAGE:  [u,jd,depth]=ecomvel(cdf,i,j,lev,[tind]);
0005 %
0006 % cdf = cdf file e.g. 'ecomsi.cdf'
0007 % (i,j) = station location in grid coordinates
0008 % lev = sigma level from which to extract time series (0 is surface)
0009 % w = complex vector of velocity with dimension time
0010 % jd = time vector  (digital Julian day)
0011 % depth = depth of time series (meters)
0012 
0013 % Rich Signell (rsignell@usgs.gov)
0014 
0015 % 3-27-96 fixed bug that was returning a slightly bad depth value
0016 %         due to averaging the sigma levels twice
0017 
0018 ncid=mexcdf('open',cdf,'nowrite');
0019 [nam,nz]=mexcdf('diminq',ncid,'zpos');
0020 [nam,nt]=mexcdf('diminq',ncid,'time');
0021 if(exist('iind')==1),
0022   istart=min(iind);
0023   icount=min(max(iind),nt)-istart+1;
0024 else
0025   istart=0;
0026   icount=nt;
0027 end
0028 sigma=mexcdf('varget',ncid,'sigma',0,nz);
0029 depth=mexcdf('varget',ncid,'depth',[j-1 i-1],[1 1]);
0030 depth=.5*(sigma(lev)+sigma(lev+1))*depth;
0031 
0032 base_date=[0 0 0 0 0 0];
0033 base_date(1:3)=mexcdf('attget',ncid,'global','base_date');
0034 t=mexcdf('varget',ncid,'time',istart,icount);
0035 jd0=julian(base_date);
0036 jd=jd0+t;
0037 
0038 ang=mexcdf('varget',ncid,'ang',[j-1 i-1],[1 1],1);
0039 u=mexcdf('varget',ncid,'u',[istart lev-1 j-1 i-1],[icount 1 1 2],1);
0040 v=mexcdf('varget',ncid,'v',[istart lev-1 j-1 i-1],[icount 1 2 1],1);
0041 
0042 w=mean(u)+sqrt(-1)*mean(v);       %average u and v to center of grid
0043 w=w(:);
0044 w=w.*exp(sqrt(-1)*ang);   %rotate vectors from grid coordinates to E,N