[R] Volume of polygon

Jeff Newmiller jdnewm|| @end|ng |rom dcn@d@v|@@c@@u@
Wed Dec 6 06:38:56 CET 2023


A raster is just a matrix of elevations. Each element of that matrix has a horizontal area. If you subtract that elevation from a reference elevation and zero out all negative values then you are left with a bunch of rectangular parallelopipeds of varying height. Add those heights up and multiply by the scalar area of the parallelopipeds and you have the volume. Repeat this for various elevations and form a set of estimates of volume vs elevation. You can use approxfun or splinefun to create a function from those points for more general use.

There are some subtleties such as converting between elevation and depth, and whether you are concerned about disjoint depressions (they would not drain)... but OP really should have read the posting guide and asked this question on R-sig-geo mailing list.

On December 5, 2023 9:09:12 PM PST, Bert Gunter <bgunter.4567 using gmail.com> wrote:
>The volume of a polygon  = 0.  Polyhedra  have volumes.
>
>This may be irrelevant, but if the lake is cylindrical == constant cross
>sectional area at all depths, then height doubles when the volume does and
>vice versa.  Otherwise you have to know how area varies with height or use
>more sensible approximations thereto.
>
>Cheers,
>Bert
>
>On Tue, Dec 5, 2023, 20:13 javad bayat <j.bayat194 using gmail.com> wrote:
>
>>  Dear all;
>> I am trying to calculate the volume of a polygon shapefile according to a
>> DEM raster. I have provided some codes at the end of this email.I dont know
>> if the codes are correct or not. Following this, I have another question
>> too.
>> I want to know if the volume of the reservoir rises or doubles, what would
>> be the elevation?
>> I would be more than happy if anyone could help me.
>> Sincerely
>>
>> "
>> library(raster)
>> library(terra)
>> library(exactextractr)
>> library(dplyr)
>> library(sf)
>> r <- raster("Base.tif")
>> p <- shapefile("p.shp")
>> r <- crop(r, p)
>> r <- mask(r, p)
>> x <- exact_extract(r, p, coverage_area = TRUE)
>>
>> x1 = as.data.frame(x[1])
>> head(x1)
>> x1 = na.omit(x1)
>>
>> x1$Height = max(x1[,1]) - x1[,1]
>>
>> x1$Vol = x1[,2] * x1[,3]
>>
>> sum(x1$Vol)
>>
>> "
>>
>> --
>> Best Regards
>> Javad Bayat
>> M.Sc. Environment Engineering
>> Alternative Mail: bayat194 using yahoo.com
>>
>>         [[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-help using 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.
>>
>
>	[[alternative HTML version deleted]]
>
>______________________________________________
>R-help using 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.

-- 
Sent from my phone. Please excuse my brevity.



More information about the R-help mailing list