python code to find lat and lon in shape file

I am tyring to find if the lon and lat values are there in the shape file. When I am reading shape file, I can see the POLYGON values POLYGON ((130.8928 -12.4831, 130.8922 -12.484,..)

When I am passing the first lon, lat values which is available as first value in polygon/shape. It's still returing the FALSE value. Please suggest what I am doing wrong here?

My python code:

import geopandas as gpd
from shapely.geometry import shape, Point, Polygon

in_shp = 'C:\\shape_file.shp'

soa_shape_map = gpd.read_file(in_shp)

# Convert to WGS84
soa_shape_map_geo = soa_shape_map.to_crs(epsg=4326)

# Extract geometry
polygon = soa_shape_map_geo['geometry'].iloc[0]

# Extract coordinates and format them
formatted_polygon = [(round(coord[0], 4), round(coord[1], 4)) for coord in polygon.exterior.coords]

# Create a Shapely Polygon object
polygon_shape = Polygon(formatted_polygon)
print(f"polygon_shape: {polygon_shape}")

def check_point_in_polygon (lon, lat):
    print(f"lon:{lon}")
    print (f"lat: {lat}")
   
    # build a shapely point from your geopoint
    point = Point(round(lon,4), round(lat,4))
    print (f"point: {point}")

    # the contains function does exactly what you want
    return polygon_shape.contains(point)

print (check_point_in_polygon(130.8928, -12.4831))
print (check_point_in_polygon(-12.4831, 130.8928))

Solution

I think there are multiple issues that might cause a False value to be returned :

wrong use of exterior to get coordinates;
rounding values to 4 digits when actually polygons in the shapefile might have 4, 6 or 8, can cause the point to fall outside of the polygon;
use of contains method for a point on the boundary will return False.

You might want to use the intersects method that will work on a boundary point also.

In the following example, the use of intersects allows for a point on the boundary to be flagged as True when it was flagged as False with contains.

The two other points are here for control :

One inside the polygon always True;
One outside the polygon always False.

Here is some working code for demonstration :

import geopandas as gpd
from shapely.geometry import Polygon, MultiPolygon, Point

# We are creating MultiPolygons with holes to mimic a shapefile with random 

# MultiPolygon 1
exterior1 = [(0.039427, 0.184263), (0.039632, 0.184263), (0.039632, 0.184528), (0.039427, 0.184528)]
interior1 = [(0.039487, 0.184313), (0.039487, 0.184478), (0.039572, 0.184478), (0.039572, 0.184313)]
multipolygon1 = MultiPolygon([Polygon(exterior1, [interior1])])

# MultiPolygon 2
exterior2 = [(0.058101, 0.194362), (0.058306, 0.194362), (0.058306, 0.194627), (0.058101, 0.194627)]
interior2 = [(0.058161, 0.194412), (0.058161, 0.194577), (0.058246, 0.194577), (0.058246, 0.194412)]
multipolygon2 = MultiPolygon([Polygon(exterior2, [interior2])])

# Create a GeoDataFrame to mimic shapefile import (already in EPSG 4326)
gdf = gpd.GeoDataFrame(
    {"id": [1, 2], 
     "geometry": [multipolygon1, multipolygon2]},
    crs="epsg:4326",
)

# Extract geometry without using exterior
polygon = gdf["geometry"].iloc[0]


# First function to work with contains (strict, not on boundaries)
def check_point_in_polygon_with_contains(lon, lat):
    return polygon.contains(Point(lon, lat))


# Second function to works with intersects
def check_point_in_polygon_with_intersects(lon, lat): 
    return Point(lon, lat).intersects(polygon)

################
# 3 test cases : 
################

# a) Point on boundary :
p1 = (exterior1[0][0], exterior1[0][1])

# b) Point slightly inside :
p2 = (exterior1[0][0] * 1.001, exterior1[0][1] * 1.001)

# c) Point slighlty outside :
p3 = (exterior1[0][0] * 0.999, exterior1[0][1] * 0.999)

# testing CONTAINS
print(check_point_in_polygon_with_contains(p1[0], p1[1])) # --> a) FALSE
print(check_point_in_polygon_with_contains(p2[0], p2[1])) # --> b) TRUE
print(check_point_in_polygon_with_contains(p3[0], p3[1])) # --> c) FALSE

# testing INTERSECTS
print(check_point_in_polygon_with_intersects(p1[0], p1[1])) # --> a) TRUE
print(check_point_in_polygon_with_intersects(p2[0], p2[1])) # --> b) TRUE
print(check_point_in_polygon_with_intersects(p3[0], p3[1])) # --> c) FALSE