diff --git a/src/coord.jl b/src/coord.jl index 1a8e406..6d8b9a2 100644 --- a/src/coord.jl +++ b/src/coord.jl @@ -88,12 +88,19 @@ end function Transformation( source_crs::GFT.CoordinateReferenceSystemFormat, target_crs::GFT.CoordinateReferenceSystemFormat; - always_xy::Bool=false, - direction::PJ_DIRECTION=PJ_FWD, - area::Ptr{PJ_AREA}=C_NULL, - ctx::Ptr{PJ_CONTEXT}=C_NULL + always_xy::Bool = false, + direction::PJ_DIRECTION = PJ_FWD, + area::Ptr{PJ_AREA} = C_NULL, + ctx::Ptr{PJ_CONTEXT} = C_NULL, ) - return Transformation(CRS(source_crs).pj, CRS(target_crs).pj; always_xy, direction, area, ctx) + return Transformation( + CRS(source_crs).pj, + CRS(target_crs).pj; + always_xy, + direction, + area, + ctx, + ) end function Transformation( diff --git a/src/crs.jl b/src/crs.jl index c0b7869..f916727 100644 --- a/src/crs.jl +++ b/src/crs.jl @@ -1,7 +1,7 @@ """ CRS(crs) -Create a CRS. `crs` can be: +Create a coordinate reference system. `crs` can be: - a proj-string, - a WKT string, - an object code (like "EPSG:4326", "urn:ogc:def:crs:EPSG::4326", "urn:ogc:def:coordinateOperation:EPSG::1671"), @@ -137,3 +137,39 @@ Base.convert(::Type{CRS}, crs::GFT.CoordinateReferenceSystemFormat) = CRS(crs) # Maybe enable later, based on https://github.com/JuliaGeo/GeoFormatTypes.jl/issues/21 # Base.convert(T::Type{<:GFT.CoordinateReferenceSystemFormat}, crs::GFT.CoordinateReferenceSystemFormat) = T(CRS(crs)) + +""" + identify(crs::CRS; auth_name = nothing) + +Returns a list of matching reference CRS and confidence values (0-100). + +# Arguments +- `crs::CRS`: Coordinate reference system +- `auth_name=nothing`: Authority name, or nothing for all authorities (e.g. "EPSG") +""" +function identify( + crs::CRS; + auth_name = nothing, +)::Vector{@NamedTuple{crs::CRS, confidence::Int32}} + + out_confidence = Ref(Ptr{Cint}(C_NULL)) + if isnothing(auth_name) + # set authority to C_NULL + auth_name = C_NULL + end + + pj_list = proj_identify(crs, auth_name, out_confidence) + list = NamedTuple{(:crs, :confidence),Tuple{CRS,Int32}}[] + + # was a match found? + if pj_list != C_NULL + n = proj_list_get_count(pj_list) + for i = 1:n + crs = CRS(proj_list_get(pj_list, i - 1)) + confidence = unsafe_load(out_confidence[], i) + push!(list, (; crs, confidence)) + end + proj_int_list_destroy(out_confidence[]) + end + return list +end diff --git a/test/libproj.jl b/test/libproj.jl index 41f6f65..5a4190c 100644 --- a/test/libproj.jl +++ b/test/libproj.jl @@ -202,7 +202,7 @@ end source_crs = GFT.EPSG("EPSG:4326") target_crs = GFT.EPSG("EPSG:32628") - trans = Proj.Transformation(source_crs, target_crs, always_xy=false) + trans = Proj.Transformation(source_crs, target_crs, always_xy = false) info = Proj.proj_pj_info(trans.pj) description1 = unsafe_string(info.definition) @@ -210,7 +210,7 @@ end source_crs = "EPSG:4326" target_crs = "EPSG:32628" - trans = Proj.Transformation(source_crs, target_crs, always_xy=false) + trans = Proj.Transformation(source_crs, target_crs, always_xy = false) info = Proj.proj_pj_info(trans.pj) description2 = unsafe_string(info.definition) @@ -477,3 +477,14 @@ end # Maybe enable later, based on https://github.com/JuliaGeo/GeoFormatTypes.jl/issues/21 # @test convert(GFT.ProjString, gftcrs) == GFT.ProjString("+proj=longlat +datum=WGS84 +no_defs +type=crs") end + +@testset "Proj.identify" begin + crs = Proj.CRS("+proj=longlat +datum=WGS84 +no_defs +type=crs") + identity = Proj.identify(crs) + @test identity[1].confidence == 70 + @test length(identity) == 4 + + identity = Proj.identify(crs; auth_name = "EPSG") + @test GFT.EPSG(identity[1].crs) == GFT.EPSG(4326) + @test identity[1].confidence == 70 +end