Jaki jest najszybszy sposób sprawdzenia, czy punkt znajduje się wewnątrz wielokąta w pytonie?

Question

Znalazłem dwie główne metody sprawdzania, czy punkt należy do wielokąta. Jedną z nich jest metoda ray tracingu wykorzystana jako here, która jest najbardziej zalecaną odpowiedzią, druga wykorzystuje matplotlib path.contains_points (co wydaje mi się nieco niejasne). Będę musiał ciągle sprawdzać wiele punktów. Czy ktokolwiek wie, czy któryś z tych dwóch jest bardziej zalecany niż drugi, czy są jeszcze lepsze trzecie opcje?

UPDATE:

Sprawdziłem dwie metody i matplotlib wygląda o wiele szybciej.

from time import time 
import numpy as np 
import matplotlib.path as mpltPath 

# regular polygon for testing 
lenpoly = 100 
polygon = [[np.sin(x)+0.5,np.cos(x)+0.5] for x in np.linspace(0,2*np.pi,lenpoly)[:-1]] 

# random points set of points to test 
N = 10000 
points = zip(np.random.random(N),np.random.random(N)) 


# Ray tracing 
def ray_tracing_method(x,y,poly): 

    n = len(poly) 
    inside = False 

    p1x,p1y = poly[0] 
    for i in range(n+1): 
     p2x,p2y = poly[i % n] 
     if y > min(p1y,p2y): 
      if y <= max(p1y,p2y): 
       if x <= max(p1x,p2x): 
        if p1y != p2y: 
         xints = (y-p1y)*(p2x-p1x)/(p2y-p1y)+p1x 
        if p1x == p2x or x <= xints: 
         inside = not inside 
     p1x,p1y = p2x,p2y 

    return inside 

start_time = time() 
inside1 = [ray_tracing_method(point[0], point[1], polygon) for point in points] 
print "Ray Tracing Elapsed time: " + str(time()-start_time) 

# Matplotlib mplPath 
start_time = time() 
path = mpltPath.Path(polygon) 
inside2 = path.contains_points(points) 
print "Matplotlib contains_points Elapsed time: " + str(time()-start_time) 

, zawierające

Ray Tracing Elapsed time: 0.441395998001 
Matplotlib contains_points Elapsed time: 0.00994491577148 

Bez różnicy względnej otrzymano stosując jeden trójkąt zamiast wieloboku 100 stron. Będę również sprawdzić zgrabna ponieważ wygląda pakiet właśnie poświęcony tego rodzaju problemów

Answer 1

można rozważyć shapely:

from shapely.geometry import Point 
from shapely.geometry.polygon import Polygon 

point = Point(0.5, 0.5) 
polygon = Polygon([(0, 0), (0, 1), (1, 1), (1, 0)]) 
print(polygon.contains(point)) 

od metod Wspominałeś Użyłem tylko drugi, path.contains_points, i działa dobrze. W każdym przypadku, w zależności od precyzji potrzebnej do testu, sugerowałbym utworzenie siatki liczbowej bool ze wszystkimi węzłami wewnątrz wielokąta, aby miała wartość True (Fałsz, jeśli nie). Jeśli masz zamiar wykonać test za dużo punktów to może być szybsze (chociaż zawiadomienie ten opiera robisz test w „piksel” tolerancja):

from matplotlib import path 
import matplotlib.pyplot as plt 
import numpy as np 

first = -3 
size = (3-first)/100 
xv,yv = np.meshgrid(np.linspace(-3,3,100),np.linspace(-3,3,100)) 
p = path.Path([(0,0), (0, 1), (1, 1), (1, 0)]) # square with legs length 1 and bottom left corner at the origin 
flags = p.contains_points(np.hstack((xv.flatten()[:,np.newaxis],yv.flatten()[:,np.newaxis]))) 
grid = np.zeros((101,101),dtype='bool') 
grid[((xv.flatten()-first)/size).astype('int'),((yv.flatten()-first)/size).astype('int')] = flags 

xi,yi = np.random.randint(-300,300,100)/100,np.random.randint(-300,300,100)/100 
vflag = grid[((xi-first)/size).astype('int'),((yi-first)/size).astype('int')] 
plt.imshow(grid.T,origin='lower',interpolation='nearest',cmap='binary') 
plt.scatter(((xi-first)/size).astype('int'),((yi-first)/size).astype('int'),c=vflag,cmap='Greens',s=90) 
plt.show() 

wyniki to:

Answer 2

Twój test jest dobry, ale to mierzy tylko niektóre szczególne sytuację: mamy jeden wielokąt z wielu wierzchołków i długi szereg punktów, aby sprawdzić je w wieloboku.

Co więcej, sądzę, że jesteś pomiarowy nie matplotlib-wewnątrz-wielokąta metodę vs Ray-metoda ale matplotlib-jakoś zoptymalizowane-iteracji vs proste-list-iteracji

Zróbmy N niezależne porównania (N par punktów i wielokątów)?

# ... your code... 
lenpoly = 100 
polygon = [[np.sin(x)+0.5,np.cos(x)+0.5] for x in np.linspace(0,2*np.pi,lenpoly)[:-1]] 

M = 10000 
start_time = time() 
# Ray tracing 
for i in range(M): 
    x,y = np.random.random(), np.random.random() 
    inside1 = ray_tracing_method(x,y, polygon) 
print "Ray Tracing Elapsed time: " + str(time()-start_time) 

# Matplotlib mplPath 
start_time = time() 
for i in range(M): 
    x,y = np.random.random(), np.random.random() 
    inside2 = path.contains_points([[x,y]]) 
print "Matplotlib contains_points Elapsed time: " + str(time()-start_time) 

Wynik:

Ray Tracing Elapsed time: 0.548588991165 
Matplotlib contains_points Elapsed time: 0.103765010834 

matplotlib jest wciąż znacznie lepiej, ale nie 100 razy lepiej. Teraz spróbujmy znacznie prostsze wielokąt ...

lenpoly = 5 
# ... same code 

wynik:

Ray Tracing Elapsed time: 0.0727779865265 
Matplotlib contains_points Elapsed time: 0.105288982391 

Answer 3

Jeśli prędkość jest co trzeba i dodatkowych zależności nie są problemem, to może znaleźć numba bardzo przydatne (teraz jest dość łatwa w instalacji, na dowolnej platformie).Proponowane podejście klasyczne: ray_tracing może być łatwo przeniesione do numba przy użyciu dekoratora numba @jit i przeniesienia wielokąta do tablicy numpy. Kod powinien wyglądać następująco:

@jit(nopython=True) 
def ray_tracing(x,y,poly): 
    n = len(poly) 
    inside = False 
    p2x = 0.0 
    p2y = 0.0 
    xints = 0.0 
    p1x,p1y = poly[0] 
    for i in range(n+1): 
     p2x,p2y = poly[i % n] 
     if y > min(p1y,p2y): 
      if y <= max(p1y,p2y): 
       if x <= max(p1x,p2x): 
        if p1y != p2y: 
        xints = (y-p1y)*(p2x-p1x)/(p2y-p1y)+p1x 
        if p1x == p2x or x <= xints: 
         inside = not inside 
     p1x,p1y = p2x,p2y 

    return inside 

Pierwsze wykonanie potrwa trochę dłużej niż wszelkie późniejsze rozmowy:

%%time 
polygon=np.array(polygon) 
inside1 = [numba_ray_tracing_method(point[0], point[1], polygon) for 
point in points] 

CPU times: user 129 ms, sys: 4.08 ms, total: 133 ms 
Wall time: 132 ms 

który po kompilacji będzie maleć do:

CPU times: user 18.7 ms, sys: 320 µs, total: 19.1 ms 
Wall time: 18.4 ms 

Jeśli potrzebujesz prędkości przy pierwszym wywołaniu funkcji, możesz następnie wstępnie skompilować kod w module przy użyciu pycc. Przechowywać funkcję w src.py jak:

from numba import jit 
from numba.pycc import CC 
cc = CC('nbspatial') 


@cc.export('ray_tracing', 'b1(f8, f8, f8[:,:])') 
@jit(nopython=True) 
def ray_tracing(x,y,poly): 
    n = len(poly) 
    inside = False 
    p2x = 0.0 
    p2y = 0.0 
    xints = 0.0 
    p1x,p1y = poly[0] 
    for i in range(n+1): 
     p2x,p2y = poly[i % n] 
     if y > min(p1y,p2y): 
      if y <= max(p1y,p2y): 
       if x <= max(p1x,p2x): 
        if p1y != p2y: 
         xints = (y-p1y)*(p2x-p1x)/(p2y-p1y)+p1x 
        if p1x == p2x or x <= xints: 
         inside = not inside 
     p1x,p1y = p2x,p2y 

    return inside 


if __name__ == "__main__": 
    cc.compile() 

zbudować i uruchomić z python src.py:

import nbspatial 

import numpy as np 
lenpoly = 100 
polygon = [[np.sin(x)+0.5,np.cos(x)+0.5] for x in 
np.linspace(0,2*np.pi,lenpoly)[:-1]] 

# random points set of points to test 
N = 10000 
# making a list instead of a generator to help debug 
points = zip(np.random.random(N),np.random.random(N)) 

polygon = np.array(polygon) 

%%time 
result = [nbspatial.ray_tracing(point[0], point[1], polygon) for point in points] 

CPU times: user 20.7 ms, sys: 64 µs, total: 20.8 ms 
Wall time: 19.9 ms 

W kodzie Numba użyłem: „B1 (F8 F8 F8 [:, :]) '

Aby skompilować z nopython=True, każdy var musi zostać zadeklarowany przed for loop.

W kodzie prebuild src linia:

@cc.export('ray_tracing', 'b1(f8, f8, f8[:,:])') 

służy do deklarowania nazwę funkcji i jej rodzaje I/O Var, wyjście logiczną b1 i dwa pływaki f8 i dwuwymiarową tablicę pływaków f8[:,:] jako dane wejściowe.