How to add variable to existing cube?

[bjoern234]

Hi!
I have a cube e.g. with axes X,Y,Time, and Bands. Bands contains "NIR" and "RED" variables. How would I add a "NDVI" variable to Bands for each time step? Is that possible or do I have to compute a new "ndvi-cube" and combine or merge them?

Thanks in advance
Bjoern

[felixcremer]

You have to compute a new ndvicube, but this is all done lazily and you don't have to save it to disk.

julia> arr = rand(300,200,90, 4)

julia> x = RangeAxis(:X, 1:300)
X                   Axis with 300 Elements from 1 to 300

julia> y = RangeAxis(:Y, 1:200)
Y                   Axis with 200 Elements from 1 to 200

julia> t = RangeAxis(:time, 1:90)
time                Axis with 90 Elements from 1 to 90


julia> b = CategoricalAxis("Bands", ["Red", "Green", "Blue", "NIR"])
Bands               Axis with 4 elements: Red Green Blue NIR 

julia> yax = YAXArray([x,y,t,b], arr)
YAXArray with the following dimensions
X                   Axis with 300 Elements from 1 to 300
Y                   Axis with 200 Elements from 1 to 200
time                Axis with 90 Elements from 1 to 90
Bands               Axis with 4 elements: Red Green Blue NIR 
Total size: 164.79 MB

julia> ndvi = map((r,nir) ->(r-nir)/(r+nir),yax[Bands="Red"], yax[Bands="NIR"])

If you would like to combine the datasets again into a cube together you would have to split the cube and concatenate it again.

julia> catax = CategoricalAxis("Bands", ["Red", "Green", "Blue", "NIR", "NDVI"])

julia> ndvicube = concatenatecubes([yax[Bands="Red"], yax[Bands="Green"], yax[Bands="Blue"], yax[Bands="NIR"], ndvi], catax)
YAXArray with the following dimensions
X                   Axis with 300 Elements from 1 to 300
Y                   Axis with 200 Elements from 1 to 200
time                Axis with 90 Elements from 1 to 90
Bands               Axis with 5 elements: Red Green Blue NIR NDVI 
Total size: 205.99 MB

[MarkusZehner]

Hi,
how would i include a single variable (say orbit no. or ascending/descending, one value, comparable to the time axis, but not distinct values) into the cube?

As in the example below, add the dir axis into the cube:

julia> caxes(cube) 
3-element Vector{CubeAxis}:
 Y                   Axis with 22500 Elements from -855560.0 to -405580.0
 X                   Axis with 15000 Elements from -253755.0 to -553735.0
 Time                Axis with 1244 Elements from 2017-01-01T16:59:45 to 2019-12-31T05:33:44

julia> ax_dir
dir                 Axis with 1244 elements: A A .. A D 

Thanks!
Markus

[MarkusZehner]

Hi, again to elaborate some stuff, maybe there is an easy fix for that, maybe its useful...

so my xarray from python, (where i added the metadata for now) looks something like this:

`<xarray.Dataset>
Dimensions: (X: 2, Y: 2, Time: 714)
Coordinates:

• X (X) float64 -3.914e+05 -3.914e+05
• Y (Y) float64 -6.672e+05 -6.672e+05
• Time (Time) datetime64[ns] 2017-01-01T16:59:45 ... 2019-12-31T05:...
  spatial_ref int64 0
  Data variables:
  layer (Time, X, Y) float32 -15.03 -14.98 nan ... -14.22 nan -14.2
  sat (Time) object 'A' 'A' 'B' 'A' 'B' 'B' ... 'B' 'A' 'B' 'A' 'A'
  orb (Time) float64 44.0 66.0 168.0 117.0 ... 168.0 117.0 44.0 66.0
  dir (Time) object 'ASCENDING' 'DESCENDING' ... 'DESCENDING'`

Trying to open this again in YAXArray for processing now fails like the following:

`julia> YAXArrays.open_dataset("cube_from_above.zarr")
ERROR: KeyError: key "vlen-utf8" not found
Stacktrace:

[1] getindex(h::Dict{String, UnionAll}, key::String)
@ Base ./dict.jl:482
[2] (::Zarr.var"#34#35")(f::Dict{String, Any})
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/Filters.jl:13
[3] iterate
@ ./generator.jl:47 [inlined]
[4] _collect(c::Vector{Any}, itr::Base.Generator{Vector{Any}, Zarr.var"#34#35"}, #unused#::Base.EltypeUnknown, isz::Base.HasShape{1})
@ Base ./array.jl:691
[5] collect_similar(cont::Vector{Any}, itr::Base.Generator{Vector{Any}, Zarr.var"#34#35"})
@ Base ./array.jl:606
[6] map(f::Function, A::Vector{Any})
@ Base ./abstractarray.jl:2294
[7] getfilters(d::Dict{String, Any})
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/Filters.jl:12
[8] Zarr.Metadata(d::Dict{String, Any})
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/metadata.jl:175
[9] Zarr.Metadata(s::String)
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/metadata.jl:166
[10] getmetadata(s::DirectoryStore)
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/Storage/Storage.jl:100
[11] ZArray(s::DirectoryStore, mode::String)
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/ZArray.jl:95
[12] zopen_noerr(s::DirectoryStore, mode::String; consolidated::Bool)
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/ZGroup.jl:39
[13] zopen_noerr(s::DirectoryStore, mode::String)
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/ZGroup.jl:37
[14] ZGroup(s::DirectoryStore, mode::String)
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/ZGroup.jl:18
[15] zopen_noerr(s::DirectoryStore, mode::String; consolidated::Bool)
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/ZGroup.jl:41
[16] #zopen#75
@ ~/.julia/packages/Zarr/2V5d9/src/ZGroup.jl:76 [inlined]
[17] zopen(s::DirectoryStore, mode::String)
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/ZGroup.jl:76
[18] zopen(s::String, mode::String; kwargs::Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
@ Zarr ~/.julia/packages/Zarr/2V5d9/src/ZGroup.jl:91
[19] zopen (repeats 2 times)
@ ~/.julia/packages/Zarr/2V5d9/src/ZGroup.jl:91 [inlined]
[20] YAXArrayBase.ZarrDataset(g::String)
@ YAXArrayBase ~/.julia/packages/YAXArrayBase/OgYsV/src/datasets/zarr.jl:7
[21] to_dataset(T::Type{YAXArrayBase.ZarrDataset}, g::String)
@ YAXArrayBase ~/.julia/packages/YAXArrayBase/OgYsV/src/datasets/datasetinterface.jl:76
[22] to_dataset(g::String; driver::Symbol)
@ YAXArrayBase ~/.julia/packages/YAXArrayBase/OgYsV/src/datasets/datasetinterface.jl:67
[23] open_dataset(g::String; driver::Symbol)
@ YAXArrays.Datasets ~/.julia/packages/YAXArrays/HzZlX/src/DatasetAPI/Datasets.jl:171
[24] open_dataset(g::String)
@ YAXArrays.Datasets ~/.julia/packages/YAXArrays/HzZlX/src/DatasetAPI/Datasets.jl:171
[25] top-level scope
@ REPL[24]:1
`

However changing the String Variables to Bool seems to work just fine, so Float and Bool do not throw errors.

`<xarray.Dataset>
Dimensions: (X: 2, Y: 2, Time: 714)
Coordinates:

• X (X) float64 -3.914e+05 -3.914e+05
• Y (Y) float64 -6.672e+05 -6.672e+05
• Time (Time) datetime64[ns] 2017-01-01T16:59:45 ... 2019-12-31T05:33:44
  Data variables:
  layer (Time, X, Y) float32 -15.03 -14.98 nan -15.35 ... -14.22 nan -14.2
  orb (Time) float64 44.0 66.0 168.0 117.0 44.0 ... 168.0 117.0 44.0 66.0
  sat_bool (Time) bool True True False True False ... True False True True
  dir_bool (Time) bool True False False True True ... False True True False`

`julia> YAXArrays.open_dataset("cube_bool.zarr")

YAXArray Dataset
Dimensions:
Y Axis with 2 Elements from -667240.0 to -667220.0
X Axis with 2 Elements from -391375.0 to -391395.0
Time Axis with 714 Elements from 2017-01-01T16:59:45 to 2019-12-31T05:33:44
Variables: layer dir_bool sat_bool orb `

Cheers!

[meggart]

Hi, thanks for the report. What you describe seems to be a Zarr error and should be fixed in Zarr.jl. Could you post the contents of the file "cube_from_above.zarr/sat/.zarray" and "cube_from_above.zarr/dir/.zarray"? We support several string data types in Zarr.jl, but there seem to be too many of them in numcodecs.

[MarkusZehner]

cube_from_above.zarr/sat/.zarray:

{
    "chunks": [
        714
    ],
    "compressor": {
        "blocksize": 0,
        "clevel": 5,
        "cname": "lz4",
        "id": "blosc",
        "shuffle": 1
    },
    "dtype": "|O",
    "fill_value": null,
    "filters": [
        {
            "id": "vlen-utf8"
        }
    ],
    "order": "C",
    "shape": [
        714
    ],
    "zarr_format": 2
}

cube_from_above.zarr/dir/.zarray:

{
    "chunks": [
        714
    ],
    "compressor": {
        "blocksize": 0,
        "clevel": 5,
        "cname": "lz4",
        "id": "blosc",
        "shuffle": 1
    },
    "dtype": "|O",
    "fill_value": null,
    "filters": [
        {
            "id": "vlen-utf8"
        }
    ],
    "order": "C",
    "shape": [
        714
    ],
    "zarr_format": 2
}

[MarkusZehner]

just an update: push!(yax.properties, "orb"=>y.orb) works fine with storing dataarrays within properties, credits to @felixcremer.

However, when saving, those properties are not kept:

julia> yax.properties
Dict{String, Any} with 2 entries:
  "sat_bin" => YAXArray with the following dimensions…
  "orb"     => YAXArray with the following dimensions…

julia> test = savecube(yax, "testcube_delete")
YAXArray with the following dimensions
Y                   Axis with 2 Elements from -667240.0 to -667220.0
X                   Axis with 2 Elements from -391375.0 to -391395.0
Time                Axis with 714 Elements from 2017-01-01T16:59:45 to 2019-12-31T05:33:44
TSdecomp            Axis with 8 elements: imf_1 imf_2 imf_3 imf_4 imf_5 imf_6 imf_7 raw 
Total size: 178.5 KB


julia> test
YAXArray with the following dimensions
Y                   Axis with 2 Elements from -667240.0 to -667220.0
X                   Axis with 2 Elements from -391375.0 to -391395.0
Time                Axis with 714 Elements from 2017-01-01T16:59:45 to 2019-12-31T05:33:44
TSdecomp            Axis with 8 elements: imf_1 imf_2 imf_3 imf_4 imf_5 imf_6 imf_7 raw 
Total size: 178.5 KB


julia> test.properties
Dict{String, Any}()

[felixcremer]

The saving of the properties will be fixed by #117. Is there anything else, that should be done here?

[lazarusA]

see https://juliadatacubes.github.io/YAXArrays.jl/dev/UserGuide/write.html#Append-to-a-Dataset for a cleaner way when it comes to append new variables.