diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml index cb17061..562fe0f 100644 --- a/.github/workflows/R-CMD-check.yaml +++ b/.github/workflows/R-CMD-check.yaml @@ -1,16 +1,13 @@ -# For help debugging build failures open an issue on the RStudio community with the 'github-actions' tag. -# https://community.rstudio.com/new-topic?category=Package%20development&tags=github-actions +# Workflow derived from https://github.com/r-lib/actions/tree/v2/examples +# Need help debugging build failures? Start at https://github.com/r-lib/actions#where-to-find-help on: push: - branches: - - main - - master + branches: [main, master] pull_request: - branches: - - main - - master -name: R-CMD-check +name: R-CMD-check.yaml + +permissions: read-all jobs: R-CMD-check: @@ -22,64 +19,33 @@ jobs: fail-fast: false matrix: config: + - {os: macos-latest, r: 'release'} - {os: windows-latest, r: 'release'} - - {os: macOS-latest, r: 'release'} - - {os: ubuntu-20.04, r: 'release', rspm: "https://packagemanager.rstudio.com/cran/__linux__/focal/latest"} - - {os: ubuntu-20.04, r: 'devel', rspm: "https://packagemanager.rstudio.com/cran/__linux__/focal/latest", http-user-agent: "R/4.1.0 (ubuntu-20.04) R (4.1.0 x86_64-pc-linux-gnu x86_64 linux-gnu) on GitHub Actions" } + - {os: ubuntu-latest, r: 'devel', http-user-agent: 'release'} + - {os: ubuntu-latest, r: 'release'} + - {os: ubuntu-latest, r: 'oldrel-1'} env: - R_REMOTES_NO_ERRORS_FROM_WARNINGS: true - RSPM: ${{ matrix.config.rspm }} GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }} + R_KEEP_PKG_SOURCE: yes steps: - - uses: actions/checkout@v2 + - uses: actions/checkout@v4 + + - uses: r-lib/actions/setup-pandoc@v2 - - uses: r-lib/actions/setup-r@v1 + - uses: r-lib/actions/setup-r@v2 with: r-version: ${{ matrix.config.r }} + http-user-agent: ${{ matrix.config.http-user-agent }} + use-public-rspm: true - - uses: r-lib/actions/setup-pandoc@v1 - - - name: Query dependencies - run: | - install.packages('remotes') - saveRDS(remotes::dev_package_deps(dependencies = TRUE), ".github/depends.Rds", version = 2) - writeLines(sprintf("R-%i.%i", getRversion()$major, getRversion()$minor), ".github/R-version") - shell: Rscript {0} - - - name: Restore R package cache - uses: actions/cache@v2 + - uses: r-lib/actions/setup-r-dependencies@v2 with: - path: ${{ env.R_LIBS_USER }} - key: ${{ runner.os }}-${{ hashFiles('.github/R-version') }}-1-${{ hashFiles('.github/depends.Rds') }} - restore-keys: ${{ runner.os }}-${{ hashFiles('.github/R-version') }}-1- - - - name: Install system dependencies - if: runner.os == 'Linux' - run: | - while read -r cmd - do - eval sudo $cmd - done < <(Rscript -e 'writeLines(c(remotes::system_requirements("ubuntu", "20.04"), remotes::system_requirements("ubuntu", "20.04", package="curl")))') - - - name: Install dependencies - run: | - remotes::install_deps(dependencies = TRUE) - remotes::install_cran("rcmdcheck") - shell: Rscript {0} - - - name: Check - env: - _R_CHECK_CRAN_INCOMING_REMOTE_: false - run: | - options(crayon.enabled = TRUE) - rcmdcheck::rcmdcheck(args = c("--no-manual", "--as-cran"), error_on = "warning", check_dir = "check") - shell: Rscript {0} + extra-packages: any::rcmdcheck + needs: check - - name: Upload check results - if: failure() - uses: actions/upload-artifact@main + - uses: r-lib/actions/check-r-package@v2 with: - name: ${{ runner.os }}-r${{ matrix.config.r }}-results - path: check + upload-snapshots: true + build_args: 'c("--no-manual","--compact-vignettes=gs+qpdf")' diff --git a/.github/workflows/pkgdown.yaml b/.github/workflows/pkgdown.yaml index e809466..bfc9f4d 100644 --- a/.github/workflows/pkgdown.yaml +++ b/.github/workflows/pkgdown.yaml @@ -1,50 +1,49 @@ +# Workflow derived from https://github.com/r-lib/actions/tree/v2/examples +# Need help debugging build failures? Start at https://github.com/r-lib/actions#where-to-find-help on: push: - branches: - - main - - master - tags: - -'*' + branches: [main, master] + pull_request: + release: + types: [published] + workflow_dispatch: -name: pkgdown +name: pkgdown.yaml + +permissions: read-all jobs: pkgdown: - runs-on: macOS-latest + runs-on: ubuntu-latest + # Only restrict concurrency for non-PR jobs + concurrency: + group: pkgdown-${{ github.event_name != 'pull_request' || github.run_id }} env: GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }} + permissions: + contents: write steps: - - uses: actions/checkout@v2 - - - uses: r-lib/actions/setup-r@v1 + - uses: actions/checkout@v4 - - uses: r-lib/actions/setup-pandoc@v1 + - uses: r-lib/actions/setup-pandoc@v2 - - name: Query dependencies - run: | - install.packages('remotes') - saveRDS(remotes::dev_package_deps(dependencies = TRUE), ".github/depends.Rds", version = 2) - writeLines(sprintf("R-%i.%i", getRversion()$major, getRversion()$minor), ".github/R-version") - shell: Rscript {0} + - uses: r-lib/actions/setup-r@v2 + with: + use-public-rspm: true - - name: Restore R package cache - uses: actions/cache@v2 + - uses: r-lib/actions/setup-r-dependencies@v2 with: - path: ${{ env.R_LIBS_USER }} - key: ${{ runner.os }}-${{ hashFiles('.github/R-version') }}-1-${{ hashFiles('.github/depends.Rds') }} - restore-keys: ${{ runner.os }}-${{ hashFiles('.github/R-version') }}-1- - - - name: Install dependencies - run: | - remotes::install_deps(dependencies = TRUE) - install.packages("pkgdown", type = "binary") - shell: Rscript {0} + extra-packages: any::pkgdown, local::. + needs: website - - name: Install package - run: R CMD INSTALL . + - name: Build site + run: pkgdown::build_site_github_pages(new_process = FALSE, install = FALSE) + shell: Rscript {0} - - name: Deploy package - run: | - git config --local user.email "actions@github.com" - git config --local user.name "GitHub Actions" - Rscript -e 'pkgdown::deploy_to_branch(new_process = FALSE)' + - name: Deploy to GitHub pages 🚀 + if: github.event_name != 'pull_request' + uses: JamesIves/github-pages-deploy-action@v4.5.0 + with: + clean: false + branch: gh-pages + folder: docs diff --git a/DESCRIPTION b/DESCRIPTION index 67466b1..029824f 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Type: Package Package: natcpp Title: Fast C++ Primitives for the 'NeuroAnatomy Toolbox' -Version: 0.1.1 +Version: 0.1.1.9000 Authors@R: person(given = "Gregory", family = "Jefferis", @@ -16,7 +16,7 @@ Description: Fast functions implemented in C++ via 'Rcpp' to support the package will automatically use routines from this package when it is available to enable large performance gains. License: GPL (>= 3) -URL: https://github.com/natverse/natcpp +URL: https://github.com/natverse/natcpp, https://natverse.org/natcpp/ BugReports: https://github.com/natverse/natcpp/issues Imports: Rcpp (>= 1.0.6) diff --git a/NAMESPACE b/NAMESPACE index d8be76d..504bfff 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -2,8 +2,11 @@ export(c_EdgeListFromSegList) export(c_ListofMatrixRows) +export(c_coords21dindex) +export(c_ijkpos) export(c_listlengths) export(c_seglengths) +export(c_sub2ind) export(c_topntail) export(c_topntail_list) export(c_total_cable) diff --git a/R/RcppExports.R b/R/RcppExports.R index d29e627..c90d42d 100644 --- a/R/RcppExports.R +++ b/R/RcppExports.R @@ -20,6 +20,43 @@ c_total_cable <- function(sl, x, y, z) { .Call(`_natcpp_c_total_cable`, sl, x, y, z) } +#' Convert physical coordinates to pixel coordinates +#' +#' @param xyz Nx3 matrix of physical coordinates +#' @param dims Integer dimensions of the 3d image array +#' @param origin Numeric: 3d coordinates of the origin +#' @param voxdims Numeric: 3 numbers describing the voxel dimensions +#' @param clamp Logical: whether or not to clamp values within the pixel +#' boundaries of the image. +#' @return Nx3 integer matrix of pixel coordinates +#' @export +c_ijkpos <- function(xyz, origin, voxdims, dims, clamp = FALSE) { + .Call(`_natcpp_c_ijkpos`, xyz, origin, voxdims, dims, clamp) +} + +#' Find 1D index given n-dimensional indices +#' @param dims Integer dimensions of the array (usually 3d) +#' @param indices Nx3 integer matrix of pixel coordinates +#' @return numeric vector of linear indices into the array +#' @export +c_sub2ind <- function(dims, indices) { + .Call(`_natcpp_c_sub2ind`, dims, indices) +} + +#' Convert physical coordinates to 1d indices into image array +#' +#' @param xyz Nx3 matrix or data.frame of physical coordinates +#' @param dims Integer dimensions of the 3d image array +#' @param origin Numeric: 3d coordinates of the origin +#' @param voxdims Numeric: 3 numbers describing the voxel dimensions +#' @param clamp Logical: whether or not to clamp values within the pixel +#' boundaries of the image. +#' @return Nx3 integer matrix of pixel coordinates +#' @export +c_coords21dindex <- function(xyz, origin, voxdims, dims, clamp = FALSE) { + .Call(`_natcpp_c_coords21dindex`, xyz, origin, voxdims, dims, clamp) +} + #' Convert a matrix into list of row vectors #' #' @details Typically this will be for 3D coordinates but there are no limits diff --git a/README.md b/README.md index edb99a9..9fb1d3d 100644 --- a/README.md +++ b/README.md @@ -7,6 +7,7 @@ [![R-CMD-check](https://github.com/natverse/natcpp/workflows/R-CMD-check/badge.svg)](https://github.com/natverse/natcpp/actions) [![Lifecycle: experimental](https://img.shields.io/badge/lifecycle-experimental-orange.svg)](https://lifecycle.r-lib.org/articles/stages.html#experimental) [![CRAN status](https://www.r-pkg.org/badges/version/natcpp)](https://CRAN.R-project.org/package=natcpp) +[![R-CMD-check](https://github.com/natverse/natcpp/actions/workflows/R-CMD-check.yaml/badge.svg)](https://github.com/natverse/natcpp/actions/workflows/R-CMD-check.yaml) The goal of *natcpp* is to provide accelerated routines through compiled C++ diff --git a/_pkgdown.yml b/_pkgdown.yml index 69f0e17..cc171ac 100644 --- a/_pkgdown.yml +++ b/_pkgdown.yml @@ -1,26 +1,33 @@ +url: https://natverse.org/natcpp/ +template: + bootstrap: 5 navbar: structure: left: - - home - intro - reference - articles - tutorials - - news - help + - news right: - - natverse + - search - github + - natverse + - lightswitch components: - home: - icon: fas fa-home fa-lg - href: index.html reference: text: Reference href: reference/index.html + search: + search: [] + news: + text: Changelog + href: news/index.html github: icon: fab fa-github fa-lg href: https://github.com/natverse/natcpp/ + aria-label: GitHub natverse: text: natverse href: https://natverse.github.io diff --git a/man/c_coords21dindex.Rd b/man/c_coords21dindex.Rd new file mode 100644 index 0000000..4c80bfe --- /dev/null +++ b/man/c_coords21dindex.Rd @@ -0,0 +1,26 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/RcppExports.R +\name{c_coords21dindex} +\alias{c_coords21dindex} +\title{Convert physical coordinates to 1d indices into image array} +\usage{ +c_coords21dindex(xyz, origin, voxdims, dims, clamp = FALSE) +} +\arguments{ +\item{xyz}{Nx3 matrix or data.frame of physical coordinates} + +\item{origin}{Numeric: 3d coordinates of the origin} + +\item{voxdims}{Numeric: 3 numbers describing the voxel dimensions} + +\item{dims}{Integer dimensions of the 3d image array} + +\item{clamp}{Logical: whether or not to clamp values within the pixel +boundaries of the image.} +} +\value{ +Nx3 integer matrix of pixel coordinates +} +\description{ +Convert physical coordinates to 1d indices into image array +} diff --git a/man/c_ijkpos.Rd b/man/c_ijkpos.Rd new file mode 100644 index 0000000..9750998 --- /dev/null +++ b/man/c_ijkpos.Rd @@ -0,0 +1,26 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/RcppExports.R +\name{c_ijkpos} +\alias{c_ijkpos} +\title{Convert physical coordinates to pixel coordinates} +\usage{ +c_ijkpos(xyz, origin, voxdims, dims, clamp = FALSE) +} +\arguments{ +\item{xyz}{Nx3 matrix of physical coordinates} + +\item{origin}{Numeric: 3d coordinates of the origin} + +\item{voxdims}{Numeric: 3 numbers describing the voxel dimensions} + +\item{dims}{Integer dimensions of the 3d image array} + +\item{clamp}{Logical: whether or not to clamp values within the pixel +boundaries of the image.} +} +\value{ +Nx3 integer matrix of pixel coordinates +} +\description{ +Convert physical coordinates to pixel coordinates +} diff --git a/man/c_sub2ind.Rd b/man/c_sub2ind.Rd new file mode 100644 index 0000000..9d9a653 --- /dev/null +++ b/man/c_sub2ind.Rd @@ -0,0 +1,19 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/RcppExports.R +\name{c_sub2ind} +\alias{c_sub2ind} +\title{Find 1D index given n-dimensional indices} +\usage{ +c_sub2ind(dims, indices) +} +\arguments{ +\item{dims}{Integer dimensions of the array (usually 3d)} + +\item{indices}{Nx3 integer matrix of pixel coordinates} +} +\value{ +numeric vector of linear indices into the array +} +\description{ +Find 1D index given n-dimensional indices +} diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp index aed385b..3415f85 100644 --- a/src/RcppExports.cpp +++ b/src/RcppExports.cpp @@ -38,6 +38,48 @@ BEGIN_RCPP return rcpp_result_gen; END_RCPP } +// c_ijkpos +IntegerMatrix c_ijkpos(SEXP xyz, NumericVector origin, NumericVector voxdims, IntegerVector dims, bool clamp); +RcppExport SEXP _natcpp_c_ijkpos(SEXP xyzSEXP, SEXP originSEXP, SEXP voxdimsSEXP, SEXP dimsSEXP, SEXP clampSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< SEXP >::type xyz(xyzSEXP); + Rcpp::traits::input_parameter< NumericVector >::type origin(originSEXP); + Rcpp::traits::input_parameter< NumericVector >::type voxdims(voxdimsSEXP); + Rcpp::traits::input_parameter< IntegerVector >::type dims(dimsSEXP); + Rcpp::traits::input_parameter< bool >::type clamp(clampSEXP); + rcpp_result_gen = Rcpp::wrap(c_ijkpos(xyz, origin, voxdims, dims, clamp)); + return rcpp_result_gen; +END_RCPP +} +// c_sub2ind +NumericVector c_sub2ind(IntegerVector dims, IntegerMatrix indices); +RcppExport SEXP _natcpp_c_sub2ind(SEXP dimsSEXP, SEXP indicesSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< IntegerVector >::type dims(dimsSEXP); + Rcpp::traits::input_parameter< IntegerMatrix >::type indices(indicesSEXP); + rcpp_result_gen = Rcpp::wrap(c_sub2ind(dims, indices)); + return rcpp_result_gen; +END_RCPP +} +// c_coords21dindex +NumericVector c_coords21dindex(SEXP xyz, NumericVector origin, NumericVector voxdims, IntegerVector dims, bool clamp); +RcppExport SEXP _natcpp_c_coords21dindex(SEXP xyzSEXP, SEXP originSEXP, SEXP voxdimsSEXP, SEXP dimsSEXP, SEXP clampSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< SEXP >::type xyz(xyzSEXP); + Rcpp::traits::input_parameter< NumericVector >::type origin(originSEXP); + Rcpp::traits::input_parameter< NumericVector >::type voxdims(voxdimsSEXP); + Rcpp::traits::input_parameter< IntegerVector >::type dims(dimsSEXP); + Rcpp::traits::input_parameter< bool >::type clamp(clampSEXP); + rcpp_result_gen = Rcpp::wrap(c_coords21dindex(xyz, origin, voxdims, dims, clamp)); + return rcpp_result_gen; +END_RCPP +} // c_ListofMatrixRows List c_ListofMatrixRows(const SEXP& object); RcppExport SEXP _natcpp_c_ListofMatrixRows(SEXP objectSEXP) { @@ -97,6 +139,9 @@ END_RCPP static const R_CallMethodDef CallEntries[] = { {"_natcpp_c_seglengths", (DL_FUNC) &_natcpp_c_seglengths, 4}, {"_natcpp_c_total_cable", (DL_FUNC) &_natcpp_c_total_cable, 4}, + {"_natcpp_c_ijkpos", (DL_FUNC) &_natcpp_c_ijkpos, 5}, + {"_natcpp_c_sub2ind", (DL_FUNC) &_natcpp_c_sub2ind, 2}, + {"_natcpp_c_coords21dindex", (DL_FUNC) &_natcpp_c_coords21dindex, 5}, {"_natcpp_c_ListofMatrixRows", (DL_FUNC) &_natcpp_c_ListofMatrixRows, 1}, {"_natcpp_c_listlengths", (DL_FUNC) &_natcpp_c_listlengths, 1}, {"_natcpp_c_topntail", (DL_FUNC) &_natcpp_c_topntail, 1}, diff --git a/src/coordinates.cpp b/src/coordinates.cpp new file mode 100644 index 0000000..42e16d2 --- /dev/null +++ b/src/coordinates.cpp @@ -0,0 +1,226 @@ +// [[Rcpp::depends(Rcpp)]] +#include +using namespace Rcpp; + +// Define a lightweight accessor template that provides (i,j) access + +// Generic accessor for NumericMatrix (no bounds checking needed) +template +struct PtAccessor { + int ncol, nrow; + const T& xyz; + + PtAccessor(const T& xyz_) : xyz(xyz_) { + nrow = xyz_.nrow(); + ncol = xyz_.ncol(); + } + inline double operator()(int i, int j) const { + return xyz(i, j); // NumericMatrix handles row & column bounds + } +}; + +// Specialization for DataFrame +template <> +struct PtAccessor { + int ncol, nrow; + std::vector cols; + + PtAccessor(const DataFrame& df) { + ncol = df.size(); + nrow = df.nrows(); + cols.reserve(ncol); + for (int j = 0; j < ncol; ++j) + cols.push_back(as(df[j])); // shallow copy of each column + } + + inline double operator()(int i, int j) const { + if (j < 0 || j >= ncol) + Rcpp::stop("Column index out of range"); + return cols[j][i]; // row bounds handled by NumericVector + } +}; + + +// Core computation, templated on accessor type +template +IntegerMatrix ijkpos_core(const PtAccessor& xyz, + const NumericVector& origin, + const NumericVector& voxdims, + const IntegerVector& dims, + bool clamp) { + int n = xyz.nrow; + int p = xyz.ncol; + + if (p != 3 || origin.size() != 3 || voxdims.size() != 3 || dims.size() != 3) + stop("All coordinate-related vectors/matrices must have length or ncol = 3"); + + IntegerMatrix ijk(n, 3); + bool warn = false; + + for (int i = 0; i < n; ++i) { + for (int j = 0; j < 3; ++j) { + double val = (xyz(i, j) - origin[j]) / voxdims[j]; + int idx = static_cast(std::round(val)) + 1; + + if (idx < 1 || idx > dims[j]) { + warn = true; + if (clamp) { + if (idx < 1) idx = 1; + else if (idx > dims[j]) idx = dims[j]; + } + } + + ijk(i, j) = idx; + } + } + + if (warn) + Rcpp::warning("pixel coordinates outside image data"); + + return ijk; +} + + +//' Convert physical coordinates to pixel coordinates +//' +//' @param xyz Nx3 matrix of physical coordinates +//' @param dims Integer dimensions of the 3d image array +//' @param origin Numeric: 3d coordinates of the origin +//' @param voxdims Numeric: 3 numbers describing the voxel dimensions +//' @param clamp Logical: whether or not to clamp values within the pixel +//' boundaries of the image. +//' @return Nx3 integer matrix of pixel coordinates +//' @export +// [[Rcpp::export]] +IntegerMatrix c_ijkpos(SEXP xyz, + NumericVector origin, + NumericVector voxdims, + IntegerVector dims, + bool clamp = false) { + // Dispatch on type without copying + if (is(xyz)) { + NumericMatrix mat = as(xyz); + return ijkpos_core(PtAccessor(mat), origin, voxdims, dims, clamp); + } else if (is(xyz)) { + DataFrame df = as(xyz); + return ijkpos_core(PtAccessor(df), origin, voxdims, dims, clamp); + } else { + stop("xyz must be a numeric matrix or a data frame with 3 numeric columns"); + } +} + + +template +NumericVector sub2ind_core(IntegerVector dims, const T& indices) { + int n = indices.nrow(); + int d = indices.ncol(); + + if (dims.size() != d) + stop("indices must have the same number of columns as dimensions in dims"); + + // Compute stride multipliers (cumprod(c(1, dims[-length(dims)]))) + std::vector stride(d); + stride[0] = 1.0; + for (int i = 1; i < d; ++i) + stride[i] = stride[i - 1] * dims[i - 1]; + + NumericVector ndx(n); + for (int r = 0; r < n; ++r) { + double idx = 1.0; // R is 1-based + for (int j = 0; j < d; ++j) { + double v = indices(r, j); + if (v < 1 || v > dims[j]) + Rcpp::warning("index out of range"); + idx += (v - 1.0) * stride[j]; + } + ndx[r] = idx; + } + + return ndx; +} + +//' Find 1D index given n-dimensional indices +//' @param dims Integer dimensions of the array (usually 3d) +//' @param indices Nx3 integer matrix of pixel coordinates +//' @return numeric vector of linear indices into the array +//' @export +// [[Rcpp::export]] +NumericVector c_sub2ind(IntegerVector dims, IntegerMatrix indices) { + return sub2ind_core(dims, indices); +} + +NumericVector c_sub2ind(IntegerVector dims, NumericMatrix indices) { + return sub2ind_core(dims, indices); +} + +template +NumericVector coords21dindex_core(const PtAccessor& xyz, + const NumericVector& origin, + const NumericVector& voxdims, + const IntegerVector& dims, + bool clamp = false) { + int n = xyz.nrow; + int d = xyz.ncol; + if (origin.size() != d || voxdims.size() != d || dims.size() != d) + Rcpp::stop("origin, voxdims, and dims must have same length as number of columns in xyz"); + + // Stride for linear indexing. Note R has column-major orientation + std::vector stride(d); + stride[0] = 1.0; + for (int j = 1; j < d; ++j) + stride[j] = stride[j-1] * dims[j-1]; + + NumericVector linidx(n); + bool warn_flag = false; + + for (int i = 0; i < n; ++i) { + double idx = 1.0; // 1-based R index + for (int j = 0; j < d; ++j) { + double pix = (xyz(i,j) - origin[j]) / voxdims[j] + 1.0; + int pixidx = static_cast(std::round(pix)); + + if (pixidx < 1 || pixidx > dims[j]) { + warn_flag = true; + if (clamp) { + if (pixidx < 1) pixidx = 1; + else if (pixidx > dims[j]) pixidx = dims[j]; + } + } + + idx += (pixidx - 1.0) * stride[j]; + } + linidx[i] = idx; + } + + if (warn_flag) + Rcpp::warning("pixel coordinates outside image data"); + + return linidx; +} + +//' Convert physical coordinates to 1d indices into image array +//' +//' @param xyz Nx3 matrix or data.frame of physical coordinates +//' @param dims Integer dimensions of the 3d image array +//' @param origin Numeric: 3d coordinates of the origin +//' @param voxdims Numeric: 3 numbers describing the voxel dimensions +//' @param clamp Logical: whether or not to clamp values within the pixel +//' boundaries of the image. +//' @return Nx3 integer matrix of pixel coordinates +//' @export +// [[Rcpp::export]] +NumericVector c_coords21dindex(SEXP xyz, + NumericVector origin, + NumericVector voxdims, + IntegerVector dims, + bool clamp = false) { + if (Rcpp::is(xyz)) { + NumericMatrix mat = as(xyz); + return coords21dindex_core(PtAccessor(mat), origin, voxdims, dims, clamp); + } else if (Rcpp::is(xyz)) { + DataFrame df = as(xyz); + return coords21dindex_core(PtAccessor(df), origin, voxdims, dims, clamp); + } else { + Rcpp::stop("xyz must be a numeric matrix or a data frame"); + } +}