2

I am calculating the volume in a raster .tif file.

My current approach is area * height where I calculate the area in the .tif file and multiply it by elevation value of .tif file.

I have used Raster Surface Volume algorithm of QGIS, Raster Volume of SAGA GIS and other techniques. But I want to use only Rasterio package.

Is there any in-built functionality of Rasterio which calculates the volume of a raster .tif file ?

Improve this question
0

1 Answer 1

Reset to default
1

Once you open a raster with rasterio basically you have a numpy array. The area of each pixel is fixed, so you can just multiply the array for the area value, and it will multiply elementwise, and then just sum all the values for the total volume.

In code:

with rasterio.open(filename) as src:
    array = src.read()
    meta = src.meta
    area = meta['transform'][0] ** 2 # Pixel size in X
volumes = array * area
total_volume = volumes.sum()
Improve this answer
6
  • The correct answer for total_volume is a value between 9500 m^3 and 10000 m^3. Unfortunately, according to your code, I am getting a value 1270057.54 which is actually wrong !!!
    – Sreeraj
    Nov 27, 2019 at 13:53
  • Are all units the same?
    – ImanolUr
    Nov 27, 2019 at 14:22
  • ImanolUr If I consider the untis, then according to your code, I think that array is in cm and area is in m^2. If I am right, I got the area = 0.00033562239.... m^2 and total_array = array.sum() = 3784183645.836426 cm = 37841836.45836426 m. So total_volume = area in m^2 * array in m = total_volume in m^3 = 0.00033562239 m^2 * 37841836.45836426 = 12700.567972543868 m^3. So, I got the volume = 12700 m^3, but it is wrong. I have to get a value between 9500 m^3 to 10000 m^3 (maximum compromisable value is around 10300 m^3. Anything above 10400 m^3 is surely incorrect).
    – Sreeraj
    Nov 28, 2019 at 8:12
  • Please note that I am using this code with input .tif file in UTM projection (and not in latlong projection). latlong projection is giving the totally incorrect values for volumes in all algorithms of QGIS, SAGA GIS, and also as per your code because of the difference in units for the latlong projection. Since, all units of input .tif file in UTM projection is in terms of meters, it is fine. But, I am getting the volume 12700 m^3, which is around 2700 m^3 more than the actual value of around 10000 m^3 as per your code, even though the logic of your code seems to be right.
    – Sreeraj
    Nov 28, 2019 at 8:45
  • What is actually your raster, and how do you know this limit of 10400m³?
    – ImanolUr
    Nov 28, 2019 at 14:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.