Writing points to vtk with WriteVTK.jl

Specific Domains Visualization
,
1

Hi, I am trying to export a set of points with corresponding value at each point using WriteVTK.jl. More precisely, I am trying to convert following codes given here.


from evtk.hl import pointsToVTK 
 import numpy as np  
 npoints = 100  
 x = np.random.rand(npoints)  
 y = np.random.rand(npoints)  
 z = np.random.rand(npoints)  
 pressure = np.random.rand(npoints)  
 temp = np.random.rand(npoints)  
 pointsToVTK("./points", x, y, z, data = {"temp" : temp, "pressure" : pressure})

I have reviewed this WriteVTK page but I could not find corresponding codes for the above python codes.
Any help is appreciated.

2

Hi, the above Python code exports an unstructured dataset (.vtu format), where each point corresponds to a “vertex” cell.

Unfortunately the equivalent in WriteVTK is a bit more verbose, as there is no direct analogous to pointsToVTK (but feel free to open an issue or PR for that!). For reference you should look at the unstructured datasets section of the README.

The following code is the equivalent to your example in Julia:

using WriteVTK

npoints = 100
x = rand(npoints)
y = rand(npoints)
z = rand(npoints)
pressure = rand(npoints)
temp = rand(npoints)
cells = [MeshCell(VTKCellTypes.VTK_VERTEX, (i, )) for i = 1:npoints]
vtk_grid("./points", x, y, z, cells) do vtk
    vtk["pressure", VTKPointData()] = pressure
    vtk["temp", VTKPointData()] = temp
end
4 Likes
3

Thank you very much. It worked!

4

Based on the solution provided by @jipolanco, I wrote this function, which works with vectors. It assumes x is always passed. Maybe someone else finds it useful.

using WriteVTK

"""
Save a set of points and associated field data to a VTK file.

# Arguments
- `save_file_path::String`: Path to save the VTK file (without `.vtu` extension)
- `x::Vector{T}`: x-coordinates of the points.
- `y::Union{Vector{T}, Nothing}`: y-coordinates of the points
    (defaults to zeros if `nothing`).
- `z::Union{Vector{T}, Nothing}`: z-coordinates of the points
    (defaults to zeros if `nothing`).
- `fields::Dict{String, Vector{T}}`: Dictionary of field data to associate with
    the points.

# Example
```julia
x = [1.0, 2.0, 3.0]
y = [4.0, 5.0, 6.0]
T = [10.0, 20.0, 30.0]
pointsToVTK("./output", x, y, fields=Dict("T" => T))
```
"""
function pointsToVTK(
    save_file_path::String,
    x::Vector{T},
    y::Union{Vector{T}, Nothing}=nothing,
    z::Union{Vector{T}, Nothing}=nothing;
    fields::Dict{String, Vector{T}}
) where {T <: Real}

    # Remove the .vtu extension if it exists
    if endswith(save_file_path, ".vtu")
        save_file_path = replace(save_file_path, ".vtu" => "")
    end

    # Default y to zeros if not provided
    if y === nothing
        y = zeros(T, length(x))
    end

    # Default z to zeros if not provided
    if z === nothing
        z = zeros(T, length(x))
    end

    # Ensure x, y, and z have the same length
    @assert length(x) == length(y) == length(z) "x, y, z lengths must match"

    npoints = length(x)

    # Create cells for each point
    cells = [MeshCell(VTKCellTypes.VTK_VERTEX, (i, )) for i = 1:npoints]

    # Write the VTK file
    vtk_grid(save_file_path, x, y, z, cells) do vtk
        for (name, values) in fields
            vtk[name, VTKPointData()] = values
        end
    end
end

I used it like this in my application:


include("control_points.jl")
XY, TT = get_control_points_and_temperatures(1.0)
x = [xy[1] for xy in XY[:]] # vector
y = [xy[2] for xy in XY[:]] # vector
T = TT[:] # vector
T2 = T .- 5 # vector
pointsToVTK("./points", x, y, fields=Dict("Temp"=>T, "Temp_shifted"=>T2))