search

JuliaDataCubes / YAXArrays.jl

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

mapCube with two input cubes and shared dimensions fails #256

Open gdkrmr opened 1 year ago

gdkrmr commented 1 year ago

I have two input cubes with a shared "Variable" dimension. In my understanding calling mapCube with both cubes should repeat over the "Variable" dimension of both cubes. This works in the second toy example but fails with an error message in the first one, even though 

julia> getAxis("variable", c) == getAxis("variable", c_stat)
true

I am not sure what is happening here.

using YAXArrays
using Zarr

c_path = "http://data.rsc4earth.de:9000/earthsystemdatacube/v3.0.2/esdc-8d-0.25deg-1x720x1440-3.0.2.zarr"
c_zarr = Zarr.zopen(c_path)
c_dataset = YAXArrays.open_dataset(c_zarr)
c_full = YAXArrays.Cube(c_dataset)
c = c_full[variable=["sensible_heat", "air_temperature_2m", "radiation_era5"]]

m2 = rand(3)
s2 = rand(3)

v_ax = getAxis("variable", c)
stat_ax = CategoricalAxis("statistic", ["mean", "std"])
c_stat = YAXArray([v_ax, stat_ax], [m2 s2])
c

getAxis("variable", c) == getAxis("variable", c_stat)

mapCube(
    (c, c_stat),
    indims=(InDims(), InDims("statistic")),
    outdims=OutDims()
) do xout, xin1, xin2
    xout[1] = (xin1[1] - xin2[1]) / xin2[2]
end
Error Message ```julia julia> getAxis("variable", c) == getAxis("variable", c_stat) true julia> mapCube( (c, c_stat), indims=(InDims(), InDims("statistic")), outdims=OutDims() ) do xout, xin1, xin2 xout[1] = (xin1[1] - xin2[1]) / xin2[2] end ERROR: TaskFailedException nested task error: DimensionMismatch: array could not be broadcast to match destination Stacktrace: [1] check_broadcast_shape @ ./broadcast.jl:553 [inlined] [2] check_broadcast_axes @ ./broadcast.jl:556 [inlined] [3] instantiate @ ./broadcast.jl:297 [inlined] [4] materialize! @ ./broadcast.jl:884 [inlined] [5] materialize! @ ./broadcast.jl:881 [inlined] [6] (::YAXArrays.DAT.var"#161#166"{CartesianIndex{4}})(iw::Vector{Float64}, x::YAXArrays.YAXTools.PickAxisArray{Float64, 5, Matrix{Float64}, (Colon(), 4), nothing}) @ YAXArrays.DAT ~/.julia/packages/YAXArrays/au5n4/src/DAT/DAT.jl:1112 [7] (::Base.var"#61#62"{YAXArrays.DAT.var"#161#166"{CartesianIndex{4}}})(#unused#::Nothing, xs::Tuple{Vector{Float64}, YAXArrays.YAXTools.PickAxisArray{Float64, 5, Matrix{Float64}, (Colon(), 4), nothing}}) @ Base ./tuple.jl:603 [8] BottomRF @ ./reduce.jl:81 [inlined] [9] _foldl_impl(op::Base.BottomRF{Base.var"#61#62"{YAXArrays.DAT.var"#161#166"{CartesianIndex{4}}}}, init::Nothing, itr::Base.Iterators.Zip{Tuple{Tuple{Array{Float32, 0}, Vector{Float64}}, Tuple{YAXArrays.YAXTools.PickAxisArray{Float32, 4, Array{Float32, 4}, (1, 2, 3, 4), nothing}, YAXArrays.YAXTools.PickAxisArray{Float64, 5, Matrix{Float64}, (Colon(), 4), nothing}}}} ) @ Base ./reduce.jl:62 [10] foldl_impl @ ./reduce.jl:48 [inlined] [11] mapfoldl_impl @ ./reduce.jl:44 [inlined] [12] #mapfoldl#288 @ ./reduce.jl:170 [inlined] [13] mapfoldl @ ./reduce.jl:170 [inlined] [14] #foldl#289 @ ./reduce.jl:193 [inlined] [15] foldl @ ./reduce.jl:193 [inlined] [16] foreach @ ./tuple.jl:603 [inlined] [17] innercode(f::Function, cI::CartesianIndex{4}, xinBC::Tuple{YAXArrays.YAXTools.PickAxisArray{Float32, 4, Array{Float32, 4}, (1, 2, 3, 4), nothing}, YAXArrays.YAXTools.PickAxisArray{Float64, 5, Matrix{Float64}, (Colon(), 4), nothing}}, xoutBC::Tuple{YAXArrays.YAXTools.PickAxisArray{Union{Missing, Float32}, 4, Array{Union{Missing, Float32}, 4}, (1, 2, 3, 4), nothing}}, filters::Tuple{Tuple{YAXArrays.DAT.AllMissing}, Tuple{YAXArrays.DAT.AllMissing}}, inwork::Tuple{Vector{Array{Float32, 0}}, Vector{Vector{Float64}}}, outwork::Tuple{Vector{Array{Union{Missing, Float32}, 0}}}, axvalcreator::YAXArrays.DAT.NoLoopAxes, offscur::NTuple{4, Int64}, addargs::Tuple{}, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}) @ YAXArrays.DAT ~/.julia/packages/YAXArrays/au5n4/src/DAT/DAT.jl:1111 [18] macro expansion @ ~/.julia/packages/YAXArrays/au5n4/src/DAT/DAT.jl:1162 [inlined] [19] (::YAXArrays.DAT.var"#415#threadsfor_fun#176"{YAXArrays.DAT.var"#415#threadsfor_fun#172#177"{var"#5#6", Tuple{YAXArrays.YAXTools.PickAxisArray{Float32, 4, Array{Float32, 4}, (1, 2, 3, 4), nothing}, YAXArrays.YAXTools.PickAxisArray{Float64, 5, Matrix{Float64}, (Colon(), 4), nothing}}, Tuple{YAXArrays.YAXTools.PickAxisArray{Union{Missing, Float32}, 4, Array{Union{Missing, Float32}, 4}, (1, 2, 3, 4), nothing}}, Tuple{Tuple{YAXArrays.DAT.AllMissing}, Tuple{YAXArrays.DAT.AllMissing}}, Tuple{Vector{Array{Float32, 0}}, Vector{Vector{Float64}}}, Tuple{Vector{Array{Union{Missing, Float32}, 0}}}, YAXArrays.DAT.NoLoopAxes, Tuple{}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, NTuple{4, Int64}, CartesianIndices{4, NTuple{4, UnitRange{Int64}}}}})(tid::Int64; onethread::Bool) @ YAXArrays.DAT ./threadingconstructs.jl:163 [20] #415#threadsfor_fun @ ./threadingconstructs.jl:130 [inlined] [21] (::Base.Threads.var"#1#2"{YAXArrays.DAT.var"#415#threadsfor_fun#176"{YAXArrays.DAT.var"#415#threadsfor_fun#172#177"{var"#5#6", Tuple{YAXArrays.YAXTools.PickAxisArray{Float32, 4, Array{Float32, 4}, (1, 2, 3, 4), nothing}, YAXArrays.YAXTools.PickAxisArray{Float64, 5, Matrix{Float64}, (Colon(), 4), nothing}}, Tuple{YAXArrays.YAXTools.PickAxisArray{Union{Missing, Float32}, 4, Array{Union{Missing, Float32}, 4}, (1, 2, 3, 4), nothing}}, Tuple{Tuple{YAXArrays.DAT.AllMissing}, Tuple{YAXArrays.DAT.AllMissing}}, Tuple{Vector{Array{Float32, 0}}, Vector{Vector{Float64}}}, Tuple{Vector{Array{Union{Missing, Float32}, 0}}}, YAXArrays.DAT.NoLoopAxes, Tuple{}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, NTuple{4, Int64}, CartesianIndices{4, NTuple{4, UnitRange{Int64}}}}}, Int64})() @ Base.Threads ./threadingconstructs.jl:108 ...and 7 more exceptions. Stacktrace: [1] threading_run(fun::YAXArrays.DAT.var"#415#threadsfor_fun#176"{YAXArrays.DAT.var"#415#threadsfor_fun#172#177"{var"#5#6", Tuple{YAXArrays.YAXTools.PickAxisArray{Float32, 4, Array{Float32, 4}, (1, 2, 3, 4), nothing}, YAXArrays.YAXTools.PickAxisArray{Float64, 5, Matrix{Float64}, (Colon(), 4), nothing}}, Tuple{YAXArrays.YAXTools.PickAxisArray{Union{Missing, Float32}, 4, Array{Union{Missing, Float32}, 4}, (1, 2, 3, 4), nothing}}, Tuple{Tuple{YAXArrays.DAT.AllMissing}, Tuple{YAXArrays.DAT.AllMissing}}, Tuple{Vector{Array{Float32, 0}}, Vector{Vector{Float64}}}, Tuple{Vector{Array{Union{Missing, Float32}, 0}}}, YAXArrays.DAT.NoLoopAxes, Tuple{}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, NTuple{4, Int64}, CartesianIndices{4, NTuple{4, UnitRange{Int64}}}}}, static::Bool) @ Base.Threads ./threadingconstructs.jl:120 [2] macro expansion @ ./threadingconstructs.jl:168 [inlined] [3] innerLoop(loopRanges::NTuple{4, UnitRange{Int64}}, f::Function, xinBC::Tuple{YAXArrays.YAXTools.PickAxisArray{Float32, 4, Array{Float32, 4}, (1, 2, 3, 4), nothing}, YAXArrays.YAXTools.PickAxisArray{Float64, 5, Matrix{Float64}, (Colon(), 4), nothing}}, xoutBC::Tuple{YAXArrays.YAXTools.PickAxisArray{Union{Missing, Float32}, 4, Array{Union{Missing, Float32}, 4}, (1, 2, 3, 4), nothing}}, filters::Tuple{Tuple{YAXArrays.DAT.AllMissing}, Tuple{YAXArrays.DAT.AllMissing}}, inwork::Tuple{Vector{Array{Float32, 0}}, Vector{Vector{Float64}}}, outwork::Tuple{Vector{Array{Union{Missing, Float32}, 0}}}, axvalcreator::YAXArrays.DAT.NoLoopAxes, addargs::Tuple{}, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}) @ YAXArrays.DAT ~/.julia/packages/YAXArrays/au5n4/src/DAT/DAT.jl:1161 [4] (::YAXArrays.DAT.var"#108#112"{YAXArrays.DAT.DATConfig{2, 1}})(r::NTuple{4, UnitRange{Int64}}) @ YAXArrays.DAT ~/.julia/packages/YAXArrays/au5n4/src/DAT/DAT.jl:709 [5] (::ProgressMeter.var"#56#59"{Distributed.RemoteChannel{Channel{Bool}}, YAXArrays.DAT.var"#108#112"{YAXArrays.DAT.DATConfig{2, 1}}})(x::NTuple{4, UnitRange{Int64}}) @ ProgressMeter ~/.julia/packages/ProgressMeter/sN2xr/src/ProgressMeter.jl:1016 [6] iterate @ ./generator.jl:47 [inlined] [7] _collect(c::DiskArrays.GridChunks{4}, itr::Base.Generator{DiskArrays.GridChunks{4}, ProgressMeter.var"#56#59"{Distributed.RemoteChannel{Channel{Bool}}, YAXArrays.DAT.var"#108#112"{YAXArrays.DAT.DATConfig{2, 1}}}}, #unused#::Base.EltypeUnknown, isz::Base.HasShape{4}) @ Base ./array.jl:802 [8] collect_similar(cont::DiskArrays.GridChunks{4}, itr::Base.Generator{DiskArrays.GridChunks{4}, ProgressMeter.var"#56#59"{Distributed.RemoteChannel{Channel{Bool}}, YAXArrays.DAT.var"#108#112"{YAXArrays.DAT.DATConfig{2, 1}}}}) @ Base ./array.jl:711 [9] map(f::Function, A::DiskArrays.GridChunks{4}) @ Base ./abstractarray.jl:3261 [10] macro expansion @ ~/.julia/packages/ProgressMeter/sN2xr/src/ProgressMeter.jl:1015 [inlined] [11] macro expansion @ ./task.jl:476 [inlined] [12] macro expansion @ ~/.julia/packages/ProgressMeter/sN2xr/src/ProgressMeter.jl:1014 [inlined] [13] macro expansion @ ./task.jl:476 [inlined] [14] progress_map(::Function, ::Vararg{Any}; mapfun::Function, progress::ProgressMeter.Progress, channel_bufflen::Int64, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}) @ ProgressMeter ~/.julia/packages/ProgressMeter/sN2xr/src/ProgressMeter.jl:1007 [15] progress_map(::Function, ::Vararg{Any}) @ ProgressMeter ~/.julia/packages/ProgressMeter/sN2xr/src/ProgressMeter.jl:999 [16] runLoop(dc::YAXArrays.DAT.DATConfig{2, 1}, showprog::Bool) @ YAXArrays.DAT ~/.julia/packages/YAXArrays/au5n4/src/DAT/DAT.jl:707 [17] mapCube(::var"#5#6", ::Tuple{YAXArray{Float32, 4, DiskArrayTools.DiskArrayStack{Float32, 4, DiskArrays.SubDiskArray{Float32, 3}, 1}, Vector{CubeAxis}}, YAXArray{Float64, 2, Matrix{Float64}, Vector{CategoricalAxis{String, S, Vector{String}} where S}}}; max_cache::Float64, indims::Tuple{InDims, InDims}, outdims::OutDims, inplace::Bool, ispar::Bool, debug::Bool, include_loopvars::Bool, showprog::Bool, irregular_loopranges::Bool, nthreads::Vector{Int64}, loopchunksize::Dict{Any, Any}, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}) @ YAXArrays.DAT ~/.julia/packages/YAXArrays/au5n4/src/DAT/DAT.jl:475 [18] top-level scope @ REPL[17]:1 ```

Am I doing something wrong? The following minimal example works:


using YAXArrays
using Dates

v_ax = CategoricalAxis("v", ["v1", "v2", "v3"])
t_ax = RangeAxis("t", Date(2000, 1, 1):Day(1):Date(2000, 1, 3))
s_ax = CategoricalAxis("s", ["s1", "s2"])

c1 = YAXArray([v_ax, t_ax], rand(3, 3))
c2 = YAXArray([v_ax, s_ax], rand(3, 2))

mapCube(
    (c1, c2),
    indims=(InDims(), InDims("s")),
    outdims=OutDims()
) do xout, xin1, xin2
    xout .= (xin1 .- xin2[1]) ./ xin2[2]
end

Trying to use Debugger.jl crashes the Julia session really hard https://github.com/JuliaDebug/Debugger.jl/issues/333