1

我想在街道（平面）上创建一个坑洼（椭圆体）。经过一些测试，我意识到布尔联合可能是使用的那个（因为我想结合两个网格但不希望街道掩盖坑洼）。
当我使用 Sphere 尝试此操作时，它工作得非常好。
但是，我使用的是 ParameticEllipsoid。在那里我经常收到错误“vtkMath::Jacobi: Error extracting eigenfunctions”，并且想要的效果只适用于布尔差异而不是联合。我想它必须与 vtk/py 及其几何对象创建有关。

虽然采用布尔差异带来了想要的结果，但它需要很多时间（我不会太介意）并且很少工作（并且只有在椭圆体不与其他物体啮合的情况下）。

有没有办法避免这种情况？我的想法是提取网格的点（NumPy 数组）并计算与它们的联合，但我无法做到。

def add_pothole():
    street_mesh = pv.Plane(i_size=100, j_size=100, i_resolution=100, j_resolution=100)

    points = street_mesh.points
    for i in range(4):
        pothole_mesh = pv.ParametricEllipsoid(10,5,3)
        # alternatively to get the half ellipsoid: (but boolean difference seemed to only work with the closed one
        # pothole_mesh = pv.ParametricEllipsoid(10,5,3, max_v=pi /2)
        pothole_mesh.translate(choice(points))
    
        # "add" the pothole to the street
        street_mesh = pothole_mesh.boolean_difference(street_mesh.triangulate())
    
    street_mesh.plot()

或者，我考虑直接定义椭圆形凹陷，使用高度作为平面内位置的函数。是否可以在不使用布尔运算的情况下做到这一点？

4

1 回答 1

1

为了解决您当前的问题，您当前的代码存在两个问题：

  1. 不像Sphere,ParametricEllipsoid由于 VTK 构造它的方式，有很多杂散点，所以最终结果在技术上并不是一个封闭的表面。这是我计划在 pyvista 中修复的一个已知问题。现在你可以通过调用clean(tolerance=1e-5)你的椭球来解决这个问题。这应该阻止这些vtkMath错误弹出。
  2. 从坑洼中减去街道将在每次调用 后反转飞机的法线boolean_difference，这将使您在飞机的交替两侧产生颠簸，这可能不是您想要的。street_mesh.flip_normals()您可以通过在循环结束时调用来解决此问题。

通过这两个更改，您的代码将运行，尽管速度很慢（我真的不知道为什么，但我也不知道布尔运算如何在幕后工作）：

from random import choice
import numpy as np
import pyvista as pv

def add_pothole():
    street_mesh = pv.Plane(i_size=100, j_size=100, i_resolution=100, j_resolution=100)
    pothole_template = pv.ParametricEllipsoid(10, 5, 3).clean(tolerance=1e-5).triangulate(

    points = street_mesh.points
    for i in range(4):
        pothole_mesh = pothole_template.copy()
        pothole_mesh.translate(choice(points))

        # "add" the pothole to the street
        street_mesh = pothole_mesh.boolean_difference(street_mesh.triangulate())
        street_mesh.flip_normals()

    street_mesh.plot()

add_pothole()

z(x, y)但是要解决您的根本问题，如果您的街道一开始就是一个表面，那么您不必首先这样做。考虑到 Perlin 噪声，您可能从一个平面开始，其中标量从 Perlin 噪声中采样，然后调用warp_by_scalar()以实现噪声调制。好吧，您可以将椭圆形的凸块添加到相同的标量，然后在最后一步进行变形（与叠加多个 Perlin 噪声样本的方式完全相同）。

为此，您必须计算z(x, y)椭球的参数化。这不是微不足道的，但也不是很困难。可能更广为人知的是，R以原点为半径的球体方程为

x^2 + y^2 + z^2= R^2.

当除以时，这是更好的（和规范的）R^2：

(x/R)^2 + (y/R)^2 + (z/R)^2 = 1.

你得到一个椭球的方法是线性变换每个轴，改变它们各自的半径（所以它们不再是 all R）。这是椭圆体的隐式方程：

(x/a)^2 + (y/b)^2 + (z/c)^2 = 1.

如果你想为椭球有一个重要的中心，你必须翻译每个坐标：

(x - x0)^2/a^2 + (y - y0)^2/b^2 + (z - z0)^2/c^2 = 1.

这是 3 维椭球的一般形式，假设它的轴与笛卡尔坐标系的轴定向。

嗯，这很棒，因为我们可以重新排列它以获得

(z - z0)^2/c^2 = 1 - (x - x0)^2/a^2 - (y - y0)^2/b^2
(z - z0)^2 = c^2 - (c/a)^2 (x - x0)^2 - (c/b)^2 (y - y0)^2
z = z0 +- sqrt(c^2 - (c/a)^2 (x - x0)^2 - (c/b)^2 (y - y0)^2)

加减部分对应于椭圆体不是z(x, y)函数的事实，因为对于每个(x, y)点，您都有两个可能的值。但是你可以从两个函数构造一个椭球：一个顶面和一个底面。

最后回到你的问题，你可以选择一个随机(x0, y0, z0)点，然后选择上面椭圆体的底面（那个z = z0 - sqrt(...)）。你唯一需要担心的是，对于(x, y)没有椭球的点，你会得到一个假想的平方根。所以你必须首先过滤掉椭球内的点。最简单的方法是计算平方根，然后丢弃 NaN：

import numpy as np
import pyvista as pv

# denser plane for sharper ellipsoid
street_mesh = pv.Plane(i_size=100, j_size=100, i_resolution=1000, j_resolution=1000)
street_mesh['Normals'] *= -1  # make normals point up for warping later
# add Perlin noise for example
freq = [0.03, 0.03, 0.03]
noise = pv.perlin_noise(5, freq, [0, 0, 0])
street_mesh['height'] = [noise.EvaluateFunction(point) for point in street_mesh.points]

# the relevant part: ellipsoid function
x0, y0, z0 = street_mesh.points[street_mesh.points.size//6, :]  # or whatever random point you want
x, y, _ = street_mesh.points.T  # two 1d arrays of coordinates
a, b, c = 10, 5, 3  # semi-axes
ellipsoid_fun = z0 - np.sqrt(c**2 - (c/a)**2 * (x - x0)**2 - (c/b)**2 *(y - y0)**2)  # RuntimeWarning
keep_inds = ~np.isnan(ellipsoid_fun)
street_mesh['height'][keep_inds] += ellipsoid_fun[keep_inds]
street_mesh.set_active_scalars('height')

# warp by 'height' and do a quick plot
street_mesh.warp_by_scalar().plot()

以上将发出关于虚平方根的警告（导致数据中的 NaN）。如果这让您感到困扰，您可以明确地将其静音。或者，我们可以按照 LBYL 并在计算平方根之前自己检查值：

arg = c**2 - (c/a)**2 * (x - x0)**2 - (c/b)**2 *(y - y0)**2
keep_inds = arg >= 0
ellipsoid = z0 - np.sqrt(arg[keep_inds])
street_mesh['height'][keep_inds] = ellipsoid

这是结果，使用增加的平面密度以获得更清晰的椭圆体： 略微凹凸不平的表面，其中心有一个尖锐的椭圆形凹陷

（我最初建议您可以使用 numpy 广播计算一次坐着的每个坑洼，但这实际上并不容易，甚至可能是不可能的，因为 NaN 过滤打破了这一点。）

于 2021-08-20T19:32:08.357 回答