Distance from Cartopy feature polygon

cartopy_feature_distance.py

import cartopy.feature as cfeature
import pyproj

from shapely.geometry import Point
from shapely.ops import transform as shapely_transform

# Some location of interest. Also used as center of projection.
wgs84_pt = Point(-72.2495, 43.886)

# Set up projections
wgs84 = pyproj.Proj(proj='latlon')
aeqd = pyproj.Proj(proj='aeqd', lon_0=wgs84_pt.x, lat_0=wgs84_pt.y)
project = pyproj.Transformer.from_proj(wgs84, aeqd, always_xy=True).transform
aeqd_point = shapely_transform(project, wgs84_pt)

# Loop through shape components, searching only for polygons,
# and generate Shapely polygons in map projection
def gen_transformed_polys(shape):
    for geom in shape.geometries():
        typestr = str(type(geom))
        if 'Polygon' in typestr:
            if 'Multi' in typestr:
                for poly in geom:
                    yield shapely_transform(project, poly)
            else:
                yield shapely_transform(project, geom)

# Test with the Natural Earth land polygons
land = cfeature.NaturalEarthFeature('physical', 'land', '50m',
                                        edgecolor='face',
                                        facecolor=cfeature.COLORS['land'])

land_distances_km = [int(one_shape.exterior.distance(aeqd_point)/1000.0) 
                for one_shape in gen_transformed_polys(land)]

print(min(land_distances_km), max(land_distances_km))

This is based on https://gis.stackexchange.com/questions/127427/transforming-shapely-polygon-and-multipolygon-objects