Title: | Multilevel B-Spline Approximation |
---|---|
Description: | Functions to interpolate irregularly and regularly spaced data using Multilevel B-spline Approximation (MBA). Functions call portions of the SINTEF Multilevel B-spline Library written by Øyvind Hjelle which implements methods developed by Lee, Wolberg and Shin (1997; <doi:10.1109/2945.620490>). |
Authors: | Andrew Finley [aut, cre], Sudipto Banerjee [aut], Øyvind Hjelle [aut], Roger Bivand [ctb] |
Maintainer: | Andrew Finley <[email protected]> |
License: | GPL (>= 2) |
Version: | 0.1-2 |
Built: | 2024-10-24 06:32:02 UTC |
Source: | CRAN |
This is a small portion of Light Detection and Ranging (LIDAR) data taken over a forested landscape in Wisconsin, USA.
data(LIDAR)
data(LIDAR)
A data frame containing 10123 rows and 3 columns corresponding to longitude, latitude, and elevation.
Data provided by: Dr. Paul V. Bolstad, Department of Forest Resources, University of Minnesota, [email protected]
The function mba.points
returns points on a surface
approximated from a bivariate scatter of points using multilevel B-splines.
mba.points(xyz, xy.est, n = 1, m = 1, h = 8, extend = TRUE, verbose = TRUE, ...)
mba.points(xyz, xy.est, n = 1, m = 1, h = 8, extend = TRUE, verbose = TRUE, ...)
xyz |
a |
xy.est |
a |
n |
initial size of the spline space in the hierarchical construction along the x axis. If the rectangular domain is a square, n = m = 1 is recommended. If the x axis is k times the length of the y axis, n = 1, m = k is recommended. The default is n = 1. |
m |
initial size of the spline space in the hierarchical construction along the y axis. If the y axis is k times the length of the x axis, m = 1, n = k is recommended. The default is m = 1. |
h |
Number of levels in the hierarchical construction. If, e.g.,
n = m = 1 and h = 8, the resulting spline surface has a coefficient
grid of size |
extend |
if FALSE, points in |
verbose |
if TRUE, warning messages are printed to the screen. |
... |
currently no additional arguments. |
List with 1 component:
xyz.est |
a |
The function mba.points
relies on the Multilevel B-spline
Approximation (MBA) algorithm. The underlying code was developed at
SINTEF Applied Mathematics by Dr. Øyvind Hjelle. Dr. Øyvind Hjelle based the
algorithm on the paper by the originators of Multilevel B-splines:
S. Lee, G. Wolberg, and S. Y. Shin. (1997) Scattered data interpolation with multilevel B-splines. IEEE Transactions on Visualization and Computer Graphics, 3(3):229–244.
For additional documentation and references see:
https://www.sintef.no/upload/IKT/9011/geometri/MBA/mba_doc/index.html.
data(LIDAR) ## Split the LIDAR dataset into training and validation sets. tr <- sample(1:nrow(LIDAR), trunc(0.5*nrow(LIDAR))) ## Look at how smoothing changes z-approximation, ## careful the number of B-spline surface coefficients ## increases at ~2^h in each direction. for(i in 1:10){ mba.pts <- mba.points(LIDAR[tr,], LIDAR[-tr,c("x","y")], h=i)$xyz.est print(sum(abs(LIDAR[-tr,"z"]-mba.pts[,"z"]))/nrow(mba.pts)) } ## Not run: ## rgl or scatterplot3d libraries can be fun. library(rgl) # Exaggerate z a bit for effect and take a smaller subset of LIDAR. LIDAR[,"z"] <- 10*LIDAR[,"z"] tr <- sample(1:nrow(LIDAR), trunc(0.99*nrow(LIDAR))) # Get the "true" surface. mba.int <- mba.surf(LIDAR[tr,], 100, 100, extend=TRUE)$xyz.est # Make nice colors for the rgl surface. zlim <- range(mba.int$z) zlen <- zlim[2] - zlim[1] + 1 colorlut <- terrain.colors(zlen) # Height color lookup table. col <- colorlut[mba.int$z - zlim[1] + 1 ] # Assign colors to heights for each point. open3d() surface3d(mba.int$x, mba.int$y, mba.int$z, color = col) # Now add the point estimates. mba.pts <- mba.points(LIDAR[tr,], LIDAR[-tr,c("x","y")])$xyz.est spheres3d(mba.pts[,"x"], mba.pts[,"y"], mba.pts[,"z"], radius=5, color="red") ## End(Not run)
data(LIDAR) ## Split the LIDAR dataset into training and validation sets. tr <- sample(1:nrow(LIDAR), trunc(0.5*nrow(LIDAR))) ## Look at how smoothing changes z-approximation, ## careful the number of B-spline surface coefficients ## increases at ~2^h in each direction. for(i in 1:10){ mba.pts <- mba.points(LIDAR[tr,], LIDAR[-tr,c("x","y")], h=i)$xyz.est print(sum(abs(LIDAR[-tr,"z"]-mba.pts[,"z"]))/nrow(mba.pts)) } ## Not run: ## rgl or scatterplot3d libraries can be fun. library(rgl) # Exaggerate z a bit for effect and take a smaller subset of LIDAR. LIDAR[,"z"] <- 10*LIDAR[,"z"] tr <- sample(1:nrow(LIDAR), trunc(0.99*nrow(LIDAR))) # Get the "true" surface. mba.int <- mba.surf(LIDAR[tr,], 100, 100, extend=TRUE)$xyz.est # Make nice colors for the rgl surface. zlim <- range(mba.int$z) zlen <- zlim[2] - zlim[1] + 1 colorlut <- terrain.colors(zlen) # Height color lookup table. col <- colorlut[mba.int$z - zlim[1] + 1 ] # Assign colors to heights for each point. open3d() surface3d(mba.int$x, mba.int$y, mba.int$z, color = col) # Now add the point estimates. mba.pts <- mba.points(LIDAR[tr,], LIDAR[-tr,c("x","y")])$xyz.est spheres3d(mba.pts[,"x"], mba.pts[,"y"], mba.pts[,"z"], radius=5, color="red") ## End(Not run)
The function mba.surf
returns a surface approximated from a
bivariate scatter of data points using multilevel B-splines.
mba.surf(xyz, no.X, no.Y, n = 1, m = 1, h = 8, extend=FALSE, sp=FALSE, ...)
mba.surf(xyz, no.X, no.Y, n = 1, m = 1, h = 8, extend=FALSE, sp=FALSE, ...)
xyz |
a |
no.X |
resolution of the approximated surface along the x axis. |
no.Y |
resolution of the approximated surface along the y axis. |
n |
initial size of the spline space in the hierarchical construction along the x axis. If the rectangular domain is a square, n = m = 1 is recommended. If the x axis is k times the length of the y axis, n = 1, m = k is recommended. The default is n = 1. |
m |
initial size of the spline space in the hierarchical construction along the y axis. If the y axis is k times the length of the x axis, m = 1, n = k is recommended. The default is m = 1. |
h |
Number of levels in the hierarchical construction. If, e.g.,
n = m = 1 and h = 8, the resulting spline surface has a coefficient
grid of size |
extend |
if FALSE, a convex hull is computed for the input points
and all matrix elements in z that have centers outside of this
polygon are set to |
sp |
if TRUE, the resulting surface is returned as a
|
... |
|
List with 8 component:
xyz.est |
a list that contains vectors x, y and the |
no.X |
|
no.Y |
|
n |
|
m |
|
h |
|
extend |
|
sp |
|
b.box |
|
If no.X != no.Y
then use sp=TRUE
for compatibility with
the image
function.
The function mba.surf
relies on the Multilevel B-spline
Approximation (MBA) algorithm. The underlying code was developed at
SINTEF Applied Mathematics by Dr. Øyvind Hjelle. Dr. Øyvind Hjelle
based the algorithm on the paper by the originators of Multilevel B-splines:
S. Lee, G. Wolberg, and S. Y. Shin. (1997) Scattered data interpolation with multilevel B-splines. IEEE Transactions on Visualization and Computer Graphics, 3(3):229–244.
For additional documentation and references see:
https://www.sintef.no/upload/IKT/9011/geometri/MBA/mba_doc/index.html.
## Not run: data(LIDAR) mba.int <- mba.surf(LIDAR, 300, 300, extend=TRUE)$xyz.est # Image plot of the surface. image(mba.int, xaxs = "r", yaxs = "r") # Perspective plot of the surface. persp(mba.int, theta = 135, phi = 30, col = "green3", scale = FALSE, ltheta = -120, shade = 0.75, expand = 10, border = NA, box = FALSE) # For a good time, I recommend using rgl. library(rgl) # Exaggerate z a bit for effect. mba.int$z <- 10*mba.int$z # Make nice colors for the rgl surface. zlim <- range(mba.int$z) zlen <- zlim[2] - zlim[1] + 1 colorlut <- terrain.colors(zlen) # Height color lookup table. col <- colorlut[mba.int$z - zlim[1] + 1 ] # Assign colors to heights for each point. open3d() surface3d(mba.int$x, mba.int$y, mba.int$z, color = col) ## End(Not run)
## Not run: data(LIDAR) mba.int <- mba.surf(LIDAR, 300, 300, extend=TRUE)$xyz.est # Image plot of the surface. image(mba.int, xaxs = "r", yaxs = "r") # Perspective plot of the surface. persp(mba.int, theta = 135, phi = 30, col = "green3", scale = FALSE, ltheta = -120, shade = 0.75, expand = 10, border = NA, box = FALSE) # For a good time, I recommend using rgl. library(rgl) # Exaggerate z a bit for effect. mba.int$z <- 10*mba.int$z # Make nice colors for the rgl surface. zlim <- range(mba.int$z) zlen <- zlim[2] - zlim[1] + 1 colorlut <- terrain.colors(zlen) # Height color lookup table. col <- colorlut[mba.int$z - zlim[1] + 1 ] # Assign colors to heights for each point. open3d() surface3d(mba.int$x, mba.int$y, mba.int$z, color = col) ## End(Not run)