объединить несколько устаревших ASCII .vtk файлов в один файл, используя python-vtk

Question

во время моих расчетов я сохраняю кучу файлов .vtk, которые связаны с одним временным шагом. Каждый файл описывает многогранник (класс C ++, который я определил), используя формат файла POLYDATA ASCII. Многогранники записываются простой функцией-членом класса многогранников.

Чтобы не определять совершенно новый класс для коллекции многогранников, который мне нужен для визуализации и путаницы в моем коде C ++, я хочу объединить несколько файлов .vtk в один файл .vtk.

Использование python-vtk вызывает у меня некоторые проблемы:

from vtk import * 

reader = vtkPolyDataReader()
reader.SetFileName("file1.vtk")

reader.Update()
polyData1 = reader.GetOutput()

reader.SetFileName('file2.vtk')

reader.Update()
polyData2 = reader.GetOutput()


# Expand the output points

points1 = polyData1.GetPoints()

points2 = polyData2.GetPoints()

insertPosition = points1.GetNumberOfPoints()

for i in xrange(points2.GetNumberOfPoints()):
    insertPoint = points2.GetPoint(i)
    points1.InsertPoint(insertPosition, 
                        insertPoint[0], insertPoint[1], insertPoint[2])
    insertPosition += 1

print points1.GetNumberOfPoints()

# Change the cell ids of every cell in the polydata2 to correspond with 
# the new points (appended point array)

increment = points1.GetNumberOfPoints();

for i in xrange(polyData2.GetNumberOfCells()):
    cell  = polyData2.GetCell(i)
    cellIds = cell.GetPointIds()
    for j in xrange(cellIds.GetNumberOfIds()):
        oldId = cellIds.GetId(j)
        cellIds.SetId(j, oldId + increment)


polyData1.Allocate(polyData1.GetNumberOfCells(), 1)

for i in xrange(polyData2.GetNumberOfCells()):
    cell = polyData2.GetCell(i)
    polyData1.InsertNextCell(cell.GetCellType(), cell.GetPointIds())


writer = vtkPolyDataWriter()
writer.SetFileName("output.vtk")
writer.SetInput(polyData1)
writer.Write()

Делая это, я получу дубликаты очков, и это нормально. Проблема заключается в том, что этот скрипт выполняется для следующих файлов .vtk:

file1:

# vtk DataFile Version 2.0
surface written 2011-12-19T15:30:18
ASCII

DATASET POLYDATA
POINTS 8 float
0.48999999999999999112 0.4000000000000000222 0.5999999999999999778
0.48999999999999999112 0.5 0.5999999999999999778
0.48999999999999999112 0.5 0.69999999999999995559
0.48999999999999999112 0.4000000000000000222 0.69999999999999995559
0.5 0.5 0.5999999999999999778
0.5 0.5 0.69999999999999995559
0.5 0.4000000000000000222 0.69999999999999995559
0.5 0.4000000000000000222 0.5999999999999999778

POLYGONS 6 30
4 0 1 2 3 
4 4 5 6 7 
4 4 1 2 5 
4 6 5 2 3 
4 7 0 1 4 
4 0 7 6 3 

CELL_DATA 6
FIELD attributes 1
zone 1 6 float
1 1 1 1 1 1

file2:

# vtk DataFile Version 2.0
surface written 2011-12-19T15:30:18
ASCII

DATASET POLYDATA
POINTS 8 float
0.58999999999999996891 0.5999999999999999778 0.5
0.58999999999999996891 0.69999999999999995559 0.5
0.58999999999999996891 0.69999999999999995559 0.5999999999999999778
0.58999999999999996891 0.5999999999999999778 0.5999999999999999778
0.5999999999999999778 0.69999999999999995559 0.5
0.5999999999999999778 0.69999999999999995559 0.5999999999999999778
0.5999999999999999778 0.5999999999999999778 0.5999999999999999778
0.5999999999999999778 0.5999999999999999778 0.5

POLYGONS 6 30
4 0 1 2 3 
4 4 5 6 7 
4 4 1 2 5 
4 6 5 2 3 
4 7 0 1 4 
4 0 7 6 3 

CELL_DATA 6
FIELD attributes 1
zone 1 6 float
1 1 1 1 1 1

Результат: точки с 8 координатами и ячейки (грани) вообще не добавляются.

Это только у меня или у Python-обёрток для vtkArray и подобных vtkObjects нет возможности перебора?

Answer 1

Я отказался от API Python-VTK для этой цели, вот небольшой быстро написанный класс, который делает эту работу, кто-то может найти это полезным:

class polyDataVtk(object):
    """Class representing the polydata vtk information stored in legacy ASCII .vtk POLYDATA files."""
    def __init__(self, fileName = None):
        self.__points = []
        self.__polygons = []
        if fileName is not None:
            self.__fileName = fileName 

    def parse(self, fileName):
        """Parse the POLYDATA information from a .vtk file and append the data to the object.
           Does not check for the file consistency."""
        file = open(fileName, 'r')

        # Use local data first. 
        points = []
        polygons = []

        line = ""
        while(True):
            line = file.readline()
            if 'POINTS' in line:
                break

        nPoints = 0

        if (line == ""):
            print "No POINTS defined in the .vtk file"
            return

        # Set the number of points
        nPoints = int(line.split()[1])

        # Append the numbers.
        for i in xrange(nPoints):
            points.append(map(lambda x : float(x), file.readline().split()))

        # Append polygons.
        line = ""

        while(True):
            line = file.readline()
            if 'POLYGONS' in line:
                break
        if (line == ""):
            print "No POLYGONS defined in the .vtk file"
            return
        # Set the number of polygons.
        nPolygons = int(line.split()[1])

        # Read the polygons.
        for i in xrange(nPolygons):
            line = file.readline()
            polygons.append(map(lambda x : int(x) + len(self.__points), line.split())[1:])

        # File parsed without a problem.
        self.__points.extend(points)
        self.__polygons.extend(polygons)
        file.close()

    def write(self,fileName=None, append=False):
        # Overwrite the filename provided to the constructor.
        if fileName is not None:
            self.__fileName = fileName

        # No fileName is provided by the constructor or the write method.
        if self.__fileName is None:
            self.__fileName = "result.vtk"

        # Append or overwrite?
        if append:
            file = open(self.__fileName, 'a')
        else:
            file = open(self.__fileName, 'w')

        file.writelines("# vtk DataFile Version 2.0\n")
        file.writelines("appended vtk files\n")
        file.writelines("ASCII\n")
        file.writelines("DATASET POLYDATA\n")

        file.writelines("POINTS %d float \n" % len(self.__points))

        for point in self.__points:
            file.writelines("%.10f %.10f %.10f\n" % (point[0], point[1], point[2]))

        size = 0
        for polygon in self.__polygons:
            size += len(polygon)
        file.writelines("POLYGONS %d %d \n" % (len(self.__polygons), 
                                               size 
                                             + len(self.__polygons)))

        for polygon in self.__polygons:
            file.writelines("%d " % len(polygon))
            for label in polygon:
                # I don't know how many labels are there in general.
                file.writelines("%d " % label)
            file.writelines("\n")

        file.close()

    def fileName(self):
        return self.__fileName

Answer 2

Вы должны рассмотреть возможность использования vtkAppendPolyData для объединения файлов.

from vtk import * 

reader = vtkPolyDataReader()
append = vtkAppendPolyData()

filenames = ['file1.vtk', 'file2.vtk']
for file in filenames:
    reader.SetFileName(file)
    reader.Update()
    polydata = vtkPolyData()
    polydata.ShallowCopy(reader.GetOutput())
    append.AddInputData(polydata)

append.Update()    

writer = vtkPolyDataWriter()
writer.SetFileName('output.vtk')
writer.SetInput(append.GetOutput())
writer.Write()

Также обратите внимание, что при повторном вызове считывателя выходной набор данных из первого входного файла будет перезаписан вторым входным файлом. Вы должны сделать ShallowCopy (как в сценарии выше), если вы хотите использовать тот же читатель.