%
% PURPOSE:  Load and plot a .grd (xyz) file from MB System.
%
% AUTHORS: (1) Rob McEwen, 22 May 07
%
% INPUTS:  ncname - a string containing the name of the .grd file, including
%                   the .grd extension.
%
%

if(exist('plg')==1)
    fprintf('loadGrd: plg is already on the stack.\n')
    return
end

if( strncmp(computer,'PCWIN',5))
    pathRoot = '//atlas';
    [status hname] = system('hostname');
    if(regexp(hname,'NORTH'))
      pathRoot = 'C:\Users\dorado1\Documents\MATLAB\Maps';
    end
elseif( strncmp(computer,'MAC',3))
    pathRoot = '/Volumes';
else
    fprintf('Unknown host %s\n',computer);
    return
end

if( exist('loc') )
    if( strncmpi('Port2008Small',loc,9) )
        ncname = 'PortugueseLedge20080424TopoUTM_NoNanSmall.grd';
        ncpath = '~/Maps/';
        decimate=0;
        if(exist('PortugueseLedge20080424TopoUTM_NoNanSmall.mat'))
            fprintf('Loading .mat version from /MATLAB\n');
            tic
            load('PortugueseLedge20080424TopoUTM_NoNanSmall.mat');
            toc
            return
        end
    elseif( strncmpi('portWall2lam50',loc,14) )
        ncname = 'portWall2lam50.mat';
        if( exist(ncname) )
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load(ncname);
	  toc
	  return
	else
	  fprintf([ncname ' not found.  See notes to create the .mat\n']);
	end
    elseif( strncmpi('Port2008',loc,8) )
        ncname = 'PortugueseLedge20080424TopoUTM_NoNan.grd';
        ncpath = '~/Maps/';
        decimate=0;
        if(exist('PortugueseLedge20080424TopoUTM_NoNan.mat'))
            fprintf('Loading .mat version from /MATLAB\n');
            tic
            load('PortugueseLedge20080424TopoUTM_NoNan.mat');
            toc
            return
        end
    elseif( strncmpi('portWall',loc,7) )
      ncname = 'portWall.mat';
        if( exist(ncname) )
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load(ncname);
	  toc
	  return
	else
	  fprintf([ncname ' not found.  See notes to create the .mat\n']);
	end
    elseif( strncmpi('portSlope',loc,10) )
      ncname = 'portSlope.mat';
      ncname = 'portSlope22V2.mat';
        if( exist(ncname) )
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load(ncname);
	  toc
	  return
	else
	  fprintf([ncname ' not found.  See notes to create the .mat\n']);
	end
      elseif( strncmpi('Portuguese',loc,3) )
        %
        % This map is offset from the 2008 versions above.
        %
        ncname = 'PortugueseLedgeTopo1m_v20140527.grd';
        ncpath = strcat(pathRoot,...
                        '/MappingAUVOps2014/PortugueseLedge2014/Figures/');
        decimate=0;
        if(exist('PortugueseLedgeTopo1m_v20140527.mat'))
            fprintf('Loading .mat version from /MATLAB\n');
            tic
            load('PortugueseLedgeTopo1m_v20140527.mat');
            toc
            return
        end
    elseif( strncmpi('MC2800',loc,5) )
        ncname = 'MontereyCanyon2800m_1m_Geo.grd';
        ncpath = strcat(pathRoot,...
                        '/OceanImaging2014/20141217AGUPoster/');
        if(exist('MontereyCanyon2800m_1m_Geo.mat'))
            fprintf('Loading .mat version from /MATLAB\n');
            tic
            load('MontereyCanyon2800m_1m_Geo.mat');
            toc
            return
        end
    elseif( strncmpi('MC1850',loc,8) )
        ncname = 'MontereyCanyon1850BendNorth_MAUV1m.grd';
        ncname = 'MontereyCanyon1850BendEast_MAUV1m.grd';
        ncpath = strcat(pathRoot,...
                        '/MappingAUVOps2014/MontereyBay/MissionPlanning/');
	ncname = 'MAUVTopo1m_CCEDeep.grd';
	ncpath = './';
        %if(exist('MontereyCanyon1850BendEast_MAUV1m.mat'))
        if(exist('MAUVTopo1m_CCEDeep.mat'))
            fprintf('Loading .mat version from /MATLAB\n');
            tic
	    %load('MontereyCanyon1850BendEast_MAUV1m.mat');
            load('MAUVTopo1m_CCEDeep.mat');
            toc
            return
        end
    elseif( strncmpi('MC1850-2',loc,8) )
        ncname = 'MontereyCanyon_2m_WB.grd';
        ncname = 'MontereyBayCanyon_10m_WB.grd';
        ncname = 'MontereyCanyon_5m_WB.grd';
        ncpath = strcat(pathRoot,...
                        '/MappingAUVOpsStuff/ALLMontereyBay/Figures/');
        if(exist('MontereyCanyon_5m_WB.mat'))
            fprintf('Loading .mat version from /MATLAB\n');
            tic
            load('MontereyCanyon_5m_WB.mat');
            toc
            return
        end
    elseif( strncmpi('wallside',loc,5) )
      %
      % This has to be ahead of 'wall',loc,4.
      %
      ncname = 'TheWall_2m.grd';
      ncpath = '/Users/rob/Maps/';
        if(exist('TheWall_2m.mat'))
            fprintf('Loading .mat version from /MATLAB\n');
            tic
            load('TheWall_2m.mat');
            toc
            return
        end
    elseif( strncmpi('wall',loc,4) )
      %
      % Way, way too big:
      % ncname = 'Topo_MAUV_SoquelCanyonTo1850Bend_2m.grd';
      % ncpath = strcat(pathRoot,...
      %                 '/MappingAUVOps2014/MontereyBay/MissionPlanning/');
      % ncpath = '/Users/rob/Maps/';
      %
      % Close, but no cigar:
      % ncname = '20081201Topoaround1000m.grd';
      % ncpath = strcat(pathRoot,...
      %                 '/MappingAUVOpsStuff/ALLMontereyBay/Figures/');
      %   if(exist('20081201Topoaround1000m.mat'))
      %       fprintf('Loading .mat version from /MATLAB\n');
      %       tic
      %       load('20081201Topoaround1000m.mat');
      %       toc
      %       return
      %   end
      %
      % Dave found it.  The real deal.  This is the original 2012 collage .grd
      % from which he made the route file.  The side-looking data exsits as a
      % .grd, but is too narrow a swath to be useful for simulation.
      % ncname = 'Topo_UpperCanyonTesting_2m.grd';
      % ncpath = strcat(pathRoot,...
      %                 '/MappingAUVOps2012/IcebergTesting/MissionPlanning/');
      %   if(exist('Topo_UpperCanyonTesting_2m.mat'))
      %       fprintf('Loading .mat version from /MATLAB\n');
      %       tic
      %       load('Topo_UpperCanyonTesting_2m.mat');
      %       toc
      %       return
      %   end
      %
      % But this is also way too big.  I used "gmt grdcut" to cut out a subregion.
      ncname = 'Topo_UpperCanyonTesting_2m_Wall.grd';
      ncpath = '/Users/rob/Maps/';
        if(exist('Topo_UpperCanyonTesting_2m_Wall.mat'))
            fprintf('Loading .mat version from /MATLAB\n');
            tic
            load('Topo_UpperCanyonTesting_2m_Wall.mat');
            toc
            return
        end
    elseif( strncmpi('upper',loc,5) )
      ncname = 'UpperMontereyCanyon10m.grd';
      ncpath = '/Users/rob/Maps/';
        if(exist('UpperMontereyCanyon10m.mat'))
            fprintf('Loading .mat version from /MATLAB\n');
            tic
            load('UpperMontereyCanyon10m.mat');
            toc
            return
        end
    elseif( strncmpi('plateau',loc,3) )
        if(exist('MontereyCanyonShallowTopo10m.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('MontereyCanyonShallowTopo10m.mat');
	  toc
	  return
	else
	  if( strcmp(computer,'PCWIN64'))
              ncpath = strcat(pathRoot,...
                              '/MappingAUVTesting/MissionPlanning/');
              %               '/MappingAUVOps2015/20150928m3/');
	    ncname = 'MontereyCanyonShallowTopo10m.grd';
	  elseif( strcmp(computer,'MACI64'))
	    ncpath = strcat('/Users/rob/Documents/MATLAB/');
	    ncname = 'MontereyCanyonShallowTopo10m.grd';
	  else
	    fprintf('Unrecognized host.  Exit.\n');
	    return
	  end
	end
    elseif( strncmpi('axial',loc,3) )
      ncname  = 'AxialSummit_MAUV_All_Topo1m_UTM.grd'
      mapName = 'AxialSummit_MAUV_All_Topo1m_UTM.mat'
      if(exist(mapName) )
	fprintf('Loading .mat version from /MATLAB\n');
	tic
	load(mapName);
	toc
	return
      else
	fprintf('%s Map not found.\n',mapName);
	return
      end
    elseif( length(loc) == 9 && strncmpi('Axial2023',loc) )
      ncname  = 'Axial2023.grd'
      mapName = 'Axial2023.mat'
      if(exist(mapName) )
	fprintf('Loading .mat version from /MATLAB\n');
	tic
	load(mapName);
	toc
	return
      else
	fprintf('%s Map not found.\n',mapName);
	return
      end
    elseif( strncmpi('axial5',loc,6) )
        if(exist('AxialEM_post2015_AUV2015Bounds_Topo5m.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('AxialEM_post2015_AUV2015Bounds_Topo5m.mat');
	  toc
	  return
	else
	  if( strcmp(computer,'PCWIN64'))
              ncpath = strcat(pathRoot,...
                              '/MappingAUVTesting/MissionPlanning/');
              %               '/MappingAUVOps2015/20150928m3/');
	    ncname = 'MontereyCanyonShallowTopo10m.grd';
	  elseif( strcmp(computer,'MACI64'))
	    ncpath = strcat('/Users/rob/Maps/');
	    ncname = 'AxialEM_post2015_AUV2015Bounds_Topo5m.grd';
	  else
	    fprintf('Unrecognized host.  Exit.\n');
	    return
	  end
	end
    elseif( strncmpi('cascadia',loc,3) )
      ncname  = 'CascadiaOps_Topo_25m_UTM.grd';
      mapName = 'CascadiaOps_Topo_25m_UTM.mat';
      if(exist(mapName) )
	fprintf('Loading %s from /MATLAB\n',mapName);
	tic
	load(mapName);
	toc
	return
      else
	fprintf('%s Map not found.\n',mapName);
	return
      end
    elseif( strncmpi('lowertrin',loc,3) )
      ncname  = 'Lower_TrinidadCanyon_Long.grd';
      mapName = 'Lower_TrinidadCanyon_Long.mat';
      if(exist(mapName) )
	fprintf('Loading %s from /MATLAB\n',mapName);
	tic
	load(mapName);
	toc
	return
      else
	fprintf('%s Map not found.\n',mapName);
	return
      end
    elseif( strncmpi('berg',loc,4) )
        if(exist('theBergSerm.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('theBergSerm.mat');
	  toc
	  return
	else
	  fprintf(['.mat version not found.  Run octree2patches to create' ...
		   ' the .mat\n']);
	end
    elseif( strcmpi('newBerg5',loc) )
        if(exist('theNewBerg5.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('theNewBerg5.mat');
	  toc
	  return
	else
	  fprintf(['.mat version not found.  Run octree2patches to create' ...
		   ' the .mat\n']);
	end
    elseif( strcmpi('newBerg',loc) )
        if(exist('theNewBerg.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('theNewBerg.mat');
	  toc
	  return
	else
	  fprintf(['.mat version not found.  Run octree2patches to create' ...
		   ' the .mat\n']);
	end
    elseif( strcmpi('CleanBerg5',loc) )
      % 5 meter resolution.
        if(exist('CleanBerg5.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('CleanBerg5.mat');
	  toc
	  return
	else
	  fprintf(['.mat version not found.  Run octree2patches to create' ...
		   ' the .mat\n']);
	end
    elseif( strcmpi('CleanBerg2Serm',loc) )
      % 2 meter resolution, in Sermilik UTM coordinates.
        if(exist('CleanBerg2Serm.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('CleanBerg2Serm.mat');
	  toc
	  return
	else
	  fprintf(['.mat version not found.  Run octree2patches to create' ...
		   ' the .mat\n']);
	end
    elseif( strcmpi('CleanBerg',loc) )
      % 1 meter resolution, UTM coord in Sermilik.  CleanBerg00 is this
      % translated to local coordinates.
        if(exist('CleanBerg.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('CleanBerg.mat');
	  toc
	  return
	else
	  fprintf(['.mat version not found.  Run octree2patches to create' ...
		   ' the .mat\n']);
	end
    elseif( strcmpi('CleanBerg00',loc) )
      % 1 meter resolution, UTM coord in local coordinates.
        if(exist('CleanBerg00.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('CleanBerg00.mat');
	  toc
	  return
	else
	  fprintf(['.mat version not found.  Run octree2patches to create' ...
		   ' the .mat\n']);
	end
    elseif( strcmpi('CleanBerg_v2R5',loc) )
      % 5 meter resolution, local coord.
        if(exist('CleanBerg_v2R5.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('CleanBerg_v2R5.mat');
	  toc
	  return
	else
	  fprintf(['.mat version not found.  Run octree2patches to create' ...
		   ' the .mat\n']);
	end
    elseif( strcmpi('CleanBerg_v2R2',loc) )
      % 2m Resolution, local coordinates.  This one contains Steve's smoothed
      % Auv data in auvE, N and U.
        if(exist('CleanBerg_v2R2.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('CleanBerg_v2R2.mat');
	  toc
	  return
	else
	  fprintf(['.mat version not found.  Run octree2patches to create' ...
		   ' the .mat\n']);
	end
    elseif( strcmpi('CleanBerg_v2R2Serm',loc) )
      %
      % Same as the above, except translated to Serm in UTM.  No Auv track.
      %
        if(exist('CleanBerg_v2R2Serm.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('CleanBerg_v2R2Serm.mat');
	  toc
	  return
	else
	  fprintf(['.mat version not found.  Run octree2patches to create' ...
		   ' the .mat\n']);
	end
    elseif( strncmpi('EMARK',loc,4) )
      ncname = 'EMARK_vents_MAUVoverEM124_Topo1m_UTM.mat';
        if( exist(ncname) )
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load(ncname);
	  toc
	  return
	else
	  fprintf([ncname ' not found.  See notes to create the .mat\n']);
	  %
	  % Dave has already converted this one to UTM.  So just read it
          % in directly.  Don't do the conversion at the end of this
          % file. 
	  %
	  % [plg.eastings plg.northings plg.z]=...
	  % grdread2('EMARK_vents_MAUVoverEM124_Topo1m_UTM.grd');
	  % plg.northings=plg.northings';
	  % save SantaMonica800mMound_Topo1m_UTM_Small1m
	end
    elseif( strncmpi('shapes',loc,5) )
        if(exist('shapes.mat'))
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load('shapes.mat');
	  toc
	  return
	else
	  fprintf(['.mat version not found.  Run octree2patches to create' ...
		   ' the .mat\n']);
	end
    elseif( strncmpi('santav2',loc,7) )
      ncname = 'SantaMonica800mMound_Topo1m_UTM_SmallV2.mat';
        if( exist(ncname) )
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load(ncname);
	  toc
	  return
	else
	  fprintf([ncname ' not found.  See notes to create the .mat\n']);
	  %
	  % Dave has already converted this one to UTM.  So just read it
          % in directly.  Don't do the conversion at the end of this
          % file. 
	  %
	  % [plg.eastings plg.northings plg.z]=...
	  % grdread2('SantaMonica800mMound_Topo1m_UTM_Small1m.grd');
	  % plg.northings=plg.northings';
	  % save SantaMonica800mMound_Topo1m_UTM_Small1m
	end
    elseif( strncmpi('santa',loc,5) )
      ncname = 'SantaMonica800mMound_Topo1m_UTM_Small1m.mat';
        if( exist(ncname) )
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load(ncname);
	  toc
	  return
	else
	  fprintf([ncname ' not found.  See notes to create the .mat\n']);
	  %
	  % Dave has already converted this one to UTM.  So just read it
          % in directly.  Don't do the conversion at the end of this
          % file. 
	  %
	  % [plg.eastings plg.northings plg.z]=...
	  % grdread2('SantaMonica800mMound_Topo1m_UTM_Small1m.grd');
	  % plg.northings=plg.northings';
	  % save SantaMonica800mMound_Topo1m_UTM_Small1m
	end
    elseif( strncmpi('sur',loc,3) )
      ncname = 'SurRidge.mat';
        if( exist(ncname) )
	  fprintf('Loading .mat version from /MATLAB\n');
	  tic
	  load(ncname);
	  toc
	  return
	else
	  fprintf([ncname ' not found.  See notes to create the .mat\n']);
	  %
	  % Dave has already converted this one to UTM.  So just read it
          % in directly.  Don't do the conversion at the end of this
          % file. 
	  %
	  % [plg.eastings plg.northings plg.z]=...
	  % grdread2('SantaMonica800mMound_Topo1m_UTM_Small1m.grd');
	  % plg.northings=plg.northings';
	  % save SantaMonica800mMound_Topo1m_UTM_Small1m
	end
    end
end
%
% Single quotes specifiy the variable name.
if(exist('ncname') == 0)
    fprintf('Please put the file name in ''ncname''\n');
    return
end
if(exist('ncpath') == 0)
    fprintf('Please put the path name in ''ncpath''\n');
    return
end
ncpathname = [ncpath,ncname];
%
% Absence of single quotes passes what the variable contains to "exist".
if(exist(ncpathname) ~=2)
    fprintf('Can''t find the file %s.  Exit.\n',ncpathname);
    return
end
fprintf('Loading Netcdf file %s.\n',ncpathname);
tic
if( 0 )
    %
    % This almost worked, but somehow the axes were reversed.  I couldn't get
    % past it, so I found a routine grdread2 on the Matlab File Exchange
    % that reads in a .grd correctly.
    %
    %ncdisp(ncname,'x');
    %ncdisp(ncname,'y');
    x = ncread(ncname,'x');
    y = ncread(ncname,'y');
    %
    % Pre-allocate memory:
    z=NaN(length(y),length(x));
    z = ncread(ncname,'z');
else
    [plg.lon plg.lat plg.z] = grdread2(ncpathname);
    plg.z = double(plg.z);
end
toc
fprintf('Finished loading bathymetry data.\n');
%
% Here I'm using the mapping convention; x is eastings, y northings.
%

%
% Convert to UTM for plotting
%
tic
fprintf('Converting to UTM.\n');
len = length(plg.lon);
indx = find( plg.lon > 180 );
plg.lon(indx) = plg.lon(indx)-360;
zone = utmzone(plg.lat(1), plg.lon(1));
zoneEnd = utmzone(plg.lat(end), plg.lon(end));
if( ~strcmp(zone,zoneEnd))
    fprintf('This data crosses UTM zone boundary. Using zone %s\n', ...
            zone);
end
mstruct = defaultm('utm');
mstruct.zone = zone;
mstruct  = defaultm(mstruct);
%
% This gets around mfwdtran's requirement that the x and y axes be
% matrices of the same dimension of z.
%
[plg.glon,plg.glat]=meshgrid(plg.lon,plg.lat); 
%
% LL2utm is ryjo's code that runs mfwdtran.  There is also a ll2utm on the
% File Exchange that transforms lat lon to utm using first principles
% without dependence on the mapping toolbox.
%
[plg.eastings,plg.northings]=LL2utm(plg.glon,plg.glat);
%rmfield(plg,'gx')
%rmfield(plg,'gy')
plg.glon = 0;
plg.glat = 0;
plg.eastings=plg.eastings(1,:); 
plg.northings=plg.northings(:,1);
toc
fprintf('Finished UTM conversion. End loadGrd.\n');
%
