CadQuery中圆锥减切半球时的异常问题求助

我需要从球壳（由半径10和8的球体差集得到）中减切4个分别沿(1,0,0)、(-1,0,0)、(0,1,0)、(0,-1,0)方向且相互交叉的圆锥。但交叉时其中一个部分无法被减切，非交叉圆锥能得到预期结果；替换为圆柱减切则一切正常。相关代码如下：

import math
import cadquery as cq
from ocp_vscode import show_object, reset_show, set_defaults

def align(obj, dir):
    l = math.sqrt(dir[0]**2+dir[1]**2+dir[2]**2); ## radial distance
    b = math.degrees(math.acos(dir[2]/l));        ## inclination angle
    c = math.degrees(math.atan2(dir[1],dir[0]));  ## azimuthal angle
    return obj.rotate((0,0,0),(1,0,0),0).rotate((0,0,0),(0,1,0),b).rotate((0,0,0),(0,0,1),c)

def cone(height, angle):
    radius = height*math.tan(math.radians(angle))
    cone = cq.Solid.makeCone(0, radius, height)
    return cone

def cylinder(height, angle):
    radius = height*math.tan(math.radians(angle))
    cylinder = cq.Solid.makeCylinder(radius, height)
    return cylinder

# BUILD
# ----------------------------------------------------------------
reset_show()
angle = 30

def piece(listCone):
    ## construct general parts
    shellOuter  = cq.Solid.makeSphere(10).cut(cq.Solid.makeSphere(8))
    pieceOuter = shellOuter
    ## generate the piece
    for element in listCone:
        ## outer shell
        obj = cone(10.01,angle)
        obj = align(obj, element)
        pieceOuter = pieceOuter.cut(obj)
    return pieceOuter

center = piece([
    (0,1,0),
    (0,-1,0),
    (1,0,0),
    (-1,0,0)
    ])
show_object(center, options=dict(alpha=0.5,color='red'))

核心是圆锥顶点重合导致的布尔运算精度问题：

• 默认创建的圆锥顶点在原点，旋转后四个圆锥的顶点都落在球壳中心位置。
• OpenCASCADE（CadQuery依赖的几何内核）处理顶点完全重合的布尔减切时，容易出现边界判定错误，导致其中一个圆锥的减切失效。
• 圆柱无此问题是因为它没有尖锐顶点，布尔运算的边界判定更稳定。

通过让圆锥顶点避开球壳中心，就能解决这个问题，有两种简单实现方式：

方式一：平移圆锥避开中心 #

在piece函数的循环中，给旋转后的圆锥加一段沿自身方向的微小平移：

obj = align(obj, element).translate(element)  # 沿圆锥方向平移单位向量距离，让顶点离开中心

方式二：修改圆锥生成逻辑 #

调整cone函数，让圆锥顶点初始就偏移一段微小距离：

def cone(height, angle):
    radius = height*math.tan(math.radians(angle))
    # 让圆锥顶点沿Z轴偏移0.01，避免在原点
    cone = cq.Solid.makeCone(0, radius, height).translate((0, 0, 0.01))
    return cone

另外可以优化align函数，用更高效的向量旋转方式替代多次旋转：

def align(obj, dir):
    dir_vec = cq.Vector(dir).normalized()
    z_axis = cq.Vector(0,0,1)
    if dir_vec.isParallelTo(z_axis):
        return obj
    rotation_axis = z_axis.cross(dir_vec)
    rotation_angle = math.degrees(z_axis.getAngle(dir_vec))
    return obj.rotate((0,0,0), rotation_axis, rotation_angle)

完整修改后的代码 #

import math
import cadquery as cq
from ocp_vscode import show_object, reset_show, set_defaults

def align(obj, dir):
    dir_vec = cq.Vector(dir).normalized()
    z_axis = cq.Vector(0,0,1)
    if dir_vec.isParallelTo(z_axis):
        return obj
    rotation_axis = z_axis.cross(dir_vec)
    rotation_angle = math.degrees(z_axis.getAngle(dir_vec))
    return obj.rotate((0,0,0), rotation_axis, rotation_angle)

def cone(height, angle):
    radius = height*math.tan(math.radians(angle))
    cone = cq.Solid.makeCone(0, radius, height)
    return cone

def cylinder(height, angle):
    radius = height*math.tan(math.radians(angle))
    cylinder = cq.Solid.makeCylinder(radius, height)
    return cylinder

# BUILD
# ----------------------------------------------------------------
reset_show()
angle = 30

def piece(listCone):
    shellOuter  = cq.Solid.makeSphere(10).cut(cq.Solid.makeSphere(8))
    pieceOuter = shellOuter
    for element in listCone:
        obj = cone(10.01,angle)
        obj = align(obj, element).translate(element)
        pieceOuter = pieceOuter.cut(obj)
    return pieceOuter

center = piece([
    (0,1,0),
    (0,-1,0),
    (1,0,0),
    (-1,0,0)
    ])
show_object(center, options=dict(alpha=0.5,color='red'))

内容的提问来源于stack exchange，提问作者Pa Dalton