Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Proj.Transformation should error if it cannot create a transformation #114

Open
tiemvanderdeure opened this issue Feb 15, 2025 · 4 comments · Fixed by #115
Open

Proj.Transformation should error if it cannot create a transformation #114

tiemvanderdeure opened this issue Feb 15, 2025 · 4 comments · Fixed by #115

Comments

@tiemvanderdeure
Copy link
Contributor

This issue came up in Rasters where Proj segfaults when attempting to create a CRS transformation: rafaqz/Rasters.jl#895

MWE:

using Proj, GeoFormatTypes
crs_string = "LOCAL_CS[\"unnamed\",UNIT[\"metre\",1,AUTHORITY[\"EPSG\",\"9001\"]],AXIS[\"Easting\",EAST],AXIS[\"Northing\",NORTH]]"
crs = WellKnownText(GeoFormatTypes.CRS(), crs_string)
transform = Proj.Transformation(crs, EPSG(4326)) # works, but show errors
transform.pj === C_NULL # true
Proj.Transformation(crs, EPSG(4326); always_xy = true) # segfaults!!

This is a CRS that is valid but cannot be transformed (according to @evetion), Proj.Transformation returns a null pointer, and if always_xy is set to true, normalize_axis_order! is called on it, which causes the segfault.

Probably Proj should just error instead of returning a null string, similar to ArchGDAL:

using ArchGDAL
gdalcrs = ArchGDAL.importCRS(crs)
target = ArchGDAL.importEPSG(4326)
ArchGDAL.createcoordtrans(gdalcrs, target) do transform end
ERROR: GDALError (CE_Failure, code 6):
        Cannot find coordinate operations from `{
  "$schema": "https://proj.org/schemas/v0.7/projjson.schema.json",
  "type": "EngineeringCRS",
  "name": "unnamed",
  "datum": {
    "name": "Unknown engineering datum"
  },
  "coordinate_system": {
    "subtype": "Cartesian",
    "axis": [
      {
        "name": "Easting",
        "abbreviation": "",
        "direction": "east",
        "unit": "metre"
      },
      {
        "name": "Northing",
        "abbreviation": "",
        "direction": "north",
        "unit": "metre"
      }
    ]
  }
}' to `EPSG:4326'

Stacktrace:
 [1] maybe_throw()
   @ GDAL C:\Users\tsh371\.julia\packages\GDAL\8oAvY\src\error.jl:42
 [2] aftercare
   @ C:\Users\tsh371\.julia\packages\GDAL\8oAvY\src\error.jl:59 [inlined]
 [3] octnewcoordinatetransformation
   @ C:\Users\tsh371\.julia\packages\GDAL\8oAvY\src\libgdal.jl:38029 [inlined]
 [4] unsafe_createcoordtrans(source::ArchGDAL.ISpatialRef, target::ArchGDAL.ISpatialRef)
   @ ArchGDAL C:\Users\tsh371\.julia\packages\ArchGDAL\ujstt\src\spatialref.jl:812
 [5] createcoordtrans(::var"#11#12", ::ArchGDAL.ISpatialRef, ::Vararg{ArchGDAL.ISpatialRef}; kwargs::@Kwargs{})
   @ ArchGDAL C:\Users\tsh371\.julia\packages\ArchGDAL\ujstt\src\context.jl:266
 [6] createcoordtrans(::Function, ::ArchGDAL.ISpatialRef, ::Vararg{ArchGDAL.ISpatialRef})
   @ ArchGDAL C:\Users\tsh371\.julia\packages\ArchGDAL\ujstt\src\context.jl:265
 [7] top-level scope
@asinghvi17
Copy link
Member

We would also need an error hinter in Rasters telling users to use Planar.

But the broader question still remains ... how can we detect if a crs is not geographically linked so we avoid erroring in the first place?

@evetion
Copy link
Member

evetion commented Feb 15, 2025

Implement crstrait ofcourse

@asinghvi17
Copy link
Member

We could do a basic parser from ProjJSON I guess...at least checking ellipsoid parameters and the nature of the CRS

@asinghvi17
Copy link
Member

Reopening because I have a thought about what the error should be - will make a PR this week

@asinghvi17 asinghvi17 reopened this Feb 18, 2025
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging a pull request may close this issue.

3 participants