zippy84 / vtkbool

A new boolean operations filter for VTK
Apache License 2.0
166 stars 40 forks source link

cylinder and sphere cannot be merged #82

Open zippy84 opened 2 months ago

zippy84 commented 2 months ago

@adelka-m Open your own issue, the next time.

I am trying to combine meshes with this code:

prod1 = vtk.vtkTrivialProducer()
prod1.SetOutput(mesh1)
prod1.Update()

booleanFilter = vtkPolyDataBooleanFilter()
booleanFilter.SetInputConnection(0, prod1.GetOutputPort())

prod2 = vtk.vtkTrivialProducer()
prod2.SetOutput(mesh2)
prod2.Update()
booleanFilter.SetInputConnection(1, prod2.GetOutputPort())

booleanFilter.SetOperModeToUnion()
booleanFilter.Update()
mesh = booleanFilter.GetOutput()
mesh.GetNumberOfPoints()

I create these meshes, so that I create a cylinder and a ball and at one of the ends of the cylinder and combine it. Every of these meshes is generated with the same code. So I really do not understand, why it works for the first two but not for the last cylinder. I show also the yellow lines from which is the cylinder object defined.

This code works for these two meshes. mesh.GetNumberOfPoints() > 0

But if I run it for these two meshes it returns an empty mesh. The only difference of the meshes are that there is a part of the ball sticking out on the other side of cylinder. mesh.GetNumberOfPoints() = 0

I even created the diameter of the meshes smaller, so there is no "part of the ball" sticking out from the other side of the mesh. But the code still does not work.

Interestingly, if I run the code for these meshes now its working (with the ball sticking out, it was not working for this pair, when the ball was sticking out at first).

zippy84 commented 2 months ago
#!/usr/bin/env python
# *-* coding: UTF-8 *-*

import numpy as np

import sys
sys.path.extend(['/home/zippy/vtkbool/build/lib/python3.12/site-packages/vtkbool'])

from vtkmodules.vtkIOLegacy import vtkPolyDataReader, vtkPolyDataWriter
from vtkmodules.vtkIOGeometry import vtkSTLReader, vtkSTLWriter
from vtkmodules.vtkFiltersCore import vtkCleanPolyData, vtkPolyDataNormals, vtkTriangleFilter
from vtkmodules.vtkCommonExecutionModel import vtkTrivialProducer
from vtkmodules.vtkFiltersSources import vtkCubeSource, vtkCylinderSource, vtkSphereSource, vtkPolyLineSource
from vtkmodules.vtkFiltersGeneral import vtkTransformPolyDataFilter
from vtkmodules.vtkCommonTransforms import vtkTransform
from vtkmodules.vtkCommonDataModel import vtkPolyData, VTK_POLYGON
from vtkmodules.vtkCommonCore import vtkObject

from vtkBool import vtkPolyDataBooleanFilter

def createCombinedMesh(point1, point2, radius=20):
    # vtkObject.GlobalWarningDisplayOff() # silence the warnings :D
    cylinder = createCylinderFromPoints(point1, point2, radius, returnFilter=False, returnTriangulated=True)
    prodCyl = vtkTrivialProducer()
    prodCyl.SetOutput(cylinder)
    sphere   = createSphereFromPoint(point2, radius + 1e-2) #tolerance factor, without it also works but it gives warnings
    # Combine the cylinder and sphere into a single mesh
    booleanFilter = vtkPolyDataBooleanFilter()
    booleanFilter.SetInputConnection(0, sphere.GetOutputPort())
    booleanFilter.SetInputConnection(1, prodCyl.GetOutputPort())
    booleanFilter.SetOperModeToUnion()
    booleanFilter.Update()
    return createOneComponentMesh(booleanFilter.GetOutput())

def createCylinderFromPoints(point1, point2, radius=20, returnFilter=False, returnTriangulated=False):
    # Create the cylinder
    cylinder = vtkCylinderSource()
    cylinder.SetRadius(radius)
    cylinder.SetResolution(50)
    direction = point2 - point1
    midpoint = (point1 + point2) / 2.0
    length = np.linalg.norm(direction)
    cylinder.SetHeight(length)
    transform = vtkTransform()
    transform.Translate(midpoint)
    # Calculate the rotation to align with the direction vector
    if np.linalg.norm(direction) != 0:
        norm_dir = direction / np.linalg.norm(direction)
        angle = np.arccos(np.dot(norm_dir, [0, 1, 0]))
        axis = np.cross([0, 1, 0], norm_dir)
        if np.linalg.norm(axis) != 0:
            axis = axis / np.linalg.norm(axis)
            transform.RotateWXYZ(np.degrees(angle), axis)
    transformFilter = vtkTransformPolyDataFilter()
    transformFilter.SetTransform(transform)
    transformFilter.SetInputConnection(cylinder.GetOutputPort())
    transformFilter.Update()
    if returnFilter:
        return transformFilter
    elif returnTriangulated:
        return triangulate(transformFilter.GetOutput())
    else:
        return transformFilter.GetOutput()

def createSphereFromPoint(point, radius=20):
    # Create the sphere
    sphere = vtkSphereSource()
    sphere.SetCenter(point)
    sphere.SetRadius(radius)
    sphere.SetPhiResolution(50)
    sphere.SetThetaResolution(50)
    sphere.Update()
    return sphere

def triangulate(inputData):
    triangleFilter = vtkTriangleFilter()
    triangleFilter.SetInputData(inputData)
    triangleFilter.Update()
    return triangleFilter.GetOutput()

def createOneComponentMesh(mesh):
    points = mesh.GetPoints()
    cells  = mesh.GetPolys()
    newMesh = vtkPolyData()
    newMesh.SetPoints(points)
    newMesh.SetPolys(cells)
    return newMesh

list = [
    [np.asarray([120.57039642,  57.58803177, 117.53710938]), np.asarray([180.38616943,  83.88054657,  87.29520416]) ],
    # [np.asarray([160.03898621,  58.91154099,  56.11118698]), np.asarray([156.74688721,  89.19313812,  80.92210388]) ]
]

combined = []
for l in list:
    com = createCombinedMesh(l[0], l[1])
    combined.append(com)
adelka-m commented 2 months ago

@adelka-m Open your own issue, the next time. Hi @zippy84, sorry, I felt it is the very similar issue, so I commented there :)

zippy84 commented 2 months ago

The boolean operation between the sphere and the cylinder will work when you change the output precision to double:

from vtkmodules.vtkCommonExecutionModel import vtkAlgorithm

# within createCylinderFromPoints
cylinder.SetOutputPointsPrecision(vtkAlgorithm.DOUBLE_PRECISION)

# within createSphereFromPoint
sphere.SetOutputPointsPrecision(vtkAlgorithm.DOUBLE_PRECISION)