JuliaDataCubes / YAXArrays.jl

Yet Another XArray-like Julia package
https://juliadatacubes.github.io/YAXArrays.jl/
Other
103 stars 18 forks source link

Batchindex with a Table can give wrong samples #194

Closed felixcremer closed 1 year ago

felixcremer commented 1 year ago

There is an axis mismatch in the batchextract function so that you can get an edge pixel of the data cube multiple times. I have this function to sample data from a data cube inside a given shapefile.

function maskcube(cube, shp)
    bbox = Extent((X=extrema(cube.Lon), Y=extrema(cube.Lat)))
    res = step(cube.Lon.values)
    shpras = rasterize(shp, to=bbox, res=res, fill=1)
    m = yaxconvert(YAXArray, DimArray(shpras, metadata=Dict("shape" => shp)))
    renameaxis!(m, :X=>"Lon")
    renameaxis!(m, :Y=>"Lat")
    @show bbox
    ct = CubeTable(Sample=m)
    ctfiltered = TableOperations.filter.(x -> !ismissing(x.Sample), ct)
    masked = cube[ctfiltered[1]]
    renameaxis!(masked, "Sample"=> CategoricalAxis("Sample", 1:length(masked.Sample.values)))

end

The result is a data cube with the right amount of samples, but the samples have all the same values. This happens because of a mismatch between the tcols and axinds in the batchextract function here: https://github.com/JuliaDataCubes/YAXArrays.jl/blob/9865671d6bd49687915b67a8022e1b325fbd9e41/src/Cubes/Cubes.jl#L253 We would need to make sure, that the axinds and tcols variables have the same length.

felixcremer commented 1 year ago

I can circumvent the problem for now, when I permute my indexing table to make sure, that the Sample entry is at the end.

felixcremer commented 1 year ago

I managed to kind of solve it when I add the following line into the if case that checks for nothingness in the axinds. The idea is to also filter the tcols, that are not in the cube.

        tcols = (;[p[1:2] for p in zip(keys(tcols), values(tcols), axinds) if !isnothing(last(p))]...)

I am going to make a PR from this next week, and to do proper testing on that as well.