get.map.lines.M3.proj {M3}

Get map lines in the model projection units

Description

Get map lines in the model projection units.

Usage

get.map.lines.M3.proj(file, database = "state", units, ...)

Arguments

file

File name of Models3-formatted file providing the model projection.

database

Geographical database to use. Choices include “state” (default), “world”, “worldHires”, “canusamex”, etc. Use “canusamex” to get the national boundaries of the Canada, the USA, and Mexico, along with the boundaries of the states. The other choices (“state”, “world”, etc.) are the names of databases included with the maps and mapdata packages. (See the map function.) Default is “state”.

units

Units for coordinates of grid rows or columns. Must be one of “m”, “km”, or “deg”. If unspecified, the default is “deg” if the file has a longitude/latitude grid, and “km” otherwise.

...

Other arguments to pass to get.proj.info.M3 function. In this case, the only relevant argument would be the earth radius to use for the projection in file.

Details

This function depends on the maps and mapdata packages to get the appropriate map boundary lines (for states, countries, etc.), ncdf4 to read the projection information from the Models3-formatted file (using a call to function get.proj.info.M3), and sf to project the boundary lines to the specified projection.

Value

Map lines for the projection described in file in either kilometers or meters (depending on value of units.km). This is a matrix, with x-coordinates in the first column and y-coordinates in the second column.

Warning

This function will only work with Lambert conic conformal or polar stereographic projections.

Author(s)

Jenise Swall

See Also

get.proj.info.M3, map, sf_project

Examples

## Find the path to the demo file.
polar.file <- system.file("extdata/surfinfo_polar.ncf", package="M3")

## Read in the terrain elevation variable.
elev <- get.M3.var(file=polar.file, var="HT")
## Make a plot.
image(elev$x.cell.ctr, elev$y.cell.ctr, elev$data[,,1],
      xlab="Projection x-coord (km)", ylab="Projection y-coord (km)",
      zlim=range(elev$data[,,1]), col=heat.colors(15))

## Superimpose national boundaries on the plot
nat.bds <- get.map.lines.M3.proj(file=polar.file, database="world")$coords
lines(nat.bds)


## Subset to a smaller geographic area in southwestern U.S.
subset.elev <- var.subset(elev, llx=-2000, urx=0, lly=-6500, ury=-4000)
## Make a plot of this subset.
image(subset.elev$x.cell.ctr, subset.elev$y.cell.ctr,
      subset.elev$data[,,1], xlab="Projection x-coord (km)",
      ylab="Projection y-coord (km)", zlim=range(subset.elev$data[,,1]),
      col=heat.colors(15))

## Superimpose Mexico, US, and Candadian national borders on the plot,
## along with state borders.
canusamex.borders <- get.map.lines.M3.proj(file=polar.file, "canusamex")$coords
lines(canusamex.borders)