从深度图计算视野

Question

我有一个深度 numpy n*n 数组。其中每个单元格的值代表该单元格的高度。我想计算视野。换句话说：

if Pixel_height > viewer_height:
    pixels_behind = -1

有没有简单的方法可以做到这一点？

编辑：要清楚，我想创建一个类似于这个绿色蒙版的位置矩阵，知道我们有一个高度矩阵 n*m enter image description here

Answer 1

我想出了这个功能，它不是很完美，但做的事情和你说的差不多。它有一个resolution 参数，表示您要考虑的圆周上有多少个“角度步长”。我不确定这是解决这个问题的“正确”方法，而且我觉得应该有更聪明的方法来做一些操作，但无论如何。

import numpy as np

def view_field(height_map, viewer_height, viewer_pos, resolution=360):
    height_map = np.asarray(height_map)
    h, w = height_map.shape
    vi, vj = viewer_pos
    # Find locations higher than viewer
    m = height_map > viewer_height
    # Find angle of each pixel relative to viewer
    ii, jj = np.ogrid[-vi:h - vi, -vj:w - vj]
    a = np.arctan2(ii, jj)
    # Distance of each pixel to viewer
    d2 = np.square(ii) + np.square(jj)
    d = np.sqrt(d2)
    # Find angle range "behind" each pixel
    pix_size = 0.5
    ad = np.arccos(d / np.sqrt(d2 + np.square(pix_size)))
    # Minimum and maximum angle encompassed by each pixel
    amin = a - ad
    amax = a + ad
    # Define angle "bins"
    ar = np.linspace(-np.pi, np.pi, resolution + 1)
    # Find the bin corresponding to each pixel
    b = np.digitize(a, ar) % resolution
    bmin = np.digitize(amin, ar) % resolution
    bmax = np.digitize(amax, ar) % resolution
    # Find the closest distance to a high pixel for each angle bin
    angdist = np.full_like(ar, np.inf)
    np.minimum.at(angdist, bmin[m], d[m])
    np.minimum.at(angdist, bmax[m], d[m])
    # Visibility is true if the pixel distance is less than the
    # visibility distance for its angle bin
    return d <= angdist[b]

接下来是这样的：

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# A wavy function for testing
def ackley(x, y):
    return (-20 * np.exp(-0.2 * np.sqrt(0.5 * (np.square(x) + np.square(y))))
            - np.exp(0.5 * (np.cos(2 * np.pi * x) + np.cos(2 * np.pi * y)))
            + np.e + 20)

# Make input data
x, y = np.ogrid[-4:4:100j, -4:4:100j]
height_map = ackley(x, y)
viewer_height = 7.5
# Viewer coordinates are in "pixel space"
viewer_pos = (50, 70)
# Compute visibility
f = view_field(height_map, viewer_height, viewer_pos)
# Plot height map and result
plt.figure(figsize=(8, 4))
ax = plt.subplot(121, projection='3d')
ax.plot_surface(np.arange(100)[:, np.newaxis], np.arange(100),
                height_map, color='b', alpha=0.6)
ax.scatter3D(*viewer_pos, viewer_height, s=5, c='r', label='Viewer')
plt.title('Height map')
plt.legend()
plt.subplot(122)
plt.imshow(f)
plt.title('Visibility map')
plt.tight_layout()

以及由此产生的情节：