Numpy y líneas de intersecciones

¿Cómo usaría numpy para calcular la intersección entre dos segmentos de línea?

En el código tengo segment1 = ((x1, y1), (x2, y2)) y segment2 = ((x1, y1), (x2, y2)). Note que el segmento 1 no es igual a segment2. Así que en mi código también he estado calculando la pendiente y la intersección con y, sería bueno si se pudiera evitar, pero no sé cómo.

He estado usando la regla de Cramer con una función que escribí en Python, pero me gustaría encontrar una forma más rápida de hacerlo.

Robado directamente de http://www.cs.mun.ca/~rod/2500/notes/numpy-arrays/numpy-arrays.html

# # line segment intersection using vectors # see Computer Graphics by FS Hill # from numpy import * def perp( a ) : b = empty_like(a) b[0] = -a[1] b[1] = a[0] return b # line segment a given by endpoints a1, a2 # line segment b given by endpoints b1, b2 # return def seg_intersect(a1,a2, b1,b2) : da = a2-a1 db = b2-b1 dp = a1-b1 dap = perp(da) denom = dot( dap, db) num = dot( dap, dp ) return (num / denom.astype(float))*db + b1 p1 = array( [0.0, 0.0] ) p2 = array( [1.0, 0.0] ) p3 = array( [4.0, -5.0] ) p4 = array( [4.0, 2.0] ) print seg_intersect( p1,p2, p3,p4) p1 = array( [2.0, 2.0] ) p2 = array( [4.0, 3.0] ) p3 = array( [6.0, 0.0] ) p4 = array( [6.0, 3.0] ) print seg_intersect( p1,p2, p3,p4) 
 import numpy as np def get_intersect(a1, a2, b1, b2): """ Returns the point of intersection of the lines passing through a2,a1 and b2,b1. a1: [x, y] a point on the first line a2: [x, y] another point on the first line b1: [x, y] a point on the second line b2: [x, y] another point on the second line """ s = np.vstack([a1,a2,b1,b2]) # s for stacked h = np.hstack((s, np.ones((4, 1)))) # h for homogeneous l1 = np.cross(h[0], h[1]) # get first line l2 = np.cross(h[2], h[3]) # get second line x, y, z = np.cross(l1, l2) # point of intersection if z == 0: # lines are parallel return (float('inf'), float('inf')) return (x/z, y/z) if __name__ == "__main__": print get_intersect((0, 1), (0, 2), (1, 10), (1, 9)) # parallel lines print get_intersect((0, 1), (0, 2), (1, 10), (2, 10)) # vertical and horizontal lines print get_intersect((0, 1), (1, 2), (0, 10), (1, 9)) # another line for fun 

Explicación

Tenga en cuenta que la ecuación de una línea es ax+by+c=0 . Entonces, si un punto está en esta línea, entonces es una solución para (a,b,c).(x,y,1)=0 ( . Es el producto puntual)

sea l1=(a1,b1,c1) , l2=(a2,b2,c2) sean dos líneas y p1=(x1,y1,1) , p2=(x2,y2,1) sean dos puntos.

Encontrando la línea que pasa por dos puntos:

sea t=p1xp2 (el producto cruzado de dos puntos) sea un vector que representa una línea.

Sabemos que p1 está en la línea t porque t.p1 = (p1xp2).p1=0 . También sabemos que p2 está en t porque t.p2 = (p1xp2).p2=0 . Entonces t debe ser la línea que pasa por p1 y p2 .

Esto significa que podemos obtener la representación vectorial de una línea tomando el producto cruzado de dos puntos en esa línea .

Encontrando el punto de intersección:

Ahora, r=l1xl2 (el producto cruzado de dos líneas) sea un vector que representa un punto

Sabemos que r encuentra en l1 porque r.l1=(l1xl2).l1=0 . También sabemos que r encuentra en l2 porque r.l2=(l1xl2).l2=0 . Entonces r debe ser el punto de intersección de las líneas l1 y l2 .

Curiosamente, podemos encontrar el punto de intersección tomando el producto cruzado de dos líneas .

Esta es una respuesta tardía, tal vez, pero fue el primer éxito cuando busqué en Google “intersecciones de línea numpy”. En mi caso, tengo dos líneas en un plano, y quería obtener rápidamente cualquier intersección entre ellas, y la solución de Hamish sería lenta, lo que requeriría un bucle nested para todos los segmentos de línea.

Aquí está cómo hacerlo sin un bucle for (es bastante rápido):

 from numpy import where, dstack, diff, meshgrid def find_intersections(A, B): # min, max and all for arrays amin = lambda x1, x2: where(x1x2, x1, x2) aall = lambda abools: dstack(abools).all(axis=2) slope = lambda line: (lambda d: d[:,1]/d[:,0])(diff(line, axis=0)) x11, x21 = meshgrid(A[:-1, 0], B[:-1, 0]) x12, x22 = meshgrid(A[1:, 0], B[1:, 0]) y11, y21 = meshgrid(A[:-1, 1], B[:-1, 1]) y12, y22 = meshgrid(A[1:, 1], B[1:, 1]) m1, m2 = meshgrid(slope(A), slope(B)) m1inv, m2inv = 1/m1, 1/m2 yi = (m1*(x21-x11-m2inv*y21) + y11)/(1 - m1*m2inv) xi = (yi - y21)*m2inv + x21 xconds = (amin(x11, x12) < xi, xi <= amax(x11, x12), amin(x21, x22) < xi, xi <= amax(x21, x22) ) yconds = (amin(y11, y12) < yi, yi <= amax(y11, y12), amin(y21, y22) < yi, yi <= amax(y21, y22) ) return xi[aall(xconds)], yi[aall(yconds)] 

Luego, para usarlo, proporcione dos líneas como argumentos, donde arg es una matriz de 2 columnas, cada fila correspondiente a un punto (x, y):

 # example from matplotlib contour plots Acs = contour(...) Bsc = contour(...) # A and B are the two lines, each is a # two column matrix A = Acs.collections[0].get_paths()[0].vertices B = Bcs.collections[0].get_paths()[0].vertices # do it x, y = find_intersections(A, B) 

que te diviertas

Esta es una versión de la respuesta de @Hamish Grubijan que también funciona para múltiples puntos en cada uno de los argumentos de entrada, es decir, a1 , a2 , b1 , b2 pueden ser matrices de filas Nx2 de puntos 2D. La función perp es reemplazada por un producto de puntos.

 T = np.array([[0, -1], [1, 0]]) def line_intersect(a1, a2, b1, b2): da = np.atleast_2d(a2 - a1) db = np.atleast_2d(b2 - b1) dp = np.atleast_2d(a1 - b1) dap = np.dot(da, T) denom = np.sum(dap * db, axis=1) num = np.sum(dap * dp, axis=1) return np.atleast_2d(num / denom).T * db + b1 

Aquí hay una (un poco forzada) de una sola línea:

 import numpy as np from scipy.interpolate import interp1d interp1d(segment1-segment2,np.arange(segment1.shape[0]))(0) 

Interpolar la diferencia (el valor predeterminado es lineal) y encontrar un 0 de la inversa.

¡Aclamaciones!

Esto es lo que uso para encontrar la intersección de líneas, funciona teniendo 2 puntos de cada línea, o solo un punto y su pendiente. Básicamente resuelvo el sistema de ecuaciones lineales.

 def line_intersect(p0, p1, m0=None, m1=None, q0=None, q1=None): ''' intersect 2 lines given 2 points and (either associated slopes or one extra point) Inputs: p0 - first point of first line [x,y] p1 - fist point of second line [x,y] m0 - slope of first line m1 - slope of second line q0 - second point of first line [x,y] q1 - second point of second line [x,y] ''' if m0 is None: if q0 is None: raise ValueError('either m0 or q0 is needed') dy = q0[1] - p0[1] dx = q0[0] - p0[0] lhs0 = [-dy, dx] rhs0 = p0[1] * dx - dy * p0[0] else: lhs0 = [-m0, 1] rhs0 = p0[1] - m0 * p0[0] if m1 is None: if q1 is None: raise ValueError('either m1 or q1 is needed') dy = q1[1] - p1[1] dx = q1[0] - p1[0] lhs1 = [-dy, dx] rhs1 = p1[1] * dx - dy * p1[0] else: lhs1 = [-m1, 1] rhs1 = p1[1] - m1 * p1[0] a = np.array([lhs0, lhs1]) b = np.array([rhs0, rhs1]) try: px = np.linalg.solve(a, b) except: px = np.array([np.nan, np.nan]) return px