function hout = worlddatamap(varargin)
%WORLDDATAMAP Map a region using data from a shapefile or data grid
%
% WORLDDATAMAP with no arguments presents a menu of region names from
% which to choose, and then displays the boundaries of the selected area
% using vector data from the default shapefile 'landareas.shp'. The
% boundaries are symbolized as lines.
%
% WORLDDATAMAP(REGION) makes a map of the specified region using vector
% data from the default shapefile, 'landareas.shp'. REGION is a string or
% a cell array of strings. The projection and limits are suitable to the
% part of the world specified by REGION. Permissible strings include
% names of continents and countries as well as 'World', 'North Pole',
% 'South Pole', and 'Pacific'. The region boundaries are symbolized as
% lines, provided vector data from the default shapefile exist within the
% region boundary.
%
% WORLDDATAMAP(REGIONONLY) adds only the specified region to the map. If
% any of the names in the REGIONONLY string or cell array of strings end
% with 'only', only the requested regions are displayed; therfore, only
% one string needs the 'only' suffix.
%
% WORLDDATAMAP(LATLIM, LONLIM) makes a map of data within a geographic
% region specified as latitude and longitude limits in degrees. LATLIM
% and LONLIM are two-element vectors of the form [southern_limit
% northern_limit] and [western_limit eastern_limit], respectively.
% Latitudes range from [-90 90] and longitudes range from [-180 180] for
% the default shapefile, 'landareas.shp'.
%
% WORLDDATAMAP(...,VECDISPLAYTYPE) where VECDISPLAYTYPE is a string,
% controls the rendering of the vector data. If omitted, VECDISPLAYTYPE
% defaults to 'line' Permissable values for VECDISPLAYTYPE are:
%
% Value Type of render
% -------- --------------
% 'line' line (default)
% 'patch' patch
%
% WORLDATADMAP(Z, REFVEC) derives the map limits from the extent of a
% regular data grid Z, with a 1-by-3 referencing vector REFVEC. The
% regular data grid is displayed using MESHM, and the region boundaries
% are sybmolized as lines, from the default shapefile, 'landareas.shp'.
%
% WORLDDATAMAP(Z, REFVEC, SURFDISPLAYTYPE) where SURFDISPLAYTYPE is a
% string, controls the rendering of the surface data. If omitted,
% SURFDISPLAYTYPE defaults to 'mesh'. Permissable values for
% SURFDISPLAYTYPE are:
%
% Value Type of render
% -------- -----------------------------------------
% 'mesh' surface in z=0 plane using MESHM (default)
% 'mesh3d' surface, set aspect ratio
% 'dem' surface, set digital elevation colormap
% 'dem3d' surface, set digital elevation colormap, aspect ratio
% 'lmesh3d' surface, apply lighting
% 'ldem3d' surface, apply lighting, digital elevation colormap
% 'line' only line boundaries defined by extent of data grid
% 'patch' only patch boundaries defined by extent of data grid
%
% WORLDDATAMAP(...,FEATURETYPE) where FEATURETYPE is a string or cell
% array of strings, controls the display of additional feature data. The
% features are obtained from the shapefiles, worldcities.shp,
% worldrivers.shp, or worldlakes.shp. If omitted, FEATURETYPE defaults to
% 'none'. Multiple features are allowed if FEATURETYPE is a cell array of
% strings. Permissable values for FEATURETYPE are:
%
% Value Type of feature
% -------- --------------
% 'none' no additional features (default)
% 'cities' draw all World cities from shapefile
% 'rivers' draw all World rivers from shapefile
% 'lakes' draw all World lakes from shapefile
% 'labels' draw labels (if LabelLat, LabelLon defined in shapefile)
%
% WORLDDATAMAP(...,DISPLAYTYPEONLY, FEATURETYPE) where DISPLAYTYPEONLY is
% either VECDISPLAYTYPE ('lineonly' or 'patchonly') or SURFDISPLAYTYPE
% ('meshonly', 'mesh3donly', etc) with the suffix 'only', prohibits the
% features defined by FEATURETYPE ('cities', 'rivers', etc) from being
% rendered.
%
% WORLDDATAMAP(...,DISPLAYTYPE, FEATURETYPEONLY) where FEATURETYPEONLY is
% a FEATURETYPE with the suffix 'only' ('citiesonly', 'riversonly', etc),
% and DISPLAYTYPE is either a VECDISPLAYTYPE ('line', 'patch') or
% SURFDISPLAYTYPE ('mesh', 'dem', etc), symbolizes only the requested
% feature.
%
% WORLDDATAMAP(FILENAME, ...) where FILENAME is a string ending in .shp
% or .dbf, reads the vector data from shapefile FILENAME instead of the
% default file, landareas.shp.
%
% WORLDDATAMAP({FILENAME, REGIONATTRIBUTE}, ...) uses the REGIONATTRIBUTE
% field name rather than the default, 'Name' to search for a specified
% REGION. If REGIONATTRIBUTE is missing, it defaults to the first match
% to the following case-insensitive list:
% 'name','cntry_name','country_name','long_name'.
%
% h = WORLDDATAMAP(...) returns the handle of the map axes.
%
% WORLDDATAMAP uses TIGHTMAP to set the axis limits tight around the map.
% If you change the projection, or need more white space around the map
% frame, use TIGHTMAP again or AXIS AUTO.
%
% Example 1
% ---------
% % World map with coastlines
% figure
% worlddatamap('World')
%
% Example 2
% ---------
% % Worldmap with land areas, major lakes and rivers, and cities and
% % populated places
% figure
% worlddatamap('World',{'lakes','rivers','cities'})
%
% Example 3
% ---------
% % Patch map of Antarctica
% figure
% worlddatamap('antarctica','patch')
%
% Example 4
% ---------
% % Map of Africa and India with major cities and populated places
% figure
% worlddatamap({'Africa','India'},'cities')
%
% Example 5
% ---------
% % Map of the geoid over South America and the central Pacific
% figure
% worlddatamap([-50 50],[160 -30])
% load geoid
% geoshow(geoid, geoidrefvec, 'DisplayType', 'texturemap');
% load coast
% geoshow(lat, long)
%
% Example 6
% ---------
% % Map of terrain elevations in Korea
% figure
% load korea
% worlddatamap(map, refvec,'dem')
%
% Example 7
% ---------
% % Map of the United States of America
% figure
% worlddatamap USA
%
% Remarks:
% WORLDDATAMAP requires MATLAB version 7.0.4 (R14SP2) (or higher) and the
% Mapping Toolbox version 2.1 (or higher).
%
% See also MESHM, SHAPEREAD, USAMAP, WORLDMAP.
% Copyright 2005 The MathWorks, Inc.
% Author: Kelly G. Luetkemeyer
% Define the external data API
api = defineExternalDataApi;
% Verify Mapping Toolbox version for ancillary data
verifyVersion(api);
% Parse and check input arguments
[worldmapArgs, dataset, attribName, displayType, featureType, latlim, lonlim, ...
region, Z, refvec, regionOnly, listSelect] = parseInputs(varargin{:});
% Verify the dataset
[worldShapeFile, attribName, usingDefaults] = verifyDataSet(api, dataset, attribName);
api.ExternalData(1).Filename = worldShapeFile;
% Change region to lat-lon limit if not using default attribute name.
if ~usingDefaults && ~isempty(region)
worldmapArgs = getLatLonLimitsFromFile(dataset, attribName, region, worldmapArgs);
end
% Display the base worldmap
if ~listSelect
h = worldmap(worldmapArgs{:});
else
% Select the region from the list dialog
try
h = worldmap(worldmapArgs{:});
catch
% Cancel was clicked.
close
return
end
end
tightmap loose
% Obtain latlim and lonlim from the projection if not specified
if ~isempty('latlim') || ~isempty('lonlim')
m = gcm;
latlim = m.maplatlimit;
lonlim = m.maplonlimit;
end
% Change [0 360] longtitude to [-180 180]
latlim = latlim(:)';
lonlim = lonlim(:)';
if any(lonlim>180)
lonlim = npi2pi(lonlim);
end
if lonlim == [0 0]
lonlim = [-180 180];
end
% Display the grid data if specified
if ~isempty(Z) && ~isempty(refvec)
addGridData(Z, refvec, displayType);
end
% Display the vector data
addWorldVectorData(api, region, regionOnly, latlim, lonlim, ...
displayType, featureType, worldShapeFile, attribName);
% Return handle if required
if nargout == 1
hout = h;
end
%--------------------------------------------------------------------------
function api = defineExternalDataApi
% Assign required version (or greater) and data filesnames
api.requiredVersion = '2.1';
api.defShapeName = 'landareas';
api.defAttribName = 'Name';
api.ExternalData(1).Filename = 'landareas.shp';
api.ExternalData(1).DrawFcn = '';
api.ExternalData(1).ValidTypes = {};
api.ExternalData(2).Filename = 'worldcities.shp';
api.ExternalData(2).DrawFcn = @drawCities;
api.ExternalData(2).ValidTypes = {'cities','citiesonly', 'all'};
api.ExternalData(3).Filename = 'worldlakes.shp';
api.ExternalData(3).DrawFcn = @drawLakes;
api.ExternalData(3).ValidTypes = {'lakes','lakesonly', 'all'};
api.ExternalData(4).Filename = 'worldrivers.shp';
api.ExternalData(4).DrawFcn = @drawRivers;
api.ExternalData(4).ValidTypes = {'rivers','riversonly', 'all'};
api.ExternalData(5).Filename = 'DEFAULTFILE';
api.ExternalData(5).DrawFcn = @drawLabels;
api.ExternalData(5).ValidTypes = {'labels','labelsonly', 'all'};
api.mappingShapeFiles = {api.ExternalData(1:4).Filename};
%--------------------------------------------------------------------------
function [worldmapArgs, dataset, attribName, displayType, featureType, latlim, lonlim, ...
region, Z, refvec, regionOnly, listSelect] = parseInputs(varargin)
% Default values for inputs
dataset = '';
attribName = '';
latlim = [];
lonlim = [];
region = [];
Z = [];
refvec = [];
regionOnly = false;
if numel(varargin) == 0
worldmapArgs = varargin;
displayType = 'line';
featureType = {'none'};
listSelect = true;
return
else
listSelect = false;
end
% Determine if resolution 'hi', 'lo', has been set.
[worldmapArgs, resolution] = findResolutionArg(varargin);
% Determine if a dataset has been set.
[worldmapArgs, dataset, attribName] = findDataSetArgs(worldmapArgs);
% Determine if type arguments have been set.
[worldmapArgs, displayType, featureType] = findTypeStrings(worldmapArgs);
switch numel(worldmapArgs)
case 0 % WOLRDDATAMAPA
listSelect = true;
case 1 % WORLDDATAMAP(REGION)
region = worldmapArgs{1};
%ensure region is a cell array of strings
if ischar(region)
region = cellstr(region);
elseif ~iscell(region)
error('%s\n\n%s\n\n%s%s%s', ...
'Function WORLDDATAMAP expected REGION to be one of these types:', ...
'string or string cell array', ...
'Instead its type was ', class(region),'.');
end
% Parse only from region, if present
if ~isempty(region)
[region,regionOnly]=removeOnlyFromRegion(region);
if regionOnly
% Remove warning from WORLDMAP for only flag
% by re-setting the region name
worldmapArgs{1} = region;
end
else
region = {'world'};
end
case 2 % WORLDDATAMAP(Z, REFVEC)
% WORLDDATAMAP(LATLIM, LONLIM)
if numel(worldmapArgs{1}) ~= 2 && numel(worldmapArgs{2}) ~= 2
% WORLDDATAMAP(Z, REFVEC)
Z = worldmapArgs{1};
refvec = worldmapArgs{2};
if isempty(displayType)
displayType = 'mesh';
end
else
% WORLDDATAMAP(LATLIM, LONLIM)
latlim = worldmapArgs{1};
lonlim = worldmapArgs{2};
end
otherwise
if isempty(displayType) && ischar(worldmapArgs{end})
typeStrings = getDisplayTypeStrings;
error('%s\n\n%s\n%s\n','Invalid DisplayType string.', ...
'Valid DisplayType strings are:',typeStrings{:});
end
error('Incorrect number of input arguments.')
end
% Set default displayType if not specified
if isempty(displayType)
displayType = 'line';
end
% Add resolution if present
if ~isempty(resolution)
worldmapArgs = {resolution worldmapArgs{:}};
end
%--------------------------------------------------------------------------
function typeStrings = getDisplayTypeStrings
typeStrings = {'patch','line','patchonly','lineonly','mesh',...
'meshonly','mesh3d', 'mesh3donly', 'dem','demonly','dem3d','dem3donly',...
'ldem3d','ldem3donly','lmesh3d','lmesh3donly'};
%--------------------------------------------------------------------------
function typeStrings = getFeatureTypeStrings
typeStrings = { 'cities','citiesonly', ...
'lakes','lakesonly','rivers','riversonly', ...
'labels', 'labelsonly', ...
'all', 'none'};
%--------------------------------------------------------------------------
function [worldmapArgs, resolution] = findResolutionArg(worldmapArgs)
% Find the resolution arg if present.
resolutionArgs = {'hi', 'lo', 'allhi'};
resolution = worldmapArgs{1};
firstArgMatches = ischar(resolution) && any(strcmpi(resolution,resolutionArgs));
if firstArgMatches
worldmapArgs(1) = [];
else
resolution = '';
end
%--------------------------------------------------------------------------
function [worldmapArgs, dataset, attribName] = findDataSetArgs(worldmapArgs)
% Find the dataset arg if present.
attribName = '';
if iscell(worldmapArgs{1})
dataset = worldmapArgs{1}{1};
if numel(worldmapArgs{1}) > 1
attribName = worldmapArgs{1}{2};
end
else
dataset = worldmapArgs{1};
end
if containsShapeExt(dataset);
worldmapArgs(1) = [];
if ~isempty(attribName) && ~ischar(attribName)
error('%s\n%s%s%s', ...
'Function WORLDDATAMAP expected REGIONATTRIBUTE to be type char.', ...
'Instead its type was ', class(attribName),'.');
end
else
dataset = '';
attribName = '';
end
%--------------------------------------------------------------------------
function shapeFilePresent = containsShapeExt(filename)
if ~ischar(filename)
shapeFilePresent = false;
else
filename = lower(filename);
shapeFilePresent = ~isempty( findstr('.shp', filename ) ) || ...
~isempty( findstr('.dbf', filename ) ) ;
end
%--------------------------------------------------------------------------
function verifyVersion(api)
requiredVersion = api.requiredVersion;
mappingShapeFiles = api.mappingShapeFiles;
for j=1:numel(mappingShapeFiles)
if ~exist(mappingShapeFiles{j},'file')
msg1 = sprintf('%s%s\n','Unable to find shapefile: ',mappingShapeFiles{j});
if ~isdeployed
mappingVersion = ver('map');
if isempty(mappingVersion) || isempty(mappingVersion.Version)
installedVersion = '0.0.0';
versionStatement = 'The Mapping Toolbox is not installed.';
else
installedVersion = mappingVersion.Version;
versionStatement = ['The installed version of the Mapping Toolbox is ', ...
installedVersion]; ...
end
if ( installedVersion(1) < requiredVersion(1) ) || ...
(( installedVersion(1) == requiredVersion(1) ) && ...
( installedVersion(3) < requiredVersion(3) ) )
msg2 = sprintf('%s%s%s%s%s\n%s\n%s', ...
'Function ',upper(mfilename), ...
' requires Mapping Toolbox version ',requiredVersion, ' or greater.', ...
versionStatement, ...
'Please install the latest version of the Mapping Toolbox.');
error([msg1 msg2]);
else
error(msg1);
end
else
msg2 = sprintf('%s%s\n%s%s%s',...
'The application is missing shapefile:',mappingShapeFiles{j}, ...
'Compile the application using the -a ', ...
mappingShapeFiles{j},' option to include the file.');
error([msg1 msg2]);
end
end
end
%--------------------------------------------------------------------------
function [worldShapeFile, attribName, usingDefaults] = verifyDataSet(...
api, dataset, attribName)
defAttribName = api.defAttribName;
defShapeFile = api.defShapeName;
% Assign default shapefile
defaultDataSetFile = [defShapeFile '.shp'];
% Assign default cases if empty
if isempty(dataset)
dataSetFile = defaultDataSetFile;
attribName = defAttribName;
usingDefaults = true;
else
dataSetFile = dataset;
if ~exist(dataSetFile, 'file')
[path, name] = fileparts(dataSetFile);
defaultDir = fullfile(fileparts(which(mfilename)), name);
if ~exist(fullfile(defaultDir,dataSetFile),'file');
msg = sprintf('%s%s','Failed to locate SHAPEFILE: ',dataSetFile);
error(msg);
else
fprintf('%s%s%s\n','Adding directory ''', defaultDir, ''' to the path.');
addpath(defaultDir);
end
end
usingDefaults = false;
if isempty(attribName)
allowableAttribNames = {'name','cntry_name','country_name','long_name'};
else
allowableAttribNames = {attribName};
end
info = shapeinfo(dataSetFile);
for j=1:numel(allowableAttribNames)
containsName = strcmpi(allowableAttribNames{j}, {info.Attributes.Name});
if any(containsName)
attribName = info.Attributes(containsName).Name;
break;
end
end
end
if isempty(attribName)
msg = sprintf('%s%s%s','REGIONATTRIBUTE for file ''',dataSetFile,''' is empty.');
error(msg);
end
worldShapeFile = dataSetFile;
%--------------------------------------------------------------------------
function worldmapArgs = getLatLonLimitsFromFile(dataset, attribName, region, worldmapArgs)
% Find the bounding box of the region
s=[];
for k=1:numel(region)
try
shape = shaperead(dataset, 'UseGeoCoords', true, 'Selector',...
{@(name) strcmpi(name,region{k}), attribName});
s = [s, shape];
catch
% If any region error's,
% then the region is not found in the file,
% but might be found by worldmap
return
end
end
switch numel(s)
case 0 % s is empty
return % Region might be found from worldmap
case 1
bbox = [s.BoundingBox];
otherwise
bbox = s(1).BoundingBox;
for num = 2:numel(s)
bbox = [ min(bbox(1,1)',[s(num).BoundingBox(1,1)]'), ...
min(bbox(1,2)', [s(num).BoundingBox(1,2)]'); ...
max(bbox(2,1)', [s(num).BoundingBox(2,1)]'), ...
max(bbox(2,2)', [s(num).BoundingBox(2,2)]') ];
end
end
lon = bbox(:,1:2:end);
lat = bbox(:,2:2:end);
lonlim = [min(lon(:)) max(lon(:))];
latlim = [min(lat(:)) max(lat(:))];
% Snap map limits to increments of INC, with a 1 degree buffer, except
% for limits that are already exact multiples of INC.
buffer = 1;
if strmatch('usa', dataset)
inc = 1; % 1 degree for usa states
else
inc = 5; % compatible with worldmap (5 degree border)
end
if mod(latlim(1),inc) ~= 0
latlim(1) = inc * floor((latlim(1) - buffer)/inc);
end
if mod(lonlim(1),inc) ~= 0
lonlim(1) = inc * floor((lonlim(1) - buffer)/inc);
end
if mod(latlim(2),inc) ~= 0
latlim(2) = inc * ceil((latlim(2) + buffer)/inc);
end
if mod(lonlim(2),inc) ~= 0
lonlim(2) = inc * ceil((lonlim(2) + buffer)/inc);
end
% Ensure that latitude limits remain within [-90 90].
if latlim(1) < -90
latlim(1) = -90;
end
if latlim(2) > 90
latlim(2) = 90;
end
if latlim(1) <= -80 && latlim(2) >=80
% changes projection
latlim = [-90 90];
end
worldmapArgs(1) = [];
worldmapArgs = {latlim, lonlim, worldmapArgs{:}};
%--------------------------------------------------------------------------
function [worldmapArgs, displayType, featureType] = findTypeStrings(worldmapArgs)
% Find the typeString if present
switch numel(worldmapArgs)
case 0
displayType = '';
featureType = {'none'};
case 1
[displayType, featureType, worldmapArgs] = findDisplayFeatureTypeAtEnd(worldmapArgs);
otherwise
% ..., DisplayType, FeatureType
% ..., FeatureType
% ..., DisplayType
[displayType, featureType, worldmapArgs] = findDisplayFeatureTypeAtEnd(worldmapArgs);
if isempty(displayType) && isempty(featureType)
displayType = '';
featureType = {'none'};
elseif isempty(displayType) && ischar(worldmapArgs{end})
displayTypeStrings = getDisplayTypeStrings;
argMatchesDisplayTypeString = any(strcmpi(worldmapArgs{end},displayTypeStrings));
if argMatchesDisplayTypeString
displayType = lower(worldmapArgs{end});
worldmapArgs(end) = [];
end
end
end
if ~iscell(featureType)
featureType = {featureType};
end
%--------------------------------------------------------------------------
function [displayType, featureType, worldmapArgs] = findDisplayFeatureTypeAtEnd(worldmapArgs)
displayType = '';
featureType = '';
argMatchesDisplayTypeString = false;
argMatchesFeatureTypeString = false;
if ischar(worldmapArgs{end})
displayTypeStrings = getDisplayTypeStrings;
featureTypeStrings = getFeatureTypeStrings;
argMatchesDisplayTypeString = any(strcmpi(worldmapArgs{end},displayTypeStrings));
argMatchesFeatureTypeString = any(strcmpi(worldmapArgs{end},featureTypeStrings));
elseif iscell(worldmapArgs{end})
featureTypeStrings = getFeatureTypeStrings;
type = worldmapArgs{end};
matchesFeatureTypeString(1:numel(type)) = false;
for k=1:numel(type)
matchesFeatureTypeString(k) = any(strcmpi(type{k},featureTypeStrings));
end
argMatchesFeatureTypeString = all(matchesFeatureTypeString);
end
if argMatchesDisplayTypeString
displayType = lower(worldmapArgs{end});
worldmapArgs(end) = [];
elseif argMatchesFeatureTypeString
featureType = lower(worldmapArgs{end});
worldmapArgs(end) = [];
end
%--------------------------------------------------------------------------
function [region,regionOnly] = removeOnlyFromRegion(region)
% detect the 'only case',
regionOnly = false;
%only case with cell array of strings
for i=1:length(region)
thisregion = region{i};
if length(thisregion) >= 5 && strcmp(thisregion(end-4:end),' only')
thisregion = thisregion(1:end-5);
regionOnly = true;
elseif length(thisregion) >= 4 && strcmp(thisregion(end-3:end),'only')
thisregion = thisregion(1:end-4);
regionOnly = true;
end
region{i} = thisregion;
end
%--------------------------------------------------------------------------
function s = loadshape(shapefile, attribName, varargin)
% LOADSHAPE(SHAPEFILE, ATTRIBNAME)
% LOADSHAPE(SHAPEFILE, ATTRIBNAME, BBOX)
args = {shapefile, 'UseGeoCoords', true};
switch numel(varargin)
case 0
case 1
bbox = varargin{1};
args = {args{:}, 'BoundingBox', bbox};
otherwise
bbox = varargin{1};
args = {args{:}, 'BoundingBox', bbox};
end
if ~isempty(attribName)
args = {args{:}, 'Attributes', {attribName}};
end
s = shaperead(args{:});
%--------------------------------------------------------------------------
function addGridData(Z, refvec, displayType)
switch displayType
case {'mesh','dem','meshonly','demonly'}
meshm(Z,refvec)
case {'dem3d','mesh3d','ldem3d','lmesh3d',...
'dem3donly','mesh3donly','ldem3donly','lmesh3donly'}
meshm(Z,refvec,size(Z),Z)
otherwise
end
% Control aspect ratio for 3d plots
switch displayType
case {'dem3d','ldem3d'}
da = daspect;
pba = pbaspect;
da(3) = 7.5*pba(3)/da(3);
daspect(da);
case {'mesh3d','lmesh3d'}
da = daspect;
pba = pbaspect;
da(3) = 2*pba(3)/da(3);
daspect(da);
otherwise
end
% Add digital elevation colormap
switch displayType
case {'dem','dem3d','ldem3d','demonly','dem3donly','ldem3donly'}
demcmap(Z)
otherwise
end
% Add lighting effects
switch displayType
case {'lmesh3d','ldem3d'}
camlight(90,5);
camlight(0,5);
lighting phong
material([0.25 0.8 0])
otherwise
end
%--------------------------------------------------------------------------
function addWorldVectorData(api, region, regionOnly, latlim, lonlim, ...
displayType, featureType, worldShapeFile, attribName)
if regionOnly
% REGIONONLY is requested with 'only'
s = readRegionOnly(worldShapeFile, attribName);
if isempty(s)
return
end
if isempty(regexp(featureType{1}, 'only', 'once'))
switch displayType
case {'line','lineonly'}
drawRegionOnlyAsLines(s, region, attribName)
case {'patch','patchonly'}
drawRegionOnlyAsPatches(s, region, attribName);
end
end
if ~isequal('none',featureType{1}) && isempty(regexp(displayType, 'only', 'once'))
for j=1:numel(region)
indx = find(strcmpi(region{j},{s.(attribName)}));
if ~isempty(indx) && numel(indx) == 1
bbox = s(indx).BoundingBox;
else
bbox = [ min(lonlim), min(latlim); max(lonlim), max(latlim)];
end
drawExternalData(api.ExternalData, featureType, attribName, bbox);
end
end
else
% Full region is specified
[s, bbox] = readRegion(worldShapeFile, latlim, lonlim, attribName);
if isempty(s)
return
end
if isempty(regexp(featureType{1}, 'only', 'once'))
switch displayType
case {'lineonly', 'line', 'mesh', 'dem', 'ldem'}
drawRegionAsLines(s, attribName)
case {'patch','patchonly'}
drawRegionAsPatches(s, attribName)
end
end
if ~isequal('none',featureType{1}) && ...
isempty(regexp(displayType, 'only', 'once'))
% Add the external data
drawExternalData(api.ExternalData, featureType, attribName, bbox);
end
end
%--------------------------------------------------------------------------
function s = readRegionOnly(filename, attribName)
s = loadshape(filename, attribName);
if ~isempty(s) && numel(s(end).Lat) < 3
% Remove degenerative point
s(end) = [];
end
%--------------------------------------------------------------------------
function [s, bbox] = readRegion(filename, latlim, lonlim, attribName)
bbox = [ min(lonlim), min(latlim); max(lonlim), max(latlim)];
if (bbox == [lonlim', latlim'])
s = loadshape(filename,attribName,bbox);
else
% May have longitude wrapping, read the whole file
s = loadshape(filename, attribName);
end
if ~isempty(s) && numel(s(end).Lat) < 3
% Remove degenerative point
s(end) = [];
end
%--------------------------------------------------------------------------
function drawRegionOnlyAsLines(s, region, attribName)
for j=1:numel(region)
indx = find(strcmpi(region{j},{s.(attribName)}));
if ~isempty(indx)
lat = [s(indx).Lat];
lon = [s(indx).Lon];
tag = s(indx(1)).(attribName);
else
lat = [s.Lat];
lon = [s.Lon];
tag = region{j};
end
plotm(lat, lon, 'k', 'Tag', tag);
end
%--------------------------------------------------------------------------
function drawRegionAsLines(s, attribName)
for i=1:numel(s)
plotm(s(i).Lat,s(i).Lon,'k','Tag',s(i).(attribName))
end
%--------------------------------------------------------------------------
function drawRegionOnlyAsPatches(s, region, attribName)
cmap = polcmap(numel(region));
for k=1:numel(region)
indx = find(strcmpi(region{k},{s.(attribName)}));
if ~isempty(indx)
lat = [s(indx).Lat];
lon = [s(indx).Lon];
tag = s(indx(1)).(attribName);
h=patchm(lat, lon, cmap(k,:));
set(h, 'Tag', tag);
else
for j=1:numel(s)
h = patchesm(s(j).Lat, s(j).Lon, 0.0, ...
'Cdata', k, 'FaceColor', 'flat');
set(h,'Tag',region{k});
end
end
end
%--------------------------------------------------------------------------
function drawRegionAsPatches(s, attribName)
alt = 0;
cmap = polcmap(numel(s));
for j=1:numel(s)
h = patchesm(s(j).Lat, s(j).Lon, alt, ...
'Cdata', j, 'FaceColor', 'flat');
set(h, 'Tag', s(j).(attribName));
end
set(gcf,'colormap',cmap);
%--------------------------------------------------------------------------
function drawExternalData(ExternalData, featureType, attribName, bbox)
for j=1:numel(ExternalData)
for k = 1:numel(featureType)
if any(strcmpi(featureType{k}, ExternalData(j).ValidTypes))
if isequal(ExternalData(j).Filename, 'DEFAULTFILE')
filename = ExternalData(1).Filename;
else
filename = ExternalData(j).Filename;
end
ExternalData(j).DrawFcn(filename, attribName, bbox);
end
end
end
%--------------------------------------------------------------------------
function drawLabels(filename, attribName, bbox)
info = shapeinfo(filename);
if any(strcmpi('LabelLat', {info.Attributes.Name})) && ...
any(strcmpi('LabelLon', {info.Attributes.Name}))
m = gcm;
s = shaperead(filename, 'UseGeoCoords', true, ...
'Attributes', {'LabelLat','LabelLon', attribName});
xv = [m.maplonlimit(1), m.maplonlimit(1), m.maplonlimit(2), m.maplonlimit(2)];
yv = [m.maplatlimit(1), m.maplatlimit(2), m.maplatlimit(2), m.maplatlimit(1)];
for k=1:numel(s)
if inpolygon(s(k).LabelLon, s(k).LabelLat, xv, yv)
h=textm(s(k).LabelLat, s(k).LabelLon, s(k).(attribName), ...
'fontsize',9, 'fontweight','bold','Tag', s(k).(attribName));
end
end
end
%--------------------------------------------------------------------------
function drawRivers(shapefile, attribName, bbox)
rivers = loadshape(shapefile, 'Name', bbox);
for i=1:length(rivers)
linem(rivers(i).Lat, rivers(i).Lon, 'b', 'Tag', rivers(i).Name);
end
%--------------------------------------------------------------------------
function drawLakes(shapefile, attribName, bbox)
lakes = loadshape(shapefile, 'Name', bbox);
for i=1:length(lakes)
h=patchm(lakes(i).Lat, lakes(i).Lon, 'b');
set(h,'Tag', lakes(i).Name);
end
%--------------------------------------------------------------------------
function drawCities(shapefile, attribName, bbox)
cities = loadshape(shapefile, 'Name', bbox);
for i=1:length(cities)
linem(cities(i).Lat,cities(i).Lon,0, '.', ...
'color' , [0.9000 0.6000 0.1000], ...
'markersize', 10, 'Tag', cities(i).Name)
textm(cities(i).Lat, cities(i).Lon, [' ' cities(i).Name], ...
'fontsize',9, 'Tag', cities(i).Name)
end
