Geomat are matrices of geographically referenced data. These are essentially georeferenced rectangular, regular grids of points. Data can be numeric (reals), integer, or logical (booleans). Objects 'geotm' are special 'geomat' matrices containing always integers and representing terrain models. Objects 'geomask' are also special 'geomat' that only contain logical values. They are mainly used to define a mask on top of a grid (which points to consider and which ones to eliminate from a calculation).
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43  geomat(x, size, xcenter, ycenter, coords = c(size = size, x = xcenter,
y = ycenter), datatype = c("numeric", "integer", "logical"), nodata = NA)
geotm(x, size, xcenter, ycenter, coords = c(size = size, x = xcenter,
y = ycenter))
geomask(x, size, xcenter, ycenter, coords = c(size = size, x = xcenter,
y = ycenter))
read.geomat(file, type = "ascii", datatype = c("numeric", "integer", "logical"),
...)
read.geotm(file, type = "ascii", ...)
read.geomask(file, type = "ascii", threshold = 0, ...)
write.geomat(x, file, type = "ascii", integers = FALSE, nodata = 9999, ...)
write.geotm(x, file, type = "ascii", nodata = 9999, ...)
write.geomask(x, file, type = "ascii", nodata = 9999, ...)
as.geomat(x, ...)
## S3 method for class 'geomat'
print(x, ...)
## S3 method for class 'geomat'
coords(x, type = "par", ...)
## S3 method for class 'geomat'
resample(x, x0 = 1, y0 = 1, step = NULL, nx = 100, ny = nx,
strict = FALSE, ...)
## S3 method for class 'geomat'
window(x, xlim, ylim, ...)
## S3 method for class 'geomat'
plot(x, y = NULL, max.xgrid = 100, nlevels = 50,
color.palette = terrain.colors, xlab = "Longitude", ylab = "Latitude",
asp = 1, ...)
## S3 method for class 'geomat'
image(x, max.xgrid = 500, col = terrain.colors(50),
add = FALSE, xlab = if (add) "" else "Longitude",
ylab = if (add) "" else "Latitude", asp = 1, ...)
## S3 method for class 'geomat'
contour(x, max.xgrid = 100, nlevels = 10, col = par("fg"),
add = FALSE, xlab = if (add) "" else "Longitude",
ylab = if (add) "" else "Latitude", asp = 1, ...)
## S3 method for class 'geomat'
persp(x, max.xgrid = 500, col = "green3",
xlab = "Longitude", ylab = "Latitude", asp = 1, theta = 10, phi = 30,
expand = 1, shade = 0.75, border = NA, box = TRUE, ...)

x 
An object (a matrix or data frame for 
size 
The size of a grid square (in decimal degrees) 
xcenter 
The position of the center of the topleft square of the grid, that is, its longitude in decimal degrees 
ycenter 
Idem, but latitude in decimal degrees 
coords 
A named vector of three numbers: 'size', 'x' and 'y' as above 
datatype 
The type of data to store in the grid, ort to read/write on the file. Can be 'numeric' (reals), 'integer', or 'logical' (booleans) 
nodata 
The number to use to represent missing data in the grid (by
default it is 
file 
The path to the file used for reading or writing data 
type 
The type of data to read/write. Currently, only \"ascii\", which
means ARC/INFO ASCII GRID format (.asc file). For 
threshold 
Value (single integer) above which all data are converted to

integers 
Do we read/write integers (saves memory and disk space used to represent the grid) 
x0 
The X origin of the new grid 
y0 
The Y origin of the new grid 
step 
The step to use for resampling ( 
nx 
The desired number of points in the X direction (longitude).

ny 
idem than nx, but in the Y direction (latitude) 
strict 
do we interpolated the grid so that we obtain exactly 
xlim 
A vector of two numbers defining the limits to use in X direction (longitude) for the window 
ylim 
A vector of two numbers defining the limits to use in Y direction (latitude) for the window 
y 
Unused argument to match 
max.xgrid 
The maximum number of points in x direction to use. If the grid that is plotted is denser, it is furst resampled to avoid drawing a graph with too much points 
nlevels 
the number of contour levles to calculate 
color.palette 
a color palette generation function 
col 
A vector of colors to use for the plot 
xlab 
The label of the X axis ( 
ylab 
The label of the Y axis ( 
asp 
The aspect ratio between 'x' and 'y'. The default value of

add 
Do we add the graph to an existing graph device, or do we plot a fresh new graph? 
theta 
angles defining the viewing direction. 
phi 

expand 
the expansion level to use for the zaxis in the perspective 
shade 
the shade at a surface facet is computed as 
border 
the color of the borders of facets. If 
box 
If 
... 
Further arguments passed to the functions (only used for the plotting method) 
An object of class, respectively 'geomat', 'geotm' or 'geomask' inheriting from 'matrix' is created. Methods either return an object of same class, or are used for their side effect of plotting a graph. Objects 'geotm' and 'geomask' also inherit from 'geomat'.
A 'geomat' object. For the print()
method, size of the grid is presented
in km.
Philippe Grosjean <phgrosjean@sciviews.org>
aurelhy
, auremask
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48  # Create a simple geomat object containing random numbers
(gm < geomat(matrix(rnorm(120), nrow = 10), 0.1, 10, 20))
# Get coordinates for this grid
coords(gm)
# Longitudes (x) and latitudes (y) for the center of all squares
coords(gm, type = "x")
coords(gm, type = "y")
# Coordinates of the center of all squares
coords(gm, type = "xy")
# Resample the grid to take one point every second points in the original grid
resample(gm, step = 2)
# Extract a window from the grid (keep only squares with centers in the window)
window(gm, xlim = c(9.5, 10.2), ylim = c(19.5, 20.6))
# Plot this grid in different ways
plot(gm)
image(gm)
contour(gm)
persp(gm, expand = 100)
# Now load real data (Morocco terrain model)
data(morocco)
morocco
image(morocco)
contour(morocco, add = TRUE)
grid()
# The mask of points inside Morocco territory was obtained like that:
#library(splancs)
#data(mbord)
#inm < inout(coords(morocco, "xy"), mbord[[1]])
#mmask < morocco
#mmask[inm] < 1
#mmask[!inm] < 0
#mmask[is.na(morocco)] < NA
#mmask < geomask(mmask, coords = coords(mmask))
data(mmask)
image(mmask)
# Get Morocco frontiers from a shapefile
# To read it from an ESRI shape
#mbord < read.geoshapes("morocco_border.shp")
data(mbord)
lines(mbord, col = "red")

Questions? Problems? Suggestions? Tweet to @rdrrHQ or email at ian@mutexlabs.com.
All documentation is copyright its authors; we didn't write any of that.