Get the most separated points

Hello,

Given a cloud of points, what would be the smartest way to determine the two points that are most widely separated?
get_max_distance

It is possible to loop through all points and calculate the distances between points so as to get the maximum, it would require n x n iterations. I wish to know if there exists a class that can do that, preferably not Python exclusive.

Thank you.

This is a classic computational geometry problem with some smart solutions. I would expect that you do not have to implement one, but you can find free, open- source implementations in Python and C++.

Thank you for the link.

vtkSphere offers to compute a bounding sphere. The problem is that it requires a hints[2] field, which is a guess of the 2 most distant points. Random pairs give different results, so getting the exact pair of points is my purpose. Iā€™ll delve further into this slowly.

Implementation in VTK seems to be a standard one (same as in three.js) that computes a bounding sphere with a simple method, which is not guaranteed to be minimal, but it is close. Hints are not required. If you donā€™t provide hints then the bounding box extreme points are used to initialize the computation.

If you only need an approximate binding sphere (that can be a few percent larger than needed) then this should work well. If you need the optimal sphere then it will require a more complex, more costly computation.

What are you planning to use this bounding sphere for?

Itā€™s intended to add 1 more option to this module.

It needs not have an exact precision. However, the ā€˜hints[2]ā€™ parameter can be passed as ā€˜nullptrā€™ in C++, while passing ā€˜Noneā€™ fails in Python. Passing arbitrary values like ā€˜[0, 0]ā€™ is meaningless. How to handle this situation in Python?

If Python Wrappers - VTK documentation does not help then you can ask on the VTK forum.

Check 2nd questionā€™s answer here:
https://www.perplexity.ai/search/declare-a-vtkidtype-2-array-on-3XtK9t28SqKQerTfUb5X8g#1

HIH

Thanks for the link.

Unfortunately, I needed to pass a ā€˜Noneā€™ parameter that would mean ā€˜nullptrā€™, but the python wrapper always wants to pass ā€˜somethingā€™ else than 0x0 to the underlying C++ function.

It would have been nice if this C++ function signature defined ā€˜hints[2] = nullptrā€™ since it is optional. I suppose it should have worked.

Anyway, Iā€™m baking a solution that should not be resource intensive.

Agreed. You can send a merge request to VTK with this simple change. It should go through review really quickly. You can also @mention David Gobbi in the merge request, asking him about the VTK null pointer syntax in Python (maybe there is a way to specify a null pointer just not documented well).

Ok, Iā€™ll do that ASAP. Just hoping they donā€™t require a test case for such trivial change.

1 Like

I came up to this function to determine the two most distant points in a cloud of points. Itā€™s decently fast, 0.06 s for 10ā¶ points, 1.01 s for 10ā· points, 11.69 s for 10āø points and bumps to 135.67 s for 10ā¹ points.

The logic seems OK to me, dropping it here for any purpose and remaining open for comments.

def findMostDistantPoints(points : vtk.vtkPoints):
    pointData = points.GetData()
    numberOfPoints = points.GetNumberOfPoints()
    
    # Sort the points on each axis, get the bounding points.
    sorter = vtk.vtkSortDataArray()

    sorter.SortArrayByComponent(pointData, 0, 0)
    xMin = pointData.GetTuple3(0)
    xMax = pointData.GetTuple3(numberOfPoints - 1)

    sorter.SortArrayByComponent(pointData, 1, 0)
    yMin = pointData.GetTuple3(0)
    yMax = pointData.GetTuple3(numberOfPoints - 1)

    sorter.SortArrayByComponent(pointData, 2, 0)
    zMin = pointData.GetTuple3(0)
    zMax = pointData.GetTuple3(numberOfPoints - 1)

    boundingPoints = vtk.vtkPoints()
    boundingPoints.InsertNextPoint(xMin)
    boundingPoints.InsertNextPoint(xMax)
    boundingPoints.InsertNextPoint(yMin)
    boundingPoints.InsertNextPoint(yMax)
    boundingPoints.InsertNextPoint(zMin)
    boundingPoints.InsertNextPoint(zMax)
    
    # Calculate the distances between each pair of points: 6 x 6 iterations.
    distances = vtk.vtkDoubleArray()
    distances.SetNumberOfComponents(7)
    numberOfBoundingPoints = boundingPoints.GetNumberOfPoints()
    for i in range(numberOfBoundingPoints):
        referencePoint = boundingPoints.GetPoint(i)
        for j in range(numberOfBoundingPoints):
            nextPoint = boundingPoints.GetPoint(j)
            distance2 = vtk.vtkMath.Distance2BetweenPoints(referencePoint, nextPoint)
            distances.InsertNextTuple([referencePoint[0], referencePoint[1], referencePoint[2], nextPoint[0], nextPoint[1], nextPoint[2], distance2])
    
    # Sort the bounding points in descending order by distance squared.
    distanceSorter = vtk.vtkSortDataArray()
    distanceSorter.SortArrayByComponent(distances, 6, 1)
    
    # Get the most distant points.
    mostDistantTuple = distances.GetTuple(0)
    p1 = [mostDistantTuple[0], mostDistantTuple[1], mostDistantTuple[2]]
    p2 = [mostDistantTuple[3], mostDistantTuple[4], mostDistantTuple[5]]
    
    # Return as vtkPoints.
    result = vtk.vtkPoints()
    result.InsertNextPoint(p1)
    result.InsertNextPoint(p2)

    return result

Computation is quite fast and the method is simple, but I donā€™t think it gives you the two farthest points (when the two farthest points are near the corners of the bounding box, in two opposite corners).

Instead of reinventing the wheel, you can use a free, open-source, restriction-free implementation like this.

Thank you for the interesting link and comment, they highlight situations termed as ā€˜degenerateā€™ which have not been accounted for.

@lassoan

vtkSphere::ComputeBoundingSphere() has now an overload in VTK that does not mandate a hint parameter. Thank you for your support.

1 Like

This implementation should overcome these limitations.

@lassoan

How does it work? Have you you implemented a known algorithm, or you invented something new? If it is something new, how it compares to known algorithms?

On test data, it seems to be accurate. Iā€™ll compare when VTK 8fcf3271ab is available in Slicerā€™s VTK. I just wanted to propose another possible solution to a question I opened.

1 Like