function string = caris2gmt(file,titlestring,nlon,nlat)
% Converts a CARIS Base to ASCII file to a grd file and makes a gmt plot
%
% Usage
%   caris2gmt(file,titlestring);
% Inputs
%   FILE - filename of CARIS Base to ASCII lat, long, depth file
%          'N' and 'W' letters must be edited out of the ascii file but
%           signs of longitude and depth will be adjusted automatically
%   TITLESTRING - Plot title
%   NLON - Number of longitude grid points (default - automatic)
%   NLAT - Number of latitude grid points (default - automatic)
% Outputs
%   None

%file = 'x1027_prelim.xyz';
%titlestring ='Hole 1027 Langseth EM122 Bathymetry';

if nargin<2
  titlestring = ' ';
end

disp('Warning - function caris2gmt must be run from terminal window')

% Load data into variable a
eval(['load ' file ' -ascii']);
if abs(file(1))>=48 & abs(file(1))<=57
  file2 = ['X' file];
else
  file2 = file;
end
eval(['a = ' file2(1:end-4) ';']);

% Correct signs and save
if (a(1,2)>0); a(:,2) = -a(:,2); end
if (a(1,3)>0); a(:,3) = -a(:,3); end;
save junk a -ascii

% Grid limits
latmax = round(max(a(:,1))*10000)/10000;
latmin = round(min(a(:,1))*10000)/10000;
lonmax = round(max(a(:,2))*10000)/10000;
lonmin = round(min(a(:,2))*10000)/10000;

% Grid increments
i = find(a(:,1)~=a(1,1),1,'first');
if nargin<3
  dlon = abs(a(2,2) - a(1,2));
  nlon = round((lonmax - lonmin)/dlon) + 1;
end
dlon = floor((lonmax - lonmin)/(nlon -1)*1e8)/1e8;
if nargin<4
  dlat = a(i,1) - a(1,1);
  nlat = round((latmax - latmin)/dlat) + 1;
end
dlat = floor((latmax - latmin)/(nlat - 1)*1e8)/1e8;
disp(['Grid size is (' int2str(nlon) ',' int2str(nlat) ')'])

% Make grd file
eval(['!xyz2grd junk -G' file(1:end-3) 'grd -I' num2str(dlon,8) '/' ...
       num2str(dlat,8) ' -R' num2str(lonmin,6) '/' num2str(lonmax,6) ...
       '/' num2str(latmin) '/' num2str(latmax) ' -:']);
!rm -f junk

% Map scale
scale1 = 7/((lonmax-lonmin)*cosd(48));
scale2 = 8.0/(latmax-latmin);
scale = min(scale1,scale2)*cosd(48);
%xpos = ((lonmax-lonmin)*cosd(48)*scale)/2+1.0;
     

% Make map
string = ['!gmtmap.sh ' file(1:end-4) ' ' num2str(lonmin,6) ' ' ...
       num2str(lonmax,6) ' ' num2str(latmin,6) ' ' num2str(latmax,6) ...
       ' ' num2str(scale,6) ' 2.0 -0.5 ''' titlestring ''''];
eval(string);