search

pyproj4 / pyproj

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

Proj4 export from a bound CRS object fails #1404

Closed phaarnes closed 2 months ago

phaarnes commented 3 months ago

Code Sample

from pyproj import Transformer, CRS
from pyproj.crs import BoundCRS

source_crs = CRS.from_user_input("EPSG:4267")
target_crs = CRS.from_user_input("EPSG:4326")
ct = Transformer.from_pipeline("EPSG:15851")

bound_crs = BoundCRS(
    source_crs = source_crs,
    target_crs = target_crs,
    transformation = ct
)

bound_crs.to_proj4()

Problem description

The example code above is failing. Pyproj throws an CRSError

  File ~\AppData\Local\miniforge3\lib\site-packages\pyproj\crs\crs.py:1295 in to_proj4
    raise CRSError("CRS cannot be converted to a PROJ string.")

CRSError: CRS cannot be converted to a PROJ string.

I suspect this issue is related to the parameters key of the transformer object being a list of length two. However, this is just a theory.

The description dictionary for this transformation object (EPSG:15851): 

{'$schema': 'https://proj.org/schemas/v0.7/projjson.schema.json',
 'type': 'Transformation',
 'name': 'NAD27 to WGS 84 (79)',
 'source_crs': {'type': 'GeographicCRS',
  'name': 'NAD27',
  'datum': {'type': 'GeodeticReferenceFrame',
   'name': 'North American Datum 1927',
   'ellipsoid': {'name': 'Clarke 1866',
    'semi_major_axis': 6378206.4,
    'semi_minor_axis': 6356583.8}},
  'coordinate_system': {'subtype': 'ellipsoidal',
   'axis': [{'name': 'Geodetic latitude',
     'abbreviation': 'Lat',
     'direction': 'north',
     'unit': 'degree'},
    {'name': 'Geodetic longitude',
     'abbreviation': 'Lon',
     'direction': 'east',
     'unit': 'degree'}]},
  'id': {'authority': 'EPSG', 'code': 4267}},
 'target_crs': {'type': 'GeographicCRS',
  'name': 'WGS 84',
  'datum_ensemble': {'name': 'World Geodetic System 1984 ensemble',
   'members': [{'name': 'World Geodetic System 1984 (Transit)',
     'id': {'authority': 'EPSG', 'code': 1166}},
    {'name': 'World Geodetic System 1984 (G730)',
     'id': {'authority': 'EPSG', 'code': 1152}},
    {'name': 'World Geodetic System 1984 (G873)',
     'id': {'authority': 'EPSG', 'code': 1153}},
    {'name': 'World Geodetic System 1984 (G1150)',
     'id': {'authority': 'EPSG', 'code': 1154}},
    {'name': 'World Geodetic System 1984 (G1674)',
     'id': {'authority': 'EPSG', 'code': 1155}},
    {'name': 'World Geodetic System 1984 (G1762)',
     'id': {'authority': 'EPSG', 'code': 1156}},
    {'name': 'World Geodetic System 1984 (G2139)',
     'id': {'authority': 'EPSG', 'code': 1309}}],
   'ellipsoid': {'name': 'WGS 84',
    'semi_major_axis': 6378137,
    'inverse_flattening': 298.257223563},
   'accuracy': '2.0',
   'id': {'authority': 'EPSG', 'code': 6326}},
  'coordinate_system': {'subtype': 'ellipsoidal',
   'axis': [{'name': 'Geodetic latitude',
     'abbreviation': 'Lat',
     'direction': 'north',
     'unit': 'degree'},
    {'name': 'Geodetic longitude',
     'abbreviation': 'Lon',
     'direction': 'east',
     'unit': 'degree'}]},
  'id': {'authority': 'EPSG', 'code': 4326}},
 'method': {'name': 'NADCON', 'id': {'authority': 'EPSG', 'code': 9613}},
 'parameters': [{'name': 'Latitude difference file',
   'value': 'conus.las',
   'id': {'authority': 'EPSG', 'code': 8657}},
  {'name': 'Longitude difference file',
   'value': 'conus.los',
   'id': {'authority': 'EPSG', 'code': 8658}}],
 'accuracy': '5.0',
 'scope': '(null/copy) Approximation for medium and low accuracy applications assuming equality between plate-fixed static and earth-fixed dynamic CRSs, ignoring static/dynamic CRS differences.',
 'area': 'United States (USA) - CONUS including EEZ - onshore and offshore - Alabama; Arizona; Arkansas; California; Colorado; Connecticut; Delaware; Florida; Georgia; Idaho; Illinois; Indiana; Iowa; Kansas; Kentucky; Louisiana; Maine; Maryland; Massachusetts; Michigan; Minnesota; Mississippi; Missouri; Montana; Nebraska; Nevada; New Hampshire; New Jersey; New Mexico; New York; North Carolina; North Dakota; Ohio; Oklahoma; Oregon; Pennsylvania; Rhode Island; South Carolina; South Dakota; Tennessee; Texas; Utah; Vermont; Virginia; Washington; West Virginia; Wisconsin; Wyoming. US Gulf of Mexico (GoM) OCS.',
 'bbox': {'south_latitude': 23.81,
  'west_longitude': -129.17,
  'north_latitude': 49.38,
  'east_longitude': -65.69},
 'id': {'authority': 'EPSG', 'code': 15851},
 'remarks': 'Parameter files are from NAD27 to NAD83 (1) (code 1241) assuming that NAD83 is equivalent to WGS 84 within the accuracy of this tfm. Uses NADCON method which expects longitudes positive west; EPSG CRS codes 4267 and 4326 have longitudes positive east.'}

Expected Output

The correct Proj4 string for this bound CRS would be: '+proj=longlat +datum=NAD27 +nadgrids=us_noaa_conus.tif +no_defs +type=crs'

Environment Information

snowman2 commented 3 months ago 
from pyproj.crs import BoundCRS
crs = BoundCRS(source_crs="EPSG:4627", target_crs="EPSG:4326", transformation="EPSG:15851")
projinfo 'BOUNDCRS[SOURCECRS[GEOGCRS["RGR92",DATUM["Reseau Geodesique de la Reunion 1992",ELLIPSOID["GRS 1980",6378137,298.257222101,LENGTHUNIT["metre",1]]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433]],CS[ellipsoidal,2],AXIS["geodetic latitude (Lat)",north,ORDER[1],ANGLEUNIT["degree",0.0174532925199433]],AXIS["geodetic longitude (Lon)",east,ORDER[2],ANGLEUNIT["degree",0.0174532925199433]],USAGE[SCOPE["Horizontal component of 3D system."],AREA["Reunion - onshore and offshore."],BBOX[-24.72,51.83,-18.28,58.24]],ID["EPSG",4627]]],TARGETCRS[GEOGCRS["WGS 84",ENSEMBLE["World Geodetic System 1984 ensemble",MEMBER["World Geodetic System 1984 (Transit)"],MEMBER["World Geodetic System 1984 (G730)"],MEMBER["World Geodetic System 1984 (G873)"],MEMBER["World Geodetic System 1984 (G1150)"],MEMBER["World Geodetic System 1984 (G1674)"],MEMBER["World Geodetic System 1984 (G1762)"],MEMBER["World Geodetic System 1984 (G2139)"],ELLIPSOID["WGS 84",6378137,298.257223563,LENGTHUNIT["metre",1]],ENSEMBLEACCURACY[2.0]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433]],CS[ellipsoidal,2],AXIS["geodetic latitude (Lat)",north,ORDER[1],ANGLEUNIT["degree",0.0174532925199433]],AXIS["geodetic longitude (Lon)",east,ORDER[2],ANGLEUNIT["degree",0.0174532925199433]],USAGE[SCOPE["Horizontal component of 3D system."],AREA["World."],BBOX[-90,-180,90,180]],ID["EPSG",4326]]],ABRIDGEDTRANSFORMATION["NAD27 to WGS 84 (79)",METHOD["NADCON",ID["EPSG",9613]],PARAMETERFILE["Latitude difference file","conus.las"],PARAMETERFILE["Longitude difference file","conus.los"],USAGE[SCOPE["(null/copy) Approximation for medium and low accuracy applications assuming equality between plate-fixed static and earth-fixed dynamic CRSs, ignoring static/dynamic CRS differences."],AREA["United States (USA) - CONUS including EEZ - onshore and offshore - Alabama; Arizona; Arkansas; California; Colorado; Connecticut; Delaware; Florida; Georgia; Idaho; Illinois; Indiana; Iowa; Kansas; Kentucky; Louisiana; Maine; Maryland; Massachusetts; Michigan; Minnesota; Mississippi; Missouri; Montana; Nebraska; Nevada; New Hampshire; New Jersey; New Mexico; New York; North Carolina; North Dakota; Ohio; Oklahoma; Oregon; Pennsylvania; Rhode Island; South Carolina; South Dakota; Tennessee; Texas; Utah; Vermont; Virginia; Washington; West Virginia; Wisconsin; Wyoming. US Gulf of Mexico (GoM) OCS."],BBOX[23.81,-129.17,49.38,-65.69]],ID["EPSG",15851],REMARK["Parameter files are from NAD27 to NAD83 (1) (code 1241) assuming that NAD83 is equivalent to WGS 84 within the accuracy of this tfm. Uses NADCON method which expects longitudes positive west; EPSG CRS codes 4267 and 4326 have longitudes positive east."]]]' -o PROJ
PROJ.4 string:
Error when exporting to PROJ string: Transformation cannot be formatted as WKT1 TOWGS84 parameters
snowman2 commented 3 months ago

I recommend raising an issue in the PROJ repository as that is where the core logic is managed: https://github.com/OSGeo/PROJ/

phaarnes commented 3 months ago

Thank you for your quick respoonse @snowman2. I have now raised an issue in the Proj repository- Proj4 export from a bound CRS object fails

snowman2 commented 2 months ago

Fix coming in PROJ 9.5: https://github.com/OSGeo/PROJ/pull/4168