我正在尝试计算 TreeLS 中的 FastPointMetrics 函数使用的 EigenVector12 指标。由于此函数是针对 TLS 数据而设计的,因此它不适用于我正在使用的 ALS 点云。该函数的确切措辞是:“3D 特征向量系数,第二特征向量的第一次加载”。然后我找到了 point_metrics() 函数并能够创建我自己的函数。
我想做的是为每个点的 20 个点的最近邻域获取第二个特征向量的第一次加载,然后取 20m 像素大小的值的平均值,最后得到值的栅格。
这是我写的函数(如有错误请指正):
EV12_metrics = function(x,y,z, th1 = 20) {
xyz <- cbind(x,y,z) # create XYZ matrix
cov_m <- cov(xyz) # compute covariance of the point locations
e = eigen(cov_m) # get eigenvalues and eigenvectors from the matrix
vect12 = e$vectors[4] # Select the 1st load of the second eigenvector
return(list(EV12 = vect12)) # return with a list of values
}
M = point_metrics(laz, ~EV12_metrics(X,Y,Z), k = 20) # apply the function to the point cloud
laz = add_attribute(laz, M$EV12, "EV12") # Add the attribute to the LAS points
gridmets = grid_metrics(laz, mean(EV12), res = 20) # Create raster of the mean of the values at 20m resolution
这很有效,只是需要一段时间来处理。然后我找到了 point_eigenvalues() 函数,但这只给了我最大、中和最小的特征值以及其他我不需要的指标。我想知道这个 point_eigenvalues() 函数是否可以执行我需要的操作,或者是否有任何方法可以加速我上面编写的函数。谢谢!
point_metrics
的文档明确提到在这种情况下使用C++代码。它为您提供了一个示例,您可以复制并粘贴仅更改输出变量
Rcpp::sourceCpp(code = "
#include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
SEXP EV12_metrics_cpp(arma::mat A) {
arma::mat coeff;
arma::mat score;
arma::vec latent;
arma::princomp(coeff, score, latent, A);
return(Rcpp::wrap(coeff));
}")
EV12_metrics2 = function(x,y,z, th1 = 25, th2 = 6) {
xyz <- cbind(x,y,z)
e <- EV12_metrics_cpp(xyz)
vect12 = e[4]
return(list(EV12 = vect12)) # return with a list of values
}
M2 = point_metrics(laz, ~EV12_metrics2(X,Y,Z), k = 20) # apply the function to the point cloud
这快了 10 倍。
使用 github 上提供的
lidR
开发版本(预计将于 1 月 24 日发布),函数 point_eigenvalue()
有一个新参数 coeff
M3 = point_eigenvalues(las, 20, coeffs = T)
M3 = M3$coeff01
使用
library(lidR)
LASfile <- system.file("extdata", "Megaplot.laz", package="lidR")
laz = readLAS(LASfile)
EV12_metrics
:36秒EV12_metrics2
:3秒point_eigenvalues
:1.2秒