pkern_snap: Snap a set of points to a grid

In deankoch/pkern: Fast Kriging with Product Kernels

Snap a set of points to a grid

Description

Maps the input points in from to the closest grid points in the extension of g covering the bounding box of from (ie. the lattice of which g is a sub-grid). In cases of duplicate mappings, the function returns the first matches only.

Usage

pkern_snap(from, g = NULL, crop_from = FALSE, crop_g = FALSE, quiet = FALSE)

Arguments

from	matrix, data frame, or points object from sp or sf, the source points
g	any grid object accepted or returned by pkern_grid, the destination grid
crop_from	logical, indicating to omit points not overlying g.
crop_g	logical, indicating to trim g to the extent of from.

Details

from can be a geometry collection from packages sf or sp, or a matrix or list of y and x coordinates. When from is a matrix, its first two columns should be the y and x coordinates (in that order), and the (optional) third column should be the data. When from is a list, the function expects (two or three) vectors of equal length, ordered as above.

When from is a geometry collection with a coordinates reference system (CRS) string, points are first transformed to the CRS of g. If one or both of from and g are missing a CRS definition, the function assumes the same one is shared in both.

g can be a raster geometry object (such as SpatRaster), in which case the function behaves like terra::rasterize. It can also be a matrix (supplying dimensions) or a list containing either gdim orgres, from which an appropriately spaced set of grid lines is derived, centered under the bounding box of the points.

crop_from and crop_g control the extent of the output grid. If both are FALSE (the default) the function returns the smallest regular grid containing both g and the snapped from points. If crop_from=TRUE and crop_g=FALSE the output grid will match g exactly. If crop_from=FALSE and crop_g=TRUE the output grid will include all snapped points, and possibly omit some or all of g. And if both are TRUE, the output grid encloses the intersection of the points with the bounding box of g.

Value

list of form returned by pkern_grid, defining a grid containing the snapped points. These are assigned the corresponding data value in from, or if from has no data, an integer mapping to the points in from. Un-mapped grid points are set to NA.

Examples

# functions to scale arbitrary inverval to (1, 2,... 100) and make color palettes
num_to_cent = function(x) 1L + floor(99*( x-min(x) ) / diff(range(x)))
my_pal = function(x) hcl.colors(x, 'Spectral', rev=T)
my_col = function(x) my_pal(1e2)[ num_to_cent(x) ]

# create a grid object
gdim = c(40, 30)
g = pkern_grid(list(gdim=gdim, gres=1.1))

# randomly position points within bounding box of g
n_pts = 10
from = lapply(g$gyx, function(yx) runif(n_pts, min(yx), max(yx)) )

# translate away from g (no overlap is required)
from[['y']] = from[['y']] + 5
from[['x']] = from[['x']] + 15

# add example data values and plot
from[['z']] = rnorm(length(from[['y']]))
pkern_plot(g, reset=FALSE)
graphics::points(from[c('x', 'y')], pch=16, col=my_col(from[['z']]))
graphics::points(from[c('x', 'y')])

# snap only the points overlying the input grid
g_snap = pkern_snap(from, g, crop_from=TRUE)
pkern_plot(g_snap, col_grid='black', reset=FALSE, leg=FALSE)
graphics::points(from[c('x', 'y')], pch=16, col=my_col(from[['z']]))
graphics::points(from[c('x', 'y')])

# snap points with  and plot (default settings)
g_snap = pkern_snap(from, g, crop_from=FALSE, crop_g=FALSE)
pkern_plot(g_snap, col_grid='black', reset=FALSE)
graphics::points(from[c('x', 'y')], pch=16, col=my_col(from[['z']]))
graphics::points(from[c('x', 'y')])

# find smallest subgrid enclosing all snapped grid points
g_snap = pkern_snap(from, g, crop_g=TRUE)
pkern_plot(g_snap, col_grid='black', reset=FALSE)
graphics::points(from[c('x', 'y')], pch=16, col=my_col(from[['z']]))
graphics::points(from[c('x', 'y')])

# create a new grid of different resolution enclosing all input points
g_snap = pkern_snap(from, g=list(gres=c(0.5, 0.5)))
pkern_plot(g_snap, reset=FALSE, col_grid='black')
graphics::points(from[c('x', 'y')], pch=16, col=my_col(from[['z']]))
graphics::points(from[c('x', 'y')])

if( requireNamespace('sf') ) {

# a different example, snapping mis-aligned subgrid
g_pts = pkern_grid(list(gdim=c(15, 8), gres=1.7), vals=FALSE)
g_pts[['gyx']][['y']] = g_pts[['gyx']][['y']] + 5
g_pts[['gyx']][['x']] = g_pts[['gyx']][['x']] + 5
n_pts = prod(g_pts$gdim)
from = pkern_coords(g_pts, out='list')

# convert to sf
eg_sfc = sf::st_geometry(pkern_coords(g_pts, out='sf'))
pkern_plot(g, reset=FALSE)
plot(eg_sfc, add=TRUE)

# generate example data and plot
eg_sf = sf::st_sf(data.frame(z=rnorm(n_pts)), geometry=eg_sfc)
pkern_plot(g, reset=FALSE)
plot(eg_sf, pch=16, add=TRUE, pal=my_pal)
plot(eg_sfc, add=TRUE)

# snap points
g_snap = pkern_snap(from=eg_sf, g)
pkern_plot(g_snap, reset=FALSE, col_grid='black')
plot(eg_sf, pch=16, add=TRUE, pal=my_pal)
plot(eg_sfc, add=TRUE)

# snapping points without data produces the mapping (non-NA values index "from")
g_snap = pkern_snap(from=eg_sfc, g)
pkern_plot(g_snap, ij=TRUE, reset=FALSE, col_grid='black')
plot(eg_sfc, add=TRUE)

# with crop_g=TRUE)
g_snap = pkern_snap(from=eg_sfc, g, crop_g=TRUE)
pkern_plot(g_snap, reset=FALSE, col_grid='black')
plot(eg_sfc, add=TRUE)

# test with sp class
eg_sp = as(eg_sf,'Spatial')
g_snap = pkern_snap(from=eg_sp, g)
pkern_plot(g_snap, reset=FALSE, col_grid='black')
plot(eg_sf, pch=16, add=TRUE, pal=my_pal)
plot(eg_sfc, add=TRUE)

}

deankoch/pkern documentation built on Oct. 26, 2023, 8:54 p.m.