軌道（軌道）の転回点/ピボット点を計算する

Question

私は、x/y座標の軌跡で転回点を決定するアルゴリズムを考え出しています。次の図は、私が意味するものを示しています：緑は開始点を示し、赤は軌道の最終点を示します（全体の軌跡は約1500点で構成されています）。

次の図では、

明らかに、真の転換点は常に議論の余地があると持っているいずれかを指定しますが、ポイントの間に存在することを角度に依存します：アルゴリズムが返すことができるグローバル）ターニングポイント。さらに、グローバルスケール（黒丸で何をしようとしたか）に転換点を定義できますが、高解像度のローカルスケールで定義することもできます。私はグローバルな（全体的な）方向の変化に興味がありますが、グローバルなソリューションとローカルなソリューションを区別するためのさまざまなアプローチについての議論が大好きです。残念ながら、後続の点以降の点の間の

見てどのように距離/角度の変更

間

• 計算距離
以降の点の間の計算角度を：私がこれまで試したどのような

これは私に頑強な結果をもたらしません。私はおそらく、複数の点に沿って曲率を計算しているかもしれませんが、それは単なるアイデアです。 私はここで私を助けるかもしれないアルゴリズムやアイデアを本当に感謝しています。コードは任意のプログラミング言語にすることができますが、matlabまたはpythonが優先されます。

EDITがここに生のデータです（ケース誰かにのは、それを再生したい）：

• mat file
• text file（xがyが2行目に座標、第1の座標）

Answer 1

パスを簡略化するためにRamer-Douglas-Peucker (RDP) algorithmを使用できます。次に、単純化されたパスの各セグメントに沿った方向の変化を計算できます。最大の方向変化に対応するポイントは、転換点と呼ばれます。

PythonのRDPアルゴリズムの実装は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() 

2つのパラメータは、上記で使用した：

RDPアルゴリズムは を簡略化パス が元のパスから逸脱することができる最大距離を表す一つのパラメータ、 toleranceをとる

1. 。 toleranceが大きいほど、パスは単純化されます。
2. もう1つのパラメータは、何がターニングポイントとみなされるのかを定義するmin_angleです。 （私は、単純化されたパス上のベクトルの入力と出力の間の角度がmin_angleより大きい元のパス上の任意の点に転換点を取っています）。

Answer 2

あなたが取ったアプローチは有望だと思うが、あなたのデータは大きくオーバーサンプリングされている。 xとyの座標を最初にフィルタリングすることができます。たとえば、幅の広いガウス関数を使用し、次にダウンサンプルします。

MATLABでは、x = conv(x, normpdf(-10 : 10, 0, 5))、次にx = x(1 : 5 : end)を使用できます。トラッキングしているオブジェクトの本来の永続性とポイント間の平均距離に応じて、これらの数値を微調整する必要があります。

次に、私が想像しているスカラー製品に基づいて、これまでに試みたのと同じアプローチを使用して、方向の変化を非常に確実に検出することができます。

Answer 3

非常に興味深い質問です。ここで私の解決策は、可変解像度を可能にします。微調整は主に絞ることを意図しているので簡単ではないかもしれませんが

凸包を計算してセットとして保存します。ほとんどのk点を通過し、点が元の順序を失わないように、凸包にない点を取り除く。

ここでの目的は、凸包がフィルタとして機能し、極端な点だけを残している「重要でない点」をすべて削除することです。もちろん、k値が大きすぎると、実際に必要なものではなく、実際の凸包に近すぎるものになります。

これは小さなk（少なくとも4）で始まり、次に求めるものが得られるまで増やしてください。おそらく角度がある量d以下の3点ごとに中点だけを含めるべきです。d。これは、すべてのターンが少なくともd度（以下のコードでは実装されていない）であることを保証します。しかし、これはおそらくk値を増やすのと同じように、情報の損失を避けるために段階的に行われるべきです。もう1つの可能な改善は、除去された点で実際に再実行し、両方の凸包にない点のみを削除することです。

次のコードかなりうまくいくように見えますが、効率とノイズ除去のために改善を使用することができます。いつ停止すべきかを決めるのはむしろ面白味がないので、コードは実際には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個の赤い点は、フィルターの後ろに残るものです。

Answer 4

もう一つのアイデアは、すべての点で左と右の周囲を検討することです。これは、各点の前後にN点の線形回帰を作成することによって行うことができる。ポイント間の交差角度がある閾値を下回っている場合は、角があります。

現在、線形回帰のポイントのキューを維持し、古いポイントを新しいポイントに置き換えることによって、これは効率的に実行できます。これは、実行中の平均と同様です。

最後に、隣接するコーナーを1つのコーナーにマージする必要があります。例えば。最も強いコーナー特性を有する点を選択する。

Answer 5

私はMatlabの経験がほとんどないので、以下のnumpy/scipyコードを与えます。

カーブが十分に滑らかであれば、あなたの転換点は最高のものであることがわかります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() 

一部がオーダーである：あなたが

smoothing_radiusを確認したい点の数を適用する平滑コンボリューションの半径

• turning_pointsですあなたのデータに曲率を計算する前に
• cluster_radiusは、曲率の高いポイントから他のポイントを候補として考慮しない転回ポイントとして選択されたポイントまでの距離です。

あなたは少しのパラメータをいじっする必要がありますが、私はこのようなものだ：

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

おそらくない良い十分な完全に自動化された検出のために、それはかなりですあなたが望むものに近い。