# Extracting GPS points inside the region boundary [closed]

I have large set of GPS coordinates of human mobility in different parts of China, but I am only interested in points in Beijing City. These points are in excel csv file that I read into pandas dataframe using python 3.

I want to filter all GPS points inside Beijing city and discard all other points in the dataframe. One approach that I have in my mind is to define a rectangular bounding box and note down the following coordinates (latitude, longitudes) of bounding box.

``````Bottom Left  (x1,y1) = 39.77750000, 116.17944444
Top Left     (x1,y2) = 40.04722222, 116.17944444
Bottom Right (x2,y1) = 39.77750000, 116.58888889
Top Right    (x2,y2) = 40.04722222, 116.58888889
``````

Then for each point, compare its latitude and longitude values with all the vertices using multiple if else conditions, but I think their should be more optimal approaches to solve this in simple way. Is there any easier way or supporting library to solve this problem?

Example image is attached for clear understanding for the problem

The input dataset in given in the link https://drive.google.com/file/d/1JjvS7igTmrtLA4E5Rs5D6tsdAXqzpYqX/view

Use within:

Returns a Series of dtype('bool') with value True for each geometry that is within other.

``````from shapely import geometry
import geopandas as gpd
import pandas as pd

file = r'C:\folder\user141.csv'
crs = {'init': 'epsg:4326'}

#Create a geoseries holding the single polygon. Coordinates in counter-clockwise order
pointList = [(39.77750000, 116.17944444),(39.77750000, 116.58888889),(40.04722222, 116.58888889),(40.04722222, 116.17944444)]
poly = geometry.Polygon(pointList)
spoly = gpd.GeoSeries([poly],crs=crs)

#Create geodataframe of points
geometry = [geometry.Point(xy) for xy in zip(dfcsv.lat, dfcsv.long)]
dfpoints = gpd.GeoDataFrame(dfcsv, crs=crs, geometry=geometry)

#Create a subset dataframe of points within the polygon
subset = dfpoints[dfpoints.within(spoly.geometry.iloc[0])]
print('Number of points within polygon: ', subset.shape[0])
``````

Outputs:

``````Number of points within polygon:  58772
``````
• Works great! Two things: Commented Feb 16, 2022 at 23:06
• In `geometry = [geometry.Point(xy) for xy in zip(dfcsv.lat, dfcsv.long)]` I think you mean `dfcsv.lng` (`dfcsv.long` throws an `AttributeError`). And FWIW, I was able to run it w/o the `crs` parameter. Commented Feb 16, 2022 at 23:13