mapping - para - coordenadas x y latitud longitud
Conversión de longitud / latitud a coordenadas cartesianas (7)
Tengo algunos puntos de coordenadas centrados en la Tierra como latitud y longitud ( WGS-84 ).
¿Cómo puedo convertirlos en coordenadas cartesianas (x, y, z) con el origen en el centro de la tierra?
¿Por qué implementar algo que ya se ha implementado y probado?
C #, por NetTopologySuite , tiene NetTopologySuite que es el puerto .NET de JTS Topology Suite.
Específicamente, usted tiene un defecto severo en su cálculo. La tierra no es una esfera perfecta, y la aproximación del radio de la tierra podría no cortarla para mediciones precisas.
Si en algunos casos es aceptable usar funciones homebrew, GIS es un buen ejemplo de un campo en el que se prefiere usar una biblioteca confiable, probada en pruebas.
Aquí está la respuesta que encontré:
Solo para completar la definición, en el sistema de coordenadas cartesianas:
- el eje x atraviesa largo, lat (0,0), de modo que la longitud 0 se encuentra con el ecuador;
- el eje y pasa (0,90);
- y el eje z pasa por los polos.
La conversión es:
x = R * cos(lat) * cos(lon)
y = R * cos(lat) * sin(lon)
z = R *sin(lat)
Donde R es el radio aproximado de la tierra (por ejemplo, 6371KM).
Si sus funciones trigonométricas esperan radianes (lo que probablemente hagan), primero deberá convertir su longitud y latitud a radianes. Obviamente, necesitas una representación decimal, no grados / minutos / segundos (ver, por ejemplo, aquí sobre la conversión).
La fórmula para la conversión de vuelta:
lat = asin(z / R)
lon = atan2(y, x)
asin es, por supuesto, arco senoidal. lee sobre atan2 en wikipedia . No te olvides de convertir de radianes a grados.
Esta página proporciona el código c # para esto (tenga en cuenta que es muy diferente de las fórmulas), y también algunas explicaciones y un buen diagrama de por qué esto es correcto,
El software proj.4 proporciona un programa de línea de comandos que puede hacer la conversión, por ejemplo
LAT=40
LON=-110
echo $LON $LAT | cs2cs +proj=latlong +datum=WGS84 +to +proj=geocent +datum=WGS84
También proporciona una API de C En particular, la función pj_geodetic_to_geocentric
hará la conversión sin tener que configurar primero un objeto de proyección.
Recientemente hice algo similar a esto utilizando la "Fórmula Haversine" en los datos WGS-84, que es un derivado de la "Ley de Haversines" con resultados muy satisfactorios.
Sí, WGS-84 supone que la Tierra es un elipsoide, pero creo que solo obtiene un error promedio de 0.5% usando un enfoque como la "Fórmula Haversine", que puede ser una cantidad aceptable de error en su caso. Siempre tendrá un poco de error a menos que esté hablando de una distancia de unos pocos pies e incluso entonces hay una curvatura teórica de la Tierra ... Si necesita un enfoque más rígido compatible con WGS-84 revise la "Fórmula Vincenty".
Entiendo de dónde viene starblue , pero la buena ingeniería de software a menudo se trata de compensaciones, por lo que todo depende de la precisión que requiera para lo que está haciendo. Por ejemplo, el resultado calculado a partir de la "Fórmula de distancia de Manhattan" frente al resultado de la "Fórmula de distancia" puede ser mejor para ciertas situaciones, ya que es menos costoso desde el punto de vista informático. Piensa "¿qué punto está más cerca?" escenarios donde no necesita una medición de distancia precisa.
En cuanto a la "Fórmula Haversine", es fácil de implementar y es agradable porque utiliza "Trigonometría esférica" en lugar de un enfoque basado en la Ley de Coseno basada en la trigonometría bidimensional, por lo que obtienes un buen equilibrio de precisión sobre la complejidad
Un caballero con el nombre de Chris Veness tiene un gran sitio web en http://www.movable-type.co.uk/scripts/latlong.html que explica algunos de los conceptos que le interesan y demuestra varias implementaciones programáticas; esto también debería responder a tu pregunta de conversión x / y.
Si le interesa obtener coordenadas basadas en un elipsoide en lugar de una esfera, eche un vistazo a http://en.wikipedia.org/wiki/Geodetic_system#From_geodetic_to_ECEF - proporciona las fórmulas y las constantes WGS84 que necesita para la conversión .
Las fórmulas allí también tienen en cuenta la altitud relativa a la superficie del elipsoide de referencia (útil si se obtienen datos de altitud de un dispositivo GPS).
Teoría para convertir GPS(WGS84)
a coordenadas cartesianas https://en.wikipedia.org/wiki/Geographic_coordinate_conversion#From_geodetic_to_ECEF_coordinates
Lo siguiente es lo que estoy usando:
- La longitud en GPS (WGS84) y las coordenadas cartesianas son las mismas.
- La latitud necesita ser convertida por los parámetros del elipsoide WGS 84 semieje mayor es 6378137 m, y
- Reciprocal de aplanamiento es 298.257223563.
Adjunté un código VB que escribí:
Imports System.Math
''Input GPSLatitude is WGS84 Latitude,h is altitude above the WGS 84 ellipsoid
Public Function GetSphericalLatitude(ByVal GPSLatitude As Double, ByVal h As Double) As Double
Dim A As Double = 6378137 ''semi-major axis
Dim f As Double = 1 / 298.257223563 ''1/f Reciprocal of flattening
Dim e2 As Double = f * (2 - f)
Dim Rc As Double = A / (Sqrt(1 - e2 * (Sin(GPSLatitude * PI / 180) ^ 2)))
Dim p As Double = (Rc + h) * Cos(GPSLatitude * PI / 180)
Dim z As Double = (Rc * (1 - e2) + h) * Sin(GPSLatitude * PI / 180)
Dim r As Double = Sqrt(p ^ 2 + z ^ 2)
Dim SphericalLatitude As Double = Asin(z / r) * 180 / PI
Return SphericalLatitude
End Function
Tenga en cuenta que la h
es la altitud sobre el WGS 84 ellipsoid
.
Por lo general, el GPS
nos dará H
de la altura MSL
superior. La altura de MSL
debe convertirse a la altura h
encima del WGS 84 ellipsoid
utilizando el modelo geopotencial EGM96
( Lemoine et al, 1998 ).
Esto se hace interpolando una cuadrícula del archivo de altura del geoide con una resolución espacial de 15 minutos de arco.
O si tiene algún nivel de GPS
profesional tiene Altitud H
( msl, altura sobre el nivel medio del mar ) y UNDULATION
, la relación entre el geoid
y el ellipsoid (m)
del dato de salida elegido de la tabla interna. puedes obtener h = H(msl) + undulation
Para XYZ por coordenadas cartesianas:
x = R * cos(lat) * cos(lon)
y = R * cos(lat) * sin(lon)
z = R *sin(lat)
Coordinate[] coordinates = new Coordinate[3];
coordinates[0] = new Coordinate(102, 26);
coordinates[1] = new Coordinate(103, 25.12);
coordinates[2] = new Coordinate(104, 16.11);
CoordinateSequence coordinateSequence = new CoordinateArraySequence(coordinates);
Geometry geo = new LineString(coordinateSequence, geometryFactory);
CoordinateReferenceSystem wgs84 = DefaultGeographicCRS.WGS84;
CoordinateReferenceSystem cartesinaCrs = DefaultGeocentricCRS.CARTESIAN;
MathTransform mathTransform = CRS.findMathTransform(wgs84, cartesinaCrs, true);
Geometry geo1 = JTS.transform(geo, mathTransform);