diff --git a/docs/src/grids/datasets.md b/docs/src/grids/datasets.md index b01e9ebb..ee1495af 100644 --- a/docs/src/grids/datasets.md +++ b/docs/src/grids/datasets.md @@ -14,8 +14,8 @@ vtk_grid("fields", x, y, z) do vtk vtk["Pressure"] = rand(Nx - 1, Ny - 1, Nz - 1) # scalar field attached to cells vtk["Velocity"] = rand(3, Nx, Ny, Nz) # vector field attached to points vtk["VelocityGradients"] = rand(3, 3, Nx, Ny, Nz) # 3×3 tensor field attached to points - vtk["date"] = "31/10/2021" # metadata ("field data" in VTK) - vtk["time"] = 0.42 # metadata ("field data" in VTK) + vtk["Date"] = "31/10/2021" # metadata ("field data" in VTK) + vtk["TimeValue"] = 0.42 # metadata ("field data" in VTK) end ``` @@ -25,6 +25,17 @@ input. In particular, note that the `Pressure` field is attached to cells instead of points, since it has dimensions ``(N_x - 1) × (N_y - 1) × (N_z - 1)``, which is the number of cells in structured files (see [Cells in structured formats](@ref)). For more control, see [Dataset types](@ref) below. +!!! note "Time values" + + As a sidenote, the `TimeValue` field used above [is special](https://docs.vtk.org/en/latest/design_documents/IOXMLTimeInFieldData.html#field-data-as-time-meta-data-in-vtk-xml-file-formats), + as it is interpreted by VTK and ParaView as the time associated to the dataset. + This can be convenient when writing time series data (for example, results from a simulation at different timesteps). + + Other alternatives for including time information are to use [ParaView collections](@ref) (`.pvd` format), or to generate + [JSON `.series` files](https://gitlab.kitware.com/paraview/paraview/-/blob/v5.5.0/Documentation/release/ParaView-5.5.0.md#json-based-new-meta-file-format-for-series-added) + (the latter can't be currently done by WriteVTK, but the format is very simple). + + ## Passing tuples of arrays Note that, in the example above, the `Velocity` vector field is passed as a single ``3 × N_x × N_y × N_z`` array, which is the layout ultimately used in VTK formats. @@ -124,8 +135,8 @@ vtk_grid("fields_explicit", x, y, z) do vtk vtk["Pressure", VTKCellData()] = rand(Nx - 1, Ny - 1, Nz - 1) vtk["Velocity", VTKPointData()] = rand(3, Nx, Ny, Nz) vtk["VelocityGradients", VTKPointData()] = rand(3, 3, Nx, Ny, Nz) - vtk["date", VTKFieldData()] = "31/10/2021" - vtk["time", VTKFieldData()] = 0.42 + vtk["Date", VTKFieldData()] = "31/10/2021" + vtk["TimeValue", VTKFieldData()] = 0.42 end ``` diff --git a/src/WriteVTK.jl b/src/WriteVTK.jl index a41ca6cf..c58ad54e 100644 --- a/src/WriteVTK.jl +++ b/src/WriteVTK.jl @@ -107,6 +107,16 @@ end DatasetFile(dtype, xdoc::XMLDocument, fname::AbstractString, args...; kwargs...) = DatasetFile(xdoc, add_extension(fname, dtype), xml_name(dtype), args...; kwargs...) +function data_format(vtk::DatasetFile) + if vtk.appended + :appended + elseif vtk.ascii + :ascii + else + :inline + end +end + function show(io::IO, vtk::DatasetFile) open_str = isopen(vtk) ? "open" : "closed" print(io, "VTK file '$(vtk.path)' ($(vtk.grid_type) file, $open_str)") diff --git a/src/gridtypes/pvtk_grid.jl b/src/gridtypes/pvtk_grid.jl index f80ed481..031890a0 100644 --- a/src/gridtypes/pvtk_grid.jl +++ b/src/gridtypes/pvtk_grid.jl @@ -14,6 +14,11 @@ struct PVTKFile <: VTKFile path::String end +# This is just to make a PVTKFile work like a DatasetFile. +# Only used when writing VTKFieldData to a PVTK file. +# Note that data is always written as text (ASCII). +data_format(::PVTKFile) = :ascii + # Returns true if the arguments do *not* contain any cell vectors. _pvtk_is_structured(x::AbstractVector{<:AbstractMeshCell}, args...) = Val(false) _pvtk_is_structured(x, args...) = _pvtk_is_structured(args...) @@ -130,11 +135,44 @@ end # Add point and cell data as usual function Base.setindex!( - pvtk::PVTKFile, data, name::AbstractString, args..., + pvtk::PVTKFile, data, name::AbstractString, loc::AbstractFieldData, + ) + pvtk.vtk[name, loc] = data +end + +# In the case of field data, also add it to the pvtk file. +# Otherwise field data (typically the time / "TimeValue") is not seen by ParaView. +function Base.setindex!( + pvtk::PVTKFile, data, name::AbstractString, loc::VTKFieldData, + ) + pvtk.vtk[name, loc] = data + if pvtk.pvtkargs.ismain + add_field_data(pvtk, data, name, loc) + end + data +end + +# Used in add_field_data. +# We need to find the PUnstructuredGrid / PStructuredGrid / ... node. +function find_base_xml_node_to_add_field(pvtk::PVTKFile, loc) + xroot = root(pvtk.xdoc) + pgrid_type = "P" * pvtk.vtk.grid_type + find_element(xroot, pgrid_type) +end + +# If `loc` was not passed, try to guess which kind of data was passed. +function Base.setindex!( + pvtk::PVTKFile, data, name::AbstractString, ) - pvtk.vtk[name, args...] = data + loc = guess_data_location(data, pvtk.vtk) :: AbstractFieldData + pvtk[name, loc] = data end +# Write XML attribute. +# Example: +# +# pvtk[VTKCellData()] = Dict("HigherOrderDegrees" => "HigherOrderDegreesDataset") +# function Base.setindex!( pvtk::PVTKFile, attributes, loc::AbstractFieldData, ) diff --git a/src/write_data.jl b/src/write_data.jl index 5b2fe9d3..10a63c2a 100644 --- a/src/write_data.jl +++ b/src/write_data.jl @@ -110,17 +110,23 @@ function write_array(io, data::Tuple) n end -function add_data_ascii(xml, x::Union{Number,String}) +function add_data_ascii(xml, x::Union{Number,AbstractString}) add_text(xml, " ") add_text(xml, string(x)) end -add_data_ascii(xml, x::AbstractArray) = map(v -> add_data_ascii(xml, v), x) +function add_data_ascii(xml, x::AbstractArray) + for v ∈ x + add_data_ascii(xml, v) + end + nothing +end function add_data_ascii(xml, data::Tuple) for i in eachindex(data...), x in data add_data_ascii(xml, x[i]) end + nothing end function set_num_components(xDA, vtk, data, loc) @@ -173,9 +179,10 @@ function data_to_xml(vtk, xParent, data, name, set_attribute(xDA, "ComponentName$(i - 1)", n) # 0-based end end - if vtk.appended + fmt = data_format(vtk) + if fmt === :appended data_to_xml_appended(vtk, xDA, data) - elseif vtk.ascii + elseif fmt === :ascii data_to_xml_ascii(vtk, xDA, data) else data_to_xml_inline(vtk, xDA, data) @@ -209,7 +216,7 @@ containing the size of the data array in bytes. """ function data_to_xml_appended(vtk::DatasetFile, xDA::XMLElement, data) - @assert vtk.appended + @assert data_format(vtk) === :appended buf = vtk.buf # append buffer compress = vtk.compression_level > 0 @@ -258,7 +265,7 @@ end Add inline, base64-encoded data to VTK XML file. """ function data_to_xml_inline(vtk::DatasetFile, xDA::XMLElement, data) - @assert !vtk.appended && !vtk.ascii + @assert data_format(vtk) === :inline compress = vtk.compression_level > 0 # DataArray node @@ -305,8 +312,8 @@ end Add inline data to VTK XML file in ASCII format. """ -function data_to_xml_ascii(vtk::DatasetFile, xDA::XMLElement, data) - @assert !vtk.appended && vtk.ascii +function data_to_xml_ascii(vtk::VTKFile, xDA::XMLElement, data) + @assert data_format(vtk) === :ascii set_attribute(xDA, "format", "ascii") add_text(xDA, "\n") add_data_ascii(xDA, data) @@ -320,18 +327,10 @@ end Add either point or cell data to VTK file. """ -function add_field_data(vtk::DatasetFile, data, name::AbstractString, +function add_field_data(vtk::VTKFile, data, name::AbstractString, loc::AbstractFieldData; component_names::Union{AbstractVector, Nothing}=nothing) - # Find Piece node. - xroot = root(vtk.xdoc) - xGrid = find_element(xroot, vtk.grid_type) - - xbase = if loc === VTKFieldData() - xGrid - else - find_element(xGrid, "Piece") - end + xbase = find_base_xml_node_to_add_field(vtk, loc) # Find or create "nodetype" (PointData, CellData or FieldData) node. nodetype = node_type(loc) @@ -344,6 +343,18 @@ function add_field_data(vtk::DatasetFile, data, name::AbstractString, xDA end +function find_base_xml_node_to_add_field(vtk::DatasetFile, loc) + # Find Piece node. + xroot = root(vtk.xdoc) + xGrid = find_element(xroot, vtk.grid_type) + xbase = if loc === VTKFieldData() + xGrid + else + find_element(xGrid, "Piece") + end + xbase +end + vtk_point_data(args...; kwargs...) = add_field_data(args..., VTKPointData(); kwargs...) vtk_cell_data(args...; kwargs...) = add_field_data(args..., VTKCellData(); kwargs...) vtk_field_data(args...; kwargs...) = add_field_data(args..., VTKFieldData(); kwargs...) diff --git a/test/checksums.sha1 b/test/checksums.sha1 index 650ddcde..ebbf13a4 100644 --- a/test/checksums.sha1 +++ b/test/checksums.sha1 @@ -55,28 +55,28 @@ b20bf14b2531aaf323b891e090d0615d0e943ef6 collection_03.vtr 8005b27b071ed50765d0d859a9e4e84a52c3963c arraydata_2D.vti d5215c84b3c1630bedd98994963560ff14200dc5 arraydata_3D.vti aa49366e66d8930fa639a8361d153a95baca361e arrays.vti -2626bef24152fcfe28c810ca68adf473c8fad90c simulation.pvtu -4a0be2eaf5c7ffdb60970f0396a11ef2ab3affca simulation/simulation_1.vtu +3edeb95b6d2f955c8ee5539cb0a64be1f53f0bf1 simulation.pvtu +1348c734f1b79881f6b74461aaee059fa2f01c09 simulation/simulation_1.vtu ee3ea76e543c6e3a9c1c28d80596b68c5b22d7ac higherorderdegrees.pvtu 9a524fb558d6e8e2a42c6bf704da88c2a835b296 higherorderdegrees/higherorderdegrees_1.vtu -e23c8f0db331e1db635aa5d98b9f1375132e7ffd pimage.pvti -e581d4824d67655584f573743f42788fb26e162c pimage/pimage_1.vti -150ce17cb8ddbb974344954a3193a3373ee1c3fc pimage/pimage_2.vti -e389b7948af0bd9cd9d178dcc6d0fb86493738f3 pimage/pimage_3.vti -669fdb06ab07c289cb5c01dd511da380fb06bdfd pimage/pimage_4.vti -1ea223999c053e38ab84726c0a4bfdd0468435ef pimage/pimage_5.vti -21dd9d15e51720f9d5822bf663c718974094d779 pimage/pimage_6.vti -d35ceff83c6c43b43eefe07007cbba8752cead15 prectilinear.pvtr -b7edbc7bb8907b472551d456f60cbc7f0c915b53 prectilinear/prectilinear_1.vtr -d26dd7cf0570b1fbad8714e16e784f09bb6affed prectilinear/prectilinear_2.vtr -e36a533248051994e7da804bff1b168d2d866bba prectilinear/prectilinear_3.vtr -4b73039f56d1a92fe958a98c86dde64b3fb86b63 prectilinear/prectilinear_4.vtr -87a07f8911949b002501da5bdd5daf497a65d09a prectilinear/prectilinear_5.vtr -4fd1d4de3867230a0636bef0260417ed17b7c17f prectilinear/prectilinear_6.vtr -fc3a26b061585c093520bbb720714d180ca22335 pstructured.pvts -584520780f94c58a9919cb216d116733f8466e0f pstructured/pstructured_1.vts -0ff33138288a96f8f6f4377f2e061748aa702d97 pstructured/pstructured_2.vts -e365dc765a17a288640e04f726cf8a7421559f39 pstructured/pstructured_3.vts -0daf5cccfaec20cdcfae94294513ef211954c9e1 pstructured/pstructured_4.vts -a2b4e0e63dda63eee19466b08d5d4ae8d34dcfd5 pstructured/pstructured_5.vts -0adbea4b9bebd2b95569398666028c19f3d8cb4a pstructured/pstructured_6.vts +03f54866860e984b2537a118def1c857f855d6ab pimage.pvti +3badee9575f3feb593109644d6124416731085da pimage/pimage_1.vti +a898a771d8e0962628118d73a216a3e4660a5295 pimage/pimage_2.vti +436cb1c2024e5795664a1355c8d19b7f336c2184 pimage/pimage_3.vti +4badae7ec804773110e662e105531e5e3d4fd288 pimage/pimage_4.vti +99843bc0f87c436d9ed385b40707ee1e0327f709 pimage/pimage_5.vti +4e37ea3120610664425712254b88d88f78d4e45d pimage/pimage_6.vti +9c0d3e7ee95c6e7bfc2f0554972d31fb9d832172 prectilinear.pvtr +e7ff55d7c8bb2c1739a2dc73c0c1fdeff2f7d536 prectilinear/prectilinear_1.vtr +6b0b2f8e5d381ec4f9644905677be8f29f432af2 prectilinear/prectilinear_2.vtr +68f065b66a0ff9ba80bc18776912f248ebc83390 prectilinear/prectilinear_3.vtr +b8b456b99f45f45e2e0fcc0de56963185b3f24ec prectilinear/prectilinear_4.vtr +a6fd8a0367efc31e7616dd7046f048c424aa8ae0 prectilinear/prectilinear_5.vtr +b1ef0b3a9e3c25781c31e0d2ceb552c78dab23f8 prectilinear/prectilinear_6.vtr +d9c0b2367bf7ba42534a134b3ef4dace753560a9 pstructured.pvts +92ff18a7c1c631800de3adc3b0bd4631e16be845 pstructured/pstructured_1.vts +9bb94bd6b2374156e60df995340e3b395f2ae6b9 pstructured/pstructured_2.vts +ad4b841fdf277c6138bbd763f4db2a359c46e400 pstructured/pstructured_3.vts +651905ba5698b503f58d40fc1051267fefa4b49e pstructured/pstructured_4.vts +bba1504744a03b3f9594fcfe63435b7311e030fe pstructured/pstructured_5.vts +2e8acea4c25a7ef64c25c05ff3f65b3890a52e72 pstructured/pstructured_6.vts diff --git a/test/pvtk_grid.jl b/test/pvtk_grid.jl index 660619b3..d3e60e24 100644 --- a/test/pvtk_grid.jl +++ b/test/pvtk_grid.jl @@ -18,6 +18,7 @@ function pvtk_unstructured() pvtk["Processor"] = [1,2] pvtk["Temperature", VTKPointData()] = y pvtk["Id", VTKCellData()] = [2, 1] + pvtk["TimeValue"] = 4.2 end @test isfile(vtufile) @test vtufile ∈ outfiles @@ -73,6 +74,7 @@ function pvtk_imagedata() ) do vtk vtk["point_data"] = point_data vtk["process_id"] = processid + vtk["TimeValue"] = 1.2 end end end @@ -95,6 +97,7 @@ function pvtk_rectilinear() ) do vtk vtk["point_data"] = point_data vtk["process_id"] = processid + vtk["TimeValue"] = 3.4 end end end @@ -124,6 +127,7 @@ function pvtk_structured() ) do vtk vtk["point_data"] = point_data vtk["process_id"] = processid + vtk["TimeValue"] = 3.5 end end end