保持Numpy插值中的尖角

cgfeq70w  于 2023-02-12  发布在  其他
关注(0)|答案(4)|浏览(157)

我正在使用gboffigboffi的优秀代码(包含在下面),用numpy的linspaceinterp插入一个形状。
这样做效果很好,但是有时会漏掉角,从而导致不希望的软化形状。

我想用一个Angular 阈值参数来保持形状的尖角。如果下一个线段的Angular 足够尖锐,有没有办法保持形状插值的尖角?谢谢!

from matplotlib import pyplot as plt
import numpy as np

x = np.array([815.9, 693.2, 570.4, 462.4, 354.4, 469.9, 585.4, 700.6, 815.9])
y = np.array([529.9, 637.9, 746, 623.2, 500.5, 326.9, 153.3, 341.6, 529.9])

fig, ax = plt.subplots(1)
ax.set_aspect('equal')
ax.scatter(x, y, s=40, zorder=3, alpha=0.3)

# compute the distances, ds, between points
dx, dy = x[+1:]-x[:-1],  y[+1:]-y[:-1]
ds = np.array((0, *np.sqrt(dx*dx+dy*dy)))

# compute the total distance from the 1st point, measured on the curve
s = np.cumsum(ds)

# interpolate
xinter = np.interp(np.linspace(0,s[-1], 30), s, x)
yinter = np.interp(np.linspace(0,s[-1], 30), s, y)

# plot the interpolated points
ax.plot(xinter, yinter)
ax.scatter(xinter, yinter, s=5, zorder=4)
plt.show()
wswtfjt7

wswtfjt71#

∮找到尖锐的Angular ∮
为了找到Angular 大于某个截止值的所有点,可以计算归一化差向量的点积的反余弦。

# cosines of the angles are the normalized dotproduct of the difference vectors, ignore first and last point
cosines = np.array( [1, *((dx[:-1]*dx[1:] +  dy[:-1]*dy[1:]) / ds[1:-1] /ds[2:]),1])

另一个公式是叉积和点积的atan2,它避免了除法。

from matplotlib import pyplot as plt
import numpy as np

x = np.array([815.9, 693.2, 570.4, 462.4, 354.4, 469.9, 585.4, 700.6, 815.9])
y = np.array([529.9, 637.9, 746, 623.2, 500.5, 326.9, 153.3, 341.6, 529.9])

fig, ax = plt.subplots(1)
ax.set_aspect('equal')
ax.scatter(x, y, s=40, zorder=3, alpha=0.3)

# compute the distances, ds, between points
dx, dy = x[+1:] - x[:-1], y[+1:] - y[:-1]
ds = np.append(0, np.sqrt(dx * dx + dy * dy))

# compute the total distance from the 1st point, measured on the curve
s = np.cumsum(ds)

# angle: atan2 of cross product and dot product
cut_off_angle = 10  # angles larger than this will be considered sharp
angles_rad = np.arctan2(dx[:-1] * dy[1:] - dy[:-1] * dx[1:],
                        dx[:-1] * dx[1:] + dy[:-1] * dy[1:])
# convert to degrees, and pad with zeros for first and last point
angles = np.pad(np.degrees(angles_rad), 1)

# interpolates
s_long = np.sort(np.append(np.linspace(0, s[-1], 30),
                           s[np.abs(angles) > cut_off_angle]))
xinter = np.interp(s_long, s, x)
yinter = np.interp(s_long, s, y)

# plot the interpolated points
ax.plot(xinter, yinter)
ax.scatter(xinter, yinter, s=5, zorder=4)
plt.show()

回答原始问题

一个简单的扩展是将s的点插入到用于中间点(np.linspace(0,s[-1], 30))的数组中。

from matplotlib import pyplot as plt
import numpy as np

x = np.array([100, 200, 200, 100, 100])
y = np.array([100, 100, 200, 200, 100])

fig, ax = plt.subplots(1)
ax.set_aspect('equal')
ax.scatter(x, y, s=40, zorder=3, alpha=0.3)

# compute the distances, ds, between points
dx, dy = x[+1:] - x[:-1], y[+1:] - y[:-1]
ds = np.array((0, *np.sqrt(dx * dx + dy * dy)))

# compute the total distance from the 1st point, measured on the curve
s = np.cumsum(ds)

# interpolates
s_long = np.sort(np.append(np.linspace(0, s[-1], 30)[1:-1], s))
xinter = np.interp(s_long, s, x)
yinter = np.interp(s_long, s, y)

# plot the interpolated points
ax.plot(xinter, yinter)
ax.scatter(xinter, yinter, s=5, zorder=4)
plt.show()

mlmc2os5

mlmc2os52#

你应该先确定角;计算曲线每段之间的Angular ,如果Angular 大于某个阈值,则可以将该点视为角点。然后,可以在角点处分割曲线,计算点之间的距离,并单独对每段进行插值。
我就是这么做的。

from matplotlib import pyplot as plt
import numpy as np

x = np.array([815.9, 693.2, 570.4, 462.4, 354.4, 469.9, 585.4, 700.6, 815.9])
y = np.array([529.9, 637.9, 746, 623.2, 500.5, 326.9, 153.3, 341.6, 529.9])

fig, ax = plt.subplots(1)
ax.set_aspect('equal')
ax.scatter(x, y, s=40, zorder=3, alpha=0.3)

# Angle threshold
threshold = np.pi / 6

# Identify the corners
corners = [0, len(x)-1]
for i in range(len(x) - 2):
    dx1, dy1 = x[i+1]-x[i], y[i+1]-y[i]
    dx2, dy2 = x[i+2]-x[i+1], y[i+2]-y[i+1]
    angle = np.arctan2(dy2, dx2) - np.arctan2(dy1, dx1)
    if abs(angle) > threshold:
        corners.append(i+1)
corners.sort()

# Split the curve into segments
segments = [(start, i) for start, i in zip([0] + corners, corners)]

# Interpolate each segment
xinter_list, yinter_list = [], []
for start, end in segments:
    xseg, yseg = x[start:end+1], y[start:end+1]
    
    # Compute the distances
    dx, dy = xseg[+1:]-xseg[:-1],  yseg[+1:]-yseg[:-1]
    ds = np.array((0, *np.sqrt(dx*dx+dy*dy)))

    # Compute the total distance from the 1st point
    s = np.cumsum(ds)
    
    # Interpolate  
    xinter_list.append(np.interp(np.linspace(0,s[-1], 30), s, xseg))
    yinter_list.append(np.interp(np.linspace(0,s[-1], 30), s, yseg))

# Plot the interpolated points
xinter = np.concatenate(xinter_list)
yinter = np.concatenate(yinter_list)
ax.plot(xinter, yinter)
ax.scatter(xinter, yinter, s=5, zorder=4)
plt.show()

whhtz7ly

whhtz7ly3#

假设您有一个“corner”列表,例如:

corners = 0, 12,43,48

(note corners包括第一点和最后一点的索引,即,在该示例中x的长度应该是49)。
接下来计算s,然后计算角位置列表

c_pos = [s[c] for c in corners]

你最终可以计算出你想要计算插值的点的列表

ds = (s[-1]-s[0])/npoints # e.g., npoints = 200
arrays = [np.linspace(c_pos0, c_pos1, (c_pos1-c_pos0/ds), endpoint=c_pos1==s[-1])
              for c_pos0, c_pos1 in zip(c_pos[0:], c_pos[1:])
s_new = np.concatenate(arrays)

xinterp = np.interp(s_new, s, x)
...

几句话

np.linspace(c_pos0, c_pos1, (c_pos1-c_pos0/ds), endpoint=c_pos1==s[-1])

这意味着:列表arrays中的数组是“开区间”[c0, c0+Δ,..., c1),除了最后一个,因为我们想要画一条闭合曲线。
这是回答你问题的容易部分,难的是找到角落...

mf98qq94

mf98qq944#

如果xinter和yinter步长不是边长的分割线,就不可能制作尖角。在这种情况下,你需要步长25,因为你的边长是100。
试试这个:

# interpolate
xinter = np.interp(np.linspace(0,s[-1], 25), s, x)
yinter = np.interp(np.linspace(0,s[-1], 25), s, y)

相关问题