# Test validity of CRS for a single coordinate in python

Does anyone know if there is a simple way to check if a CRS (given by EPSG code) is valid for a single coordinate (lon, lat)?

• What do you mean by "check validity"?
– Erik
Aug 24, 2022 at 14:41
• Sorry wasn't very clear of me. I mean check if the coordinate lies inside or outside the valid lon/lat ranges for the CRS. Aug 24, 2022 at 15:14
• Look up the geographic bounding box and do a point in polygon test. Aug 24, 2022 at 16:10

### Using `pyproj` and `shapely`

``````from shapely.geometry import box, Point
from pyproj import CRS

def within_crs_bounds(epsg_code, lon, lat):
""" Returns `True` if point (`lon`, `lat`) is within the
bounds of the given projection (`epsg_code`).
Otherwise returns `False`. """

crs = CRS.from_user_input(epsg_code)

bounding_box = box(*crs.area_of_use.bounds)

return Point(lon, lat).within(bounding_box)
``````

### Using PyQGIS

``````from qgis.core import QgsCoordinateReferenceSystem, QgsPointXY

def within_crs_bounds(epsg_code, lon, lat):
""" Returns `True` if point (`lon`, `lat`) is within the
bounds of the given projection (`epsg_code`).
Otherwise returns `False`. """

spatial_reference = QgsCoordinateReferenceSystem.fromEpsgId(epsg_code)

bounding_box = spatial_reference.bounds()

return bounding_box.contains(QgsPointXY(lon, lat))
``````

### Using GDAL/OGR

``````from osgeo import osr, ogr

def within_crs_bounds(epsg_code, lon, lat):
""" Returns `True` if point (`lon`, `lat`) is within the
bounds of the given projection (`epsg_code`).
Otherwise returns `False`. """

srs = osr.SpatialReference()
srs.ImportFromEPSG(epsg_code)

area = srs.GetAreaOfUse()

point = ogr.Geometry(ogr.wkbPoint)

min_x = area.west_lon_degree
min_y = area.south_lat_degree
max_x = area.east_lon_degree
max_y = area.north_lat_degree

ring = ogr.Geometry(ogr.wkbLinearRing)

box = ogr.Geometry(ogr.wkbPolygon)

return point.Within(box)
``````

Since creating geometries in GDAL leads to very verbose code, I would probably also check the spatial relation myself:

``````from osgeo import osr

def within_crs_bounds(epsg_code, lon, lat):
""" Returns `True` if point (`lon`, `lat`) is within the
bounds of the given projection (`epsg_code`).
Otherwise returns `False`. """

srs = osr.SpatialReference()
srs.ImportFromEPSG(epsg_code)

area = srs.GetAreaOfUse()

return (
area. west_lon_degree <= lon <= area. east_lon_degree
and area.south_lat_degree <= lat <= area.north_lat_degree
)
``````

### Test Code

``````epsg_code = 2056  # Projection: CH1903+ / LV95

# Aarhus (Denmark)
assert not within_crs_bounds(epsg_code, lon=10.2107, lat=56.1572)

# Lausanne (Switzerland)
assert within_crs_bounds(epsg_code, lon=6.6335, lat=46.519833)
``````
• Amazing! Thank you for such a thorough answer Thomas, I've learnt a lot from it :) Aug 25, 2022 at 12:48

In the end I wrote the following function which seems to do the job:

``````def verify_crs_bounds(lon, lat, epsg_code):
crs = pyproj.CRS.from_user_input(epsg_code)
if (crs.area_of_use.west <= lon <= crs.area_of_use.east) and (
crs.area_of_use.south <= lat <= crs.area_of_use.north
):
result = True
else:
print(f"EPSG:{epsg_code} is invalid as the coordinate given is not within the permitted area of use.")
result = False
return result
``````