Added latitude band for the UTMZ point type

README.md

@@ -52,9 +52,9 @@ point_enu = ENU(point_lla, origin_lla, wgs84)
 
 Similarly, we could convert to UTM/UPS coordinates, and two types are provided
 for this - `UTM` stores 3D coordinates `x`, `y`, and `z` in an unspecified zone,
-while `UTMZ` includes the `zone` number and `hemisphere` bool (where `true` =
-northern, `false` = southern). To get the canonical zone for your coordinates,
-simply use:
+while `UTMZ` includes the `zone` number, `hemisphere` bool (where `true` =
+northern, `false` = southern) and `latletter` char with the latitude band letter. 
+To get the canonical zone for your coordinates, simply use:
 ```julia
 x_utmz = UTMZ(x_lla, wgs84)
 ```
@@ -317,11 +317,12 @@ universal polar-stereographic (UPS) coordinates (where the zone is `0`).
 ##### `UTMZ{T}` - universal transverse-Mercator + zone
 
 In addition to the easting `x`, northing `y` and height `z`, the global `UTMZ` type
-also encodes the UTM `zone` and `hemisphere`, where `zone` is a `UInt8` and
-`hemisphere` is a `Bool` for compact storage. The northern hemisphere is
-denoted as `true`, and the southern as `false`. Zone `0` corresponds to the UPS
+also encodes the UTM `zone` , `hemisphere` and `latletter`, where `zone` is a `UInt8`,
+`hemisphere` is a `Bool` for compact storage and `latletter` a char. The northern hemisphere 
+is denoted as `true`, and the southern as `false`. Zone `0` corresponds to the UPS
 projection about the corresponding pole, otherwise `zone` is an integer between
-`1` and `60`.
+`1` and `60`. The latitude band letter`latletter` char with the latitude band is a char
+that can go from `C` to `X`.
 
 ##### `ENU{T}` - east-north-up
 

src/Geodesy.jl

@@ -49,7 +49,7 @@ export
     ITRF, GDA94,
 
     # UTM helpers
-    utm_zone
+    utm_zone, utm_lat_letter
 
 include("ellipsoids.jl")
 include("datums.jl")

src/conversion.jl

@@ -86,5 +86,5 @@ UTM(enu::ENU, zone::Integer, hemisphere::Bool, origin, datum) = UTMfromENU(origi
 ### UTMZ <-> UTM coordinates ###
 ###############################
 
-UTMZ(utm::UTM, zone::Integer, hemisphere::Bool, datum) = UTMZfromUTM(zone, hemisphere, datum)(utm)
-UTM(utmz::UTMZ, zone::Integer, hemisphere::Bool, datum) = UTMfromUTMZ(zone, hemisphere, datum)(utmz)
+UTMZ(utm::UTM, zone::Integer, hemisphere::Bool, latletter::Char, datum) = UTMZfromUTM(zone, hemisphere, latletter, datum)(utm)
+UTM(utmz::UTMZ, zone::Integer, hemisphere::Bool, latletter::Char, datum) = UTMfromUTMZ(zone, hemisphere, latletter, datum)(utmz)

src/points.jl

@@ -118,10 +118,11 @@ struct UTMZ{T <: Number}
     z::T
     zone::UInt8
     isnorth::Bool # which hemisphere
+    latletter::Char # which latitude band
 end
-UTMZ(x::T, y::T, zone::Integer, isnorth::Bool) where {T} = UTMZ(x, y, zero(T), UInt8(zone), isnorth)
-UTMZ(x, y, z, zone::Integer, isnorth::Bool) = UTMZ(promote(x, y, z)..., UInt8(zone), isnorth)
-UTMZ(utm::UTM, zone::Integer, isnorth::Bool) = UTMZ(utm.x, utm.y, utm.z, UInt8(zone), isnorth)
+UTMZ(x::T, y::T, zone::Integer, isnorth::Bool, latletter::Char) where {T} = UTMZ(x, y, zero(T), UInt8(zone), isnorth, latletter)
+UTMZ(x, y, z, zone::Integer, isnorth::Bool, latletter::Char) = UTMZ(promote(x, y, z)..., UInt8(zone), isnorth, latletter)
+UTMZ(utm::UTM, zone::Integer, isnorth::Bool, latletter::Char) = UTMZ(utm.x, utm.y, utm.z, UInt8(zone), isnorth, latletter)
 UTM(utmz::UTMZ) = UTM(utmz.x, utmz.y, utmz.z)
-Base.isapprox(utm1::UTMZ, utm2::UTMZ; atol = 1e-6, kwargs...) = (utm1.zone == utm2.zone) & (utm1.isnorth == utm2.isnorth) & isapprox(utm1.x, utm2.x; atol = atol, kwargs...) & isapprox(utm1.y, utm2.y; atol = atol, kwargs...) & isapprox(utm1.z, utm2.z; atol = atol, kwargs...) # atol in metres (1μm)
-Base.show(io::IO, utm::UTMZ) = print(io, "UTMZ($(utm.x), $(utm.y), $(utm.z), zone=$(utm.zone == 0 ? "polar" : Int(utm.zone)) ($(utm.isnorth ? "north" : "south")))")
+Base.isapprox(utm1::UTMZ, utm2::UTMZ; atol = 1e-6, kwargs...) = (utm1.zone == utm2.zone) & (utm1.isnorth == utm2.isnorth) & (utm1.latletter == utm2.latletter) & isapprox(utm1.x, utm2.x; atol = atol, kwargs...) & isapprox(utm1.y, utm2.y; atol = atol, kwargs...) & isapprox(utm1.z, utm2.z; atol = atol, kwargs...) # atol in metres (1μm)
+Base.show(io::IO, utm::UTMZ) = print(io, "UTMZ($(utm.x), $(utm.y), $(utm.z), zone=$(utm.zone == 0 ? "polar" : Int(utm.zone)) $(utm.latletter) ($(utm.isnorth ? "north" : "south")))") 

src/transformations.jl

@@ -476,14 +476,15 @@ Base.show(io::IO, trans::UTMZfromLLA) = print(io, "UTMZfromLLA($(trans.datum))")
 
 function (trans::UTMZfromLLA)(lla::LLA)
     (zone, isnorth) = utm_zone(lla)
+    latletter = utm_lat_letter(lla.lat)
     if zone == 0
         # Do polar steriographic projection
         k0 = 0.994
         (x,y,γ,k) = polarst_fwd(isnorth, k0, trans.tm, lla.lat, lla.lon)
         x = x + 2e6
         y = y + 2e6
 
-        return UTMZ(x, y, lla.alt, 0, isnorth)
+        return UTMZ(x, y, lla.alt, 0, isnorth, latletter)
     else
         lon_ref = Float64(utm_meridian(zone))
         k0 = 0.9996 # Horizontal scaling factor
@@ -492,7 +493,7 @@ function (trans::UTMZfromLLA)(lla::LLA)
         x = 5e5 + x # Convention has 500km offset for easting
         y = (isnorth ? 0.0 : 1e7) + y # Northing offset for southern hemisphere
         # also, k = k * k0
-        return UTMZ(x, y, lla.alt, zone, isnorth) # Note: scaling not applied to vertical dimension
+        return UTMZ(x, y, lla.alt, zone, isnorth, latletter) # Note: scaling not applied to vertical dimension
     end
 end
 
@@ -511,9 +512,10 @@ Transformation to append the UTM/UPS zone and hemisphere to a raw `UTM` data poi
 struct UTMZfromUTM{Datum} <: Transformation
     zone::Int
     isnorth::Bool
+    latletter::Char
     datum::Datum
 end
-UTMZfromUTM(zone::Integer, h, d::D) where {D} = UTMZfromUTM{D}(UInt8(zone), h, d)
+UTMZfromUTM(zone::Integer, h, l, d::D) where {D} = UTMZfromUTM{D}(UInt8(zone), h, l, d)
 Base.show(io::IO, trans::UTMZfromUTM) = print(io, "UTMZfromUTM(zone=$(trans.zone == 0 ? "polar" : trans.zone) ($(trans.isnorth ? "north" : "south")), $(trans.datum))")
 #Base.show(io::IO, trans::UTMZfromUTM) = print(io, "UTMZfromUTM(zone=$(trans.zone == 0 ? "polar" : trans.zone) ($(trans.isnorth ? "north" : "south")))")
 
@@ -526,15 +528,16 @@ automatically convert the data to the specified zone if necessary.
 struct UTMfromUTMZ{Datum} <: Transformation
     zone::Int
     isnorth::Bool
+    latletter::Char
     datum::Datum
 end
-UTMfromUTMZ(zone::Integer, h, d::D) where {D} = UTMfromUTMZ{D}(UInt8(zone), h, d)
+UTMfromUTMZ(zone::Integer, h, l, d::D) where {D} = UTMfromUTMZ{D}(UInt8(zone), h, l, d)
 Base.show(io::IO, trans::UTMfromUTMZ) = print(io, "UTMfromUTMZ(zone=$(trans.zone == 0 ? "polar" : trans.zone) ($(trans.isnorth ? "north" : "south")), $(trans.datum))")
 #Base.show(io::IO, trans::UTMfromUTMZ) = print(io, "UTMfromUTMZ(zone=$(trans.zone == 0 ? "polar" : trans.zone) ($(trans.isnorth ? "north" : "south")))")
 
-(trans::UTMZfromUTM)(utm::UTM) = UTMZ(utm.x, utm.y, utm.z, trans.zone, trans.isnorth)
+(trans::UTMZfromUTM)(utm::UTM) = UTMZ(utm.x, utm.y, utm.z, trans.zone, trans.isnorth, trans.latletter)
 function (trans::UTMfromUTMZ)(utm::UTMZ)
-    if trans.zone == utm.zone && trans.isnorth == utm.isnorth
+    if trans.zone == utm.zone && trans.isnorth == utm.isnorth && trans.latletter == utm.latletter
         UTM(utm.x, utm.y, utm.z)
     else
         # Should this be an error or an automatic transformation to the correct zone?
@@ -543,8 +546,8 @@ function (trans::UTMfromUTMZ)(utm::UTMZ)
     end
 end
 
-Base.inv(trans::UTMfromUTMZ) = UTMZfromUTM(trans.zone, trans.isnorth, trans.datum)
-Base.inv(trans::UTMZfromUTM) = UTMfromUTMZ(trans.zone, trans.isnorth, trans.datum)
+Base.inv(trans::UTMfromUTMZ) = UTMZfromUTM(trans.zone, trans.isnorth, trans.latletter, trans.datum)
+Base.inv(trans::UTMZfromUTM) = UTMfromUTMZ(trans.zone, trans.isnorth, trans.latletter, trans.datum)
 
 ###################
 ## ECEF <-> UTMZ ##

src/utm.jl

@@ -54,3 +54,60 @@ function utm_meridian(zone::Integer)
     end
     return 6 * zone - 183
 end
+
+"""
+    utm_lat_letter(lat)
+
+Calculates the latitude bands from the latitude poistion.
+This is part of the military grid reference system (MGRS).
+A letter coming before "N" in the alphabet is in the southern 
+hemisphere and any letter "N" or after is in the northern hemisphere.
+"""
+function utm_lat_letter(lat::T) where T
+
+    lat_letter::Char = ' '
+    if lat <-72 
+        lat_letter='C'
+    elseif lat <-64
+        lat_letter='D'
+    elseif lat <-56
+        lat_letter='E'
+    elseif lat <-48
+        lat_letter='F'
+    elseif lat <-40
+        lat_letter='G'
+    elseif lat <-32
+        lat_letter='H'
+    elseif lat <-24
+        lat_letter='J'
+    elseif lat <-16
+        lat_letter='K'
+    elseif lat <-8
+        lat_letter='L'
+    elseif lat <0
+        lat_letter='M'
+    elseif lat <8
+        lat_letter='N'
+    elseif lat <16
+        lat_letter='P'
+    elseif lat <24
+        lat_letter='Q'
+    elseif lat <32
+        lat_letter='R'
+    elseif lat <40
+        lat_letter='S'
+    elseif lat <48
+        lat_letter='T'
+    elseif lat <56
+        lat_letter='U'
+    elseif lat <64
+        lat_letter='V'
+    elseif lat <72
+        lat_letter='W'
+    else
+        lat_letter='X'
+    end
+
+    return lat_letter
+
+end

test/conversion.jl

@@ -22,7 +22,7 @@
         @test ECEF(enu, lla_ref, wgs84) ≈ ecef
 
         # LLA <-> UTM
-        (z, h) = (19, true)
+        (z, h, l) = (19, true, 'T')
         utm = UTM(327412.48528248386, 4.692686244318043e6, 0.0)
         @test UTM(lla, z, h, wgs84) ≈ utm
         @test LLA(utm, z, h, wgs84) ≈ lla
@@ -36,7 +36,7 @@
         @test ENU(utm, z, h, lla_ref, wgs84) ≈ enu
 
         # LLA <-> UTMZ
-        utmz = UTMZ(327412.48528248386, 4.692686244318043e6, 0, z, h)
+        utmz = UTMZ(327412.48528248386, 4.692686244318043e6, 0, z, h, l)
         @test UTMZ(lla, wgs84) ≈ utmz
         @test LLA(utmz, wgs84) ≈ lla
 
@@ -49,8 +49,8 @@
         @test ENU(utmz, lla_ref, wgs84) ≈ enu
 
         # UTM <-> UTMZ
-        @test UTMZ(utm, z, h, wgs84) ≈ utmz
-        @test UTM(utmz, z, h, wgs84) ≈ utm
+        @test UTMZ(utm, z, h, l, wgs84) ≈ utmz
+        @test UTM(utmz, z, h, l, wgs84) ≈ utm
 
         # Identity
         @test ECEF(ecef, wgs84) == ecef

test/points.jl

@@ -19,9 +19,9 @@
     utm = UTM(10., 10., 0.)
     @test UTM(10., 10.) == utm
 
-    utmz = UTMZ(10., 10., 0., 1, true)
-    @test UTMZ(10., 10., 1, true) == utmz
-    @test UTMZ(utm, 1, true) == utmz
+    utmz = UTMZ(10., 10., 0., 1, true, 'D')
+    @test UTMZ(10., 10., 1, true, 'D') == utmz
+    @test UTMZ(utm, 1, true, 'D') == utmz
 
     @test UTM(utmz) == utm
 end # @testset

test/transformations.jl

@@ -30,19 +30,19 @@
     utmz = utmz_lla(lla)
     utmz2 = utmz_lla(lla2)
     # Test also the poles go to UPS
-    @test utmz_lla(LLA(89.0, 89.0, 89.0)) ≈ UTMZ(2.111009610242531e6, 1.9980623200455606e6, 89.0, 0, true)
-    @test lla_utmz(UTMZ(2.111009610242531e6, 1.9980623200455606e6, 89.0, 0, true)) ≈ LLA(89.0, 89.0, 89.0)
+    @test utmz_lla(LLA(89.0, 89.0, 89.0)) ≈ UTMZ(2.111009610242531e6, 1.9980623200455606e6, 89.0, 0, true, 'X')
+    @test lla_utmz(UTMZ(2.111009610242531e6, 1.9980623200455606e6, 89.0, 0, true, 'X')) ≈ LLA(89.0, 89.0, 89.0)
 
-    @test utmz ≈ UTMZ(239579.67583179142,7.342551466042985e6,1374.7804632852078, UInt8(47), false)
+    @test utmz ≈ UTMZ(239579.67583179142,7.342551466042985e6,1374.7804632852078, UInt8(47), false, 'J')
     @test lla ≈ lla_utmz(utmz)
 
-    utmz_utm = UTMZfromUTM(47, false, wgs84)
-    utm_utmz = UTMfromUTMZ(47, false, wgs84)
+    utmz_utm = UTMZfromUTM(47, false, 'J', wgs84)
+    utm_utmz = UTMfromUTMZ(47, false, 'J', wgs84)
 
     @test utmz == utmz_utm(utm)
     @test utm == utm_utmz(utmz)
     # Sometimes this has to do something non-trivial - make sure it does
-    utm_utmz2 = UTMfromUTMZ(46, false, wgs84)
+    utm_utmz2 = UTMfromUTMZ(46, false,'X', wgs84)
     @test utm_utmz2(utmz) ≈ UTM(850009.7418418773, 7.340641097689279e6, 1374.7804632852078)
 
     # ENU coordinates