--- a/python/moving.py	Mon Feb 08 12:24:26 2016 -0500
+++ b/python/moving.py	Thu Feb 11 11:56:50 2016 -0500
@@ -12,7 +12,7 @@
 
 try:
     from shapely.geometry import Polygon, Point as shapelyPoint
-    from shapely.prepared import prep
+    from shapely.prepared import prep, PreparedGeometry
     shapelyAvailable = True
 except ImportError:
     print('Shapely library could not be loaded')
@@ -349,8 +349,11 @@
 
 if shapelyAvailable:
     def pointsInPolygon(points, polygon):
-        '''Optimized tests of a series of points within (Shapely) polygon '''
-        prepared_polygon = prep(polygon)
+        '''Optimized tests of a series of points within (Shapely) polygon (not prepared)'''
+        if type(polygon) == PreparedGeometry:
+            prepared_polygon = polygon
+        else:
+            prepared_polygon = prep(polygon)
         return filter(prepared_polygon.contains, points)
 
 # Functions for coordinate transformation
@@ -933,32 +936,50 @@
                            self.positions[1][::step]])
 
     if shapelyAvailable:
-        def getTrajectoryInPolygon(self, polygon):
-            '''Returns the trajectory built with the set of points inside the (shapely) polygon'''
+        def getTrajectoryInPolygon(self, polygon, t2 = None):
+            '''Returns the trajectory built with the set of points inside the (shapely) polygon
+            The polygon could be a prepared polygon (faster) from prepared.prep
+
+            t2 is another trajectory (could be velocities) 
+            which is filtered based on the first (self) trajectory'''
             traj = Trajectory()
-            points = [p.asShapely() for p in self]
-            for p in pointsInPolygon(points, polygon):
-                    traj.addPositionXY(p.x, p.y)
-            return traj
+            inPolygon = []
+            for x, y in zip(self.positions[0], self.positions[1]):
+                inPolygon.append(polygon.contains(shapelyPoint(x, y)))
+                if inPolygon[-1]:
+                    traj.addPositionXY(x, y)
+            traj2 = Trajectory()
+            if t2 is not None:
+                for inp, x, y in zip(inPolygon, t2.positions[0], t2.positions[1]):
+                    if inp:
+                        traj2.addPositionXY(x, y)
+            return traj, traj2
 
         def proportionInPolygon(self, polygon, minProportion = 0.5):
-            pointsIn = pointsInPolygon([p.asShapely() for p in self], polygon)
+            inPolygon = [polygon.contains(shapelyPoint(x, y)) for x, y in zip(self.positions[0], self.positions[1])]
             lengthThreshold = float(self.length())*minProportion
-            return len(pointsIn) >= lengthThreshold
+            return sum(inPolygon) >= lengthThreshold
     else:
-        def getTrajectoryInPolygon(self, polygon):
+        def getTrajectoryInPolygon(self, polygon, t2 = None):
             '''Returns the trajectory built with the set of points inside the polygon
             (array of Nx2 coordinates of the polygon vertices)'''
             traj = Trajectory()
+            inPolygon = []
             for p in self:
-                if p.inPolygon(polygon):
+                inPolygon.append(p.inPolygon(polygon))
+                if inPolygon[-1]:
                     traj.addPosition(p)
-            return traj
+            traj2 = Trajectory()
+            if t2 is not None:
+                for inp, x, y in zip(inPolygon, t2.positions[0], t2.positions[1]):
+                    if inp:
+                        traj2.addPositionXY(p.x, p.y)
+            return traj, traj2
 
         def proportionInPolygon(self, polygon, minProportion = 0.5):
-            pointsInPolygon = [p.inPolygon(polygon) for p in self]
+            inPolygon = [p.inPolygon(polygon) for p in self]
             lengthThreshold = float(self.length())*minProportion
-            return len(pointsInPolygon) >= lengthThreshold
+            return sum(inPolygon) >= lengthThreshold
 
     @staticmethod
     def lcss(t1, t2, lcss):

--- a/python/tests/moving_shapely.txt	Mon Feb 08 12:24:26 2016 -0500
+++ b/python/tests/moving_shapely.txt	Thu Feb 11 11:56:50 2016 -0500
@@ -1,9 +1,28 @@
 >>> from moving import *
 >>> from shapely.geometry import Polygon
+>>> from shapely.prepared import prep
 
 >>> t1 = Trajectory([[0.5,1.5,2.5],[0.5,3.5,6.5]])
->>> t1.getTrajectoryInPolygon(Polygon([[0,0],[4,0],[4,3],[0,3]]))
+>>> poly = Polygon([[0,0],[4,0],[4,3],[0,3]])
+>>> sub1, sub2 = t1.getTrajectoryInPolygon(poly)
+>>> sub1
+(0.500000,0.500000)
+>>> sub1, sub2 = t1.getTrajectoryInPolygon(Polygon([[10,10],[14,10],[14,13],[10,13]]))
+>>> sub1.length()
+0
+>>> sub1, sub2 = t1.getTrajectoryInPolygon(prep(poly))
+>>> sub1
 (0.500000,0.500000)
->>> t1.getTrajectoryInPolygon(Polygon([[10,10],[14,10],[14,13],[10,13]])).length()
-0
+>>> t2 = t1.differentiate(True)
+>>> sub1, sub2 = t1.getTrajectoryInPolygon(prep(poly), t2)
+>>> sub1.length() == sub2.length()
+True
+>>> sub1
+(0.500000,0.500000)
+>>> sub2
+(1.000000,3.000000)
 
+>>> t1.proportionInPolygon(poly, 0.5)
+False
+>>> t1.proportionInPolygon(poly, 0.3)
+True