oce (version 0.9-18)

mapPlot: Plot a map

Description

Plot coordinates as a map, using one of the subset of projections provided by the rgdal package. The projection information specified with the mapPlot call is stored so that can be retrieved by related functions, making it easy to add more items so the map, including points, lines, text, images and contours.

Usage

mapPlot(longitude, latitude, longitudelim, latitudelim,
        grid=TRUE, bg, fill, border=NULL, col=NA, type='l', axes=TRUE, drawBox=TRUE,
        showHemi=TRUE, polarCircle=0,
        lonlabel=NULL, latlabel=NULL, sides=NULL,
        projection="+proj=moll", tissot=FALSE, trim=TRUE,
        debug=getOption("oceDebug"), ...)

Arguments

longitude

either a vector of longitudes of points to be plotted, or something (an oce object, a list, or a data frame) from which both longitude and latitude may be inferred (in which case the latitude argument is ignored). If

latitude

latitudes of points to be plotted (ignored if the first argument contains both latitude and longitude).

longitudelim

optional limit of longitudes to plot

latitudelim

optional limit of latitudes to plot

grid

either a number (or pair of numbers) indicating the spacing of longitude and latitude lines, in degrees, or a logical value (or pair of values) indicating whether to draw an auto-scaled grid, or whether to skip the grid drawing. In the case o

bg

background colour for plot (ignored at present).

fill

included only for compatibility with the previous version of mapPlot; will be removed in next version.

border

colour of coastlines and international borders. The default, NULL, means to use par("fg"); see polygon.

col

colour with which to fill coastline elements. The default, NA, is not to fill; see polygon.

type

value to indicate type of plot, as with par("plot").

axes

logical value indicating whether to draw longitude and latitude values in the lower and left margin, respectively. This may not work well for some projections or scales.

drawBox

logical value indicating whether to draw a box around the plot. This is helpful for many projections at sub-global scale.

showHemi

logical value indicating whether to show the hemisphere in axis tick labels.

polarCircle

a number indicating the number of degrees of latitude extending from the poles, within which zones are not drawn.

lonlabel,latlabel,sides

optional vectors of longitude and latitude to label on the indicated sides of plot, passed to plot.coastline. Using these arguments permits reasonably simple customization. If they are

projection

optional indication of projection, a character string that in the format used by the rgdal package (and in much of computer cartography). For example, projection="+proj=merc" specifies a Mercator projection. See

trim

a logical value indicating whether to trim islands or lakes containing only points that are off-scale of the current plot box. This solves the problem of Antarctica overfilling the entire domain, for an Arctic-centred stereographic projection

tissot

logical indicating whether to use mapTissot to plot Tissot indicatrices, i.e. ellipses at grid intersection points, which indicate map distortion.

debug

a flag that turns on debugging. Set to 1 to get a moderate amount of debugging information, or to 2 to get more.

...

optional arguments passed to some plotting functions. This can be useful in many ways, e.g. Example 5 shows how to use xlim etc to reproduce a scale exactly between two plots.

Available Projections

Map projections are provided by the rgdal package, but not all projections in that package are available. The available list is given in the table below. The cartographic community has set up a naming scheme in a coded scheme, e.g. projection="+proj=aea" selects the Albers equal area projection.

The allowed projections include those PROJ.4 projections provided by rgdal that have inverses, minus a few that cause problems: alsk overdraws coastlineWorld, and is a niche projection for Alaska; calcofi is not a real projection, but rather a coordinate system; gs48 overdraws coastlineWorld, and is a niche projection for the USA; gs50 overdraws coastlineWorld, and is a niche projection for the USA; gstmerc overdraws coastlineWorld; isea causes segmentation faults on OSX systems; krovak overdraws coastlineWorld, and is a niche projection for the Czech Republic; labrd returns NaN for most of the world, and is a niche projection for Madagascar; lee_os overdraws coastlineWorld; and nzmg overdraws coastlineWorld. The information in the table is reformatted from the output of the unix command proj -lP, where proj is provided by version 4.9.0 of the PROJ.4 system. Most of the arguments listed have default values. In addition, most projections can handle arguments lon_0 and lat_0, for shifting the reference point, although in some cases shifting the longitude can yield poor filling of coastlines.

Further details of the projections and the controlling arguments are provided at several websites, because PROJ.4 has been incorporated into rgdal and other R packages, plus many other software systems; a good starting point for learning is [6]. See Examples for suggested projections for some common applications, and [8] for a gallery indicating how to use every projection.

lll{ Projection Code Arguments Albers equal area aea lat_1, lat_2 Azimuthal equidistant aeqd lat_0, guam Aitoff aitoff - Mod. stererographics of Alaska alsk - Bipolar conic of western hemisphere bipc - Bonne Werner bonne lat_1 Cassini cass - Central cylindrical cc - Equal area cylindrical cea lat_ts Collignon collg - Craster parabolic Putnins P4 crast - Eckert I eck1 - Eckert II eck2 - Eckert III eck3 - Eckert IV eck4 - Eckert V eck5 - Eckert VI eck6 - Equidistant cylindrical plate (Caree) eqc lat_ts, lat_0 Equidistant conic eqdc lat_1, lat_2 Euler euler lat_1, lat_2 Extended transverse Mercator etmerc lat_ts, lat_0 Fahey fahey - Foucaut fouc - Foucaut sinusoidal fouc_s - Gall stereographic gall - Geostationary satellite view geos h General sinusoidal series gn_sinu m, n Gnomonic gnom - Goode homolosine goode - Hatano asymmetrical equal area hatano - HEALPix healpix - rHEALPix rhealpix north_square, south_square Interrupted Goode homolosine igh - Int'l map of the world polyconic imw_p lat_1, lat_2, lon_1 Kavraisky V kav5 - Kavraisky VII kav7 - Lambert azimuthal equal area laea - Lat/long lonlat - Lat/long latlon - Lambert conformal conic lcc lat_1, lat_2, lat_0 Lambert conformal conic alternative lcca lat_0 Lambert equal area conic leac lat_1, south Loximuthal loxim Space oblique for Landsat lsat lsat, path McBryde-Thomas flat-polar sine, no. 1 mbt_s McBryde-Thomas flat-polar sine, no. 2 mbt_fps McBryde-Thomas flat-polar parabolic mbtfpp McBryde-Thomas flat-polar quartic mbtfpq McBryde-Thomas flat-polar sinusoidal mbtfps Mercator merc lat_ts Miller oblated stereographic mil_os Miller cylindrical mill Mollweide moll Murdoch I murd1 lat_1, lat_2 Murdoch II murd2 lat_1, lat_2 murdoch III murd3 lat_1, lat_2 Natural earth natearth Nell nell Nell-Hammer nell_h Near-sided perspective nsper h New Zealand map grid nzmg General oblique transformation ob_tran o_proj, o_lat_p, o_lon_p, o_alpha, o_lon_c o_lat_c, o_lon_1, o_lat_1, o_lon_2, o_lat_2 Oblique cylindrical equal area ocea lat_1, lat_2, lon_1, lon_2 Oblated equal area oea n, m, theta Oblique Mercator omerc alpha, gamma, no_off, lonc, lon_1, lat_1, lon_2, lat_2 Orthographic ortho - Perspective conic pconic lat_1, lat_2 Polyconic American poly - Putnins P1 putp1 - Putnins P2 putp2 - Putnins P3 putp3 - Putnins P3' putp3p - Putnins P4' putp4p - Putnins P5 putp5 - Putnins P5' putp5p - Putnins P6 putp6 - Putnins P6' putp6p - Quartic authalic qua_aut - Quadrilateralized spherical cube qsc - Robinson robin - Roussilhe stereographic rouss - Sinusoidal aka Sanson-Flamsteed sinu - Swiss. oblique Mercator somerc - Stereographic stere lat_ts Oblique stereographic alternative sterea - Gauss-schreiber transverse Mercator gstmerc lat_0, lon_0, k_0 Transverse cylindrical equal area tcea - Tissot tissot lat_1, lat_2 Transverse Mercator tmerc - Two point equidistant tpeqd lat_1, lon_1, lat_2, lon_2 Tilted perspective tpers tilt, azi, h Universal polar stereographic ups south Urmaev flat-polar sinusoidal urmfps n Universal transverse Mercator utm zone, south van der Grinten I vandg - Vitkovsky I vitk1 lat_1, lat_2 Wagner I Kavraisky VI wag1 - Wagner II wag2 - Wagner III wag3 lat_ts Wagner IV wag4 - Wagner V wag5 - Wagner VI wag6 - Werenskiold I weren - Winkel I wink1 lat_ts Winkel Tripel wintri lat_ts }

Available ellipse formulations

In the PROJ.4 system of specifying projections, the following ellipse models are available: MERIT, SGS85, GRS80, IAU76, airy, APL4.9, NWL9D, mod_airy, andrae, aust_SA, GRS67, bessel, bess_nam, clrk66, clrk80, clrk80ign, CPM, delmbr, engelis, evrst30, evrst48, evrst56, evrst69, evrstSS, fschr60, fschr60m, fschr68, helmert, hough, intl, krass, kaula, lerch, mprts, new_intl, plessis, SEasia, walbeck, WGS60, WGS66, WGS72, WGS84, and sphere (the default). For example, use projection="+proj=aea +ellps=WGS84" for an Albers Equal Area projection using the most recent of the World Geodetic System model. It is unlikely that changing the ellipse will have a visible effect on plotted material at the plot scale appropriate to most oceanographic applications.

Available datum formulations

In the PROJ.4 system of specifying projections, the following datum formulations are available: WGS84, GGRS87, Greek_Geodetic_Reference_System_1987, NAD83, North_American_Datum_1983, NAD27, North_American_Datum_1927, potsdam, Potsdam, carthage, Carthage, hermannskogel, Hermannskogel, ire65, Ireland, nzgd49, New, OSGB36, and Airy. It is unlikely that changing the datum will have a visible effect on plotted material at the plot scale appropriate to most oceanographic applications.

Choosing a projection

The use of the PROJ.4 scheme is greatly encouraged. The best choice of projection depends on the application. Readers may find projection="+proj=moll" useful for world-wide plots, ortho for hemispheres viewed from the equator, stere for polar views, lcc for wide meridional ranges in mid latitudes, and merc in limited-area cases where angle preservation is important.

Issues

Map projection is a complicated matter that is addressed here in a limited and pragmatic way. For example, mapPlot tries to draw axes along a box containing the map, instead of trying to find spots along the ``edge'' of the map at which to put longitude and latitude labels. This design choice greatly simplifies the coding effort, freeing up time to work on issues regarded as more pressing. Chief among those issues are (a) the occurrence of horizontal lines in maps that have prime meridians (b) inaccurate filling of land regions that (again) occur with shifted meridians and (c) inaccurate filling of Antarctica in some projections. Generally, issues are tackled first for commonly used projections, such as those used in the examples.

Details

Creates a map using the indicated projection. As noted in the information on the projection argument, projections are specified in the notation used by project() in the rgdal package; see Available Projections for a list of possibilities. Once a projection is set, other map* functions may be used to add to the map. Further details on map projections are provided by [1,11], an exhaustive treatment that includes many illustrations, an overview of the history of the topic, and some notes on the strengths and weaknesses of the various formulations. See especially pages 2 through 7, which define terms and provide recommendations. Reference [2] is also useful, especially regarding datum shifts; [3] and [4] are less detailed and perhaps better for novices. See [8] for a gallery of projections.

References

1. Snyder, John P., 1987. Map Projections: A Working Manual. USGS Professional Paper: 1395 (available at pubs.usgs.gov/pp/1395/report.pdf).

2. Natural Resources Canada http://www.nrcan.gc.ca/earth-sciences/geography/topographic-information/maps/9805

3. Wikipedia page http://en.wikipedia.org/wiki/List_of_map_projections

4. Radical Cartography website http://www.radicalcartography.net/?projectionref

5. The PROJ.4 website is http://trac.osgeo.org/proj, and it is the place to start to learn about the code.

6. PROJ.4 projection details are at http://www.remotesensing.org/geotiff/proj_list/.

7. A gallery of map plots is provided at http://dankelley.github.io/r/2015/04/03/oce-proj.html.

8. A fascinating historical perspective is provided by Snyder, J. P. (1993). Two thousand years of map projections. University of Chicago Press.

See Also

Points may be added to a map with mapPoints, lines with mapLines, text with mapText, polygons with mapPolygon, images with mapImage, and scale bars with mapScalebar. Points on a map may be determined with mouse clicks using mapLocator. Great circle paths can be calculated with geodGc. See [8] for a demonstration of the available map projections (with graphs).

Examples

library(oce)
data(coastlineWorld)

# Example 1.
# Mollweide ([1] page 54) is an equal-area projection that works well
# for whole-globe views, below shown in a Pacific-focus view.
# Note that filling is not employed when the prime meridian
# is shifted, because this causes a problem with Antarctica
par(mfrow=c(2,1), mar=c(3, 3, 1, 1))
mapPlot(coastlineWorld, projection="+proj=moll", fill='gray')
mtext("Mollweide", adj=1)
cl180 <- coastlineCut(coastlineWorld, lon_0=-180)
mapPlot(cl180, projection="+proj=moll +lon_0=-180")
mtext("Mollweide", adj=1)
par(mfrow=c(1,1))

# Example 2.
# Orthographic projections resemble a globe, making them attractive for
# non-technical use, but they are neither conformal nor equal-area, so they
# are somewhat limited for serious use on large scales.  See Section 20 of
# [1]. Note that filling is not employed because it causes a problem with 
# Antarctica.
par(mar=c(3, 3, 1, 1))
mapPlot(coastlineWorld, projection="+proj=ortho +lon_0=-180")
mtext("Orthographic", adj=1)

# Example 3.
# The Lambert conformal conic projection is an equal-area projection
# recommended by [1], page 95, for regions of large east-west extent
# away from the equator, here illustrated for the USA and Canada.
par(mar=c(3, 3, 1, 1))
mapPlot(coastlineWorld, longitudelim=c(-130,-55), latitudelim=c(35,60),
        projection="+proj=lcc +lat_0=30 +lat_1=60 +lon_0=-100", fill='gray')
mtext("Lambert conformal", adj=1)

# Example 4.
# The stereographic projection [1], page 120, is conformal, used
# below for an Arctic view with a Canadian focus.  Note the trick of going
# past the pole: the second latitudelim value is 180 minus the first, and the
# second longitudelim is 180 plus the first; this uses image points "over"
# the pole.
par(mar=c(3, 3, 1, 1))
mapPlot(coastlineWorld, longitudelim=c(-130,50), latitudelim=c(70,110),
        proj="+proj=stere +lat_0=90 +lon_0=-135", fill='gray')
mtext("Stereographic", adj=1)

# Example 5.
# Spinning globe: create PNG files that can be assembled into a movie
png("539B-%03d.png")
lons <- seq(360, 0, -15)
ilon <- seq_along(lons)
par(mar=rep(0, 4))
for (i in ilon) {
    p <- paste("+proj=ortho +lat_0=30 +lon_0=", lons[i], sep="")
    if (i == 1) {
        mapPlot(coastlineWorld, projection=p, col="blue", lwd=1.4)
        xlim <- par("usr")[1:2]
        ylim <- par("usr")[3:4]
    } else {
        mapPlot(coastlineWorld, projection=p, col="blue", lwd=1.4,
                xlim=xlim, ylim=ylim, xaxs="i", yaxs="i")
    }
}
dev.off()