Segfault when calling StartAppend() in vtkUnsignedDistance()

[fraseyboo]

So I’m trying to use the Append() function in vtkSignedDistance / vtkUnsignedDistance which allows for polydata to be incrementally added which is supposedly supported from the documentation, when initially calling StartAppend() the code immediately segfaults and displays:

Segmentation fault: 11

I initially thought this was due to the point locator used being static and unable to handle the insertion of new data so I tried a variety of locators which should work but the issue persists.

I’ve produced a minimum working example which describes the issue, the code runs through the singular polydata example fine but breaks on the example with multiple polydata.

Are there some additional steps I’m missing here?

import vtk
import numpy as np


def generate_random_polydata(total_points=1000):

    points = np.random.rand(total_points,3)
    vtk_points = vtk.vtkPoints()
    
    for point in points:
        vtk_points.InsertNextPoint(point)
    vtk_polydata = vtk.vtkPolyData()
    vtk_polydata.SetPoints(vtk_points)
    vtk_polydata.Modified()

    return vtk_polydata


def generate_distance_field(bounds=[-1.5,1.5,-1.5,1.5,-1.5,1.5], voxel_spacing=[0.1,0.1,0.1], total_polydata=4):

    signed_distance = vtk.vtkSignedDistance()

    dimensions = [
        int((bounds[1] - bounds[0]) / voxel_spacing[0]),
        int((bounds[3] - bounds[2]) / voxel_spacing[1]),
        int((bounds[5] - bounds[4]) / voxel_spacing[2])
    ]


    array_size = [
        (bounds[1] - bounds[0]),
        (bounds[3] - bounds[2]),
        (bounds[5] - bounds[4])
    ]

    signed_distance = vtk.vtkUnsignedDistance() # signed distance functions require normals

    # signed_distance = vtk.vtkSignedDistance()

    signed_distance.SetBounds(bounds)
    
    signed_distance.SetDimensions(dimensions)
    signed_distance.SetRadius(np.linalg.norm(array_size))


    # These parameters do not change the error

    # signed_distance.AdjustBoundsOff()
    # signed_distance.CappingOff()


    # Changing the locator to an incremental one does not change the error

    # locator = signed_distance.GetLocator()

    # new_locator = vtk.vtkNonMergingPointLocator()
    # new_locator = vtk.vtkPointLocator()
    # new_locator = vtk.vtkIncrementalOctreePointLocator()
    # new_locator = vtk.vtkNonMergingPointLocator()
    # new_locator = vtk.vtkMergePoints()

    # signed_distance.SetLocator(new_locator)



    signed_distance.Modified()

    locator = signed_distance.GetLocator()

    # print(locator)

    if total_polydata > 1:

        print('Case with multiple polydata')

        signed_distance.StartAppend() #required to call before Append(), but just segfaults (11)
        print('Code breaks in the line above here')

        for _ in range(total_polydata):
            test_polydata = generate_random_polydata()
            signed_distance.Append(test_polydata)
        signed_distance.EndAppend()


    else:

        print('Case with single polydata')

        polydata = generate_random_polydata()
        signed_distance.SetInputData(polydata)


    signed_distance.Update()
    distance_field = signed_distance.GetOutput()

    print('Distances successfully generated\n')

    # locator = signed_distance.GetLocator()

    # print(locator)

    return distance_field

Many thanks.

[Jens_Munk_Hansen]

The vtkSignedDistance is not so user friendly. One need to set the radius, the dimensions and bounds. Then call Update(), StartAppend() and Append(polydata) (one or more time). Finally, you call EndAppend.

We should fix this error in the source code and also provide an example where Append is used.

[fraseyboo]

I can confirm that calling Update() before StartAppend() seems to be sufficient for the code to work. Hopefully some modification can be made to StartAppend() in the future to check that the parameters set in Update() are applied without it needed to be called, or a section in the documentation to clarify the behaviour.

If a new example were to be made it’d be nice to see the relative performance of all the point locators in the computation.

[fraseyboo]

For completeness on this issue it is also worth mentioning that the distance array produced with this method is not directly equivalent to running the function on pre-combined data. It appears that the distance computation is run sequentially and will override older data even if the older data had smaller values.

You can see the effect in this screenshot, the distances from the Polydata on the right overrode the older data from the left polydata.

The more ‘typical’ behaviour would be to only replace distances in which the distance is smaller than the previous value, as seen in this screenshot. AFAIK there isn’t any functionality to do this with the current code, it’d make the function far more useful though.

[Jens_Munk_Hansen]

I would guess, the typicaly behavior is to always update an average value and sometimes with a weighting function

[fraseyboo]

As far as I can tell the algorithm just does a direct override of all values in the image that are within the user-specified radius of the polydata, both behaviours have their usage (e.g the current behaviour would be useful for transient objects where overriding previous data is important, the second behaviour is more useful for observations of static objects like LIDAR scanning).

I was initially a little confused with the behaviour as I had the radius parameter set to a large value which meant that only the most recent object was shown affecting the field, with some tuning I got the first screenshot but ideally I’d prefer the second.

Looking at the source code I believe the behaviour could be modified at line 110 of the source code, behaviour could be chosen from the following:

• A direct override if a points within the radius are found (current mode of operation).
• A conditional override if the new value is smaller than the current one.
• A weighted average of the new value and the older value ((n1 * v1) + (n2 * v2) / (n1 + n2)), however this would require another array to store the value of n1.

[cory.quammen (Cory Quammen (Kitware))]

@will.schroeder do you recall if the behavior of overwriting the previously computed distance field was intended in this filter, or was it an oversight that the distance to the appended object isn’t compared to the existing distance field value value?

[will.schroeder]

As Fraseyboo suggests, the filter should support set operations (union, intersect, etc). This is consistent with other classes in VTK (e.g., see vtkGaussianSplatter and AccumulationMode). The class implementation seems to be in a prototype state, and needs polish.