function [lon, lat, vardata, jd, pn, status] = prof_xsect ( pointer_file, varname, depth_i, opt_plot )
% PROF_XSECT:  computes cross section of profiles.
% USAGE:  [lon, lat, var, jd, pn, status] = prof_xsect ( pointer_file, varname, depth_i )
%
% PARAMETERS:
% INPUT:
%   pointer_file:  
%      list of profile netcdf files
%   varname:  
%      name of epic variable to use.  Should be something like 'T_20',
%      'P_4023', etc.
%   depth_i:  
%      depth at which to take cross section
%   opt_plot:
%      Optional, if present and equal to 1, put up a gridded plot of
%      the data.
%
% OUTPUT:
%   lon, lat:  
%      coordinates of profiles
%   var:       
%      data interpolated at 'depth_i'
%   jd:     
%      julian day
%   pn: profile name
%   status:  
%      -1 if the routine fails for some reason
%
%

if ( nargin < 4 )
    opt_plot = 0;
end

status = -1;
lon = -1;
lat = -1;
vardata = -1;
jd = -1;

afid = fopen ( pointer_file, 'r' );
if ( afid == -1 )
    fprintf ( 2, 'Could not open pointer file ''%s''.\n', pointer_file );
    help prof_xsect;
    return;
end


mexcdf ( 'setopts', 0 );


%
% process each epic file
epic_file = fgetl ( afid );
data_count = 0;
while ( epic_file ~= -1 )

    ncid = mexcdf ( 'open', epic_file, 'nowrite' );
    if ( ncid == -1 )
	fprintf ( 2, 'Could not open %s.\n', epic_file );
	return;
    end
    %
    % Get the lat, lon, time data
    [lat, status] = mexcdf('VARGET1', ncid, 'lat', [0]);
    if ( status == -1 )
	fprintf ( 2, 'Could not get latitude from %s.\n', epic_file );
	return;
    end
    [lon, status] = mexcdf('VARGET1', ncid, 'lon', [0]);
    if ( status == -1 )
	fprintf ( 2, 'Could not get longitude from %s.\n', epic_file );
	return;
    end
    [time, status] = mexcdf('VARGET1', ncid, 'time', [0]);
    if ( status == -1 )
	fprintf ( 2, 'Could not get time from %s.\n', epic_file );
	return;
    end
    [time2, status] = mexcdf('VARGET1', ncid, 'time2', [0]);
    if ( status == -1 )
	fprintf ( 2, 'Could not get time2 from %s.\n', epic_file );
	return;
    end
    time = time + time2/8.64e7;


    %
    % get the indicated variable along with the depth column
    [vardata, status] = mexcdf('VARGET', ncid, varname, [0 0 0 0], [-1 -1 -1 -1]);
    if ( status == -1 )
	fprintf ( 2, 'Could not get %s from %s.\n', varname, epic_file );
	return;
    end

    [depth_data, status] = mexcdf('VARGET', ncid, 'depth', [0 0 0 0], [-1 -1 -1 -1]);
    if ( status == -1 )
	fprintf ( 2, 'Could not get depth from %s.\n', varname, epic_file );
	return;
    end


    %
    % We only use this file if the input depth is between the minimum and
    % maximum depths in the file.
    if ( (min(depth_data) <= depth_i) & (max(depth_data) >= depth_i ) )
     
        data_count = data_count+1;
        pn{data_count}=epic_file;
	keeper_lon(data_count) = lon;
	keeper_lat(data_count) = lat;
	keeper_time(data_count) = time;
	keeper_var(data_count) = interp1(depth_data, vardata, depth_i);


    end

    epic_file = fgetl ( afid );

end

fclose ( afid );

if ( opt_plot )
    min_lon = min(keeper_lon(:));
    max_lon = max(keeper_lon(:));
    min_lat = min(keeper_lat(:));
    max_lat = max(keeper_lat(:));
    
    [LON,LAT] = meshgrid ( [min_lon:(max_lon-min_lon)/10:max_lon], ...
    		       [min_lat:(max_lat-min_lat)/10:max_lat] );
    
    [lon_i, lat_i, var_i] = griddata ( keeper_lon, keeper_lat, keeper_var, LON, LAT );
    pslice(lon_i, lat_i, var_i);
    hold on
    plot (keeper_lon, keeper_lat,'k*');
end


lon = keeper_lon;
lat = keeper_lat;
vardata = keeper_var;
jd = keeper_time;
status = 1;

return;
