Como representar un punto en la superficie terrestre: Parte III - Transformación de Coordenadas con Cambio de Datum

En los dos meses anteriores vimos la diferencia entre una Coordenada Geográfica  y una Coordenada UTM, y hemos hablado sobre la importancia del Datum. Ahora hablaremos de cómo transformar una Coordenada Geográfica a una Coordenada UTM y viceversa.

Primero veamos un flujo de la transformación dicha:

Dónde se ve los pasos para transformar de una Coordenada UTM con un Datum Inicial hacia una Coordenada UTM con un Datum Final, y como paso previo se ve la transformación de la  Coordenada UTM Inicial hacia una Coordenada Geográfica Inicial, luego a una Coordenada Geocéntrica y luego a una Coordenada Geográfica Final hacia la Coordenada UTM Final.

Previamente, debemos conocer algunas “constantes” o símbolos usados en estas fórmulas:

FE: es el Falso Este, por lo general tiene un valor de 500000

FN: es el Falso Norte, para el Hemisferio Norte tiene un valor de 0, para el Hemisferio Sur tiene un valor de 10000000

ϕ : Latitud de la Coordenada Geográfica, por lo general se usa en su notación decimal (si tenemos una Latitud de 72º30’00’’, debemos escribirlo a 72.5)

λ  : Longitud de la Coordenada Geográfica, por lo general se usa en su notación decimal (si tenemos una Longitud de 12º30’00’’, debemos escribirlo a 12.5)

a  : Semieje mayor del Elipsoide de Revolución

b  : Semieje menor del Elipsoide de Revolución

1 / f: inversa del achatamiento del Elipsoide, donde: 1/ f= a / (a – b) o sino,  f= (a – b)/a

e² : primera excentricidad del Elipsoide, donde: e² = 2.f – f²

e'²: segunda excentricidad del Elipsoide, donde: e'² = e²  /(1 – f)²

ko: factor de Escala

λo : Longitud de Origen

ϕo: Latitud de Origen

Para hacer estas transformaciones, recurriremos a algunas fórmulas descritas en este documento:

http://sites.google.com/site/recursoswebgis/file-cabinet/guid7.pdf

y cuando es indicado, se usaran las fórmulas de este otro documento:

http://www.hydrometronics.com/downloads/Ellipsoidal%20Orthographic%20Projection.pdf  , página 10

Paso 1: Fórmulas para convertir de a Coordenadas UTM a  Coordenadas Geográficas:

 

Fuente: http://sites.google.com/site/recursoswebgis/file-cabinet/guid7.pdf, páginas 29 y 30

Ejemplo de uso:

Vamos a transformar la Coordenada UTM:

Este (E): 523,456.78

Norte (N): 8’123,456.78

Datum: PSAD56 - Mean

Zona UTM: 18S

 La transformaremos a una Coordenada Geográfica. Luego usando las formulas descritas arriba seria entonces:

Parámetros del Elipsoide: Para el Datum PSAD56 – Mean, este usa el Elipsoide “Internacional 1924”, luego sus parámetros son:

Semieje Menor “a”: 6’378,388

Inversa del Achatamiento “1/f”: 297

e²  = 2.f – f² , reemplazando:  e² = 2(1/297) – (1/297)² = 0.00672267002233332

e'² = e²  /(1 – f)² , reemplazando: e'² = 0.00672267002233332/(1 – (1/297))² = 0.00676817019722425

ko: 0.9996

FE: 500,000

FN: 10’000,000

λo : Por estar en la Zona UTM 18S, le corresponde un λo = -75°, que en radianes son:  -1.30899693899575

ϕo : 0°, que en radianes es: 0

 e₁  = [1 –  (1 –  0.00672267002233332)^0.5]/[1 + (1 –  0.00672267002233332)^0.5] = 0.00168634064080944

Mo = 6’378,388[(1 –   0.00672267002233332/4 –  3(0.00672267002233332)²/64 –  5(0.00672267002233332)³/256)0  –  (3(0.00672267002233332)/8 + 3(0.00672267002233332)²/32 + 45(0.00672267002233332)³/1024)sin2(0) + (15(0.00672267002233332)²/256 + 45(0.00672267002233332)³/1024)sin4(0)  –  (35(0.00672267002233332)³/3072)sin6(0)]

Mo = 0

M₁  = 0  + (8’123,456.78 –  10’000,000)/ 0.9996 = -1877294.13765506

µ₁  = -1877294.13765506/[6’378,388 (1 –0.00672267002233332/4 –  3(0.00672267002233332)²/64 –  5(0.00672267002233332)³/256)]
µ₁  = -0.294817210575242

ϕ₁  = -0.294817210575242 + (3(0.00168634064080944)/2 –  27(0.00168634064080944)³/32)sin2(-0.294817210575242) + (21(0.00168634064080944)²/16 – 55(0.00168634064080944)⁴/32)sin4(-0.294817210575242) + (151(0.00168634064080944)³/96)sin6(-0.294817210575242) + (1097(0.00168634064080944)⁴/512)sin8(-0.294817210575242)
ϕ₁  = -0.296227218504302

 T₁  = tan² (-0.296227218504302) = 0.0931514678527443

 ν₁  = 6’378,388/(1 –   (0.00672267002233332)sin²(-0.296227218504302))^0.5 = 6380215.75817981

ρ₁  = 6’378,388(1 –   0.00672267002233332)/(1 –  (0.00672267002233332)sin²(-0.296227218504302))^1.5

ρ₁ = 6340956.17534071

C₁  = (0.00676817019722425)cos²(-0.296227218504302) = 0.00619142945535154

D = (523,456.78 –  500,000)/((6380215.75817981)(0.9996))  = 0.00367795813746607

Luego la Latitud (ϕ) será:

ϕ = -0.296227218504302  –  ((6380215.75817981)tan(-0.296227218504302)/6340956.17534071)[( 0.00367795813746607)²/2 –  (5 + 3(0.0931514678527443) + 10(0.00619142945535154) –  4(0.00619142945535154)²  –  9(0.00676817019722425))(0.00367795813746607)⁴/24 + (61 + 90(0.0931514678527443) + 298(0.00619142945535154) + 45(0.0931514678527443)²  –  252(0.00676817019722425)  –  3(0.00619142945535154)²)(0.00367795813746607)⁶/720]

ϕ = -0.296225141408883, está en radianes, luego lo pasamos en grados sexagesimales y tenemos:

ϕ = -16.97245039  = -16° 58' 20.8213982''

 

Luego la Longitud (λ) será:

λ = -1.30899693899575 + [0.00367795813746607 –  (1 + 2(0.0931514678527443) + 0.00619142945535154)( 0.00367795813746607)³/6 + (5 –  2(0.00619142945535154) + 28(0.0931514678527443)  –  3(0.00619142945535154)² + 8(0.00676817019722425) + 24(0.0931514678527443)²)( 0.00367795813746607)⁵/120]/cos(-0.296227218504302)

λ = -1.30515150124643, está en radianes, luego lo pasamos en grados sexagesimales y tenemos:

λ = -74.7796726465838 = -74° 46' 46.8215277''

 

Paso 2: Fórmulas para convertir de a Coordenadas Geográficas a  Coordenadas Geocéntricas:

 

 

Fuente: http://sites.google.com/site/recursoswebgis/file-cabinet/guid7.pdf, página 54

En las formulas se ve la variable “h”, que es la altura elipsoidal, en esta presente nota, esta tendrá un valor de cero (0)

Ejemplo de uso:

Vamos a transformar la Coordenada Geográfica:

λ  =  -74° 46' 46.8215277'' = -1.30515150124643 (en radianes)

ϕ =  -16° 58' 20.8213982'' = -0.296225141408883 (en radianes)

 

Datum: PSAD56 - Mean

 La transformaremos a una Coordenada Geocéntrica. Luego usando las formulas descritas arriba seria entonces:

Parámetros del Elipsoide: Para el Datum PSAD56 – Mean, este usa el Elipsoide “Internacional 1924”, luego sus parámetros son:

Semieje Menor “a”: 6’378,388

Inversa del Achatamiento “1/f”: 297

e²  = 2.f – f² , reemplazando:  e² = 2.(1/297) – (1/297)² = 0.00672267002233332

e'² = e²  /(1 – f)² , reemplazando: e'² = 0.00672267002233332  /(1 – (1/297))² = 0.00676817019722425

 ν = 6’378,388/(1 –  (0.00672267002233332)sin²(-0.296225141408883))^0.5 = 6380215.73329148

Luego los valores X, Y, Z de la Coordenada Geocéntrica son:

X = (6380215.73329148 + 0)cos(-0.296225141408883)cos(-1.30515150124643) = 1602053.21335958

Y = (6380215.73329148 + 0)cos(-0.296225141408883)sin(-1.30515150124643) = -5888278.09383977

Z = ((1 – 0.00672267002233332)(6380215.73329148) + 0)sin(-0.296225141408883)  = -1849939.86741314

 

Paso 3: Fórmulas para convertir de a Coordenadas Geocéntricas a  Coordenadas Geográficas:

 

 

En este paso se hace el cambio de Datum. Para nuestro caso, nuestro Datum Inicial es el PSAD56 – Mean y el Datum Final es el WGS84.

 

Ejemplo de uso:

Vamos a transformar la Coordenada Geocéntrica que se encuentran en el Datum Inicial PSAD56 – Mean:

X  = 1602053.21335958

Y = -5888278.09383977

Z  = -1849939.86741314

 

La transformaremos a una Coordenada Geográfica en el Datum Final WGS84. Luego usando las formulas descritas arriba seria entonces:

Parámetros del Elipsoide: Para el Datum WGS84, este usa el Elipsoide “WGS 84”, luego sus parámetros son:

Semieje Menor “a”: 6’378,137

Inversa del Achatamiento “1/f”: 298.2572236

Semieje Menor “b”:  b = a – a.f = 6’378,137 - 6’378,137(1/298.2572236) = 6356752.31

e²  = 2.f – f² , reemplazando:  e² = 2.(1/298.2572236) – (1/298.2572236)² = 0.00669437999020854

e'² = e²  /(1 – f)² , reemplazando: e'² = 0.00669437999020854/(1 – (1/298.2572236))² =0.00673949674234457

ϕ =  -16° 58' 20.8213982'' = -0.296225141408883 (en radianes)

 

Los “shift” de transformación del Datum PSAD56 – Mean son:

Shift X_PSAD56: -288

Shift Y_PSAD56: 175

Shift Z_PSAD56: -376

 

Los “shift” de transformación del Datum WGS84 son:

Shift X _WGS84: 0

Shift Y _WGS84: 0

Shift Z _WGS84: 0

 

Luego, la transformación seria:

X_WGS84  = 1602053.21335958 + (-288) – 0 = 1601765.21335958

Y _WGS84 = -5888278.09383977 + 175 – 0 = -5888103.09383977

Z _WGS84  = -1849939.86741314 + (-376) – 0 = -1850315.86741314

 

P = (1601765.21335958² + (-5888103.09383977)²)^1/2 = 6102082.41852028

θ = atan(((-1850315.86741314)(6378137))/((6102082.41852028)(6356752.31))) = -0.295348589994358

Luego la Latitud en el Datum WGS84 (ϕ) será:

ϕ = atan ((-1850315.86741314 + 0.00673949674234457(6356752.31)sin³(-0.295348589994358))/( 6102082.41852028 – 0.00669437999020854(6378137)cos³(-0.295348589994358)))

ϕ = -0.296285086067313, está en radianes, luego lo pasamos en grados sexagesimales y tenemos:

ϕ = -16.9758849643274 = -16° 58' 33.1858716''

 

Luego la Longitud en el Datum WGS84 (λ) será:

λ = atan(-5888103.09383977/-5888103.09383977)  = -1.3051895136674, está en radianes, luego lo pasamos en grados sexagesimales y tenemos:

λ = -74.7818505978744 = -74° 46' 54.6621523

ν = 6’378,137 /(1 –  (0.00669437999020854)sin²(-0.296285086067313))^0.5  = 6379957.68227534

La altura sobre el Elipsoide es:

h = (6102082.41852028/cos(-0.296285086067313)) - 6379957.68227534 = 119.19

 

Paso 4: Fórmulas para convertir de a Coordenadas Geográficas  a  Coordenadas UTM:

Fuente: http://sites.google.com/site/recursoswebgis/file-cabinet/guid7.pdf, página 29

Ejemplo de uso:

Vamos a transformar la Coordenada Geográfica:

Latitud (ϕ) = -16° 58' 33.1858716'' = -16.9758849643274 = -0.296285086067313 (radianes)

Longitud (λ) = -74° 46' 54.6621523 = -1.3051895136674 = -1.3051895136674 (radianes)

Datum: WGS84

La transformaremos a una Coordenada UTM en la Zona UTM 18S. Luego usando las formulas descritas arriba seria entonces:

Parámetros del Elipsoide: Para el Datum WGS84, este usa el Elipsoide “WGS 84”, luego sus parámetros son:

Semieje Menor “a”: 6’378,137

Inversa del Achatamiento “1/f”: 298.2572236

Semieje Menor “b”:  b = a – a.f = 6’378,137 - 6’378,137(1/298.2572236) = 6356752.31

e²  = 2.f – f² , reemplazando:  e² = 2.(1/298.2572236) – (1/298.2572236)² = 0.00669437999020854

e'² = e²  /(1 – f)² , reemplazando: e'² = 0.00669437999020854/(1 – (1/298.2572236))² =0.00673949674234457

ko: 0.9996

FE: 500000

FN: 10000000

λo : Por estar en la Zona UTM 18S, le corresponde un λo = -75°, que en radianes son:  -1.30899693899575

ϕo : 0°, que en radianes es: 0

A = (-1.3051895136674  –  (-1.30899693899575)).cos(-0.296285086067313) = 0.00364152715114343

T = tan²(-0.296285086067313) = 0.093190086144421

ν = 6378137/(1 –  (0.00669437999020854)sin²(-0.296285086067313))^0.5 = 6379957.68227558

C = (0.00669437999020854).cos²(-0.296285086067313)/(1 –  (0.00669437999020854)) = 0.00616498157801096

M = (6378137)[(1 –  (0.00669437999020854)/4 –  3(0.00669437999020854)²/64 –  5(0.00669437999020854)³/256)(-0.296285086067313)  –  (3(0.00669437999020854)/8 + 3(0.00669437999020854)²/32 + 45(0.00669437999020854)³/1024)sin2(-0.296285086067313)+ (15(0.00669437999020854)²/256 + 45(0.00669437999020854)³/1024)sin4(-0.296285086067313)  –  (35(0.00669437999020854)³/3072).sin6(-0.296285086067313)]

M = -1877638.36860933

Mo = (6378137)[(1 –  (0.00669437999020854)/4 –  3(0.00669437999020854)²/64 –  5(0.00669437999020854)³/256)(0)  –  (3(0.00669437999020854)/8 + 3(0.00669437999020854)²/32 + 45(0.00669437999020854)³/1024)sin2(0)+ (15(0.00669437999020854)²/256 + 45(0.00669437999020854)³/1024)sin4(0)  –  (35(0.00669437999020854)³/3072)sin6(0)]

Mo = 0

Luego, las Coordenadas UTM serian:

Este (E)  =  500000 + (0.9996)(6379957.68227558)[0.00364152715114343 + (1 –  0.093190086144421 + 0.00616498157801096)(0.00364152715114343)³/6 + (5 –  18(0.093190086144421) + (0.093190086144421) ²  + 72(0.00616498157801096) –  58(0.00673949674234457))(0.00364152715114343)⁵/120]

Este (E)   = 523223.543

Norte (N) = 10000000 + (0.9996){-1877638.36860933 –  0  + (6379957.68227558)tan(-0.296285086067313)[( 0.00364152715114343)²/2 + (5 –  0.093190086144421 + 9(0.00616498157801096) + 4(0.00616498157801096)²)(0.00364152715114343)⁴/24 +(61 –  58(0.093190086144421) + (0.093190086144421)²  + 600(0.00616498157801096) – 330(0.00673949674234457))(0.00364152715114343)⁶/720]}

Norte (N) = 8123099.778

Estas son las mismas formulas usadas en la Hoja de Cálculo que se encuentra en este post:

http://www.ingeografos.com.pe/2012/10/transformacion-de-coordenadas-utm-tm.html

Esperando que sea de ayuda estos tres post sobre las Coordenadas Geográficas, Coordenadas UTM, Zonas UTM, Datum y Transformaciones de Coordenadas, me despido hasta el otro mes. Gracias y saludos