-
Notifications
You must be signed in to change notification settings - Fork 0
/
nb_is_inside_sm.py
52 lines (39 loc) · 1.5 KB
/
nb_is_inside_sm.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
from numba import jit, njit
import numba
import numpy as np
@jit(nopython=True)
def is_inside_sm(polygon, point):
length = len(polygon)-1
dy2 = point[1] - polygon[0][1]
intersections = 0
ii = 0
jj = 1
while ii < length:
dy = dy2
dy2 = point[1] - polygon[jj][1]
# consider only lines which are not completely above/bellow/right from the point
if dy*dy2 <= 0.0 and (point[0] >= polygon[ii][0] or point[0] >= polygon[jj][0]):
# non-horizontal line
if dy < 0 or dy2 < 0:
F = dy*(polygon[jj][0] - polygon[ii][0]) / \
(dy-dy2) + polygon[ii][0]
# if line is left from the point - the ray moving towards left, will intersect it
if point[0] > F:
intersections += 1
elif point[0] == F: # point on line
return 2
# point on upper peak (dy2=dx2=0) or horizontal line (dy=dy2=0 and dx*dx2<=0)
elif dy2 == 0 and (point[0] == polygon[jj][0] or (dy == 0 and (point[0]-polygon[ii][0])*(point[0]-polygon[jj][0]) <= 0)):
return 2
ii = jj
jj += 1
# print 'intersections =', intersections
return intersections & 1
# px = point_df['Easting'].values.get()
@njit(parallel=True)
def get_spatial_selection_cpu(points, polygon):
ln = len(points)
D = np.empty(ln, dtype=numba.boolean)
for i in numba.prange(ln):
D[i] = is_inside_sm(polygon, points[i])
return D