# How to determine if a point, specified by its latitude and longitude, falls within a given area using the shapely library in Python ?

Published: September 16, 2023

Tags: Python; Shapely;

Using the Shapely library in Python, it is possible to check if a point with a latitude and longitude is within an area or not. First, you need to create a shapely Geometry object for your area of interest (polygon). Then you can create a shapely Point object with your latitude and longitude. Finally, you can check if this point is within the area defined by the polygon using the contains or within methods, examples:

## Installing Shapely

Shapely can be installed, for example, along with geopandas as it is a dependency:

````pip install geopandas`
```

## Check if a point is within a given polygon

Please import the required Python libraries:

````from shapely.geometry import Point`
`from shapely.geometry.polygon import Polygon`

`import numpy as np`
```

### Define a shapely point

Define a target point, for example

````target_lat = 43.7580`
`target_lon = -65.3208`

`point = Point(target_lon,target_lat) # create point`
```

### Define a shapely polygon

Define an area using a shapely polygon

````lons_vect = np.array([-83.57072,-47.520203,-47.604248,-86.768524])`
`lats_vect = np.array([38.963406,44.391663,49.442635,43.54688,])`

`lons_lats_vect = np.column_stack((lons_vect, lats_vect)) # Reshape coordinates`

`polygon = Polygon(lons_lats_vect) # create polygon`
```

### Using shapely contains() method

Now you can use this method to check if a point with a latitude and longitude is within an area of interest. It is important to note that the coordinates should be provided in the same coordinate system as used for your polygon or else it will produce erroneous results.:

````polygon.contains(point)`
```

returns here

````True`
```

### Using shapely within() method

To check if a point is in the polygon

````point.within(polygon)`
```

returns here

````True`
```

### Data visualization

````import matplotlib.pyplot as plt`
`import matplotlib as mpl`
`import cartopy.crs as ccrs`
`import matplotlib.patches as mpatches`

`plt.figure(figsize=(16,9))`

`proj = ccrs.PlateCarree()`

`ease_extent = [-180., 180., 90., -90.]`

`ax = plt.axes(projection=proj)`

`ax.set_extent(ease_extent, crs=proj)`

`ax.gridlines(color='gray', linestyle='--')`

`ax.coastlines()`

`lat_corners = lats_vect`
`lon_corners = lons_vect`

`poly_corners = np.zeros((len(lat_corners), 2), np.float64)`
`poly_corners[:,0] = lon_corners`
`poly_corners[:,1] = lat_corners`

`poly = mpatches.Polygon(poly_corners, closed=True, fill=False, lw=1, fc=None, transform=ccrs.PlateCarree())`
`ax.add_patch(poly)`

`plt.scatter( target_lon, target_lat,`
`         color='red', linewidth=2, marker='o', s=20,`
`         transform=ccrs.PlateCarree(),`
`         )`

`plt.scatter( lons_vect, lats_vect,`
`         color='blue', linewidth=2, marker='o', s=20,`
`         transform=ccrs.PlateCarree(),`
`         )`

`plt.tight_layout()`

`plt.show()`

`plt.close()`
```

## Transforming Shapely Polygon and point in another coordinate system

Now let's consider the following example

````from shapely.geometry import Point`
`from shapely.geometry.polygon import Polygon`

`import numpy as np`

`target_lat =  47.18241500854492`
`target_lon = -115.90896606445312`

`lons_vect = np.array([-104.50501251, -140.47285461, -143.56848145, -104.57888031])`
`lats_vect = np.array([43.46247482, 38.01778412, 42.60084152, 48.49499893])`

`lons_lats_vect = np.column_stack((lons_vect, lats_vect)) # Reshape coordinates`

`polygon = Polygon(lons_lats_vect) # create polygon`

`point = Point(target_lon,target_lat) # create point`

`print(polygon.contains(point)) # check if polygon contains point`
`print(point.within(polygon)) # check if a point is in the polygon`
```

it will returns

````False`
`False`
```

Let's plot the data

````import matplotlib.pyplot as plt`
`import matplotlib as mpl`
`import cartopy.crs as ccrs`
`import matplotlib.patches as mpatches`

`plt.figure(figsize=(16,9))`

`proj = ccrs.PlateCarree()`

`ease_extent = [-180., 180., 90., -90.]`

`ax = plt.axes(projection=proj)`

`ax.set_extent(ease_extent, crs=proj)`

`ax.gridlines(color='gray', linestyle='--')`

`ax.coastlines()`

`lat_corners = lats_vect`
`lon_corners = lons_vect`

`poly_corners = np.zeros((len(lat_corners), 2), np.float64)`
`poly_corners[:,0] = lon_corners`
`poly_corners[:,1] = lat_corners`

`poly = mpatches.Polygon(poly_corners, closed=True, fill=False, lw=1, fc=None, transform=ccrs.PlateCarree())`
`ax.add_patch(poly)`

`plt.scatter( target_lon, target_lat,`
`         color='red', linewidth=2, marker='o', s=20,`
`         transform=ccrs.PlateCarree(),`
`         )`

`plt.scatter( lons_vect, lats_vect,`
`         color='blue', linewidth=2, marker='o', s=20,`
`         transform=ccrs.PlateCarree(),`
`         )`

`plt.tight_layout()`

`plt.title('How to check if a point with a latitude and longitude is in a given area or not using shapely in python ?' ,fontsize=14)`

`plt.savefig("shapely_point_within_02.png", bbox_inches='tight', dpi=100)`

`plt.show()`

`plt.close()`
```

However if we replace PlateCarree

````poly = mpatches.Polygon(poly_corners, closed=True, fill=False, lw=1, fc=None, transform=ccrs.PlateCarree())`
```

by Geodetic

````poly = mpatches.Polygon(poly_corners, closed=True, fill=False, lw=1, fc=None, transform=ccrs.Geodetic())`
```

we can see that the point is within the polygon !

Please take caution when selecting a coordinate system, as it can greatly impact your results.

### From EPSG:4326 to EPSG:32618

````from shapely.ops import transform`

`import pyproj`

`project = pyproj.Transformer.from_proj(`
`    pyproj.Proj("EPSG:4326"), `
`    pyproj.Proj('EPSG:32618'),`
`    always_xy=True`
`)`

`polygon = transform(project.transform, polygon)`
`point = transform(project.transform, point)`

`print(polygon.contains(point)) # check if polygon contains point`
`print(point.within(polygon)) # check if a point is in the polygon`
```

returns

````True`
`True`
```

## References

geopandas.GeoSeries.contains geopandas.org
geopandas.GeoSeries.within geopandas.org
shapely pypi.org
geopandas geopandas.org