8

对于那些想要将简单的 3D numpy 数组(连同轴)导出到 .vtk(或 .vtr)文件以进行后处理并在 Paraview 或 Mayavi 中显示的人来说,有一个名为PyEVTK的小模块可以做到这一点。该模块支持结构化和非结构化数据等。不幸的是,即使代码在基于 unix 的系统中运行良好,我也无法在任何 Windows 安装上使其工作(不断崩溃),这只会让事情变得复杂。我联系了开发商,但他的建议没有奏效

因此我的问题是:如何使用from vtk.util import numpy_support函数将 3D 数组(函数本身不支持 3D 数组)导出到 .vtk 文件?有没有一种简单的方法可以在不创建 vtkDatasets 等的情况下做到这一点?

非常感谢!

4

5 回答 5

10

它已经永远了,我完全忘记了问这个问题,但我最终弄明白了。我在我的博客 (PyScience) 中写了一篇关于它的文章,提供了关于如何在 NumPy 和 VTK 之间进行转换的教程。有兴趣的不妨看看:

pyscience.wordpress.com/2014/09/06/numpy-to-vtk-converting-your-numpy-arrays-to-vtk-arrays-and-files/

于 2014-09-12T00:04:33.357 回答
5

这不是您问题的直接答案,但如果您有tvtk(如果您有 mayavi,您应该有),您可以使用它将数据写入 vtk 格式。(见:http ://code.enthought.com/projects/files/ETS3_API/enthought.tvtk.misc.html )

它不使用PyEVTK,并且它支持广泛的数据源(不仅仅是结构化和非结构化网格),因此它可能会在其他东西不适用的地方工作。

举个简单的例子(Mayavi 的mlab界面可以让这变得不那么冗长,特别是如果你已经在使用它的话。):

import numpy as np
from enthought.tvtk.api import tvtk, write_data

data = np.random.random((10,10,10))

grid = tvtk.ImageData(spacing=(10, 5, -10), origin=(100, 350, 200), 
                      dimensions=data.shape)
grid.point_data.scalars = np.ravel(order='F')
grid.point_data.scalars.name = 'Test Data'

# Writes legacy ".vtk" format if filename ends with "vtk", otherwise
# this will write data using the newer xml-based format.
write_data(grid, 'test.vtk')

以及输出文件的一部分:

# vtk DataFile Version 3.0
vtk output
ASCII
DATASET STRUCTURED_POINTS
DIMENSIONS 10 10 10
SPACING 10 5 -10
ORIGIN 100 350 200
POINT_DATA 1000
SCALARS Test%20Data double
LOOKUP_TABLE default
0.598189 0.228948 0.346975 0.948916 0.0109774 0.30281 0.643976 0.17398 0.374673 
0.295613 0.664072 0.307974 0.802966 0.836823 0.827732 0.895217 0.104437 0.292796 
0.604939 0.96141 0.0837524 0.498616 0.608173 0.446545 0.364019 0.222914 0.514992 
...
...
于 2013-04-02T02:03:09.217 回答
2

Mayavi 的 TVTK 有一种编写 vtk 文件的优美方式。这是我根据@Joe 和 tvtk 文档为自己编写的测试示例。它比 evtk 的优势是同时支持 ascii 和 html。希望它能帮助其他人。

from tvtk.api import tvtk, write_data
import numpy as np

#data = np.random.random((3, 3, 3))
#
#i = tvtk.ImageData(spacing=(1, 1, 1), origin=(0, 0, 0))
#i.point_data.scalars = data.ravel()
#i.point_data.scalars.name = 'scalars'
#i.dimensions = data.shape
#
#w = tvtk.XMLImageDataWriter(input=i, file_name='spoints3d.vti')
#w.write()

points = np.array([[0,0,0], [1,0,0], [1,1,0], [0,1,0]], 'f')
(n1, n2)  = points.shape
poly_edge = np.array([[0,1,2,3]])

print n1, n2
## Scalar Data
#temperature = np.array([10., 20., 30., 40.])
#pressure = np.random.rand(n1)
#
## Vector Data
#velocity = np.random.rand(n1,n2)
#force     = np.random.rand(n1,n2)
#
##Tensor Data with 
comp = 5
stress = np.random.rand(n1,comp)
#
#print stress.shape
## The TVTK dataset.
mesh = tvtk.PolyData(points=points, polys=poly_edge)
#
## Data 0 # scalar data
#mesh.point_data.scalars = temperature
#mesh.point_data.scalars.name = 'Temperature'
#
## Data 1 # additional scalar data
#mesh.point_data.add_array(pressure)
#mesh.point_data.get_array(1).name = 'Pressure'
#mesh.update()
#
## Data 2 # Vector data
#mesh.point_data.vectors = velocity
#mesh.point_data.vectors.name = 'Velocity'
#mesh.update()
#
## Data 3 additional vector data
#mesh.point_data.add_array( force)
#mesh.point_data.get_array(3).name = 'Force'
#mesh.update()

mesh.point_data.tensors = stress
mesh.point_data.tensors.name = 'Stress'

# Data 4 additional tensor Data
#mesh.point_data.add_array(stress)
#mesh.point_data.get_array(4).name = 'Stress'
#mesh.update()

write_data(mesh, 'polydata.vtk')

# XML format 
# Method 1
#write_data(mesh, 'polydata')

# Method 2
#w = tvtk.XMLPolyDataWriter(input=mesh, file_name='polydata.vtk')
#w.write()
于 2015-06-05T05:06:46.140 回答
2

我知道这有点晚了,我很喜欢你的教程@somada141。这也应该有效。

def numpy2VTK(img, spacing=[1.0, 1.0, 1.0]):
 # evolved from code from Stou S.,
 # on http://www.siafoo.net/snippet/314
 # This function, as the name suggests, converts numpy array to VTK
 importer = vtk.vtkImageImport()

 img_data = img.astype('uint8')
 img_string = img_data.tostring()  # type short
 dim = img.shape

 importer.CopyImportVoidPointer(img_string, len(img_string))
 importer.SetDataScalarType(VTK_UNSIGNED_CHAR)
 importer.SetNumberOfScalarComponents(1)

 extent = importer.GetDataExtent()
 importer.SetDataExtent(extent[0], extent[0] + dim[2] - 1,
                       extent[2], extent[2] + dim[1] - 1,
                       extent[4], extent[4] + dim[0] - 1)
 importer.SetWholeExtent(extent[0], extent[0] + dim[2] - 1,
                        extent[2], extent[2] + dim[1] - 1,
                        extent[4], extent[4] + dim[0] - 1)

 importer.SetDataSpacing(spacing[0], spacing[1], spacing[2])
 importer.SetDataOrigin(0, 0, 0)


 return importer

希望能帮助到你!

于 2016-10-07T15:38:31.837 回答
1

这是从此处SimpleITK获取功能的版本:load_itk

import SimpleITK as sitk
import numpy as np


if len(sys.argv)<3:
    print('Wrong number of arguments.', file=sys.stderr)
    print('Usage: ' + __file__ + ' input_sitk_file' + ' output_sitk_file', file=sys.stderr)
    sys.exit(1)


def quick_read(filename):
    # Read image information without reading the bulk data.
    file_reader = sitk.ImageFileReader()
    file_reader.SetFileName(filename)
    file_reader.ReadImageInformation()
    print('image size: {0}\nimage spacing: {1}'.format(file_reader.GetSize(), file_reader.GetSpacing()))
    # Some files have a rich meta-data dictionary (e.g. DICOM)
    for key in file_reader.GetMetaDataKeys():
        print(key + ': ' + file_reader.GetMetaData(key))



def load_itk(filename):
    # Reads the image using SimpleITK
    itkimage = sitk.ReadImage(filename)
    # Convert the image to a  numpy array first and then shuffle the dimensions to get axis in the order z,y,x
    data = sitk.GetArrayFromImage(itkimage)
    # Read the origin of the ct_scan, will be used to convert the coordinates from world to voxel and vice versa.
    origin = np.array(list(reversed(itkimage.GetOrigin())))
    # Read the spacing along each dimension
    spacing = np.array(list(reversed(itkimage.GetSpacing())))
    return data, origin, spacing


def convert(data, output_filename):
    image = sitk.GetImageFromArray(data)
    writer = sitk.ImageFileWriter()
    writer.SetFileName(output_filename)
    writer.Execute(image)


def wait():
    print('Press Enter to load & convert or exit using Ctrl+C')
    input()


quick_read(sys.argv[1])
print('-'*20)
wait()
data, origin, spacing = load_itk(sys.argv[1])
convert(sys.argv[2])
于 2019-04-30T08:32:06.033 回答