3

I need to display the trajectories in space of 50 particles using vtk and paraview. Currently my data is pos(x,y,t,n) where n is the label of the n-th particle. I have saved all my data in a vtk file which is organized as:

# vtk DataFile Version 3.0
VTK from Matlab
BINARY

DATASET POLYDATA
POINTS 199250 double
[PROPERLY ORGANIZED BINARY COORDINATES OF THE 199250 POINTS]

The above file can be imported correctly into ParaView, and at a first glance I can apply the following programmagle filter

pdi = self.GetPolyDataInput()
pdo =  self.GetPolyDataOutput()
newPoints = vtk.vtkPoints()
numPoints = pdi.GetNumberOfPoints()

for i in range(0, numPoints):
    coord = pdi.GetPoint(i)
    x, y, z = coord[:3]
    x = x * 0.03
    y = y * 0.03
    z = z * 2
    newPoints.InsertPoint(i, x, y, z)

pdo.SetPoints(newPoints)
aPolyLine = vtk.vtkPolyLine()
aPolyLine.GetPointIds().SetNumberOfIds(numPoints)

for i in range(0,numPoints):
 aPolyLine.GetPointIds().SetId(i, i)

pdo.Allocate(1, 1)
pdo.InsertNextCell(aPolyLine.GetCellType(), aPolyLine.GetPointIds())

to create the following graphic

enter image description here

The filter basically creates a line connecting every other point, but in my case it connects the last point of one particle trajectory to the first one of the following trajectory.

I would like to create a different polydata object for every particle, so I could even label or color them separately (and eliminate the "all-connected" issue). My first guess for a similar Filter would be

pdi = self.GetPolyDataInput()
pdo =  self.GetPolyDataOutput()
nParts = 50
newPoints = vtk.vtkPoints()
numPoints = pdi.GetNumberOfPoints()
numPointsPart = int(numPoints/nParts)

for i in range(0, numPoints):
    coord = pdi.GetPoint(i)
    x, y, z = coord[:3]
    x = x * 0.03
    y = y * 0.03
    z = z * 2
    newPoints.InsertPoint(i, x, y, z)

pdo.SetPoints(newPoints)
pdo.Allocate(nParts, 1)

for i in range(0,nParts):
  aPolyLine = vtk.vtkPolyLine()
  aPolyLine.GetPointIds().SetNumberOfIds(numPointsPart)
  indStart=int(i*numPointsPart)
  indEnd=int((i+1)*numPointsPart)
  for i in range(indStart,indEnd):
    aPolyLine.GetPointIds().SetId(i, i)
  pdo.InsertNextCell(aPolyLine.GetCellType(), aPolyLine.GetPointIds())

However, when I run this, ParaView just closes due to segfault. Bonus: the code is written correctly: if I set nParts to 1, I get the same graph as above.

EDIT : the dataset can be found here

1 Answers1

3

The problem is in the iteration logic:

for i in range(indStart, indEnd):
    aPolyLine.GetPointIds().SetId(i, i)

Must be changed to the following:

for i in range(indStart, indEnd):
    aPolyLine.GetPointIds().SetId(i - indStart, i)

The original code was calling SetId with an incorrect index. Also, note the shadowed variable i. Although that's not what's causing the issue, if you use i after the second for loop, the result may not be what you'd be expecting.

Utkarsh
  • 1,492
  • 1
  • 11
  • 19
  • Thanks Utkarsh, the shadowed i was a transcription typo and it's right in my code (tuh! I feel embarrassed!)... Anyway, I think I'm not really getting the SetId working. Could you please provide me some example of how I have to set the indexes? – Juan Sebastian Totero Sep 21 '14 at 09:02
  • 1
    It's an idlist http://www.vtk.org/doc/nightly/html/classvtkIdList.html. For 'SetId' The first value is index, the second is the index for the point index for the point that forms that cell. – Utkarsh Sep 22 '14 at 13:30
  • 1
    oh I understand now. So if I want to set my the second point of my polyline to be the tenth point in the points list I have to set SetId(1,9) right? – Juan Sebastian Totero Sep 22 '14 at 16:24