PyProj fails to find CRS

ldesousa commented 2 years ago

Code Sample, a copy-pastable example if possible

import psycopg2
import geopandas as gpd

con = psycopg2.connect(database="data", user="postgres", password="postgres", host="localhost")
sql = "SELECT *, public.ST_Transform(wkb_geometry, 54052) as geom FROM public.waypoint"
points = gpd.read_postgis(sql, con, geom_col='geom')

Problem description

My system is Ubuntu 22.04. The use case is to read into a data frame a set of points from around the globe and plot it to map. The points are in a PostGIS database. Proj fails with an exception when reading the points to a GeoPandas data frame. It claims the CRS does not exist, but it does.

Expected Output

A GeoPandas data frame.

Environment Information

$ pyproj -v
pyproj info:
    pyproj: 3.3.1
      PROJ: 8.2.0
  data dir: /home/duque004/git/wrb-class/env/lib/python3.10/site-packages/pyproj/proj_dir/share/proj
user_data_dir: /home/duque004/.local/share/proj
PROJ DATA (recommended version): 1.8
PROJ Database: 1.2
EPSG Database: v10.038 [2021-10-21]
ESRI Database: ArcMap 12.8 [2021-05-06]
IGNF Database: 3.1.0 [2019-05-24]

System:
    python: 3.10.4 (main, Jun 29 2022, 12:14:53) [GCC 11.2.0]
executable: /home/duque004/git/wrb-class/env/bin/python3
   machine: Linux-5.15.0-46-generic-x86_64-with-glibc2.35

Python deps:
   certifi: 2022.06.15
    Cython: None
setuptools: 59.6.0
       pip: 22.0.2

Installation method

pip

Full log

This is the full environment:

$ pip freeze
attrs==22.1.0
certifi==2022.6.15
click==8.1.3
click-plugins==1.1.1
cligj==0.7.2
Fiona==1.8.21
geopandas==0.11.1
munch==2.5.0
numpy==1.23.2
packaging==21.3
pandas==1.4.4
psycopg2==2.9.3
pyparsing==3.0.9
pyproj==3.3.1
python-dateutil==2.8.2
pytz==2022.2.1
Shapely==1.8.4
six==1.16.0

The log from the interactive command line with the exception:

$ python3
Python 3.10.4 (main, Jun 29 2022, 12:14:53) [GCC 11.2.0] on linux
Type "help", "copyright", "credits" or "license" for more information.
>>> import psycopg2
>>> import geopandas as gpd
>>>
>>> con = psycopg2.connect(database="data", user="postgres", password="postgres",
...     host="localhost")
>>>
>>> sql = "SELECT *, public.ST_Transform(wkb_geometry, 54052) as geom FROM public.waypoint"
>>> points = gpd.read_postgis(sql, con, geom_col='geom')
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/home/duque004/git/wrb-class/env/lib/python3.10/site-packages/geopandas/io/sql.py", line 163, in _read_postgis
    return _df_to_geodf(df, geom_col=geom_col, crs=crs)
  File "/home/duque004/git/wrb-class/env/lib/python3.10/site-packages/geopandas/io/sql.py", line 93, in _df_to_geodf
    return GeoDataFrame(df, crs=crs, geometry=geom_col)
  File "/home/duque004/git/wrb-class/env/lib/python3.10/site-packages/geopandas/geodataframe.py", line 189, in __init__
    self.set_geometry(geometry, inplace=True, crs=crs)
  File "/home/duque004/git/wrb-class/env/lib/python3.10/site-packages/geopandas/geodataframe.py", line 347, in set_geometry
    level = _ensure_geometry(level, crs=crs)
  File "/home/duque004/git/wrb-class/env/lib/python3.10/site-packages/geopandas/geodataframe.py", line 59, in _ensure_geometry
    out = from_shapely(np.asarray(data), crs=crs)
  File "/home/duque004/git/wrb-class/env/lib/python3.10/site-packages/geopandas/array.py", line 150, in from_shapely
    return GeometryArray(vectorized.from_shapely(data), crs=crs)
  File "/home/duque004/git/wrb-class/env/lib/python3.10/site-packages/geopandas/array.py", line 285, in __init__
    self.crs = crs
  File "/home/duque004/git/wrb-class/env/lib/python3.10/site-packages/geopandas/array.py", line 335, in crs
    self._crs = None if not value else CRS.from_user_input(value)
  File "/home/duque004/git/wrb-class/env/lib/python3.10/site-packages/pyproj/crs/crs.py", line 488, in from_user_input
    return cls(value, **kwargs)
  File "/home/duque004/git/wrb-class/env/lib/python3.10/site-packages/pyproj/crs/crs.py", line 335, in __init__
    self._local.crs = _CRS(self.srs)
  File "pyproj/_crs.pyx", line 2352, in pyproj._crs._CRS.__init__
pyproj.exceptions.CRSError: Invalid projection: epsg:54052: (Internal Proj Error: proj_create: crs not found)

First check, the CRS definition in PostGIS:

$ psql -d data
psql (14.5 (Ubuntu 14.5-0ubuntu0.22.04.1))
Type "help" for help.

data=# select * from public.spatial_ref_sys where srid = 54052;

 srid  | auth_name | auth_srid |                                                                                                                                                                                                                                srtext                                                                                                                                                                                                                                |                             proj4text
-------+-----------+-----------+----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------+-------------------------------------------------------------------
 54052 | ESRI      |     54052 | PROJCS["World_Goode_Homolosine_Land",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["Degree",0.0174532925199433]],PROJECTION["Goode_Homolosine"],PARAMETER["central_meridian",0],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["ESRI","54052"]] | +proj=goode +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs
(1 row)

Then check the same in the Proj database:

$ sqlite3 ./lib/python3.10/site-packages/pyproj/proj_dir/share/proj/proj.db
SQLite version 3.37.2 2022-01-06 13:25:41
Enter ".help" for usage hints.
sqlite> select * from projected_crs where code = 54052;
ESRI|54052|World_Goode_Homolosine_Land||||EPSG|4326|||PROJCS["World_Goode_Homolosine_Land",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Goode_Homolosine"],PARAMETER["False_Easting",0.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",0.0],PARAMETER["Option",1.0],UNIT["Meter",1.0]]|0

Looks like PyProj (or GeoPandas) is inadvertently changing the CRS authority, but this might not explain the exception. The srid field in PostGIS and code in the Proj database are both unique (not enforced, but the codes do not repeat). The programme should be able to fetch the correct CRS definition with the code alone.

Thank you for reading.

snowman2 commented 2 years ago

This is something that geopandas needs to address. If you use CRS.from_authority and pass in the auth_name of ESRI and the auth_code 54052, it should work.

snowman2 commented 2 years ago

Closing as this is a geopandas issue. Please report it here: https://github.com/geopandas/geopandas/

ldesousa commented 2 years ago

For reference, this bug is now reported as issue #2545 at GeoPandas.

pyproj4 / pyproj