Convertire una coordinata geografica dal sistema UTM al sistema EPSG:4326 WGS 84 e viceversa
Da Gambas-it.org - Wikipedia.
Mostriamo un possibile codice per convertire una coordinata UTM - GPS (WGS84) UTM - nelle corrispondenti coordinate geografiche in formato decimale del sistema EPSG:4326 WGS 84:
' Procedura derivata da un'elaborazione di © 2005 Jonathan Stott
Private Const RifEllMAGGIORE As Float = 6378137
Private Const RifEllMINORE As Float = 6356752.314
Private Const UTM_F0 As Float = 0.9996
Private wgs84_ecc As Float = ((RifEllMAGGIORE ^ 2) - (RifEllMINORE ^ 2)) / (RifEllMAGGIORE ^ 2)
Public Sub Main()
' Imposta la coordinata Est, la coordinata Nord e il numero del Fuso del punto:
Dim coord_est As Float = 291954.0
Dim coord_nord As Float = 4640624.0
Dim zonaNumero As Byte = 33
' Il numero della lettera della Fascia servirà soltanto per la stampa finale:
Dim zonaLettera As String = "T"
Dim latlon As Float[]
latlon = UTMRefToLatLng(coord_est, coord_nord, zonaNumero)
' Fuso Fascia x-Est y-Nord
Print "\nUTM: "; zonaNumero; zonaLettera;; coord_est; " E "; coord_nord; " N"
Print "\nLatitudine: \e[31m"; Format(latlon[0], "0.000000"), "\e[0m Longitudine: \e[34m"; Format(latlon[1], "0.000000")
End
' Converte una coordinata da UTM nelle corrispondenti coordinate geografiche della Latitudine e della Longitude
Private Function UTMRefToLatLng(est As Float, nord As Float, zonaNumero As Byte) As Float[]
Dim ePrimeSquared As Float = wgs84_ecc / (1.0 - wgs84_ecc)
Dim e1 As Float = (1 - Sqr(1 - wgs84_ecc)) / (1 + Sqr(1 - wgs84_ecc))
Dim x As Float = est - 500000.0
Dim longitudeOrigin As Float = (zonaNumero - 1.0) * 6.0 - 180.0 + 3.0
Dim m As Float = nord / UTM_F0
Dim mu As Float = m / (RifEllMAGGIORE * (1.0 - wgs84_ecc / 4.0 - 3.0 * wgs84_ecc * wgs84_ecc / 64.0 - 5.0 * (wgs84_ecc ^ 3) / 256.0))
Dim phi1Rad As Float = mu + (3.0 * e1 / 2.0 - 27.0 * (e1 ^ 3) / 32.0) * Sin(2.0 * mu) +
(21.0 * e1 * e1 / 16.0 - 55.0 * (e1 ^ 4) / 32.0) * Sin(4.0 * mu) +
(151.0 * (e1 ^ 3) / 96.0) * Sin(6.0 * mu)
Dim n As Float = RifEllMAGGIORE / Sqr(1.0 - wgs84_ecc * Sin(phi1Rad) * Sin(phi1Rad))
Dim t As Float = Tan(phi1Rad) * Tan(phi1Rad)
Dim c As Float = ePrimeSquared * Cos(phi1Rad) * Cos(phi1Rad)
Dim r As Float = RifEllMAGGIORE * (1.0 - wgs84_ecc) / (CFloat(1.0 - wgs84_ecc * Sin(phi1Rad) * Sin(phi1Rad)) ^ 1.5)
Dim d As Float = x / (n * UTM_F0)
Dim latitudine As Float = (phi1Rad - (n * Tan(phi1Rad) / r) * (d * d / 2.0 - (5.0 + (3.0 * t) + (10.0 * c) -
(4.0 * c * c) - (9.0 * ePrimeSquared)) * (d ^ 4) / 24.0 +
(61.0 + (90.0 * t) + (298.0 * c) + (45.0 * t * t) - (252.0 * ePrimeSquared) -
(3.0 * c * c)) * (d ^ 6) / 720.0)) * (180.0 / Pi)
Dim longitudine As Float = longitudeOrigin + ((d - (1.0 + 2.0 * t + c) * (d ^ 3) / 6.0 + (5.0 - (2.0 * c) +
(28.0 * t) - (3.0 * c * c) + (8.0 * ePrimeSquared) + (24.0 * t * t)) *
(d ^ 5) / 120.0) / Cos(phi1Rad)) * (180.0 / Pi)
Return [latitudine, longitudine]
End
Convertire una coordinata geografica dal formato decimale del sistema EPSG:4326 WGS 84 al sistema UTM
Mostriamo un possibile codice per convertire una coordinata geografica (Latitudine/Longitudine) dal formato decimale del sistema EPSG:4326 WGS 84 al sistema UTM [GPS (WGS84) UTM].
Il codice è tarato per i Fusi e le Fasce UTM che comprendono il territorio della Repubblica Italiana.
' Procedura derivata da un'elaborazione di © 2005 Jonathan Stott
Private Const RifEllMAGGIORE As Float = 6378137
Private Const RifEllMINORE As Float = 6356752.314
Private Const UTM_F0 As Float = 0.9996
Private wgs84_ecc As Float = ((RifEllMAGGIORE ^ 2) - (RifEllMINORE ^ 2)) / (RifEllMAGGIORE ^ 2)
Public Sub Main()
' Imposta la Latitudine e la Longitudine del punto:
Dim latitudine As Float = 41.89018
Dim longitudine As Float = 12.49230
Dim ii As Integer[]
ii = LatLngToUTMRef(latitudine, longitudine)
Print "Longitudine: "; Longitudine, "Latitudine: "; Latitudine
' Fuso Fascia x-Est y-Nord
Print "\nUTM: \e[31m"; ii[0]; Chr(ii[1]);; ii[2]; " E "; ii[3]; " N\e[0m"
End
' Converte una coordinata da valori latitude/longitude in un corrispondente riferimento UTM
Private Function LatLngToUTMRef(lat As Float, lon As Float) As Integer[]
Dim latitudeRad As Float = lat * (Pi / 180)
Dim longitudeRad As Float = lon * (Pi / 180.0)
Dim longitudeZone As Byte = Floor((lon + 180.0) / 6.0) + 1
Dim longitudeOrigin As Float = (longitudeZone - 1) * 6 - 180 + 3
Dim longitudeOriginRad As Float = longitudeOrigin * (Pi / 180.0)
Dim ePrimeSquared As Float = wgs84_ecc / (1 - wgs84_ecc)
Dim n As Float = RifEllMAGGIORE / Sqr(1 - wgs84_ecc * Sin(latitudeRad) * Sin(latitudeRad))
Dim t As Float = Tan(latitudeRad) * Tan(latitudeRad)
Dim c As Float = ePrimeSquared * Cos(latitudeRad) * Cos(latitudeRad)
Dim A As Float = Cos(latitudeRad) * (longitudeRad - longitudeOriginRad)
Dim M As Float = RifEllMAGGIORE * ((1 - wgs84_ecc / 4 - 3 * wgs84_ecc * wgs84_ecc / 64 - 5 *
wgs84_ecc * wgs84_ecc * wgs84_ecc / 256) *
latitudeRad - (3 * wgs84_ecc / 8 + 3 * wgs84_ecc * wgs84_ecc / 32 + 45 *
wgs84_ecc * wgs84_ecc * wgs84_ecc / 1024) * Sin(2 * latitudeRad) +
(15 * wgs84_ecc * wgs84_ecc / 256 + 45 * wgs84_ecc * wgs84_ecc * wgs84_ecc / 1024) * Sin(4 *
latitudeRad) - (35 * wgs84_ecc * wgs84_ecc * wgs84_ecc / 3072) * Sin(6 * latitudeRad))
Dim UTMEasting As Float = UTM_F0 * n * (A + (1 - t + c) * (A ^ 3) / 6 + (5 - 18 * t + t * t + 72 * c - 58 * ePrimeSquared) *
(A ^ 5) / 120) + 500000.0
Dim UTMNorthing As Float = UTM_F0 * (M + n * Tan(latitudeRad) * (A * A / 2 + (5 - t + (9 * c) + (4 * c * c)) *
(A ^ 4) / 24 + (61 - (58 * t) + (t * t) + (600 * c) - (330 * ePrimeSquared)) * (A ^ 6) / 720))
Dim UTMzone As Byte
If (48 > lat) And (lat >= 40) Then UTMzone = Asc("T")
If (40 > lat) And (lat >= 32) Then UTMzone = Asc("S")
Return [longitudeZone, UTMzone, Round(UTMEasting), Round(UTMNorthing)]
End
