Search code examples
pythonshapefileesrimatplotlib-basemappyshp

Given latitude/longitude say if the coordinate is within continental US or not


I want to check if a particular latitude/longitude is within continental US or not. I don't want to use Online APIs and I'm using Python.

I downloaded this shapefile

from shapely.geometry import MultiPoint, Point, Polygon
import shapefile    
sf = shapefile.Reader("cb_2015_us_nation_20m")
shapes = sf.shapes()
fields = sf.fields
records = sf.records()
points = shapes[0].points
poly = Polygon(points)
lon = -112
lat = 48
point = Point(-112, 48)
poly.contains(point)
#should return True because it is in continental US but returns False

The sample lon, lat is within US boundary but poly.contains returns False. I'm not sure what the problem is and how to solve the issue so that I can test if a point is within continental US.


Solution

  • I ended up checking if lat/lon was in every state instead of check in continental U.S., if a point is in one of the states, then it is in continental U.S..

    from shapely.geometry import MultiPoint, Point, Polygon
    import shapefile
    #return a polygon for each state in a dictionary
    def get_us_border_polygon():
    
        sf = shapefile.Reader("./data/states/cb_2015_us_state_20m")
        shapes = sf.shapes()
        #shapes[i].points
        fields = sf.fields
        records = sf.records()
        state_polygons = {}
        for i, record in enumerate(records):
            state = record[5]
            points = shapes[i].points
            poly = Polygon(points)
            state_polygons[state] = poly
    
        return state_polygons
    
    #us border
    state_polygons = get_us_border_polygon()   
    #check if in one of the states then True, else False
    def in_us(lat, lon):
        p = Point(lon, lat)
        for state, poly in state_polygons.iteritems():
            if poly.contains(p):
                return state
        return None