Mercurial > hg > nsaunier > traffic-intelligence

comparison python/moving.py @ 372:349eb1e09f45

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Cleaned the methods/functions indicating if a point is in a polygon In general, shapely should be used, especially for lots of points: from shapely.geometry import Polygon, Point poly = Polygon(array([[0,0],[0,1],[1,1],[1,0]])) p = Point(0.5,0.5) poly.contains(p) -> returns True poly.contains(Point(-1,-1)) -> returns False You can convert a moving.Point to a shapely point: p = moving.Point(1,2) p.asShapely() returns the equivalent shapely point If you have several points to test, use moving.pointsInPolygon(points, polygon) where points are moving.Point and polygon is a shapely polygon.
author Nicolas Saunier <nicolas.saunier@polymtl.ca>
date Tue, 16 Jul 2013 17:00:17 -0400
parents 027e254f0b53
children 2aed569f39e7
comparison
equal deleted inserted replaced
371:924e38c9f70e 372:349eb1e09f45
5 import cvutils 5 import cvutils
6 6
7 from math import sqrt 7 from math import sqrt
8 from numpy import median 8 from numpy import median
9 9
10 #from shapely.geometry import Polygon 10 try:
11 from shapely.geometry import Polygon, Point as shapelyPoint
12 from shapely.prepared import prep
13 shapelyAvailable = True
14 except ImportError:
15 print('Shapely library could not be loaded')
16 shapelyAvailable = False
11 17
12 __metaclass__ = type 18 __metaclass__ = type
13 19
14 class Interval(object): 20 class Interval(object):
15 '''Generic interval: a subset of real numbers (not iterable)''' 21 '''Generic interval: a subset of real numbers (not iterable)'''
195 return (self.x, self.y) 201 return (self.x, self.y)
196 202
197 def asint(self): 203 def asint(self):
198 return Point(int(self.x), int(self.y)) 204 return Point(int(self.x), int(self.y))
199 205
206 if shapelyAvailable:
207 def asShapely(self):
208 return shapelyPoint(self.x, self.y)
209
200 def project(self, homography): 210 def project(self, homography):
201 from numpy.core.multiarray import array 211 from numpy.core.multiarray import array
202 projected = cvutils.projectArray(homography, array([[self.x], [self.y]])) 212 projected = cvutils.projectArray(homography, array([[self.x], [self.y]]))
203 return Point(projected[0], projected[1]) 213 return Point(projected[0], projected[1])
204 214
205 def inPolygon(self, poly): 215 def inPolygonNoShapely(self, polygon):
206 '''Returns if the point x, y is inside the polygon. 216 '''Indicates if the point x, y is inside the polygon
207 The polygon is defined by the ordered list of points in poly 217 (array of Nx2 coordinates of the polygon vertices)
208 218
209 taken from http://www.ariel.com.au/a/python-point-int-poly.html 219 taken from http://www.ariel.com.au/a/python-point-int-poly.html
210 220
211 Use points_inside_poly from matplotlib.nxutils''' 221 Use Polygon.contains if Shapely is installed'''
212 222
213 n = len(poly); 223 n = polygon.shape[0];
214 counter = 0; 224 counter = 0;
215 225
216 p1 = poly[0]; 226 p1 = polygon[0,:];
217 for i in range(n+1): 227 for i in range(n+1):
218 p2 = poly[i % n]; 228 p2 = polygon[i % n,:];
219 if self.y > min(p1.y,p2.y): 229 if self.y > min(p1[1],p2[1]):
220 if self.y <= max(p1.y,p2.y): 230 if self.y <= max(p1[1],p2[1]):
221 if self.x <= max(p1.x,p2.x): 231 if self.x <= max(p1[0],p2[0]):
222 if p1.y != p2.y: 232 if p1[1] != p2[1]:
223 xinters = (self.y-p1.y)*(p2.x-p1.x)/(p2.y-p1.y)+p1.x; 233 xinters = (self.y-p1[1])*(p2[0]-p1[0])/(p2[1]-p1[1])+p1[0];
224 if p1.x == p2.x or self.x <= xinters: 234 if p1[0] == p2[0] or self.x <= xinters:
225 counter+=1; 235 counter+=1;
226 p1=p2 236 p1=p2
227 return (counter%2 == 1); 237 return (counter%2 == 1);
228 238
229 @staticmethod 239 @staticmethod
242 252
243 @staticmethod 253 @staticmethod
244 def plotAll(points, **kwargs): 254 def plotAll(points, **kwargs):
245 from matplotlib.pyplot import scatter 255 from matplotlib.pyplot import scatter
246 scatter([p.x for p in points],[p.y for p in points], **kwargs) 256 scatter([p.x for p in points],[p.y for p in points], **kwargs)
257
258 if shapelyAvailable:
259 def pointsInPolygon(points, polygon):
260 '''Optimized tests of a series of points within (Shapely) polygon '''
261 prepared_polygon = prep(polygon)
262 return filter(prepared_polygon.contains, points)
263
247 264
248 class NormAngle(object): 265 class NormAngle(object):
249 '''Alternate encoding of a point, by its norm and orientation''' 266 '''Alternate encoding of a point, by its norm and orientation'''
250 267
251 def __init__(self, norm, angle): 268 def __init__(self, norm, angle):
523 return Trajectory([self.positions[0][inter.first:inter.last], 540 return Trajectory([self.positions[0][inter.first:inter.last],
524 self.positions[1][inter.first:inter.last]]) 541 self.positions[1][inter.first:inter.last]])
525 else: 542 else:
526 return None 543 return None
527 544
528 def getTrajectoryInPolygon(self, polygon): 545 def getTrajectoryInPolygonNoShapely(self, polygon):
529 '''Returns the set of points inside the polygon 546 '''Returns the trajectory built with the set of points inside the polygon
530 (array of Nx2 coordinates of the polygon vertices)''' 547 (array of Nx2 coordinates of the polygon vertices)'''
531 import matplotlib.nxutils as nx
532 traj = Trajectory() 548 traj = Trajectory()
533 result = nx.points_inside_poly(self.asArray().T, polygon) 549 for p in self:
534 for i in xrange(self.length()): 550 if p.inPolygonNoShapely(polygon):
535 if result[i]: 551 traj.addPosition(p)
536 traj.addPositionXY(self.positions[0][i], self.positions[1][i])
537 return traj 552 return traj
538 553
539 # version 2: use shapely polygon contains 554 if shapelyAvailable:
540 555 def getTrajectoryInPolygon(self, polygon):
556 '''Returns the trajectory built with the set of points inside the (shapely) polygon'''
557 traj = Trajectory()
558 points = [p.asShapely() for p in self]
559 for p in pointsInPolygon(points, polygon):
560 traj.addPositionXY(p.x, p.y)
561 return traj
562
541 @staticmethod 563 @staticmethod
542 def lcss(t1, t2, lcss): 564 def lcss(t1, t2, lcss):
543 return lcss.compute(t1, t2) 565 return lcss.compute(t1, t2)
544 566
545 class CurvilinearTrajectory(Trajectory): 567 class CurvilinearTrajectory(Trajectory):