How to convert 3D numpy array to vtk and save the .vtk file?

Support
,
1

Hi,all.Sorry for my little knowledge about VTK.
Here is the thing I want to do:
I have a scalar volume (3D numpy array)
1646504309(1)
I want to add some information describing the image and then save it as a .vtk file as follows:
image

image
image969×229 22.3 KB

What should I do to achieve this conversion?
Thank you in advance for your attention and help!

2

I don’t know about adding information to a vtk file, but for converting a numpy array to a vtk array format there is vtk.util.numpy_support.numpy_to_vtk; see here

1 Like
3

yes you are right.

import vtk.util.numpy_support as numpy_support

def numpyToVTK(data):
    data_type = vtk.VTK_FLOAT
    shape = data.shape

    flat_data_array = data.flatten()
    vtk_data = numpy_support.numpy_to_vtk(num_array=flat_data_array, deep=True, array_type=data_type)
    
    img = vtk.vtkImageData()
    img.GetPointData().SetScalars(vtk_data)
    img.SetDimensions(shape[0], shape[1], shape[2])
    return img

with the above function I did a quick test below:

test_arr = np.arange(36).reshape(3,4,3)
test_arr

--> array([[
        [ 0,  1,  2],
        [ 3,  4,  5],
        [ 6,  7,  8],
        [ 9, 10, 11]],

       [[12, 13, 14],
        [15, 16, 17],
        [18, 19, 20],
        [21, 22, 23]],

       [[24, 25, 26],
        [27, 28, 29],
        [30, 31, 32],
        [33, 34, 35]]])

img = numpyToVTK(data)

Then verify the vtkImageData is right or not

test = vtk.util.numpy_support.vtk_to_numpy(img.GetPointData().GetScalars())
test = test.reshape(3,4,3)
(test == test_arr).all()

-->True

However,I’m like you and don’t know how to save the vtk file to the hard drive

1 Like
4

If you have the array associated to a storable mrml node then you can write it to file as shown here.

To create a volume node out of your numpy array,

  1. Add a volume node: volume_node = slicer.mrmlScene.AddNewNodeByClass("vktMRMLScalarVolumeNode")
  2. Set the values of volume_node's underlying vtkImageData using slicer.util.updateVolumeFromArray described here.
2 Likes
5

Thanks, Ebrahim.
Thank you for telling me the way to create the volume node.
However,if I want to save the numpy array as .vtk file,actually it’s quite simple as follows:

f = open('test.vtk','w') # change your vtk file name
f.write('# vtk DataFile Version 2.0\n')
f.write('test\n')
f.write('ASCII\n')
f.write('DATASET STRUCTURED_POINTS\n')
f.write('DIMENSIONS 3 3 4\n') # change your dimension
f.write('SPACING 0.100000 0.100000 0.100000\n')
f.write('ORIGIN 0 0 0\n')
f.write('POINT_DATA 36\n') # change the number of point data
f.write('SCALARS volume_scalars float 1\n')
f.write('LOOKUP_TABLE default\n')
f.write(numpy_array) # change your numpy array here but need some processing
f.close()
6

do you maybe have an idea how to do this for a 2D numpy array?

7

Just found out:

8

Thanks Ebrahim,
In the begining,I just use python open() to write the numpy data,but later I found it is not efficient,taking too long to save large volume.
Your way is better so I have changed the Solution.I have followed your steps to save vtk file on disk but there’s some modification.

  1. Create a 3D scalar image just as an example
import numpy as np

numpy_arr = np.arange(1,61)
numpy_arr = numpy_arr.reshape(5,4,3)
numpy_arr.shape
--> (5, 4, 3)
  1. Create a new empty volume
nodeName = "MyNewVolume"
imageSize = [5, 4, 3]
voxelType=vtk.VTK_UNSIGNED_CHAR
imageOrigin = [0.0, 0.0, 0.0]
imageSpacing = [0.1, 0.1, 0.1]
imageDirections = [[1,0,0], [0,1,0], [0,0,1]]
fillVoxelValue = 0

# Create an empty image volume, filled with fillVoxelValue
imageData = vtk.vtkImageData()
imageData.SetDimensions(imageSize)
imageData.AllocateScalars(voxelType, 1)
imageData.GetPointData().GetScalars().Fill(fillVoxelValue)

# Create volume node
volumeNode = slicer.mrmlScene.AddNewNodeByClass("vtkMRMLScalarVolumeNode", nodeName)
volumeNode.SetOrigin(imageOrigin)
volumeNode.SetSpacing(imageSpacing)
volumeNode.SetIJKToRASDirections(imageDirections)
volumeNode.SetAndObserveImageData(imageData)
volumeNode.CreateDefaultDisplayNodes()
  1. Update voxels of my volume node from a numpy array
slicer.util.updateVolumeFromArray(volumeNode,numpy_arr)
  1. Save volumeNode to file
myStorageNode = volumeNode.CreateDefaultStorageNode()
myStorageNode.SetFileName("test_data.vtk")
myStorageNode.WriteData(volumeNode)

It is pretty okay when I open this test_data.vtk to check.

image
image790×239 8.25 KB

Also,I have loaded this file in Slicer as follows,

image
image999×819 25.2 KB

Then,I get the array back from this volume in order to double check,

volumeNode = getNode('MyNewVolume')
# volumeNode
a = arrayFromVolume(volumeNode)
a.shape
--> (5, 4, 3)

(a==numpy_arr).all()
--> True

So, I am pretty sure the data is saved correctly as vtk file.
Thank you very much for your help again,Ebrahim.

1 Like