From 29baebfba0d7ebb9e2d388e9aa14231ff1aa979e Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Bj=C3=B6rn=20Reineking?= <bjoern@mcbr2017alien.local> Date: Fri, 5 Oct 2018 17:57:30 +0200 Subject: [PATCH] Add functions used in Ilastik --- NAMESPACE | 6 +++ R/RcppExports.R | 56 +++++++++++++++++++++ R/difference_of_gaussians.R | 13 +++++ man/difference_of_gaussians.Rd | 19 +++++++ man/gaussian_gradient_magnitude.Rd | 19 +++++++ man/hessian_of_gaussian_eigenvalues.Rd | 19 +++++++ man/laplacian_of_gaussian.Rd | 19 +++++++ man/structure_tensor_eigenvalues.Rd | 19 +++++++ man/structure_tensor_eigenvalues_2.Rd | 21 ++++++++ src/RcppExports.cpp | 66 +++++++++++++++++++++++++ src/gaussian_gradient_magnitude.cpp | 31 ++++++++++++ src/gaussian_smoothing.cpp | 5 +- src/hessian_of_gaussian_eigenvalues.cpp | 36 ++++++++++++++ src/laplacian_of_gaussian.cpp | 28 +++++++++++ src/structure_tensor_eigenvalues.cpp | 50 +++++++++++++++++++ 15 files changed, 404 insertions(+), 3 deletions(-) create mode 100644 R/difference_of_gaussians.R create mode 100644 man/difference_of_gaussians.Rd create mode 100644 man/gaussian_gradient_magnitude.Rd create mode 100644 man/hessian_of_gaussian_eigenvalues.Rd create mode 100644 man/laplacian_of_gaussian.Rd create mode 100644 man/structure_tensor_eigenvalues.Rd create mode 100644 man/structure_tensor_eigenvalues_2.Rd create mode 100644 src/gaussian_gradient_magnitude.cpp create mode 100644 src/hessian_of_gaussian_eigenvalues.cpp create mode 100644 src/laplacian_of_gaussian.cpp create mode 100644 src/structure_tensor_eigenvalues.cpp diff --git a/NAMESPACE b/NAMESPACE index 2485ee9..d3182ab 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,5 +1,11 @@ # Generated by roxygen2: do not edit by hand +export(difference_of_gaussians) +export(gaussian_gradient_magnitude) export(gaussian_smoothing) +export(hessian_of_gaussian_eigenvalues) +export(laplacian_of_gaussian) +export(structure_tensor_eigenvalues) +export(structure_tensor_eigenvalues_2) importFrom(Rcpp,sourceCpp) useDynLib(vigratools, .registration = TRUE) diff --git a/R/RcppExports.R b/R/RcppExports.R index 5ffb959..6e5c777 100644 --- a/R/RcppExports.R +++ b/R/RcppExports.R @@ -1,6 +1,17 @@ # Generated by using Rcpp::compileAttributes() -> do not edit by hand # Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393 +#' Gaussian gradient magnitude +#' +#' Magnitude of Gaussian gradient. +#' @param x Numeric matrix +#' @param sigma Scale of the Gaussian kernel +#' @return Numeric matrix; smoothed version of x +#' @export +gaussian_gradient_magnitude <- function(x, sigma) { + .Call(`_vigratools_gaussian_gradient_magnitude`, x, sigma) +} + #' Gaussian smoothing #' #' Perform isotropic Gaussian convolution. The function uses BORDER_TREATMENT_REFLECT. @@ -12,3 +23,48 @@ gaussian_smoothing <- function(x, sigma) { .Call(`_vigratools_gaussian_smoothing`, x, sigma) } +#' Hessian of Gaussian eigenvalues +#' +#' Calculate eigenvalues of Hessian of Gaussian +#' @param x Numeric matrix +#' @param sigma Scale of smoothing +#' @return List of 2 matrices with 1st and 2nd eigenvalue +#' @export +hessian_of_gaussian_eigenvalues <- function(x, sigma) { + .Call(`_vigratools_hessian_of_gaussian_eigenvalues`, x, sigma) +} + +#' Laplacian of Gaussian +#' +#' Calculates laplacian of Gaussian. +#' @param x Numeric matrix +#' @param sigma Scale of the Gaussian kernel +#' @return Numeric matrix; smoothed version of x +#' @export +laplacian_of_gaussian <- function(x, sigma) { + .Call(`_vigratools_laplacian_of_gaussian`, x, sigma) +} + +#' Structure tensor eigenvalues (2 argument) +#' +#' Calculate eigenvalues of structure tensor +#' @param x Numeric matrix +#' @param inner_sigma Scale of inner smoothing +#' @param outer_sigma Scale of outer smoothing +#' @return List of 2 matrices with 1st and 2nd eigenvalue +#' @export +structure_tensor_eigenvalues_2 <- function(x, inner_sigma, outer_sigma) { + .Call(`_vigratools_structure_tensor_eigenvalues_2`, x, inner_sigma, outer_sigma) +} + +#' Structure tensor eigenvalues +#' +#' Calculate eigenvalues of structure tensor +#' @param x Numeric matrix +#' @param sigma Scale of inner smoothing; outer smoothing is sigma/2. +#' @return List of 2 matrices with 1st and 2nd eigenvalue +#' @export +structure_tensor_eigenvalues <- function(x, sigma) { + .Call(`_vigratools_structure_tensor_eigenvalues`, x, sigma) +} + diff --git a/R/difference_of_gaussians.R b/R/difference_of_gaussians.R new file mode 100644 index 0000000..ee7746a --- /dev/null +++ b/R/difference_of_gaussians.R @@ -0,0 +1,13 @@ +#' Difference of Gaussians +#' +#' Compute difference of Gaussians +#' @param x Input matrix +#' @param sigma Scale of smoothing (smaller gaussian is 0.66 times sigma) +#' @return Matrix containing difference of original image smoothed at two scales. +#' @export +#' +difference_of_gaussians <- function(x, sigma) { + smoothed_1 <- gaussian_smoothing(x, sigma) + smoothed_2 <- gaussian_smoothing(x, 0.66*sigma) + smoothed_1 - smoothed_2 +} diff --git a/man/difference_of_gaussians.Rd b/man/difference_of_gaussians.Rd new file mode 100644 index 0000000..a1066fb --- /dev/null +++ b/man/difference_of_gaussians.Rd @@ -0,0 +1,19 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/difference_of_gaussians.R +\name{difference_of_gaussians} +\alias{difference_of_gaussians} +\title{Difference of Gaussians} +\usage{ +difference_of_gaussians(x, sigma) +} +\arguments{ +\item{x}{Input matrix} + +\item{sigma}{Scale of smoothing (smaller gaussian is 0.66 times sigma)} +} +\value{ +Matrix containing difference of original image smoothed at two scales. +} +\description{ +Compute difference of Gaussians +} diff --git a/man/gaussian_gradient_magnitude.Rd b/man/gaussian_gradient_magnitude.Rd new file mode 100644 index 0000000..5c22da4 --- /dev/null +++ b/man/gaussian_gradient_magnitude.Rd @@ -0,0 +1,19 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/RcppExports.R +\name{gaussian_gradient_magnitude} +\alias{gaussian_gradient_magnitude} +\title{Gaussian gradient magnitude} +\usage{ +gaussian_gradient_magnitude(x, sigma) +} +\arguments{ +\item{x}{Numeric matrix} + +\item{sigma}{Scale of the Gaussian kernel} +} +\value{ +Numeric matrix; smoothed version of x +} +\description{ +Magnitude of Gaussian gradient. +} diff --git a/man/hessian_of_gaussian_eigenvalues.Rd b/man/hessian_of_gaussian_eigenvalues.Rd new file mode 100644 index 0000000..6c7d075 --- /dev/null +++ b/man/hessian_of_gaussian_eigenvalues.Rd @@ -0,0 +1,19 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/RcppExports.R +\name{hessian_of_gaussian_eigenvalues} +\alias{hessian_of_gaussian_eigenvalues} +\title{Hessian of Gaussian eigenvalues} +\usage{ +hessian_of_gaussian_eigenvalues(x, sigma) +} +\arguments{ +\item{x}{Numeric matrix} + +\item{sigma}{Scale of smoothing} +} +\value{ +List of 2 matrices with 1st and 2nd eigenvalue +} +\description{ +Calculate eigenvalues of Hessian of Gaussian +} diff --git a/man/laplacian_of_gaussian.Rd b/man/laplacian_of_gaussian.Rd new file mode 100644 index 0000000..3a7c47d --- /dev/null +++ b/man/laplacian_of_gaussian.Rd @@ -0,0 +1,19 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/RcppExports.R +\name{laplacian_of_gaussian} +\alias{laplacian_of_gaussian} +\title{Laplacian of Gaussian} +\usage{ +laplacian_of_gaussian(x, sigma) +} +\arguments{ +\item{x}{Numeric matrix} + +\item{sigma}{Scale of the Gaussian kernel} +} +\value{ +Numeric matrix; smoothed version of x +} +\description{ +Calculates laplacian of Gaussian. +} diff --git a/man/structure_tensor_eigenvalues.Rd b/man/structure_tensor_eigenvalues.Rd new file mode 100644 index 0000000..5f43fda --- /dev/null +++ b/man/structure_tensor_eigenvalues.Rd @@ -0,0 +1,19 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/RcppExports.R +\name{structure_tensor_eigenvalues} +\alias{structure_tensor_eigenvalues} +\title{Structure tensor eigenvalues} +\usage{ +structure_tensor_eigenvalues(x, sigma) +} +\arguments{ +\item{x}{Numeric matrix} + +\item{sigma}{Scale of inner smoothing; outer smoothing is sigma/2.} +} +\value{ +List of 2 matrices with 1st and 2nd eigenvalue +} +\description{ +Calculate eigenvalues of structure tensor +} diff --git a/man/structure_tensor_eigenvalues_2.Rd b/man/structure_tensor_eigenvalues_2.Rd new file mode 100644 index 0000000..94712cc --- /dev/null +++ b/man/structure_tensor_eigenvalues_2.Rd @@ -0,0 +1,21 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/RcppExports.R +\name{structure_tensor_eigenvalues_2} +\alias{structure_tensor_eigenvalues_2} +\title{Structure tensor eigenvalues (2 argument)} +\usage{ +structure_tensor_eigenvalues_2(x, inner_sigma, outer_sigma) +} +\arguments{ +\item{x}{Numeric matrix} + +\item{inner_sigma}{Scale of inner smoothing} + +\item{outer_sigma}{Scale of outer smoothing} +} +\value{ +List of 2 matrices with 1st and 2nd eigenvalue +} +\description{ +Calculate eigenvalues of structure tensor +} diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp index 9301353..0d71229 100644 --- a/src/RcppExports.cpp +++ b/src/RcppExports.cpp @@ -5,6 +5,18 @@ using namespace Rcpp; +// gaussian_gradient_magnitude +NumericMatrix gaussian_gradient_magnitude(NumericMatrix x, double sigma); +RcppExport SEXP _vigratools_gaussian_gradient_magnitude(SEXP xSEXP, SEXP sigmaSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< NumericMatrix >::type x(xSEXP); + Rcpp::traits::input_parameter< double >::type sigma(sigmaSEXP); + rcpp_result_gen = Rcpp::wrap(gaussian_gradient_magnitude(x, sigma)); + return rcpp_result_gen; +END_RCPP +} // gaussian_smoothing NumericMatrix gaussian_smoothing(NumericMatrix x, double sigma); RcppExport SEXP _vigratools_gaussian_smoothing(SEXP xSEXP, SEXP sigmaSEXP) { @@ -17,9 +29,63 @@ BEGIN_RCPP return rcpp_result_gen; END_RCPP } +// hessian_of_gaussian_eigenvalues +List hessian_of_gaussian_eigenvalues(NumericMatrix x, double sigma); +RcppExport SEXP _vigratools_hessian_of_gaussian_eigenvalues(SEXP xSEXP, SEXP sigmaSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< NumericMatrix >::type x(xSEXP); + Rcpp::traits::input_parameter< double >::type sigma(sigmaSEXP); + rcpp_result_gen = Rcpp::wrap(hessian_of_gaussian_eigenvalues(x, sigma)); + return rcpp_result_gen; +END_RCPP +} +// laplacian_of_gaussian +NumericMatrix laplacian_of_gaussian(NumericMatrix x, double sigma); +RcppExport SEXP _vigratools_laplacian_of_gaussian(SEXP xSEXP, SEXP sigmaSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< NumericMatrix >::type x(xSEXP); + Rcpp::traits::input_parameter< double >::type sigma(sigmaSEXP); + rcpp_result_gen = Rcpp::wrap(laplacian_of_gaussian(x, sigma)); + return rcpp_result_gen; +END_RCPP +} +// structure_tensor_eigenvalues_2 +List structure_tensor_eigenvalues_2(NumericMatrix x, double inner_sigma, double outer_sigma); +RcppExport SEXP _vigratools_structure_tensor_eigenvalues_2(SEXP xSEXP, SEXP inner_sigmaSEXP, SEXP outer_sigmaSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< NumericMatrix >::type x(xSEXP); + Rcpp::traits::input_parameter< double >::type inner_sigma(inner_sigmaSEXP); + Rcpp::traits::input_parameter< double >::type outer_sigma(outer_sigmaSEXP); + rcpp_result_gen = Rcpp::wrap(structure_tensor_eigenvalues_2(x, inner_sigma, outer_sigma)); + return rcpp_result_gen; +END_RCPP +} +// structure_tensor_eigenvalues +List structure_tensor_eigenvalues(NumericMatrix x, double sigma); +RcppExport SEXP _vigratools_structure_tensor_eigenvalues(SEXP xSEXP, SEXP sigmaSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< NumericMatrix >::type x(xSEXP); + Rcpp::traits::input_parameter< double >::type sigma(sigmaSEXP); + rcpp_result_gen = Rcpp::wrap(structure_tensor_eigenvalues(x, sigma)); + return rcpp_result_gen; +END_RCPP +} static const R_CallMethodDef CallEntries[] = { + {"_vigratools_gaussian_gradient_magnitude", (DL_FUNC) &_vigratools_gaussian_gradient_magnitude, 2}, {"_vigratools_gaussian_smoothing", (DL_FUNC) &_vigratools_gaussian_smoothing, 2}, + {"_vigratools_hessian_of_gaussian_eigenvalues", (DL_FUNC) &_vigratools_hessian_of_gaussian_eigenvalues, 2}, + {"_vigratools_laplacian_of_gaussian", (DL_FUNC) &_vigratools_laplacian_of_gaussian, 2}, + {"_vigratools_structure_tensor_eigenvalues_2", (DL_FUNC) &_vigratools_structure_tensor_eigenvalues_2, 3}, + {"_vigratools_structure_tensor_eigenvalues", (DL_FUNC) &_vigratools_structure_tensor_eigenvalues, 2}, {NULL, NULL, 0} }; diff --git a/src/gaussian_gradient_magnitude.cpp b/src/gaussian_gradient_magnitude.cpp new file mode 100644 index 0000000..73f9a19 --- /dev/null +++ b/src/gaussian_gradient_magnitude.cpp @@ -0,0 +1,31 @@ +//[[Rcpp::depends(vigRa)]] +#include <Rcpp.h> +#include <vigra/multi_array.hxx> +#include <vigra/multi_convolution.hxx> +#include <vigra/numerictraits.hxx> +//#include <vigra/multi_blockwise.hxx> +using namespace Rcpp; +using namespace vigra; + + +//' Gaussian gradient magnitude +//' +//' Magnitude of Gaussian gradient. +//' @param x Numeric matrix +//' @param sigma Scale of the Gaussian kernel +//' @return Numeric matrix; smoothed version of x +//' @export +// [[Rcpp::export]] +NumericMatrix gaussian_gradient_magnitude(NumericMatrix x, double sigma) { + const int nrows = x.nrow(); + const int ncols = x.ncol(); + Shape2 shape(nrows, ncols); + MultiArrayView<2, double> input_image(shape, x.begin()); + NumericMatrix result(nrows, ncols); + MultiArrayView<2, double> output_image(shape, result.begin()); + // gaussianGradientMagnitude(input_image, output_image, sigma); + //std::copy(output_image.begin(), output_image.end(), result.begin()); + return(result); +} + + diff --git a/src/gaussian_smoothing.cpp b/src/gaussian_smoothing.cpp index df3ea6c..a1fe864 100644 --- a/src/gaussian_smoothing.cpp +++ b/src/gaussian_smoothing.cpp @@ -2,9 +2,6 @@ #include <Rcpp.h> #include <vigra/multi_array.hxx> #include <vigra/convolution.hxx> -// #include <vigra/tensorutilities.hxx> -// #include <vigra/mathutil.hxx> -// #include <vigra/multi_tensorutilities.hxx> using namespace Rcpp; using namespace vigra; @@ -27,3 +24,5 @@ NumericMatrix gaussian_smoothing(NumericMatrix x, double sigma) { gaussianSmoothing(input_image, output_image, sigma); return(result); } + + diff --git a/src/hessian_of_gaussian_eigenvalues.cpp b/src/hessian_of_gaussian_eigenvalues.cpp new file mode 100644 index 0000000..f947b51 --- /dev/null +++ b/src/hessian_of_gaussian_eigenvalues.cpp @@ -0,0 +1,36 @@ +//[[Rcpp::depends(vigRa)]] +#include <Rcpp.h> +#include <vigra/multi_array.hxx> +#include <vigra/convolution.hxx> +#include <vigra/tensorutilities.hxx> +#include <vigra/mathutil.hxx> +#include <vigra/multi_tensorutilities.hxx> +using namespace Rcpp; +using namespace vigra; + +//' Hessian of Gaussian eigenvalues +//' +//' Calculate eigenvalues of Hessian of Gaussian +//' @param x Numeric matrix +//' @param sigma Scale of smoothing +//' @return List of 2 matrices with 1st and 2nd eigenvalue +//' @export +// [[Rcpp::export]] +List hessian_of_gaussian_eigenvalues(NumericMatrix x, double sigma) { + const int nrows = x.nrow(); + const int ncols = x.ncol(); + Shape2 shape(nrows, ncols); + MultiArrayView<2, double> input_image(shape, x.begin()); + MultiArray<2, TinyVector<float, 3> > hessian(shape); + MultiArray<2, TinyVector<float, 2> > eigen(shape); + hessianMatrixOfGaussian(input_image, hessian, sigma); + tensorEigenvaluesMultiArray(hessian, eigen); + MultiArrayView<2, float, StridedArrayTag> first = eigen.bindElementChannel(0); + MultiArrayView<2, float, StridedArrayTag> second = eigen.bindElementChannel(1); + NumericMatrix first_result(nrows, ncols); + NumericMatrix second_result(nrows, ncols); + std::copy(first.begin(), first.end(), first_result.begin()); + std::copy(second.begin(), second.end(), second_result.begin()); + return Rcpp::List::create(Rcpp::Named("ei1") = first_result, + Rcpp::Named("ei2") = second_result); +} diff --git a/src/laplacian_of_gaussian.cpp b/src/laplacian_of_gaussian.cpp new file mode 100644 index 0000000..f1e8600 --- /dev/null +++ b/src/laplacian_of_gaussian.cpp @@ -0,0 +1,28 @@ +//[[Rcpp::depends(vigRa)]] +#include <Rcpp.h> +#include <vigra/multi_array.hxx> +#include <vigra/convolution.hxx> +using namespace Rcpp; +using namespace vigra; + + +//' Laplacian of Gaussian +//' +//' Calculates laplacian of Gaussian. +//' @param x Numeric matrix +//' @param sigma Scale of the Gaussian kernel +//' @return Numeric matrix; smoothed version of x +//' @export +// [[Rcpp::export]] +NumericMatrix laplacian_of_gaussian(NumericMatrix x, double sigma) { + const int nrows = x.nrow(); + const int ncols = x.ncol(); + Shape2 shape(nrows, ncols); + MultiArrayView<2, double> input_image(shape, x.begin()); + NumericMatrix result(nrows, ncols); + MultiArrayView<2, double> output_image(shape, result.begin()); + laplacianOfGaussian(input_image, output_image, sigma); + return(result); +} + + diff --git a/src/structure_tensor_eigenvalues.cpp b/src/structure_tensor_eigenvalues.cpp new file mode 100644 index 0000000..841d636 --- /dev/null +++ b/src/structure_tensor_eigenvalues.cpp @@ -0,0 +1,50 @@ +//[[Rcpp::depends(vigRa)]] +#include <Rcpp.h> +#include <vigra/multi_array.hxx> +#include <vigra/convolution.hxx> +#include <vigra/tensorutilities.hxx> +#include <vigra/mathutil.hxx> +#include <vigra/multi_tensorutilities.hxx> +using namespace Rcpp; +using namespace vigra; + +//' Structure tensor eigenvalues (2 argument) +//' +//' Calculate eigenvalues of structure tensor +//' @param x Numeric matrix +//' @param inner_sigma Scale of inner smoothing +//' @param outer_sigma Scale of outer smoothing +//' @return List of 2 matrices with 1st and 2nd eigenvalue +//' @export +// [[Rcpp::export]] +List structure_tensor_eigenvalues_2(NumericMatrix x, double inner_sigma, double outer_sigma) { + const int nrows = x.nrow(); + const int ncols = x.ncol(); + Shape2 shape(nrows, ncols); + MultiArrayView<2, double> input_image(shape, x.begin()); + MultiArray<2, TinyVector<float, 3> > tensor(shape); + MultiArray<2, TinyVector<float, 2> > eigen(shape); + structureTensor(input_image, tensor, inner_sigma, outer_sigma); + tensorEigenvaluesMultiArray(tensor, eigen); + MultiArrayView<2, float, StridedArrayTag> first = eigen.bindElementChannel(0); + MultiArrayView<2, float, StridedArrayTag> second = eigen.bindElementChannel(1); + NumericMatrix first_result(nrows, ncols); + NumericMatrix second_result(nrows, ncols); + std::copy(first.begin(), first.end(), first_result.begin()); + std::copy(second.begin(), second.end(), second_result.begin()); + return Rcpp::List::create(Rcpp::Named("ei1") = first_result, + Rcpp::Named("ei2") = second_result); +} + + +//' Structure tensor eigenvalues +//' +//' Calculate eigenvalues of structure tensor +//' @param x Numeric matrix +//' @param sigma Scale of inner smoothing; outer smoothing is sigma/2. +//' @return List of 2 matrices with 1st and 2nd eigenvalue +//' @export +// [[Rcpp::export]] +List structure_tensor_eigenvalues(NumericMatrix x, double sigma) { + return structure_tensor_eigenvalues_2(x, sigma, sigma/2.0); +} -- GitLab