计算轨迹（路径）中的转折点/枢轴点

Posted

技术标签:

【中文标题】计算轨迹（路径）中的转折点/枢轴点【英文标题】：calculate turning points / pivot points in trajectory (path) 【发布时间】：2013-01-15 21:47:45 【问题描述】：

我正在尝试提出一种算法来确定 x/y 坐标轨迹中的转折点。下图说明了我的意思：绿色表示轨迹的起点，红色表示轨迹的终点（整个轨迹由约 1500 个点组成）：

在下图中，我手动添加了算法可能返回的可能（全局）转折点：

显然，真正的转折点总是值得商榷的，并且取决于人们指定的必须位于两点之间的角度。此外，可以在全局范围内定义转折点（我试图用黑色圆圈做的事情），但也可以在高分辨率局部范围内定义。我对全球（整体）方向的变化很感兴趣，但我很想看到关于人们用来区分全球和本地解决方案的不同方法的讨论。

到目前为止我已经尝试过：

计算后续点之间的距离 计算后续点之间的角度 看看后续点之间的距离/角度如何变化

不幸的是，这并没有给我任何可靠的结果。我可能也计算了沿多个点的曲率，但这只是一个想法。 我真的很感激任何可能对我有帮助的算法/想法。代码可以是任何编程语言，首选matlab或python。

编辑这是原始数据（以防有人想玩）：

mat file text file（x坐标在前，y坐标在第二行）

【问题讨论】：

非常有趣的问题，但我不确定这个论坛是否适合提问。我看到了许多主观方式来定义轨迹中的转折点，例如你看到它的规模。当你仔细观察时，我可以看到许多不同的转折点。继续进行的方法可能是对每个点两侧的点进行某种平滑（或仅使用 n 个点绘制一条线），并决定这两条直线之间的角度。然后，尽管有矫直算法，但您将“只有”两个参数（n 和最小角度）。也许这有帮助？ @Alex 我知道这个问题的主观性。我仍然认为这可能是一个普遍感兴趣的问题，我希望看到人们讨论可以用来区分本地转折点和全球转折点的不同方法。 【参考方案1】：

您可以使用Ramer-Douglas-Peucker (RDP) algorithm 来简化路径。然后，您可以计算沿简化路径的每一段的方向变化。方向变化最大的点可以称为转折点：

RDP 算法的 Python 实现可以在 on github 找到。

import matplotlib.pyplot as plt
import numpy as np
import os
import rdp

def angle(dir):
    """
    Returns the angles between vectors.

    Parameters:
    dir is a 2D-array of shape (N,M) representing N vectors in M-dimensional space.

    The return value is a 1D-array of values of shape (N-1,), with each value
    between 0 and pi.

    0 implies the vectors point in the same direction
    pi/2 implies the vectors are orthogonal
    pi implies the vectors point in opposite directions
    """
    dir2 = dir[1:]
    dir1 = dir[:-1]
    return np.arccos((dir1*dir2).sum(axis=1)/(
        np.sqrt((dir1**2).sum(axis=1)*(dir2**2).sum(axis=1))))

tolerance = 70
min_angle = np.pi*0.22
filename = os.path.expanduser('~/tmp/bla.data')
points = np.genfromtxt(filename).T
print(len(points))
x, y = points.T

# Use the Ramer-Douglas-Peucker algorithm to simplify the path
# http://en.wikipedia.org/wiki/Ramer-Douglas-Peucker_algorithm
# Python implementation: https://github.com/sebleier/RDP/
simplified = np.array(rdp.rdp(points.tolist(), tolerance))

print(len(simplified))
sx, sy = simplified.T

# compute the direction vectors on the simplified curve
directions = np.diff(simplified, axis=0)
theta = angle(directions)
# Select the index of the points with the greatest theta
# Large theta is associated with greatest change in direction.
idx = np.where(theta>min_angle)[0]+1

fig = plt.figure()
ax =fig.add_subplot(111)

ax.plot(x, y, 'b-', label='original path')
ax.plot(sx, sy, 'g--', label='simplified path')
ax.plot(sx[idx], sy[idx], 'ro', markersize = 10, label='turning points')
ax.invert_yaxis()
plt.legend(loc='best')
plt.show()

上面使用了两个参数：

    RDP 算法采用一个参数，tolerance，它 表示简化路径的最大距离 可能会偏离原来的路径。 tolerance 越大，简化路径越粗略。 另一个参数是min_angle，它定义了什么是转折点。 （我取一个转折点为原路径上的任意一点，简化路径上的进出向量之间的夹角大于min_angle）。

【讨论】：

这似乎是我试图通过凸包方法实现的目标，因为我有点想到了快速船体。 +1，并且必须记住这一点。我唯一的问题是它可能应该基于被视为转弯而不是点数的最小角度。 @MC：好主意。谢谢你。 （已进行更改。）【参考方案2】：

我将在下面给出 numpy/scipy 代码，因为我几乎没有 Matlab 经验。

如果你的曲线足够平滑，你可以将你的转折点识别为最高的curvature。将点索引号作为曲线参数，加上central differences scheme，可以用下面的代码计算曲率

import numpy as np
import matplotlib.pyplot as plt
import scipy.ndimage

def first_derivative(x) :
    return x[2:] - x[0:-2]

def second_derivative(x) :
    return x[2:] - 2 * x[1:-1] + x[:-2]

def curvature(x, y) :
    x_1 = first_derivative(x)
    x_2 = second_derivative(x)
    y_1 = first_derivative(y)
    y_2 = second_derivative(y)
    return np.abs(x_1 * y_2 - y_1 * x_2) / np.sqrt((x_1**2 + y_1**2)**3)

您可能希望先平滑曲线，然后计算曲率，然后确定最高曲率点。下面的函数就是这样做的：

def plot_turning_points(x, y, turning_points=10, smoothing_radius=3,
                        cluster_radius=10) :
    if smoothing_radius :
        weights = np.ones(2 * smoothing_radius + 1)
        new_x = scipy.ndimage.convolve1d(x, weights, mode='constant', cval=0.0)
        new_x = new_x[smoothing_radius:-smoothing_radius] / np.sum(weights)
        new_y = scipy.ndimage.convolve1d(y, weights, mode='constant', cval=0.0)
        new_y = new_y[smoothing_radius:-smoothing_radius] / np.sum(weights)
    else :
        new_x, new_y = x, y
    k = curvature(new_x, new_y)
    turn_point_idx = np.argsort(k)[::-1]
    t_points = []
    while len(t_points) < turning_points and len(turn_point_idx) > 0:
        t_points += [turn_point_idx[0]]
        idx = np.abs(turn_point_idx - turn_point_idx[0]) > cluster_radius
        turn_point_idx = turn_point_idx[idx]
    t_points = np.array(t_points)
    t_points += smoothing_radius + 1
    plt.plot(x,y, 'k-')
    plt.plot(new_x, new_y, 'r-')
    plt.plot(x[t_points], y[t_points], 'o')
    plt.show()

一些解释是为了：

turning_points 是您要识别的点数 smoothing_radius 是在计算曲率之前应用于数据的平滑卷积的半径 cluster_radius 是与被选为转折点的高曲率点的距离，在该转折点不应将其他点视为候选点。

您可能需要稍微调整一下参数，但我得到的是这样的：

>>> x, y = np.genfromtxt('bla.data')
>>> plot_turning_points(x, y, turning_points=20, smoothing_radius=15,
...                     cluster_radius=75)

对于全自动检测可能不够好，但它非常接近您想要的。

【讨论】：

【参考方案3】：

一个非常有趣的问题。这是我的解决方案，它允许可变分辨率。虽然，微调它可能并不简单，因为它主要是为了缩小范围

每k个点，计算凸包并将其存储为一个集合。遍历最多 k 个点，并删除所有不在凸包中的点，这样这些点就不会失去原来的顺序。

这里的目的是凸包将充当过滤器，删除所有“不重要的点”，只留下极值点。当然，如果 k 值太高，你最终会得到与实际凸包太接近的东西，而不是你真正想要的。

这应该从一个小的 k 开始，至少 4，然后增加它直到你得到你想要的。您还应该只包括角度低于一定量的每 3 个点的中点，d。这将确保所有转弯都至少为 d 度（未在下面的代码中实现）。然而，这可能应该逐步完成以避免信息丢失，就像增加 k 值一样。另一个可能的改进是实际重新运行已删除的点，并且仅删除不在两个凸包中的点，尽管这需要更高的最小 k 值，至少为 8。

以下代码似乎运行良好，但仍可以使用改进来提高效率和消除噪音。它在确定何时停止时也相当不雅，因此代码实际上只在 k=4 到 k=14 左右有效（就目前而言）。 

def convex_filter(points,k):
    new_points = []
    for pts in (points[i:i + k] for i in xrange(0, len(points), k)):
        hull = set(convex_hull(pts))
        for point in pts:
            if point in hull:
                new_points.append(point)
    return new_points

# How the points are obtained is a minor point, but they need to be in the right order.
x_coords = [float(x) for x in x.split()]
y_coords = [float(y) for y in y.split()]
points = zip(x_coords,y_coords)

k = 10

prev_length = 0
new_points = points

# Filter using the convex hull until no more points are removed
while len(new_points) != prev_length:
    prev_length = len(new_points)
    new_points = convex_filter(new_points,k)

这是上面 k=14 代码的屏幕截图。 61 个红点是过滤后剩下的。

【讨论】：

【参考方案4】：

另一个想法是检查每个点的左右环境。这可以通过在每个点之前和之后创建 N 个点的线性回归来完成。如果点之间的相交角度低于某个阈值，那么您就有了一个角。

这可以通过保持当前处于线性回归中的点的队列并用新点替换旧点来有效地完成，类似于运行平均值。

您最终必须将相邻的角合并为一个角。例如。选择具有最强角属性的点。

【讨论】：

【参考方案5】：

您采用的方法听起来很有希望，但您的数据严重过度采样。您可以先过滤 x 和 y 坐标，例如使用宽高斯然后下采样。

在 MATLAB 中，您可以使用 x = conv(x, normpdf(-10 : 10, 0, 5))，然后使用 x = x(1 : 5 : end)。您必须根据所跟踪对象的内在持久性和点之间的平均距离来调整这些数字。

然后，您将能够非常可靠地检测方向变化，使用您之前尝试过的相同方法，基于标量积，我想。

【讨论】：

以上是关于计算轨迹（路径）中的转折点/枢轴点的主要内容，如果未能解决你的问题，请参考以下文章

CCF - 201604-1 - 折点计数

CCF——折点计数（2016-4）

CCF CSP 201604-1 折点计数

CCF 201604-1折点计数 (水题，暴力)

在PS中,路径锚点分为哪些类型?各有啥特点?

CSP201604-1：折点计数