Преобразование VTK UnstructuredGrid в StructuredGrid

Причиной моей проблемы является то, что у меня есть 3D-структура, сохраненная в файле.vtk, которой мне нужно манипулировать (расширять, разрушать и т. Д.). Следующие фрагменты кода предназначены для последовательного запуска, т. Е. Если вы запускаете их один за другим, проблем быть не должно (кроме тех, о которых я упоминал!).

Я очень новичок в ВТК, поэтому прошу прощения за любые очень простые ошибки!

проблема

Моя проблема проистекает из проблемы с SimpleITK, когда он не может прочитать UnstructuredGrid или PolyData:

In [1]: import SimpleITK as sitk
In [2]: img_vtk = sitk.ReadImage(file_vtk)
Traceback (most recent call last):

  File "<ipython-input-52-435ce999db50>", line 1, in <module>
    img_vtk = sitk.ReadImage(file_vtk)

  File "/usr/local/lib/python3.5/dist-packages/SimpleITK/SimpleITK.py", line 8614, in ReadImage
    return _SimpleITK.ReadImage(*args)

RuntimeError: Exception thrown in SimpleITK ReadImage: /tmp/SimpleITK/Code/IO/src/sitkImageReaderBase.cxx:97:
sitk::ERROR: Unable to determine ImageIO reader for "/data/ROMPA_MRIandSeg/09S/Analysis/1_model/clip_dilate.vtk"

Однако SimpleITK может читать StructuredGrid, поэтому я попытался решить эту проблему путем чтения с использованием VTK и конвертации.

import vtk
reader = vtk.vtkGenericDataObjectReader() # Using generic to allow it to match either Unstructured or PolyData
reader.SetFileName(file_vtk)
reader.Update()
output = reader.GetOutput()

Однако с этого момента каждый метод, который я пробовал, кажется, потерпел неудачу.

Предлагаемые решения

Преобразование в NumPy, затем преобразование в ситк изображения

Я попытался преобразовать его в массив numpy (), а затем интерполировать регулярную сетку с фиктивной переменной 1, чтобы указать значения в структуре.

from vtk.utils import numpy_support
import scipy.interpolate
import numpy as np

nparray = numpy_support.vtk_to_numpy(output.GetPointData().GetArray(0))

output_bounds = output.GetBounds()
x_grid = range(math.floor(output_bounds[0]),math.ceil(output_bounds[1]),1)
y_grid = range(math.floor(output_bounds[2]),math.ceil(output_bounds[3]),1)
z_grid = range(math.floor(output_bounds[4]),math.ceil(output_bounds[5]),1)
grid = list()
for x in x_grid:
   for y in y_grid:
      for z in z_grid:
         grid.append((x,y,z))
dummy = np.array([1 for i in range(nparray.shape[0])])
npgrid = scipy.interpolate.griddata(nparray,dummy,grid,fill_value=0)

npgrid.reshape(len(x_grid),len(y_grid),len(z_grid))
img = sitk.GetImageFromArray(npgrid)
sitk.WriteImage(img,file_out)

Однако, когда я загружаю это в ParaView, для вывода отображается ограничивающая рамка, но контур вывода пуст.

С помощью ShepardMethod

Я попытался интерполировать с помощью встроенного ShepardMethodпосле преобразования UnstructuredGrid в PolyData (как я в основном видел, как ShepardMethod применяется к PolyData):

bounds = output.GetBounds()
spacings = [1.0,1.0,1.0] # arbitrary spacing
dimensions = [0,0,0]
for i,spacing in enumerate(spacings):
    dimensions[i] = int(math.ceil((bounds[i*2 + 1]-bounds[i*2])/spacing))

vtkPoints = vtk.vtkPoints()
for i in range(0,nparray.shape[0]):
   x=nparray[i,0]
   y=nparray[i,1]
   z=nparray[i,2]
   p=[x,y,z]
   vtkPoints.InsertNextPoint(p)
poly = vtk.vtkPolyData()
poly.SetPoints(vtkPoints)    

shepard = vtk.vtkShepardMethod()
shepard.SetInputData(poly)
shepard.SetSampleDimensions(dimensions)
shepard.SetModelBounds(output.GetBounds())
shepard.Update()
shepard_data = shepard.GetOutput().GetPointData().GetArray(0)

shepard_numpy = numpy_support.vtk_to_numpy(shepard_data)
shepard_numpy = shepard_numpy.reshape(dimensions[0],dimensions[1],dimensions[2])
shepard_img = sitk.GetImageFromArray(shepard_numpy)
sitk.WriteImage(shepard_img,file_out)

Как и в приведенном выше описании, это обеспечило ограничивающий прямоугольник в ParaView. Применение контура обеспечило структуру из двух треугольников, т.е. почти ничего не было написано успешно. В качестве альтернативы я попытался записать вывод напрямую, используя VTK.

shepard_data = shepard.GetOutput()
shepard_grid = vtk.vtkImageToStructuredGrid()
shepard_grid.SetInputData(shepard_data)
shepard_grid.Update()

writer = vtk.vtkStructuredGridWriter()
writer.SetFileName(file_out)
writer.SetInputData(shepard_grid.GetOutput())
writer.Write()

Это дало тот же результат, что и раньше.

С помощью ProbeFilter

Я попробовал выше, используя ProbeFilter вместо этого (с обоими преобразованиями в numpy и писать прямо). К сожалению, результат был такой же, как указано выше.

mesh = vtk.vtkStructuredGrid()
mesh.SetDimensions(dimensions)

probe = vtk.vtkProbeFilter()
probe.SetInputData(mesh)
probe.SetSourceData(output)
probe.Update()
probe_out = probe.GetOutput()

writer = vtk.vtkStructuredGridWriter()
writer.SetFileName(file_out)
writer.SetInputData(probe.GetOutput())
writer.Write()

probe_data = probe.GetOutput().GetPointData().GetArray(0)
probe_numpy = numpy_support.vtk_to_numpy(probe_data)
probe_numpy = probe_numpy.reshape(dimensions[0],dimensions[1],dimensions[2])
probe_img = sitk.GetImageFromArray(probe_numpy)

sitk.WriteImage (probe_img, file_out)

Однако, это, казалось, не дало жизнеспособного результата (vtkStructuredGridWriter произвел пустой файл, и probe_numpy был пуст).

Изменение вывода ParaView

Мои исходные данные поступают из файла structdGrid .vtk, который я открываю с помощью ParaView, а затем обрезаю для удаления структур, которые не требуются в сетке. Сохранение выходных данных сохраняет неструктурированную сетку, и я не смог выяснить, могу ли я это изменить, и, во-первых, избежать этого беспорядка!

Источник

1 ответ

Просто используйте фильтр Resample With Dataset в ParaView.

  • Открыть ParaView
  • Откройте файл структуры StructuredGrid с геометрией, которую вы хотите иметь
  • Откройте файл UnstructuredGrid
  • Добавить фильтр "Resample with dataset"
  • Выберите структурированные данные в качестве источника ввода
  • Применять
Источник
Другие вопросы по тегам