Hi:
Some preliminary comments first
1- You might get a better response posting to R-sig-Geo
2- Probably better to switch to the `terra` package instead of the older
`raster`
3- Furthermore, there is the package `exactextractr` that ensures that
edge pixels are taken into account, weighted by the actual, partial
coverage of each pixel by the polygon.
4- Can you post: your R version, and version of the relevant packages,
and information about the polygons and DEM (i.e. `str(p)` and `str(r)` )
On 20/11/2023 14:30, javad bayat wrote:
Dear all;
I am trying to calculate volume under each polygon of a shapefile according
to a DEM.
when I run the code, it gives me an error as follows.
"
Error in h(simpleError(msg, call)) :
error in evaluating the argument 'x' in selecting a method for function
'addAttrToGeom': sp supports Z dimension only for POINT and MULTIPOINT.
use `st_zm(...)` to coerce to XY dimensions
"
I want to have a table that contains one column corresponding to the
polygon rank and another that contains the volume on that polyon.
I would be more than happy if anyone could help me.
I provided codes at the end of this email.
Sincerely
##########################################################################################
library(raster);
library(sf)
# Load the DEM raster and shapefile
r <- raster("E:/Base1.tif")
p <- read_sf(dsn = "E:/Sites.shp", layer = " Sites")
Now to the question:
In `read_sf()`, when the source is a shapefile, you do not need the
'layer' option. It doesn't hurt, but in this case I see that you have an
extra space before the layer name. Could that be causing the problem?
# Extract the values of the DEM raster for each polygon
values <- extract(r, p)
Error in h(simpleError(msg, call)) :
error in evaluating the argument 'x' in selecting a method for function
'addAttrToGeom': sp supports Z dimension only for POINT and MULTIPOINT.
use `st_zm(...)` to coerce to XY dimensions
# Calculate the volume of each polygon
volumes <- sapply(values, function(x) rasterVolume(x, r))
What is this "rasterVolume" function that you call here?
# Print the results
for (i in 1:length(volumes)) {
cat(sprintf("Volume under polygon %d: %f\n", i, volumes[i]))
}
##########################################################################################
--
Micha Silver
Ben Gurion Univ.
Sde Boker, Remote Sensing Lab
cell: +972-523-665918
______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
