-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathinpoly2_mat.m
143 lines (99 loc) · 4.06 KB
/
inpoly2_mat.m
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
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
function [stat,bnds] = ...
inpoly2_mat(vert,node,edge,fTOL,lbar)
%INPOLY2_MAT the local m-code version of the crossing-number
%test. Loop over edges; do a binary-search for the first ve-
%rtex that intersects with the edge y-range; do crossing-nu-
%mber comparisons; break when the local y-range is exceeded.
% Darren Engwirda : 2017 --
% Email : [email protected]
% Last updated : 19/12/2020
feps = fTOL * lbar ^ +1 ;
veps = fTOL * lbar ^ +1 ;
nvrt = size (vert,1) ;
nnod = size (node,1) ;
nedg = size (edge,1) ;
stat = false(nvrt,1) ;
bnds = false(nvrt,1) ;
%----------------------------------- loop over polygon edges
for epos = +1 : size(edge,1)
inod = edge(epos,1) ;
jnod = edge(epos,2) ;
%------------------------------- calc. edge bounding-box
yone = node(inod,2) ;
ytwo = node(jnod,2) ;
xone = node(inod,1) ;
xtwo = node(jnod,1) ;
xmin = min(xone,xtwo) ;
xmax = max(xone,xtwo) ;
xmin = xmin - veps;
xmax = xmax + veps;
ymin = yone - veps;
ymax = ytwo + veps;
ydel = ytwo - yone;
xdel = xtwo - xone;
edel = abs(xdel) + ydel ;
%------------------------------- find top VERT(:,2)<YONE
ilow = +1; iupp = nvrt;
while (ilow < iupp - 1) % binary search
imid = ilow ...
+ floor((iupp-ilow) / 2);
if (vert(imid,2) < ymin)
ilow = imid ;
else
iupp = imid ;
end
end
if (vert(ilow,2) >= ymin)
ilow = ilow - 1 ;
end
%------------------------------- calc. edge-intersection
for jpos = ilow+1 : nvrt
if (bnds(jpos)), continue ; end
xpos = vert(jpos,1) ;
ypos = vert(jpos,2) ;
if (ypos <= ymax)
if (xpos >= xmin)
if (xpos <= xmax)
%------------------- compute crossing number
mul1 = ...
ydel * (xpos - xone) ;
mul2 = ...
xdel * (ypos - yone) ;
if ((feps * edel) >= ...
abs(mul2 - mul1) )
%------------------- BNDS -- approx. on edge
bnds(jpos)= true ;
stat(jpos)= true ;
elseif (ypos == yone ...
&& xpos == xone )
%------------------- BNDS -- match about ONE
bnds(jpos)= true ;
stat(jpos)= true ;
elseif (ypos == ytwo ...
&& xpos == xtwo )
%------------------- BNDS -- match about TWO
bnds(jpos)= true ;
stat(jpos)= true ;
elseif (mul1 < mul2)
if (ypos >= yone ...
&& ypos < ytwo)
%------------------- advance crossing number
stat(jpos) = ...
~stat (jpos) ;
end
end
end
else
if (ypos >= yone ...
&& ypos < ytwo)
%------------------- advance crossing number
stat(jpos) = ...
~stat (jpos) ;
end
end
else
break ; % done -- due to the sort
end
end
end
end