计算轨迹(路径)中的转折点/枢轴点

Question

我正在尝试提出一种算法来确定x/y坐标轨迹中的转折点.下图说明了我的意思:绿色表示起点,红色表示轨迹的最终点(整个轨迹由~1500点组成):

在下图中,我手动添加了算法可能返回的可能(全局)转折点:

显然,真正的转折点总是有争议的,并且将取决于一个人指定必须位于点之间的角度.此外,转折点可以在全球范围内定义(我试图用黑色圆圈做的),但也可以在高分辨率的局部尺度上定义.我对全球(整体)方向变化感兴趣,但我很乐意看到人们用来梳理全局与本地解决方案的不同方法的讨论.

到目前为止我尝试过的:

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

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

编辑这里是原始数据(如果有人想要玩它):

• 垫文件
• 文本文件(首先是x坐标,第二行是y坐标)

Answer 1

您可以使用Ramer-Douglas-Peucker(RDP)算法来简化路径.然后,您可以计算沿简化路径的每个段的方向变化.对应于方向最大变化的点可以称为转折点:

可以在github上找到RDP算法的Python实现.

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()

上面使用了两个参数:

1. RDP算法采用一个参数,即一个参数,tolerance表示简化路径可以偏离原始路径的最大距离.越大tolerance,越简化的路径.
2. 另一个参数是min_angle定义什么被视为转折点.(我正在转向原点路径上的任何一点,它在简化路径上的进入和退出向量之间的角度大于min_angle).

Answer 2

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

如果曲线足够平滑,则可以将转折点识别为曲率最高的转折点.将点索引号作为曲线参数和中心差分方案,可以使用以下代码计算曲率

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)

可能不足以进行全自动检测,但它非常接近你想要的.