pyproj4 / pyproj

Python interface to PROJ (cartographic projections and coordinate transformations library)
https://pyproj4.github.io/pyproj
MIT License
1.04k stars 211 forks source link

Transformer creation broken with empty area of interest #1326

Closed tpwrules closed 2 weeks ago

tpwrules commented 1 year ago

Code Sample, a copy-pastable example if possible

from pyproj import CRS, Transformer
from pyproj.aoi import AreaOfInterest
from pyproj.database import query_utm_crs_info

def broken(lat, lon):
    aoi = AreaOfInterest(
        west_lon_degree=lon,
        south_lat_degree=lat,
        east_lon_degree=lon,
        north_lat_degree=lat,
    )
    utm_crs_list = query_utm_crs_info(
        datum_name="WGS 84",
        area_of_interest=aoi,
    )
    utm_crs = CRS.from_epsg(utm_crs_list[0].code)

    return Transformer.from_crs(CRS("WGS 84"), utm_crs, area_of_interest=aoi)

def working(lat, lon):
    aoi = AreaOfInterest(
        west_lon_degree=lon,
        south_lat_degree=lat,
        east_lon_degree=lon+1e-8,
        north_lat_degree=lat+1e-8,
    )
    utm_crs_list = query_utm_crs_info(
        datum_name="WGS 84",
        area_of_interest=aoi,
    )
    utm_crs = CRS.from_epsg(utm_crs_list[0].code)

    return Transformer.from_crs(CRS("WGS 84"), utm_crs, area_of_interest=aoi)

lat, lon = 42.032974, -93.581543

working(lat, lon)

# worked in 3.4.1, broken in 3.5.0 and 3.6.0
# throws pyproj.exceptions.ProjError: Error creating Transformer from CRS.
broken(lat, lon)

Problem description

An area of interest with 0 area (i.e. east/west and/or north/south coordinates the same) causes some internal exception when passed to Transformer.from_crs. Enlarging the area even negligibly results in the operation working properly.

A zero-size area should still technically be valid. Both variants do work without crashing in v3.4.1.

Expected Output

A 0 area of interest should work and produce a Transformer object.

Environment Information

pyproj info:
    pyproj: 3.5.0
      PROJ: 9.2.0
  data dir: /nix/store/l21lg0k87m6z743923ld2clp76yykskk-proj-9.2.0/share/proj
user_data_dir: /home/tpwatson/.local/share/proj
PROJ DATA (recommended version): 1.13
PROJ Database: 1.2
EPSG Database: v10.082 [2023-02-06]
ESRI Database: ArcGIS Pro 3.1 [2023-19-01]
IGNF Database: 3.1.0 [2019-05-24]

System:
    python: 3.10.12 (main, Jun  6 2023, 22:43:10) [GCC 12.3.0]
executable: /nix/store/46j792vz3np5f32kg36ifa2s2qsyppq7-python3-3.10.12-env/bin/python3.10
   machine: Linux-4.15.0-169-generic-x86_64-with-glibc2.37

Python deps:
   certifi: 2022.12.7
    Cython: None
setuptools: None
       pip: None

Installation method

nixpkgs built from source

snowman2 commented 2 weeks ago

It is likely due to PROJ as I can use pyproj 3.4.1 with PROJ 9.4.1 and it fails. However, if I use PROJ 9.1.0, it succeeds. As mentioned in the issue, once I change to PROJ 9.2+ it fails.

snowman2 commented 2 weeks ago

Command:

PROJ_DEBUG=3 cs2cs "WGS 84" "EPSG:32615" --bbox -93.581543,42.032974,-93.581543,42.032974
PROJ 9.1.1: It works
...
No operation found matching criteria
Rel. 9.4.1, June 1st, 2024
<cs2cs>: 
cannot initialize transformation
cause: (null)
program abnormally terminated
snowman2 commented 2 weeks ago

This command works for both versions:

PROJ_DEBUG=3 cs2cs "WGS 84" "EPSG:32615" --bbox -93.581543,42.032974,-93.581542999,42.032974001
snowman2 commented 2 weeks ago

Forwarded: https://github.com/OSGeo/PROJ/issues/4236