Newsletters - MATLAB Digest

Maps Made Easy

by Walter Stumpf
Send email to Walter Stumpf c/o Russ Minkwitz

MATLAB's extensive collection of numerical methods and programming capabilities have led to applications that range from aerospace engineering to zoology. Many of these applications produce data that have a geographic context. The Mapping Toolbox makes it easy to create maps to display and analyze your data. This article will show you some new functions that make creating these maps even easier.

Before widespread use of computers, we probably would have taken a paper map and hand drawn the data on top. Copies of this base map might have been on hand, and each new data set would have been redrawn. We would have given little thought to the base map, though it would have been smart to work with a cartographer to ensure the map's suitability.

For the cartographer who created the map, there would have been many things to think about. Should the projection preserve area, shapes or distances? Map projections inevitably introduce some kind of distortion as the round surface of a sphere is flattened to a flat map. What class of projection is appropriate to the area shown, and which of an infinite number of projections should be used? Where should the lines of zero distortion be placed? How detailed should the coastline data be? From what source do you get the coastline and country boundary data for the base map? The Mapping Toolbox's recently introduced worldmap function takes the place of the cartographer in making these decisions. To illustrate the use of this tool, we will use an example in which we enter data about earthquake activity in the Korean peninsula. The first step is creating a base map. You can do this by specifying a part of the world either by the name of a region or country, or by latitude and longitude limits. In addition, by calling the worldmap function with no input arguments you can choose the region from a list. worldmap

Here is the result when you select Asia: Click to see enlarged view (31 k)
Let's refine the base map further. We can enter the name of the region as a string and the worldmap function will match all country names that begin with the provided string. close worldmap korea	Click to see enlarged view (28 k)

Notice how the appearance of the map changed as we asked for a smaller region. Worldmap picked a different projection, switched to higher resolution atlas data, and added cities and country labels as the risk of clutter was reduced.

We can also override some of the default choices. For example, we can specifiy all of the high-resolution atlas data, even the smallest islands, for countries with names beginning with "Korea", are displayed as lines only.

close all
worldmap('allhi','korea','lineonly')

Now we can easily add our earthquake data to this base map. The Mapping Toolbox always adds displayed data to an existing map, rather than replacing it as MATLAB does. Mapping Toolbox display commands are just like those in MATLAB, except that they use geographic spherical coordinates instead of Cartesian. Here we load our earthquake data and display it on the base map as markers with the area proportional to the square of the magnitude. We also change the marker face color to white to clarify close events, and move the markers to below the atlas data.

load koreaEQdata
whos

Name	Size	Bytes	Class
CAPTION	4x68	544	char array
kdata	803x3	19272	double array

Grand total is 2681 elements using 19816 bytes

h = scatterm(kdata(:,1),kdata(:,2),kdata(:,3).^2,'r','filled');
set(h,'markerfacecolor','w')
zdatam(h,-10)

Click to see enlarged view (30 k)

Another way to specify the region is by latitude and longitude limits. If you don't know the limits of a small region in advance, it's convenient to make a larger map, and then define the limits by clicking on the map. For example, to make a base map of the southern Japanese island of Kyushu, you can use the map of the Korean Peninsula and then use inputm to define the limits:

[latlim,lonlim]=inputm(2)

latlim =
34.4886
30.7035


lonlim =
129.0171
133.5265

figure
worldmap(latlim,lonlim,'lineonly')
h = scatterm(kdata(:,1),kdata(:,2),kdata(:,3).^2,'r','filled');
set(h,'markerfacecolor','w')
zdatam(h,-10)

Click to see enlarged view (41 k)
The matrix calling form for worldmap can do a lot of sophisticated graphic tweaking for you. The optional argument can request that matrix data be displayed as a lighted three-dimensional digital elevation map. In that case, worldmap takes control of a disparate set of Handle Graphics properties, which would require many lines of code and considerable experience to control. Here is an example of one-kilometer resolution data for Kyushu. This data was imported from the GTOPO30 dataset, high-resolution elevation data available over the internet. We used the Mapping Toolbox's demdataui GUI which is a user interface for extracting digital elevation map data from a number of external data files. worldmap(s.map,s.maplegend,'ldem3d') hidem(gca)	Click to see enlarged view (36 k)
Click to see enlarged view (23 k)	For base maps of the United States, there is the usamap function, which is similar to worldmap. It can display Alaska and Hawaii as inset maps at differing or equal scales. Each inset is in a different axes, similar to the MATLAB subplot command. Use the axes command to specify which axes will receive the results of the next plotting command. close all usamap allequal

If these functions do not produce the base map you want, you can use the lower level map projection and atlas commands to retain complete control over the base map.

axesm('mapprojection','ortho','origin',[40 120],...
'mlinelocation',15,'plinelocation',15,...
'mlinelimit',[-75 75],'mlineexception',-180:90:90,...
'grid','on','frame','on')
displaym(worldlo('POpatch'))
polcmap
hidem(gca)
set(gcf,'color','w')
tightmap loose

Click to see enlarged view (29 k)

These simple but powerful base map commands, combined with other Mapping Toolbox display commands make it easy to create geographic displays with your data. They reduce the amount of time spent preparing base maps, and let you spend more time thinking about what the results really mean.