From f15059d9869fb141a8b6e4c1428597761152a639 Mon Sep 17 00:00:00 2001 From: Elias Carvalho Date: Thu, 7 Nov 2024 16:41:15 -0300 Subject: [PATCH 1/3] Add utility functions to easily extract values from metadata --- src/GeoTIFF.jl | 1 + src/userutils.jl | 57 ++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 58 insertions(+) create mode 100644 src/userutils.jl diff --git a/src/GeoTIFF.jl b/src/GeoTIFF.jl index 813b97e..756eac9 100644 --- a/src/GeoTIFF.jl +++ b/src/GeoTIFF.jl @@ -13,6 +13,7 @@ import TiffImages include("geokeys.jl") include("metadata.jl") +include("userutils.jl") include("image.jl") include("load.jl") include("save.jl") diff --git a/src/userutils.jl b/src/userutils.jl new file mode 100644 index 0000000..68f8b99 --- /dev/null +++ b/src/userutils.jl @@ -0,0 +1,57 @@ +# ----------------------------------------------------------------- +# Licensed under the MIT License. See LICENSE in the project root. +# ----------------------------------------------------------------- + +function geokey(metadata::Metadata, id::GeoKeyID) + geokeys = metadata.geokeydirectory.geokeys + i = findfirst(geokey -> geokey.id == id, geokeys) + isnothing(i) ? nothing : geokeys[i] +end + +function geokeyvalue(metadata::Metadata, id::GeoKeyID) + geokey = geokey(metadata, id) + isnothing(geokey) ? nothing : geokey.value +end + +rastertype(metadata::Metadata) = geokeyvalue(metadata, GTRasterTypeGeoKey) + +modeltype(metadata::Metadata) = geokeyvalue(metadata, GTModelTypeGeoKey) + +function epsgcode(metadata::Metadata) + mt = modeltype(metadata) + isnothing(mt) && return nothing + if mt == Projected2D + geokeyvalue(metadata, ProjectedCRSGeoKey) + elseif mt == Geographic2D || mt == Geocentric3D + geokeyvalue(metadata, GeodeticCRSGeoKey) + else + # Undefined or UserDefined + nothing + end +end + +function affineparams(metadata::Metadata) + pixelscale = metadata.modelpixelscale + tiepoint = metadata.modeltiepoint + transformation = metadata.modeltransformation + if !isnothing(pixelscale) && !isnothing(tiepoint) + sx = pixelscale.x + sy = pixelscale.y + sz = pixelscale.z + (; i, j, k, x, y, z) = tiepoint + tx = x - i / sx + ty = y + j / sy + tz = z - k / sz + A = [ + sx 0 0 + 0 -sy 0 + 0 0 sz + ] + b = [tx, ty, tz] + A, b + elseif !isnothing(transformation) + transformation.A, transformation.b + else + nothing + end +end From 3b8ec9374b9460394da334265bc580fb5ba6b5c1 Mon Sep 17 00:00:00 2001 From: Elias Carvalho Date: Thu, 7 Nov 2024 17:12:21 -0300 Subject: [PATCH 2/3] Add 'affineparams2D' --- src/userutils.jl | 25 ++++++++++++++++++++++++- 1 file changed, 24 insertions(+), 1 deletion(-) diff --git a/src/userutils.jl b/src/userutils.jl index 68f8b99..e98c6d4 100644 --- a/src/userutils.jl +++ b/src/userutils.jl @@ -30,7 +30,30 @@ function epsgcode(metadata::Metadata) end end -function affineparams(metadata::Metadata) +function affineparams2D(metadata::Metadata) + pixelscale = metadata.modelpixelscale + tiepoint = metadata.modeltiepoint + transformation = metadata.modeltransformation + if !isnothing(pixelscale) && !isnothing(tiepoint) + sx = pixelscale.x + sy = pixelscale.y + (; i, j, x, y) = tiepoint + tx = x - i / sx + ty = y + j / sy + A = [ + sx 0 + 0 -sy + ] + b = [tx, ty] + A, b + elseif !isnothing(transformation) + transformation.A[1:2, 1:2], transformation.b[1:2] + else + nothing + end +end + +function affineparams3D(metadata::Metadata) pixelscale = metadata.modelpixelscale tiepoint = metadata.modeltiepoint transformation = metadata.modeltransformation From bb4643d646230a18433094260d1f5df63afe01f3 Mon Sep 17 00:00:00 2001 From: Elias Carvalho Date: Thu, 7 Nov 2024 17:18:50 -0300 Subject: [PATCH 3/3] Fix code --- src/metadata.jl | 8 ++++---- src/userutils.jl | 4 ++-- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/src/metadata.jl b/src/metadata.jl index 94dfb91..a52e23e 100644 --- a/src/metadata.jl +++ b/src/metadata.jl @@ -31,10 +31,10 @@ GeoKeyDirectory(; version=1, revision=1, minor=1, geokeys=GeoKey[]) = function GeoKeyDirectory(params::Vector{UInt16}) nkeys = params[4] geokeys = map((1:nkeys) * 4) do i - id = GeoKeyID(geokey[1 + i]) - tag = geokey[2 + i] - count = geokey[3 + i] - value = geokey[4 + i] + id = GeoKeyID(params[1 + i]) + tag = params[2 + i] + count = params[3 + i] + value = params[4 + i] GeoKey(id, tag, count, value) end diff --git a/src/userutils.jl b/src/userutils.jl index e98c6d4..b4b6a31 100644 --- a/src/userutils.jl +++ b/src/userutils.jl @@ -9,8 +9,8 @@ function geokey(metadata::Metadata, id::GeoKeyID) end function geokeyvalue(metadata::Metadata, id::GeoKeyID) - geokey = geokey(metadata, id) - isnothing(geokey) ? nothing : geokey.value + gk = geokey(metadata, id) + isnothing(gk) ? nothing : gk.value end rastertype(metadata::Metadata) = geokeyvalue(metadata, GTRasterTypeGeoKey)