如何计算Python中2个轮廓之间的点数答案

【问题标题】：How can I count number of points between 2 contours in Python如何计算Python中2个轮廓之间的点数
【发布时间】：2017-10-30 22:50:29
【问题描述】：

我正在使用matplotlib.pyplot 插入我的数据并创建轮廓。遵循this 答案/示例（关于如何计算轮廓内的面积），我能够得到轮廓线的顶点。有没有办法使用该信息，即一条线的顶点，来计算两个给定轮廓之间有多少点？这些点将不同于用于推导轮廓的数据。

【问题讨论】：

标签： python matplotlib contour

【解决方案1】：

通常，您不想通过逆向工程来获取一些数据。相反，您可以对稍后用于绘制轮廓的数组进行插值，并找出哪些点位于特定值的区域中。

下面将找到-0.8 和-0.4 级别之间的所有点，打印它们并在图上以红色显示。

import numpy as np; np.random.seed(1)
import matplotlib.mlab as mlab
import matplotlib.pyplot as plt
from scipy.interpolate import Rbf

X, Y = np.meshgrid(np.arange(-3.0, 3.0, 0.1), np.arange(-2.4, 1.0, 0.1))
Z1 = mlab.bivariate_normal(X, Y, 1.0, 1.0, 0.0, 0.0)
Z2 = mlab.bivariate_normal(X, Y, 1.5, 0.5, 1, 1)
Z = 10.0 * (Z2 - Z1)

points = np.random.randn(15,2)/1.2
levels = [-1.2, -0.8,-0.4,-0.2]

# interpolate points
f = Rbf(X.flatten(), Y.flatten(), Z.flatten()) 
zi = f(points[:,0], points[:,1])
# add interpolated points to array with columns x,y,z
points3d = np.zeros((points.shape[0],3))
points3d[:,:2] = points
points3d[:,2] = zi
# masking condition for points between levels
filt = (zi>levels[1]) & (zi <levels[2])
# print points between the second and third level
print(points3d[filt,:])

### plotting
fig, ax = plt.subplots()

CS = ax.contour(X, Y, Z, levels=levels)
ax.clabel(CS, inline=1, fontsize=10)
#plot points between the second and third level in red:
ax.scatter(points[:,0], points[:,1], c=filt.astype(float), cmap="bwr" )

plt.show()

【讨论】：

感谢您的回复。我认为这将帮助我更好地表达我的目标：我不需要对情节进行逆向工程。相反，我想使用从已知 X、Y、Z 的数组/数据集获得的 2 个轮廓，并希望将其用于仅已知 X 和 Y 的 other 数组。由于我的 Z 值来自耗时的模型运行，我只想知道它们落在轮廓的哪一侧（而不是它们的确切位置），我希望我可以避免为新数组再次运行模型而是使用现有轮廓的顶点。这更有意义吗？
当然，但上面的内容不正是这样做的吗？
看起来您需要知道示例中 Z 的值才能使用掩码并将其过滤掉。我只有一个数组有 Z，其余的没有。如果我在这里遗漏了一些明显的东西，请原谅我？
从本质上讲，我的任务类似于通过“光线投射”完成的“多边形中的点”问题。除了代替多边形外，我有两个轮廓。我可以将轮廓转换为多边形，但想知道是否有更直接的方法
你需要X,Y,Z（我想这些是来自你的密集模型的数据）。你需要知道两个轮廓级别（例如z=-0.8,z=-0.4），最后你需要一些points，你想知道它们是否位于轮廓内（你只有x，y的点，z=zi是在脚本中计算）。

【解决方案2】：

我不确定我是否理解您要检查哪些点，但是，如果您有线顶点（两个点）并且想要检查第三个点是否介于两者之间，您可以采用简单的 (效率不高）接近并计算三者形成的三角形的面积。如果面积为 0，则该点落在同一条线上。此外，您可以计算点之间的距离，并查看该点是在线之间还是线外（延长线）。

希望这会有所帮助！

【讨论】：

是的，我说的是检查“第三点”是否介于两者之间，除了我有一组大约 20-40K 点要检查/计数。您的评论确实帮助我考虑了更多，但我很难概念化如何为一条线（而不仅仅是一组 x,y 顶点）和许多点执行此操作。谢谢！

【解决方案3】：

Ampere's law of magnetic fields 可以在这里提供帮助 - 尽管计算成本可能很高。该定律表明，沿闭合回路的磁场的路径积分与回路内的电流成正比。

假设你有一个轮廓 C 和一个点 P (x0,y0)。想象一条位于 P 处垂直于页面（进入页面的电流）的无限导线承载一些电流。使用安培定律，我们可以证明导线在点 P (x,y) 处产生的磁场与 (x0,y0) 到 (x,y) 的距离成反比，并且与以点 P 为中心的圆相切通过点 (x,y)。因此，如果导线位于轮廓之外，则路径积分为零。

import numpy as np
import pylab as plt

# generating a mesh and values on it
delta = 0.1
x = np.arange(-3.1*2, 3.1*2, delta)
y = np.arange(-3.1*2, 3.1*2, delta)
X, Y = np.meshgrid(x, y)
Z = np.sqrt(X**2 + Y**2)

# generating the contours with some levels
levels = [1.0]
plt.figure(figsize=(10,10))
cs = plt.contour(X,Y,Z,levels=levels)

# finding vertices on a particular level
contours = cs.collections[0]
vertices_level = contours.get_paths()[0].vertices # In this example the 
shape of vertices_level is (161,2)
# converting points into two lists; one per dimension. This step can be optimized
lX, lY = list(zip(*vertices_level))
          
# computing Ampere's Law rhs 
def AmpereLaw(x0,y0,lX,lY):
    S = 0
    for ii in range(len(lX)-1):
        dx = lX[ii+1] - lX[ii]
        dy = lY[ii+1] - lY[ii]
        ds = (1/((lX[ii]-x0)**2+(lY[ii]-y0)**2))*(-(lY[ii]-y0)*dx+(lX[ii]-x0)*dy)
        if -1000 < ds < 1000: #to avoid very lare numbers when denominator is small
            S = S + ds
    return(S)
    
# we know point (0,0) is inside the contour
AmpereLaw(0,0,lX,lY)
# result: -6.271376740062852

# we know point (0,0) is inside the contour
AmpereLaw(-2,0,lX,lY)
# result:  0.00013279920934375876

您可以使用此结果找到一个轮廓内但另一轮廓外的点。

【讨论】：