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
import folium
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.19.0 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 : 1.0.1 pandas : 2.2.2 fiona : 1.9.6 numpy : 1.26.4 shapely : 2.0.4 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: >
In [ ]:
m = my_traj.explore(color="blue", style_kwds={"weight": 4}, name="Trajectory")
intersections.explore(m=m, color="red", style_kwds={"weight": 4}, name="Intersection")
folium.TileLayer("OpenStreetMap").add_to(m)
folium.LayerControl().add_to(m)
m
Out[ ]:
Make this Notebook Trusted to load map: File -> Trust Notebook
Clipping a TrajectoryCollection¶
Alternatively, using TrajectoryCollection:
In [ ]:
clipped = tc.clip(polygon)
clipped
Out[ ]:
TrajectoryCollection with 2 trajectories
In [ ]:
clipped.plot()
Out[ ]:
<Axes: >
In [ ]:
clipped.explore(
column="trajectory_id",
cmap="cool",
tiles="CartoDB positron",
style_kwds={"weight": 4},
)
Out[ ]:
Make this Notebook Trusted to load map: File -> Trust Notebook
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.explore(color="blue", style_kwds={"weight": 4})
Out[ ]:
Make this Notebook Trusted to load map: File -> Trust Notebook
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.90468) | 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.9048) | 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.explore(
column="trajectory_id",
cmap="autumn",
tiles="CartoDB positron",
style_kwds={"weight": 4},
)
Out[ ]:
Make this Notebook Trusted to load map: File -> Trust Notebook
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.90468) | 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.9048) | 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 [ ]: