Grid in JLD2 files produced by JLD2OutputWriter is not sliced

[ali-ramadhan]

Right now JLD2 files always include the full grid (with halos) in file["grid"] which can make it hard to make plots of the data, especially if it was sliced (requiring users to know how to cut out the halos).

Not sure if this is a bug as file["grid"] might be referring to the grid the simulation was run on in which case it makes sense to include all the coordinates including halos, but this might still surprise users wishing to make plots.

Thanks to @vchuravy for helping uncover this issue!

using Oceananigans, Oceananigans.OutputWriters, JLD2
grid = RegularCartesianGrid(size=(8, 8, 8), extent=(1, 2, 3))
model = IncompressibleModel(grid=grid)
simulation = Simulation(model, Δt=1, stop_iteration=1)
simulation.output_writers[:test] = JLD2OutputWriter(model, merge(model.velocities, model.tracers), schedule=IterationInterval(1), prefix="test", field_slicer=FieldSlicer(k=grid.Nz, with_halos=false))
run!(simulation)
f = jldopen("test.jld2")

so

julia> size(f["timeseries/u/1"])
(8, 8, 1)

but

julia> f["grid/zC"]
10-element Array{Float64,1}:
 -3.1875
 -2.8125
 -2.4375
 -2.0625
 -1.6875
 -1.3125
 -0.9375
 -0.5625
 -0.1875
  0.1875

[glwagner]

That was the intent of file["grid"], but we can support alternative behavior.

I think a more complete solution would abstract the concept of fields in output, and then we could use znodes(outputted_field) to obtain the z-coordinates.

[sandreza]

is this still an issue?

[glwagner]

We should close this, but I'd like to leave an explanatory note in case others want to weigh in.

This issue is really about "vision" than something missing with output writers. This issue really proposes that we think about JLD2 output differently (perhaps, in a way that's similar to how NetCDF views output). In that mode output would be somewhat independent of Oceananigans; if we output "sliced" data, then we could also output a "sliced grid" that represents the portion of the grid on which the data lives. One might then be able to load data and a grid "simply" (ie as Arrays) and make sense of it without using Oceananigans.

But since this issue was opened we have implemented Field.indices to move towards more full support for "windowed" or "sliced" fields (not only when writing output), which sends us in a slightly different direction. Now a single output writer can be used to output fields that are differently sliced or windowed. We also support loading sliced fields with FieldTimeSeries (which is then able to correctly locate sliced fields on a grid), and doing further computations with those fields (we even envision supporting computations between fields that are sliced differently by computing index "intersections").

Since we do save index information now when we slice a field, it is actually possible to manually slice a grid after data has been outputted. So savvy users can achieve this functionality with themselves if they (for some reason) are allergic to FieldTimeSeries.

[glwagner]

We don't want to slice the grid in JLD2 output.