sábado, 24 de agosto de 2013

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 – f2

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

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:

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


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

 
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

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

e'² = e2  /(1 – f)2 , reemplazando: e'² = 0.00672267002233332/(1 – (1/297))2 = 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

 e1  = [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)2/64 –  5(0.00672267002233332)3/256)0    (3(0.00672267002233332)/8 + 3(0.00672267002233332)2/32 + 45(0.00672267002233332)3/1024)sin2(0) + (15(0.00672267002233332)2/256 + 45(0.00672267002233332)3/1024)sin4(0)    (35(0.00672267002233332)3/3072)sin6(0)]
Mo = 0

M1  = 0  + (8’123,456.78  10’000,000)/ 0.9996 = -1877294.13765506
µ1  = -1877294.13765506/[6’378,388 (1 –0.00672267002233332/4 –  3(0.00672267002233332)2/64 –  5(0.00672267002233332)3/256)]
µ1  = -0.294817210575242

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

 T1  = tan2 (-0.296227218504302) = 0.0931514678527443

 ν1  = 6’378,388/(1 –   (0.00672267002233332)sin2(-0.296227218504302))0.5 = 6380215.75817981

ρ1  = 6’378,388(1 –   0.00672267002233332)/(1 –  (0.00672267002233332)sin2(-0.296227218504302))1.5
ρ1 = 6340956.17534071

C1  = (0.00676817019722425)cos2(-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/2 –  (5 + 3(0.0931514678527443) + 10(0.00619142945535154) –  4(0.00619142945535154)2    9(0.00676817019722425))(0.00367795813746607)4/24 + (61 + 90(0.0931514678527443) + 298(0.00619142945535154) + 45(0.0931514678527443)2    252(0.00676817019722425)    3(0.00619142945535154)2)(0.00367795813746607)6/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)3/6 + (5 –  2(0.00619142945535154) + 28(0.0931514678527443)    3(0.00619142945535154)2 + 8(0.00676817019722425) + 24(0.0931514678527443)2)( 0.00367795813746607)5/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:

 

 


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

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

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

 ν = 6’378,388/(1 –  (0.00672267002233332)sin2(-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

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

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

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

 

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

Shift XPSAD56: -288

Shift YPSAD56: 175

Shift ZPSAD56: -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:

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

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

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

 
P = (1601765.213359582 + (-5888103.09383977)2)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)sin3(-0.295348589994358))/( 6102082.418520280.00669437999020854(6378137)cos3(-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)sin2(-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:


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

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

e'² = e2  /(1 – f)2 , reemplazando: e'² = 0.00669437999020854/(1 – (1/298.2572236))2 =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 = tan2(-0.296285086067313) = 0.093190086144421

ν = 6378137/(1 –  (0.00669437999020854)sin2(-0.296285086067313))0.5 = 6379957.68227558

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

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

Mo = (6378137)[(1 –  (0.00669437999020854)/4 –  3(0.00669437999020854)2/64 –  5(0.00669437999020854)3/256)(0)    (3(0.00669437999020854)/8 + 3(0.00669437999020854)2/32 + 45(0.00669437999020854)3/1024)sin2(0)+ (15(0.00669437999020854)2/256 + 45(0.00669437999020854)3/1024)sin4(0)    (35(0.00669437999020854)3/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)3/6 + (5 –  18(0.093190086144421) + (0.093190086144421) 2  + 72(0.00616498157801096) –  58(0.00673949674234457))(0.00364152715114343)5/120]
Este (E)   = 523223.543

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

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


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
ª