我正在尝试使用skimage.measure.marching_cubes_lewiner
来解决一些 isosurface f(x,y,z)=0
。在我的情况下f
是强非线性的,并且当坐标以对数间距给出时最好映射。因为行进立方体需要一个规则网格来构建体素,所以我正在处理与我的原始坐标X,Y,Z
相对应的坐标网格网格log10
,因此我的等值面由 等价地给出f(10**X,10**Y,10**Z)=0
。一切都会好起来的,如果不是因为假设我正在使用X,Y,Z
in [-1.5,2]^3
(相当于x,y,z
in [0.03,100.]^3
),由 给出的解的顶点坐标skimage.measure.marching_cubes_lewiner
不在这个立方体中。
在回答了关于 SO 的另一个相关问题之后,我认为这可能是因为该算法可能在考虑单一体积的情况下工作,因此我需要spacing
在我的调用中设置正确的输入参数skimage.measure.marching_cubes_lewiner
。以这种方式,假设我将我的函数映射到每个坐标f
的点网格上N
,以便我通过numpy.diff([-1.5,2])/N
每个坐标增加指数,因此我调用:
import numpy as np
from skimage import measure as msr
def f(x,y,z):
val = ... # some lengthy code to define my implicit function
return val
# Define ranges of my coordinates
xRange = [0.03,100.]
yRange = [0.03,100.]
zRange = [0.03,100.]
XRange = np.log10(xRange)
YRange = np.log10(yRange)
ZRange = np.log10(zRange)
# Create regular grid
N = 50 # number of points per coordinate
X,Y,Z = np.mesh[XRange[0]:XRange[1]:N*1j,
YRange[0]:YRange[1]:N*1j,
ZRange[0]:ZRange[1]:N*1j]
F = f(10**X,10**Y,10**Z)
sol,_,_,_ = skimage.measure.marching_cubes_lewiner(F,0.0,spacing(np.diff(XRange)/N,np.diff(YRange)/N,np.diff(ZRange)/N))
然而,出乎意料的是,解点的坐标通常出现在[0,Vx]*[0,Vy]*[0,Vz]
,和Vx>XRange[-1]
中。我不知道为什么会发生这种情况,也不知道如何正确地将等值面解决方案的坐标重新调整为问题的真实单位。Vy>YRange[-1]
Vz>ZRange[-1]