VTK/Writing VTK files using python

From KitwarePublic
Jump to navigationJump to search

Before starting, you are strongly advised to read the VTK file formats documentation here.

Introduction

The purpose of this wiki is to show how one can write VTK files using python. Why using python ? Because it's beautiful & powerful, of course ! ;-) Obviously, you can use your favourite programming language to write such files. In the following, the data you want to write to a VTK file are supposed to arise from scipy/numpy arrays. Please browse scipy.org to get more information if needed.

Before writing VTK files, you have to know that there are two kinds of VTK files: "legacy" VTK files and "XML" VTK files. The first kind is quite obsolete now, but still widely used. Although it is easier to write than the XML kind (understand: with low-level instructions like "printf" in C or "print" in python), this wiki will focus on the last kind, the XML format.

Datasets format

Note: Well, I have initiated this wiki and never found time to write some materials to it :-(

Meanwhile, Gaël wrote a really cool doc on VTK Data Structures, in the Mayavi Users Guide.

These pages show how one can easily write snipped code to create theses datasets, using TVTK. And believe me, this is far easier to write them using TVTK than "pure" VTK ;-)

But if one of you, wikiers, think that this link is irrelevant (not "pure" VTK), no problem, you can remove it, I won't blame you ;-)

Cheers,



As you should already be aware if you have read the VTK file formats documentation, VTK can handle several kinds of datasets, for points data and cells data. The primary are:

  • StructuredPoints (aka ImageData, see later);
  • RectilinearGrid;
  • StructuredGrid;
  • UnstructuredGrid;
  • PolyData.

Briefly:

  • StructuredPoints is the "simplest" format: you only have to set the mesh dimensions (nx, ny, nz), mesh origin (x0, y0, z0) and cell dimensions (dx, dy, dz). Thus, points data are regularly and uniformly spaced.

StructuredPoints dataset

  • RectilinearGrid is also a regularly spaced points data but spacing can be not uniform. Thus, you have to specify the nodes coordinates along the three axes, Ox, Oy & Oz.

RectilinearGrid dataset

  • StructuredGrid is a not regular and not uniform points data. You have to specify the coordinates for all mesh nodes.

StructuredGrid dataset

  • UnstructuredGrid is a not regular and not uniform points data, like the StructuredGrid, but can handle all cell types (see VTK file formats documentation for cell types available).

File:Unstructured grid.png

  • PolyData can be used for any polygonal data. Thus, you have to specify the coordinates for the points and the polygons defined by the points index.

PolyData dataset

"legacy"

First, the "legacy" VTK format has only one file extension for all kind of datasets format: ".vtk", which is not the case for the "XML" as you will see later. Writing "legacy" VTK files is well documented in the VTK file formats documentation and is quite straightforward. You have just to set a few parameters, such as:

  • file format (ASCII or binary);
  • data format (char, int, float, double, etc...)
  • specific information related to the kind of the dataset you want to use.

Note: The following issue is frequently asked on the vtkusers mailing-list, so please read it carefully before complaining about errors while reading your VTK file after you wrote it. If you want to write your data in binary format in your "legacy" VTK file (XML files are not concerned), you have to write it in BigEndian, and not in LittleEndian, even if you use an intel or amd box which are known to be in LittleEndian. So, you have to convert your data in BigEndian before writing them in your VTK file. This is explicitly mentioned in the VTK file formats documentation.

Using python, you can write "legacy" VTK files with your own methods or use pyvtk for VTK version 2.0 format.

"XML Format"

The new XML format is somewhat more complex, but ultimately more flexible. Different types of data have different file extensions. Recently, a simple open-source Python class for writing unstructured data in the VTK XML format has been released. Although this class is probably too specialized for general-purpose use, it is an easily understood example that you can use to create your own classes. It utilized the built-in Python tool xml.dom.minidom to simplify the creation of XML. It also generates a .pvd file, which allows a series of VTK files to be read into Paraview simultaneously.

EVTK

EVTK (Export to VTK) package allows exporting data to binary VTK files. EVTK depends only on common libraries, so it is easy to install in different systems. EVTK is released under a BSD license.

EVTK has been released in two languages: Python and Java.

PyEVTK [1] is composed of a set of pure Python files and it only requires Numpy as additional requirement. EVTK provides low and high level interfaces. While the low level interface can be used to export data that is stored in any type of container, the high level functions make easy to export data stored in Numpy arrays.

JEVTK [2] is a newer version of EVTK that is only composed of pure Java files and it does not requiere any additional library. It provides almost the same interface as PyEVTK with low and high-level interfaces.

Examples

The following examples show how to use the high level functions to export data stored in different grid format (see above for grid definitions).

Image

To export data stored in a regular Cartesian grid using PyEVTK:

from evtk.hl import imageToVTK 
import numpy as np
# Dimensions
nx, ny, nz = 6, 6, 2
ncells = nx * ny * nz
npoints = (nx + 1) * (ny + 1) * (nz + 1)
# Variables
pressure = np.random.rand(ncells).reshape( (nx, ny, nz), order = 'C')
temp = np.random.rand(npoints).reshape( (nx + 1, ny + 1, nz + 1))
imageToVTK("./image", cellData = {"pressure" : pressure}, pointData = {"temp" : temp} )

A similar example using JEVTK:

var nc = 10; var ncells = new int[] {nc, nc, nc}; var nnpoints = (nc + 1) * (nc + 1) * (nc + 1); var nncells = nc * nc * nc;

var origin = new double[]{0.0, 0.0, 0.0}; var spacing = new double[]{1.0, 1.0, 1.0};

var cellData = makeCellData(); var temp = new double[nncells];

// Fill with some random data var rnd = new Random(); for (int i = 0; i < nncells; i++) {

   temp[i] = rnd.nextDouble();

}

cellData.addData("temp", temp);

var pointData = makePointData();

imageToVTK("image", ncells, origin, spacing, cellData, pointData);


Structured grid

To export data stored in a logically structured grid using PyEVTK:

from evtk.hl import gridToVTK 
import numpy as np
import random as rnd
# Dimensions
nx, ny, nz = 6, 6, 2
lx, ly, lz = 1.0, 1.0, 1.0
dx, dy, dz = lx/nx, ly/ny, lz/nz
ncells = nx * ny * nz
npoints = (nx + 1) * (ny + 1) * (nz + 1)
# Coordinates
X = np.arange(0, lx + 0.1*dx, dx, dtype='float64')
Y = np.arange(0, ly + 0.1*dy, dy, dtype='float64')
Z = np.arange(0, lz + 0.1*dz, dz, dtype='float64')
x = np.zeros((nx + 1, ny + 1, nz + 1))
y = np.zeros((nx + 1, ny + 1, nz + 1))
z = np.zeros((nx + 1, ny + 1, nz + 1))
# We add some random fluctuation to make the grid more interesting
for k in range(nz + 1):
for j in range(ny + 1): for i in range(nx + 1):
x[i,j,k] = X[i] + (0.5 - rnd.random()) * 0.2 * dx
y[i,j,k] = Y[j] + (0.5 - rnd.random()) * 0.2 * dy
z[i,j,k] = Z[k] + (0.5 - rnd.random()) * 0.2 * dz
# Variables
pressure = np.random.rand(ncells).reshape( (nx, ny, nz))
temp = np.random.rand(npoints).reshape( (nx + 1, ny + 1, nz + 1))
gridToVTK("./structured", x, y, z, cellData = {"pressure" : pressure}, pointData = {"temp" : temp})

Rectilinear grid

To export data stored in a rectilinear grid (AKA irregular Cartesian grid) using PyEVTK:

from evtk.hl import gridToVTK 
import numpy as np
# Dimensions
nx, ny, nz = 6, 6, 2
lx, ly, lz = 1.0, 1.0, 1.0
dx, dy, dz = lx/nx, ly/ny, lz/nz
ncells = nx * ny * nz
npoints = (nx + 1) * (ny + 1) * (nz + 1)
# Coordinates
x = np.arange(0, lx + 0.1*dx, dx, dtype='float64')
y = np.arange(0, ly + 0.1*dy, dy, dtype='float64')
z = np.arange(0, lz + 0.1*dz, dz, dtype='float64')
# Variables
pressure = np.random.rand(ncells).reshape( (nx, ny, nz))
temp = np.random.rand(npoints).reshape( (nx + 1, ny + 1, nz + 1))
gridToVTK("./rectilinear", x, y, z, cellData = {"pressure" : pressure}, pointData = {"temp" : temp})

Points

To export a set points with two scalar variables (temp and pressure) at each point using PyEVTK:

<syntaxhighlight lang="python">

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})

</syntaxhighlight>


A similar example using JEVTK:

<syntaxhighlight lang="java"> // Creates some points with random coordinates between 0 and 1 var npoints = 25; var x = new double[npoints]; var y = new double[npoints]; var z = new double[npoints];

var rnd = new Random(); for (int i = 0; i < npoints; i++) {

   x[i] = rnd.nextDouble();
   y[i] = rnd.nextDouble();
   z[i] = rnd.nextDouble();

}

// Add a variable associated to each point abd fill with some random values var temp = new double[npoints]; for (int i = 0; i < npoints; i++) {

   temp[i] = rnd.nextDouble();

}

var pointData = EVTK.makePointData(); pointData.addData("temp", temp);

EVTK.pointsToVTK("points", x, y, z, pointData); </syntaxhighlight>