Home > omviz > scrum_kslice.m

scrum_kslice

PURPOSE ^

SCRUM_KSLICE: returns horizontal slice at particular layer.

SYNOPSIS ^

function [u,x,y]=scrum_kslice(cdf,var,time,layer)

DESCRIPTION ^

SCRUM_KSLICE:  returns horizontal slice at particular layer.

 Horizontal slice is at specified s_rho layer at given time-step for a
 SCRUM 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]=scrum_kslice(cdf,var,[time],[layer])

 where 
   cdf:  file name for netCDf file (e.g. 'scrum_c.nc')
   var:  the variable to select (eg. 'salt' for salinity)
   time:  time step 
   layer:  in the NetCDF file, s_rho is dimensioned as [0 ... s_rho_length-1].
           The layer variables which are a function of s_rho run from -1 to 0,
           which is reversed from ECOM.  In order to have conformity, we will 
           have the parameter 'layer' be 1-based, and have layer 1 correspond 
           to the top layer instead of the bottom.

    
       Examples: 

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

          [elev,x,y]=scrum_kslice('scrum_c.nc','zeta',4);
              returns the elevation field from the 4th time step

          [depth,x,y]=scrum_kslice('scrum_c.nc','h');
              returns the depth field

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [u,x,y]=scrum_kslice(cdf,var,time,layer)
0002 %SCRUM_KSLICE:  returns horizontal slice at particular layer.
0003 %
0004 % Horizontal slice is at specified s_rho layer at given time-step for a
0005 % SCRUM 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]=scrum_kslice(cdf,var,[time],[layer])
0011 %
0012 % where
0013 %   cdf:  file name for netCDf file (e.g. 'scrum_c.nc')
0014 %   var:  the variable to select (eg. 'salt' for salinity)
0015 %   time:  time step
0016 %   layer:  in the NetCDF file, s_rho is dimensioned as [0 ... s_rho_length-1].
0017 %           The layer variables which are a function of s_rho run from -1 to 0,
0018 %           which is reversed from ECOM.  In order to have conformity, we will
0019 %           have the parameter 'layer' be 1-based, and have layer 1 correspond
0020 %           to the top layer instead of the bottom.
0021 %
0022 %
0023 %       Examples:
0024 %
0025 %          [s,x,y]=scrum_kslice('scrum_c.nc','salt',2,3);
0026 %              returns the salinity field from the 3rd sigma level
0027 %              at the 2nd time step.
0028 %
0029 %          [elev,x,y]=scrum_kslice('scrum_c.nc','zeta',4);
0030 %              returns the elevation field from the 4th time step
0031 %
0032 %          [depth,x,y]=scrum_kslice('scrum_c.nc','h');
0033 %              returns the depth field
0034 %
0035 if (nargin<2 | nargin>4),
0036   help scrum_kslice; return
0037 end
0038 
0039 % turn off warnings from NetCDf
0040 mexcdf('setopts',0);
0041 
0042 %
0043 % open existing file
0044 ncid=mexcdf('open',cdf,'nowrite');
0045 if(ncid==-1),
0046   disp(['file ' cdf ' not found'])
0047   return
0048 end
0049 
0050 
0051 %
0052 % Acquire the grid.
0053 %
0054 % If "lon_rho" and "lat_rho" are present, grab them.
0055 % Otherwise, get "x_rho" and "y_rho".
0056 [lon_rho_varid, rcode] = ncmex('VARID', ncid, 'lon_rho');
0057 [lat_rho_varid, rcode] = ncmex('VARID', ncid, 'lat_rho');
0058 if ( (lon_rho_varid >= 0) | (lat_rho_varid >= 0) )
0059     x=ncmex('varget',ncid,'lon_rho',[0 0],[-1 -1]);
0060     y=ncmex('varget',ncid,'lat_rho',[0 0],[-1 -1]);
0061 else
0062     x=ncmex('varget',ncid,'x_rho',[0 0],[-1 -1]);
0063     y=ncmex('varget',ncid,'y_rho',[0 0],[-1 -1]);
0064 end
0065 
0066 
0067 [xi_rho_length, eta_rho_length] = size(x);
0068 
0069 
0070 %
0071 % Determine which variables make up the grid.
0072 % This assumes that the variable being asked for is
0073 % defined something like this.
0074 %
0075 % salt ( time, depth, ydimension, xdimension );
0076 %
0077 % The x dimension must be last, and the y dimension must
0078 % be second to last.
0079 [varid, rcode] = ncmex('VARID', ncid, var);
0080 [dud, dud, ndims, vardims, dud, status] = ncmex('varinq', ncid, varid);
0081 
0082 y_dimid = vardims(ndims-1);
0083 x_dimid = vardims(ndims);
0084 
0085 [dud, y_length, status] = ncmex('diminq', ncid, y_dimid);
0086 [dud, x_length, status] = ncmex('diminq', ncid, x_dimid);
0087 
0088 %
0089 % if the lengths of the x and y dimensions do not match up with
0090 % the eta_rho and xi_rho lengths, then we must have something
0091 % that depends upon eta_u, xi_v, or similar.  Adjust the grid
0092 % accordingly.
0093 if ( x_length == (xi_rho_length-1) )
0094     x = (x(1:xi_rho_length-1,:) + x(2:xi_rho_length,:))/2;
0095     y = (y(1:xi_rho_length-1,:) + y(2:xi_rho_length,:))/2;
0096 end
0097 if ( y_length == (eta_rho_length-1) )
0098     x = (x(:,1:eta_rho_length-1) + x(:,2:eta_rho_length))/2;
0099     y = (y(:,1:eta_rho_length-1) + y(:,2:eta_rho_length))/2;
0100 end
0101 
0102 
0103 
0104 
0105 %
0106 % allow for using kslice on 2D, 3D and 4D variables
0107 switch ( nargin )
0108 
0109     %
0110     % 2D variable, such as bathymetry (h)
0111     case 2
0112         [u,ierr]=ncmex('varget',ncid,var,[0 0],[-1 -1]);
0113 
0114     %
0115     % 3D
0116     case 3
0117         [u,ierr]=ncmex('varget',ncid,var,[(time-1) 0 0],[1 -1 -1]);
0118 
0119     %
0120     % 4D
0121     case 4
0122         [s_rho_dimid, status] = ncmex ( 'dimid', ncid, 's_rho' );
0123         [dud, s_rho_length, status] = ncmex ( 'diminq', ncid, s_rho_dimid );
0124         [u,ierr]=ncmex('varget',ncid,var,[(time-1) (s_rho_length-layer) 0 0],[1 1 -1 -1]);
0125 
0126 end
0127 
0128     
0129 
0130 %
0131 % If the appropriate mask is present, use it to mask out the land.
0132 % Find the correct mask by comparing the dimension ids of the
0133 % variable in question and that of each mask variable.
0134 mask_vars = ['mask_rho'; 'mask_u  '; 'mask_v  '; 'mask_psi' ];
0135 for ind = 1:4
0136     mask_var = deblank(mask_vars(ind,:));
0137 
0138     [mask_varid, status] = ncmex ( 'varid', ncid, mask_var );
0139     if ( mask_varid ~= -1 )
0140         [dud, dud, ndims, mask_dimids, dud, status] = ncmex('varinq', ncid, mask_varid);
0141 
0142         if ( ~isempty(find(vardims==mask_dimids(1))) & ~isempty(find(vardims==mask_dimids(2))) )
0143             [mask, status] = ncmex ( 'varget', ncid, mask_var, [0 0], [-1 -1] );
0144             land = find(mask==0);
0145             u(land) = u(land)*NaN;
0146         end
0147 
0148     end
0149 
0150 end
0151 
0152 ncmex('close',ncid);
0153

Generated on Thu 01-Dec-2005 08:50:35 by m2html © 2003