Computing intersections with polygons¶
![No description has been provided for this image](https://movingpandas.github.io/movingpandas/assets/img/movingpandas.png)
Clipping and intersection functions can be used to extract trajectory segments that are located within an area of interest polygon.
In [ ]:
import pandas as pd
import geopandas as gpd
import movingpandas as mpd
import shapely as shp
import hvplot.pandas
from geopandas import GeoDataFrame, read_file
from shapely.geometry import Point, LineString, Polygon
from datetime import datetime, timedelta
from holoviews import opts
import warnings
warnings.filterwarnings('ignore')
opts.defaults(opts.Overlay(active_tools=['wheel_zoom'], frame_width=500, frame_height=400))
mpd.show_versions()
MovingPandas 0.18.1 SYSTEM INFO ----------- python : 3.10.14 | packaged by conda-forge | (main, Mar 20 2024, 12:40:08) [MSC v.1938 64 bit (AMD64)] executable : c:\Users\Agarkovam\AppData\Local\miniforge3\envs\mpd-ex\python.exe machine : Windows-10-10.0.19045-SP0 GEOS, GDAL, PROJ INFO --------------------- GEOS : None GEOS lib : None GDAL : 3.8.5 GDAL data dir: None PROJ : 9.4.0 PROJ data dir: C:\Users\Agarkovam\AppData\Local\miniforge3\envs\mpd-ex\Library\share\proj PYTHON DEPENDENCIES ------------------- geopandas : 0.14.4 pandas : 2.2.2 fiona : 1.9.6 numpy : 1.26.4 shapely : 2.0.4 rtree : 1.2.0 pyproj : 3.6.1 matplotlib : 3.8.4 mapclassify: 2.6.1 geopy : 2.4.1 holoviews : 1.17.1 hvplot : 0.8.3 geoviews : 1.9.6 stonesoup : 1.2
In [ ]:
gdf = read_file('../data/geolife_small.gpkg')
tc = mpd.TrajectoryCollection(gdf, 'trajectory_id', t='t')
Clipping a Trajectory¶
In [ ]:
xmin, xmax, ymin, ymax = 116.365035,116.3702945,39.904675,39.907728
polygon = Polygon([(xmin,ymin), (xmin,ymax), (xmax,ymax), (xmax,ymin), (xmin,ymin)])
my_traj = tc.trajectories[2]
intersections = my_traj.clip(polygon)
print("Found {} intersections".format(len(intersections)))
Found 1 intersections
In [ ]:
ax = my_traj.plot()
gpd.GeoSeries(polygon).plot(ax=ax, color='lightgray')
intersections.plot(ax=ax, color='red', linewidth=5)
Out[ ]:
<Axes: >
Clipping a TrajectoryCollection¶
Alternatively, using TrajectoryCollection:
In [ ]:
clipped = tc.clip(polygon)
clipped
Out[ ]:
TrajectoryCollection with 2 trajectories
In [ ]:
clipped.plot()
Out[ ]:
<Axes: >
Computing intersections for a Trajectory¶
In [ ]:
polygon_feature = {
"geometry": polygon,
"properties": {'field1': 'abc'}
}
In [ ]:
my_traj = tc.trajectories[2]
intersections = my_traj.intersection(polygon_feature)
intersections
Out[ ]:
TrajectoryCollection with 1 trajectories
In [ ]:
intersections.plot()
Out[ ]:
<Axes: >
In [ ]:
intersections.to_point_gdf()
Out[ ]:
id | sequence | trajectory_id | tracker | geometry | intersecting_field1 | |
---|---|---|---|---|---|---|
t | ||||||
2009-02-04 10:43:04.222205 | 3225 | 773 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36799 39.90467) | abc |
2009-02-04 10:43:05.000000 | 3226 | 774 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36798 39.90471) | abc |
2009-02-04 10:43:07.000000 | 3227 | 775 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36795 39.90480) | abc |
2009-02-04 10:43:09.000000 | 3228 | 776 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36793 39.90487) | abc |
2009-02-04 10:43:11.000000 | 3229 | 777 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36794 39.90495) | abc |
... | ... | ... | ... | ... | ... | ... |
2009-02-04 10:48:14.000000 | 3390 | 938 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36535 39.90589) | abc |
2009-02-04 10:48:15.000000 | 3391 | 939 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36525 39.90589) | abc |
2009-02-04 10:48:16.000000 | 3392 | 940 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36516 39.90589) | abc |
2009-02-04 10:48:17.000000 | 3393 | 941 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36507 39.90589) | abc |
2009-02-04 10:48:17.399995 | 3393 | 941 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36504 39.90589) | abc |
170 rows × 6 columns
Computing intersections for a TrajectoryCollection¶
In [ ]:
intersections = tc.intersection(polygon_feature)
intersections
Out[ ]:
TrajectoryCollection with 2 trajectories
In [ ]:
intersections.plot()
Out[ ]:
<Axes: >
In [ ]:
intersections.to_point_gdf()
Out[ ]:
id | sequence | trajectory_id | tracker | geometry | intersecting_field1 | |
---|---|---|---|---|---|---|
t | ||||||
2009-02-04 10:43:04.222205 | 3225 | 773 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36799 39.90467) | abc |
2009-02-04 10:43:05.000000 | 3226 | 774 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36798 39.90471) | abc |
2009-02-04 10:43:07.000000 | 3227 | 775 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36795 39.90480) | abc |
2009-02-04 10:43:09.000000 | 3228 | 776 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36793 39.90487) | abc |
2009-02-04 10:43:11.000000 | 3229 | 777 | 3_2009-02-04 10:43:04.222205 | 2 | POINT (116.36794 39.90495) | abc |
... | ... | ... | ... | ... | ... | ... |
2009-03-10 11:08:05.000000 | 4934 | 672 | 4_2009-03-10 11:01:41.826068 | 2 | POINT (116.36539 39.90588) | abc |
2009-03-10 11:08:06.000000 | 4935 | 673 | 4_2009-03-10 11:01:41.826068 | 2 | POINT (116.36529 39.90588) | abc |
2009-03-10 11:08:07.000000 | 4936 | 674 | 4_2009-03-10 11:01:41.826068 | 2 | POINT (116.36519 39.90588) | abc |
2009-03-10 11:08:08.000000 | 4937 | 675 | 4_2009-03-10 11:01:41.826068 | 2 | POINT (116.36509 39.90588) | abc |
2009-03-10 11:08:08.575754 | 4937 | 675 | 4_2009-03-10 11:01:41.826068 | 2 | POINT (116.36504 39.90588) | abc |
375 rows × 6 columns
In [ ]: