ecom_kslice

ECOM_KSLICE:  returns horizontal slice at particular layer.

 Horizontal slice is at specified sigma layer at given time-step for a
 ECOM file.

 This function can also be used to read in 2D and 3D fields such as
 bathymetry and heat_flux.  The coordinates of u are returned as x and y.

 USAGE: [u,x,y]=ecom_kslice(cdf,var,[time],[layer])

 where 
   cdf:  file name for netCDf file (e.g. 'ecom.cdf')
   var:  the variable to select (eg. 'salt' for salinity)
   time:  time step 
   layer:  sigma layer (e.g 1 for bottom layer)

    
       Examples: 

          [s,x,y]=ecom_kslice('ecom.cdf','salt',2,3);
              returns the salinity field from the 3rd sigma level
              at the 2nd time step.

          [elev,x,y]=ecom_kslice('ecom.cdf','elev',4);
              returns the elevation field from the 4th time step

          [depth,x,y]=ecom_kslice('ecom.cdf','depth');
              returns the depth field

0001 function [u,x,y]=ecom_kslice(cdf,var,time,layer)
0002 %ECOM_KSLICE:  returns horizontal slice at particular layer.
0003 %
0004 % Horizontal slice is at specified sigma layer at given time-step for a
0005 % ECOM file.
0006 %
0007 % This function can also be used to read in 2D and 3D fields such as
0008 % bathymetry and heat_flux.  The coordinates of u are returned as x and y.
0009 %
0010 % USAGE: [u,x,y]=ecom_kslice(cdf,var,[time],[layer])
0011 %
0012 % where
0013 %   cdf:  file name for netCDf file (e.g. 'ecom.cdf')
0014 %   var:  the variable to select (eg. 'salt' for salinity)
0015 %   time:  time step
0016 %   layer:  sigma layer (e.g 1 for bottom layer)
0017 %
0018 %
0019 %       Examples:
0020 %
0021 %          [s,x,y]=ecom_kslice('ecom.cdf','salt',2,3);
0022 %              returns the salinity field from the 3rd sigma level
0023 %              at the 2nd time step.
0024 %
0025 %          [elev,x,y]=ecom_kslice('ecom.cdf','elev',4);
0026 %              returns the elevation field from the 4th time step
0027 %
0028 %          [depth,x,y]=ecom_kslice('ecom.cdf','depth');
0029 %              returns the depth field
0030 %
0031 
0032 if (nargin<2 | nargin>4),
0033   help ecom_kslice; return
0034 end
0035 % turn off warnings from NetCDf
0036 mexcdf('setopts',0);
0037 %
0038 % open existing file
0039 ncid=mexcdf('open',cdf,'nowrite');
0040 if(ncid==-1),
0041   disp(['file ' cdf ' not found'])
0042   return
0043 end
0044 [name, nx]=mexcdf('diminq',ncid,'xpos');
0045 [name, ny]=mexcdf('diminq',ncid,'ypos');
0046 [name, nz]=mexcdf('diminq',ncid,'zpos');
0047 
0048 
0049 
0050 %
0051 % Acquire the grid.
0052 %
0053 % If "lon_rho" and "lat_rho" are present, grab them.
0054 % Otherwise, get "x_rho" and "y_rho".
0055 [lon_varid, rcode] = mexcdf('VARID', ncid, 'lon');
0056 [lat_varid, rcode] = mexcdf('VARID', ncid, 'lat');
0057 if ( (lon_varid >= 0) | (lat_varid >= 0) )
0058     x=ncmex('varget',ncid,'lon',[0 0],[-1 -1]);
0059     y=ncmex('varget',ncid,'lat',[0 0],[-1 -1]);
0060 else
0061     x=ncmex('varget',ncid,'x',[0 0],[-1 -1]);
0062     y=ncmex('varget',ncid,'y',[0 0],[-1 -1]);
0063 end
0064 
0065 %if(nargout==3),
0066 %  x=mexcdf('varget',ncid,'x',[0 0],[ny nx]);
0067 %  y=mexcdf('varget',ncid,'y',[0 0],[ny nx]);
0068 %end
0069 
0070 
0071 
0072 %
0073 % use depth=-99999  to mask points on land
0074 %
0075 depth=mexcdf('varget',ncid,'depth',[0 0],[ny nx]);
0076 land=find(depth==-99999);
0077 %
0078 % allow for using kslice on 2D, 3D and 4D variables
0079 %
0080 if(nargin==4),
0081  [u,ierr]=mexcdf('varget',ncid,var,[(time-1) layer-1 0 0],[1 1 ny nx],1); 
0082 elseif(nargin==3),
0083  [u,ierr]=mexcdf('varget',ncid,var,[(time-1) 0 0],[1 ny nx],1); 
0084 else
0085  [u,ierr]=mexcdf('varget',ncid,var,[0 0],[ny nx],1); 
0086 end
0087 
0088 mexcdf('close',ncid);
0089 u(land)=u(land)*NaN;
0090