Dimension values disappear again?

[kapple19]

Sorry for finding another example of dimension values disappearing (previously #868).

julia> using DimensionalData

julia> x = X(1:10)
X 1:10

julia> y = Y(2:0.1:4)
Y 2.0:0.1:4.0

julia> st1 = DimStack((; A = rand(x, y), B = rand(x, Z([1]))))
┌ 10×21×1 DimStack ┐
├──────────────────┴───────────────────────────────────────────────────────────────────────────────────────────── dims ┐  ↓ X Sampled{Int64} 1:10 ForwardOrdered Regular Points,
  → Y Sampled{Float64} 2.0:0.1:4.0 ForwardOrdered Regular Points,
  ↗ Z Sampled{Int64} [1] ForwardOrdered Irregular Points
├────────────────────────────────────────────────────────────────────────────────────────────────────────────── layers ┤  :A eltype: Float64 dims: X, Y size: 10×21
  :B eltype: Float64 dims: X, Z size: 10×1
└──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┘
julia> st2 = DimStack((; A = rand(x, y), B = rand(x, Z([2]))))
┌ 10×21×1 DimStack ┐
├──────────────────┴───────────────────────────────────────────────────────────────────────────────────────────── dims ┐  ↓ X Sampled{Int64} 1:10 ForwardOrdered Regular Points,
  → Y Sampled{Float64} 2.0:0.1:4.0 ForwardOrdered Regular Points,
  ↗ Z Sampled{Int64} [2] ForwardOrdered Irregular Points
├────────────────────────────────────────────────────────────────────────────────────────────────────────────── layers ┤  :A eltype: Float64 dims: X, Y size: 10×21
  :B eltype: Float64 dims: X, Z size: 10×1
└──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┘
julia> cat(st1, st2; dims = Z)
┌ 10×21×2 DimStack ┐
├──────────────────┴───────────────────────────────────────────────────────────────────────────────────────────── dims ┐  ↓ X Sampled{Int64} 1:10 ForwardOrdered Regular Points,
  → Y Sampled{Float64} 2.0:0.1:4.0 ForwardOrdered Regular Points,
  ↗ Z
├────────────────────────────────────────────────────────────────────────────────────────────────────────────── layers ┤  :A eltype: Float64 dims: X, Y, Z size: 10×21×2
  :B eltype: Float64 dims: X, Z size: 10×2
└──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┘
(jl_9Nf4Wv) pkg> st
Status `C:\Users\T0185914\AppData\Local\Temp\jl_9Nf4Wv\Project.toml`
  [0703355e] DimensionalData v0.29.16

[rafaqz]

Thats because your Z would be broken otherwise, with the same number repeated

[kapple19]

Do I define the values of Z wrong?
It's 1 in st1, and 2 in st2.

[rafaqz]

Oh right yes I misread, that is a bug

[tiemvanderdeure]

with 2 DimArrays this does work

julia> cat(st1.B, st2.B; dims = Z)
┌ 10×2 DimArray{Float64, 2} B ┐
├─────────────────────────────┴──────────────────────── dims ┐
  ↓ X Sampled{Int64} 1:10 ForwardOrdered Regular Points,
  → Z Sampled{Int64} [1, 2] ForwardOrdered Irregular Points
└────────────────────────────────────────────────────────────┘
  ↓ →  1          2
  1    0.242448   0.857127
  2    0.291349   0.234988
  3    0.865921   0.955676
  4    0.856839   0.410151
  ⋮               
  7    0.667107   0.91918
  8    0.0467946  0.518208
  9    0.137003   0.361783
 10    0.853074   0.208523

And even if DimArrays are misalligned this at least throws a warning:

 julia> cat(st2.B, st1.B; dims = Z)
┌ Warning: `cat` cannot concatenate `Dimension`s, falling back to `parent` type:
│ `Ordered` lookups are misaligned at 1 and 1 on dimension Z. 
│ 
│ To fix for `AbstractDimArray`, pass new lookup values as `cat(As...; dims=Z(newlookupvals))` keyword or `dims=Z()` for empty `NoLookup`.
└ @ DimensionalData ~/.julia/dev/DimensionalData/src/array/methods.jl:545
┌ 10×2 DimArray{Float64, 2} B ┐
├─────────────────────────────┴───────────────────── dims ┐
  ↓ X Sampled{Int64} 1:10 ForwardOrdered Regular Points,
  → Z
└─────────────────────────────────────────────────────────┘
  1  0.857127  0.242448
  2  0.234988  0.291349
  3  0.955676  0.865921
  4  0.410151  0.856839
  5  0.519696  0.83716
  6  0.44849   0.292439
  7  0.91918   0.667107
  8  0.518208  0.0467946
  9  0.361783  0.137003
 10  0.208523  0.853074

I guess there are some functions that the aren't hit when cat is called on DimStacks

[tiemvanderdeure]

Actually it's just because the A layer in both stacks don't have a Z dimension. So cat(st1.A, st2.A; dims = Z) returns this empty Z dim, which the stack then inherits because it is the first Z dim it encounters.

We could explicitly call _vcat_dims on the dimensions of the stack and pass that to rebuild_from_arrays. I guess it would add a tiny overhead because we concat the dimensions for each layers also, but maybe that's worth it?

[rafaqz]

Yeah that's better, we can remove the overhead by using that finished dimension to pass to cat on each array

[tiemvanderdeure]

Nice, that will even (ever so slightly) improve performance

[tiemvanderdeure]

Realized it's not so easy because in a case like this we can't really do cat at all:

x = X(1:10)
y = Y(2:0.1:4)
st1 = DimStack((; A = rand(x, y), B = rand(x, Z([1]))))
st2 = DimStack((; A = rand(x, y), B = rand(x, Z([2,3]))))
st3 = cat(st1, st2; dims = Z)

This errors because the sizes of the concatenated layers don't match.

So this is a special case that can only occur with dimensions of length 1.

[rafaqz]

Ahhh I get it you can concatenate a dimension that is in some layers but not in others, unless its of length one.

Even doing it with length one is a bit dodgy because its saying A automatically takes on the Z value of B where actually A can't really be seen to be sliced from or belonging to the Z axis when B had a Z axis in the same stack.

Probably we should error in all cases, or instead check that st1.A == st2.A and keep one, only concatenating layers with Z, or when no layers have Z?

[tiemvanderdeure]

Probably we should error in all cases

I'm not sure about that. cat(st1.A, st2.A; dims = Z) doesn't error even though there is no Z dim. It just gives a warning and a NoLookup.

It think the behaviour we get from my PR is reasonable. It makes this niche thing work, arguably the code is even a tiny bit cleaner, and the error messages are going to be the same

[rafaqz]

But... cat on an array has no ambiguity. Here we are actually combining two quite different modes if cat in one operation in the same dimension. One is adding a new dimension Z, the other growing an existing dimension Z. I think its actually undefined behaviour to do both at once with one object

[tiemvanderdeure]

So maybe it should just error?

If you think of DimStacks as similar to DimArrays of namedtuples then I'm not sure it's ambigious. Like cat(DimArray(st1), DimArray(st2), dims=Z) also has the Z dimension with values [1,2].

Currently the behaviour is a bit weird because cat(st1[(:B, :A)], st2[(:B, :A)], dims=Z) does retain Z values.

[rafaqz]

Yeah, what's weird is if there is already a Z.

Think about it like the cat on each layer will have different dimensions to cat on the stack then selecting the layers by name.

Probably making a stack out of the separate arrays method will error because Z won't match

[tiemvanderdeure]

Probably making a stack out of the separate arrays method will error because Z won't match

It doesn't error and actually it's because of something I tried to address in #879.

Basically the nolookup and the values pass comparedims. What happens is basically

using DimensionalData, DimensionalData.Lookups
import DimensionalData: format, comparedims

a = format(Z([1,2]))
b = format(Z(NoLookup(Base.OneTo(2))))
comparedims(a, b)

this doesn't error. comparedims(a, b; valtype = true) does