From 7ba182c5d7f07f79e607f218f170e3e5c7a5e1a8 Mon Sep 17 00:00:00 2001
From: rhijmans <r.hijmans@gmail.com>
Date: Fri, 26 Apr 2024 01:52:54 -0700
Subject: [PATCH] for #1487

---
 R/cells.R              | 15 ++++++++++++---
 man/cells.Rd           |  3 ++-
 src/raster_methods.cpp | 20 ++++++++++++++++----
 src/spatRaster.h       |  2 +-
 4 files changed, 31 insertions(+), 9 deletions(-)

diff --git a/R/cells.R b/R/cells.R
index bd8f83939..150f600ff 100644
--- a/R/cells.R
+++ b/R/cells.R
@@ -17,11 +17,20 @@ setMethod("cells", signature(x="SpatRaster", y="missing"),
 )
 
 setMethod("cells", signature(x="SpatRaster", y="numeric"),
-	function(x, y) {
+	function(x, y, pairs=FALSE) {
 		opt <- spatOptions()
-		v <- x@ptr$is_in_cells(y, opt)
+		v <- x@ptr$is_in_cells(y, pairs, opt)
 		x <- messages(x, "cells")
-		v <- lapply(v, function(i) i+1)
+		if (pairs) {
+			v <- lapply(v, function(i) {
+				m <- matrix(i, ncol=2)
+				m[,1] <- m[,1] + 1
+				colnames(m) <- c("cell", "value")
+				m
+			})
+		} else {
+			v <- lapply(v, function(i) i+1)
+		}
 		names(v) <- names(x)
 		v
 	}
diff --git a/man/cells.Rd b/man/cells.Rd
index 788cb12b7..a10e290e2 100644
--- a/man/cells.Rd
+++ b/man/cells.Rd
@@ -18,7 +18,7 @@ Get the cell numbers covered by a SpatVector or SpatExtent. Or that match values
 \usage{
 \S4method{cells}{SpatRaster,missing}(x, y)
 
-\S4method{cells}{SpatRaster,numeric}(x, y)
+\S4method{cells}{SpatRaster,numeric}(x, y, pairs=FALSE)
 
 \S4method{cells}{SpatRaster,SpatVector}(x, y, method="simple", weights=FALSE, exact=FALSE, 
 				touches=is.lines(y), small=TRUE)
@@ -32,6 +32,7 @@ Get the cell numbers covered by a SpatVector or SpatExtent. Or that match values
   \item{y}{SpatVector, SpatExtent, 2-column matrix representing points, numeric representing values to match, or missing}
   \item{method}{character. Method for getting cell numbers for points. The default is "simple", the alternative is "bilinear". If it is "bilinear", the four nearest cells and their weights are returned}
   \item{weights}{logical. If \code{TRUE} and \code{y} has polygons, the approximate fraction of each cell that is covered is returned as well}
+  \item{pairs}{logical. If \code{TRUE} the cell values matched area also returned}
   \item{exact}{logical. If \code{TRUE} and \code{y} has polygons, the exact fraction of each cell that is covered is returned as well}
   \item{touches}{logical. If \code{TRUE}, values for all cells touched by lines or polygons are extracted, not just those on the line render path, or whose center point is within the polygon. Not relevant for points}
   \item{small}{logical. If \code{TRUE}, values for all cells in touched polygons are extracted if none of the cells center points is within the polygon; even if \code{touches=FALSE}}
diff --git a/src/raster_methods.cpp b/src/raster_methods.cpp
index dd8d0f500..1741f1107 100644
--- a/src/raster_methods.cpp
+++ b/src/raster_methods.cpp
@@ -861,9 +861,10 @@ SpatRaster SpatRaster::is_in(std::vector<double> m, SpatOptions &opt) {
 
 
 
-std::vector<std::vector<double>> SpatRaster::is_in_cells(std::vector<double> m, SpatOptions &opt) {
+std::vector<std::vector<double>> SpatRaster::is_in_cells(std::vector<double> m, bool keepvalue, SpatOptions &opt) {
 
 	std::vector<std::vector<double>> out(nlyr());
+	std::vector<std::vector<double>> outval(nlyr());
 
 	if (m.empty()) {
 		return(out);
@@ -894,14 +895,20 @@ std::vector<std::vector<double>> SpatRaster::is_in_cells(std::vector<double> m,
 		readBlock(v, bs, i);
 		size_t cellperlayer = bs.nrows[i] * nc;
 		for (size_t j=0; j<v.size(); j++) {
-			size_t lyr = j / cellperlayer;
-			size_t cell = j % cellperlayer + bs.row[i] * nc;
 			if (std::isnan(v[j])) {
-				if (hasNAN)	out[lyr].push_back(cell);
+				if (hasNAN)	{
+					size_t cell = j % cellperlayer + bs.row[i] * nc;
+					size_t lyr = j / cellperlayer;
+					out[lyr].push_back(cell);
+					if (keepvalue) outval[lyr].push_back(NAN);
+				}
 			} else {
 				for (size_t k=0; k<m.size(); k++) {
 					if (v[j] == m[k]) {
+						size_t cell = j % cellperlayer + bs.row[i] * nc;
+						size_t lyr = j / cellperlayer;
 						out[lyr].push_back(cell);
+						if (keepvalue) outval[lyr].push_back(v[j]);
 						break;
 					}
 				}
@@ -909,6 +916,11 @@ std::vector<std::vector<double>> SpatRaster::is_in_cells(std::vector<double> m,
 		}
 	}
 	readStop();
+	if (keepvalue) {
+		for (size_t i=0; i<nlyr(); i++) {
+			out[i].insert( out[i].end(), outval[i].begin(), outval[i].end() );
+		}
+	}
 	return(out);
 }
 
diff --git a/src/spatRaster.h b/src/spatRaster.h
index 563fedab1..77dccf27d 100644
--- a/src/spatRaster.h
+++ b/src/spatRaster.h
@@ -698,7 +698,7 @@ class SpatRaster {
 		SpatRaster init(std::vector<double> values, SpatOptions &opt);
 		
 		SpatRaster is_in(std::vector<double> m, SpatOptions &opt);
-		std::vector<std::vector<double>> is_in_cells(std::vector<double> m, SpatOptions &opt);
+		std::vector<std::vector<double>> is_in_cells(std::vector<double> m, bool keepvalue, SpatOptions &opt);
 
 		std::vector<std::string> getDataType(bool unique, bool memtype);
 		std::vector<std::string> dataType();