Zonal function - TaskFailedException #438
Comments
|
Please put the file downloads in the script so its fully self contained and I can copy and paste snd run it. I will not have time to find those files. |
|
Please find the file here: using Rasters
using CairoMakie
using ArchGDAL
using DataFrames
using Statistics
# Slope export folder
slope_folder = "main_path"
# Slope VRT
slope_vrt = "index_slope.vrt"
# Load slope
slope = Raster(joinpath(slope_folder, slope_vrt), crs = EPSG(2154))
# Load RPG
dataset = ArchGDAL.read(joinpath(slope_folder, "rpg.gpkg"))
rpg = DataFrames.DataFrame(ArchGDAL.getlayer(dataset, 0))
rename!(rpg, :geom => :geometry)
rpg_slope_mean_1 = zonal(mean, slope; of=rpg, boundary=:touches)
rpg_slope_mean_2 = [zonal(mean, slope; of=rpg.geometry[i], boundary=:touches) for i in 1:size(rpg,1)]Thanks |
|
Can you put the file download in the script and run it from start to finish and include its output. So you can verify befire I try it. A full MWE makes my life so much easier |
|
Sorry, my bad... |
|
I'm not able to load your file either with Rasters.jl or ArchGDAL.jl directly, something is broken with it. Its actually really important to follow all the instructions in the template for me to be able to solve your problem, they are very specific because this happens all the time. So please:
Then we can be fairly confident its reproducible. |
|
I apologize for my misuse of the GitHub issues template. using Downloads
using Rasters
using ArchGDAL
using DataFrames
using Statistics
# Load slope rasters
println("load raster")
slope_url = "https://www.dropbox.com/s/mu7so62hqde868m/BDALTIV2_25M_FXX_0825_6525_MNT_LAMB93_IGN69_slope.tif?dl=0";
slope_file = Downloads.download(slope_url);
slope = Raster(slope_file, crs = EPSG(2154));
# Load RPG shapes
println("load vector")
rpg_url = "https://www.dropbox.com/s/1narkro2bdt12vn/rpg_1.gpkg?dl=0";
rpg_file = Downloads.download(rpg_url);
dataset = ArchGDAL.read(rpg_file);
rpg = DataFrames.DataFrame(ArchGDAL.getlayer(dataset, 0));
rename!(rpg, :geom => :geometry);
# Get zonal statistics
println("Stat zonal")
rpg_slope_mean = zonal(mean, slope; of=rpg, boundary=:touches)Note: the same script on locale files was previously failing but now works too (problem of fresh Julia ?). But, with the initial script (.vrt instead of .tif or problem size ?), the issue remains. I am still struggling to find a way to share the .vrt due to the nature of the type of file and the low internet bandwidth I have now. |
|
Thanks for the update, and no worries at all, at least it demonstrates the reproducibility problem now. |
|
So, here are a version with bigger files that fail... using Downloads
using Rasters
using ArchGDAL
using DataFrames
using Statistics
# Load slope rasters
println("load raster")
slope_url = "https://www.dropbox.com/s/xu8232h5cuym3us/slope_big.tif?dl=0";
slope_file = Downloads.download(slope_url);
slope = Raster(slope_file, crs = EPSG(2154));
# Load RPG shapes
println("load vector")
rpg_url = "https://www.dropbox.com/s/bl96xbez3ppwfg2/rpg.gpkg?dl=0";
rpg_file = Downloads.download(rpg_url);
dataset = ArchGDAL.read(rpg_file);
rpg = DataFrames.DataFrame(ArchGDAL.getlayer(dataset, 0));
rename!(rpg, :geom => :geometry);
# Get zonal statistics
println("Stat zonal")
# Version 1
rpg_slope_mean_1 = zonal(mean, slope; of=rpg, boundary=:touches)
# Version 2
rpg_slope_mean_2 = [zonal(mean, slope; of=rpg.geometry[i], boundary=:touches) for i in 1:size(rpg,1)]And here are the output |
|
For me it seems that you have missing values in your data. Could you use x->mean(skipmissing(x)) as your inner function? |
|
Well, it was I thought at the begining, but it seems not the problem... using Downloads
using Rasters
using ArchGDAL
using DataFrames
using Statistics
# Load slope rasters
println("load raster")
slope_url = "https://www.dropbox.com/s/xu8232h5cuym3us/slope_big.tif?dl=0";
slope_file = Downloads.download(slope_url);
slope = Raster(slope_file, crs = EPSG(2154));
# Load RPG shapes
println("load vector")
rpg_url = "https://www.dropbox.com/s/bl96xbez3ppwfg2/rpg.gpkg?dl=0";
rpg_file = Downloads.download(rpg_url);
dataset = ArchGDAL.read(rpg_file);
rpg = DataFrames.DataFrame(ArchGDAL.getlayer(dataset, 0));
rename!(rpg, :geom => :geometry);
# Get zonal statistics
println("Stat zonal")
# Version 1
rpg_slope_mean_1 = zonal(x->mean(skipmissing(x)), slope; of=rpg, boundary=:touches)Output: |
|
Ok I can reproduce this one! Thanks for writing that up. I'll see what I can do. |
|
Solved. Should be registered in the hour. (it was missing values coming from polygons outside the rasters extent, which I had forgotten to put in the tests) Unfortunately ArchGDAL points are super slow (yeesian/ArchGDAL.jl#369) most of the time |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
I am using
zonalto get the mean of slopes for agricultural plots (data from https://www.ign.fr/)The zonal function returns a conversion error when passing the whole shape of plots:
TaskFailedException
nested task error: MethodError: Cannot
convertan object of type Missing to an object of type Float32It succeeds for a subset of rows below an indexes under 800, but any subset above fails.
While there is no issue when applying the
zonalto features individually.The issue should not be then related to the raster or shape.
For now, I am using the second solution, but it seems quite slower.
Side question: Is it possible to pass several functions directly to the
zonal(min, max, mean, for example), or is it more efficient to reapply several times thezonal?Thanks.
The text was updated successfully, but these errors were encountered: