python - girona - qgis manual

Encuentra la distancia mínima de punto a curva complicada

Question

Tengo una curva complicada definida como un conjunto de puntos en una tabla como tal (la tabla completa está here ):

# x y
1.0577 12.0914
1.0501 11.9946
1.0465 11.9338
...

Si trazo esta tabla con los comandos:

plt.plot(x_data, y_data, c=''b'',lw=1.)
plt.scatter(x_data, y_data, marker=''o'', color=''k'', s=10, lw=0.2)

Me sale lo siguiente:

donde he añadido los puntos y segmentos rojos manualmente. Lo que necesito es una forma de calcular esos segmentos para cada uno de esos puntos, es decir: una forma de encontrar la distancia mínima desde un punto dado en este espacio 2D a la curva interpolada .

No puedo usar la distancia a los puntos de datos en sí mismos (los puntos negros que generan la curva azul) ya que no están ubicados en intervalos iguales, a veces están cerca y otras están muy alejados y esto afecta profundamente mis resultados más abajo en el línea.

Dado que esta no es una curva bien educada, no estoy realmente seguro de lo que podría hacer. He intentado interpolarlo con una línea de vida UnivariateSpline pero me da un ajuste muy pobre:

# Sort data according to x.
temp_data = zip(x_data, y_data)
temp_data.sort()
# Unpack sorted data.
x_sorted, y_sorted = zip(*temp_data)
# Generate univariate spline.
s = UnivariateSpline(x_sorted, y_sorted, k=5)
xspl = np.linspace(0.8, 1.1, 100)
yspl = s(xspl)
# Plot.
plt.scatter(xspl, yspl, marker=''o'', color=''r'', s=10, lw=0.2)

También intenté aumentar el número de puntos de interpolación, pero tuve un lío:

# Sort data according to x.
temp_data = zip(x_data, y_data)
temp_data.sort()
# Unpack sorted data.
x_sorted, y_sorted = zip(*temp_data)
t = np.linspace(0, 1, len(x_sorted))
t2 = np.linspace(0, 1, 100)
# One-dimensional linear interpolation.
x2 = np.interp(t2, t, x_sorted)
y2 = np.interp(t2, t, y_sorted)
plt.scatter(x2, y2, marker=''o'', color=''r'', s=10, lw=0.2)

Cualquier idea / puntero será muy apreciado.

Answer 1

La curva es por naturaleza paramétrica, es decir, para cada x no es necesaria una única y y viceversa. Así que no debes interpolar una función de la forma y (x) o x (y). En su lugar, debe hacer dos interpolaciones, x (t) y y (t) donde t es, por ejemplo, el índice del punto correspondiente.

Luego usa scipy.optimize.fminbound para encontrar la t óptima de tal manera que (x (t) - x0) ^ 2 + (y (t) - y0) ^ 2 sea el más pequeño, donde (x0, y0) son los puntos rojos en tu primera figura Para fminsearch, puede especificar el límite mínimo / máximo para que t sea 1 y len(x_data)

Answer 2

Puede intentar implementar un cálculo de distancia de punto a línea en pares incrementales de puntos en la curva y encontrar ese mínimo. Esto introducirá un pequeño error de la curva según se dibuja, pero debería ser muy pequeño, ya que los puntos están relativamente juntos.

http://en.wikipedia.org/wiki/Distance_from_a_point_to_a_line

Answer 3

Si está dispuesto a usar una biblioteca para esto, mire shapely : https://github.com/Toblerity/Shapely

Como un ejemplo rápido ( points.txt contiene los datos que vinculó en su pregunta):

import shapely.geometry as geom
import numpy as np
coords = np.loadtxt(''points.txt'')
line = geom.LineString(coords)
point = geom.Point(0.8, 10.5)
# Note that "line.distance(point)" would be identical
print point.distance(line)

Como ejemplo interactivo (esto también dibuja los segmentos de línea que querías):

import numpy as np
import shapely.geometry as geom
import matplotlib.pyplot as plt
class NearestPoint(object):
def __init__(self, line, ax):
self.line = line
self.ax = ax
ax.figure.canvas.mpl_connect(''button_press_event'', self)
def __call__(self, event):
x, y = event.xdata, event.ydata
point = geom.Point(x, y)
distance = self.line.distance(point)
self.draw_segment(point)
print ''Distance to line:'', distance
def draw_segment(self, point):
point_on_line = line.interpolate(line.project(point))
self.ax.plot([point.x, point_on_line.x], [point.y, point_on_line.y],
color=''red'', marker=''o'', scalex=False, scaley=False)
fig.canvas.draw()
if __name__ == ''__main__'':
coords = np.loadtxt(''points.txt'')
line = geom.LineString(coords)
fig, ax = plt.subplots()
ax.plot(*coords.T)
ax.axis(''equal'')
NearestPoint(line, ax)
plt.show()

Tenga en cuenta que he añadido ax.axis(''equal'') . shapely opera en el sistema de coordenadas en el que se encuentran los datos. Sin la gráfica de ejes iguales, la vista se distorsionará, y aunque shapely aún encontrará el punto más cercano, no se verá del todo bien en la pantalla: