function [u,y,z,jd] =islice(cdf,var,time,iindex)
%  ISLICE  Extracts a vertical slice of a scalar variable 
%          along row I from a 4D ECOM NetCDF file. 
%
%  USAGE: [u,y,z,jd]=islice(cdf,var,time,i_index)
%    cdf = NetCDF file name (string) (e.g. 'ecomsi.cdf')
%    var = variable name (string) (e.g. 'salt' for salinity)
%    time = time step (integer) (e.g. 1 for 1st time step)
%    i_index = I index along which to take vertical slice
%
%    u = 2D array of the variable "var"
%    y = distance along slice in km (assuming units in netCDF file are m)
%    z = depth in meters.
%   jd = Julian time of slice (type "help julian" for details)
%
%
%   Example:  to extract temperature along I=20
%             from the 3rd time step in the file 'ecom3d.cdf', do:
%
%          [temp,y,z]=islice('ecom3d.cdf','temp',3,20);
%
%   Tip: Visualize these arrays using EXTCONTOUR,
%       [e.g. EXTCONTOUR(y,z,temp,'label','fill');]

% Rich Signell rsignell@usgs.gov

if nargin~=4 ,
  help islice; return
end
mexcdf('setopts',0);
ncid=mexcdf('open',cdf,'nowrite');
[name, ny]=mexcdf('diminq',ncid,'ypos');
[name, nz]=mexcdf('diminq',ncid,'zpos');
u=mexcdf('varget',ncid,var,[(time-1) 0 2 (iindex-1)],[1 nz-1 ny-4 1],1);
if nargout>1,
  d=mexcdf('varget',ncid,'depth',[2 (iindex-1)],[ny-4 1],1);
  ind=find(d==-99999.);
  if(~isempty(ind)),
    d(ind)=d(ind)*0;
    u(ind,:)=u(ind,:)*NaN;
  end
  sigma=mexcdf('varget',ncid,'sigma',0,nz,1);
  sigma=0.5*(sigma(1:nz-1)+sigma(2:nz));
  h2=mexcdf('varget',ncid,'h2',[0 (iindex-1)],[ny 1],1);
  y=cumsum(h2)/1000.;
  y=y(3:ny-2);
  y=y*ones(1,nz-1);
  z=d*(sigma');
  ind=find(isnan(z));
  z(ind)=zeros(size(ind));
end
% get time stuff
if(nargout>3)
[nam,nt]=mexcdf('diminq',ncid,'time');
base_date=[0 0 0 0 0 0];
base_date(1:3)=mexcdf('attget',ncid,'global','base_date');
t=mexcdf('varget',ncid,'time',time-1,1);
jd0=julian(base_date);
jd=jd0+t;
end
mexcdf('close',ncid);
