Identificar el punto latitudinal y longitudinal directamente debajo de un objeto en el espacio

Question

Dada la azimuth and altitude / latitude and longitude del observador al satélite, ¿cómo puedo encontrar con precisión la latitud y la longitud del punto de intersección directamente bajo el satélite?

Este problema me está afectando desde hace tiempo, ¡gracias de antemano por su ayuda!

Editar: También tengo los valores L1/2 p y L1/2 c/a.

Answer 1

Suponiendo que se dispone de la distancia del observador al satélite -sin la cual el problema no tiene solución definitiva-, se trata de resolver tres subproblemas, utilizando la estrategia de calcular las coordenadas cartesianas geocéntricas (x,y,z) del satélite.

En lo siguiente, a es el semieje mayor (6.378.137,0 metros en WGS 84 ) y b es el eje semiminor (aproximadamente 6.356.752,314 245 metros). Todos los cálculos deben utilizar sistemáticamente radianes o grados, lo que prefiera su software para las funciones trigonométricas. Los cálculos de doble precisión proporcionarán una precisión submicrónica, suponiendo que dicha precisión esté presente en los datos originales.

Los ejemplos (de trabajo) utilizan Mathematica .

Convertir entre coordenadas geodésicas y cartesianas

Podemos parametrizar una sección vertical del elipsoide de la forma (x, z) = (a * cos(t), b * sin(t)) para la latitud t y luego rotarlo alrededor del eje z. La inversión de las ecuaciones nos permite recuperar (lon, lat) a partir de (x,y,z).

ellipsoidToCartesian[{lon_, lat_}, {a_,b_}] := {a Cos[lat] Cos[lon], a Cos[lat] Sin[lon], b Sin[lat]};
cartesianToEllipsoid[{x_, y_, z_}, {a_,b_}] := {ArcTan[x, y], ArcTan[Norm[{x, y}]/a, z/b]};

Obtener una base para el plano tangente en cualquier punto

Este código devuelve una lista cuyo primer elemento es el origen y cuyo segundo elemento es la base ortonormalizada {este, norte, arriba}. Calcula las derivadas de los parámetros geodésicos y las normaliza para obtener una base para el plano tangente, luego utiliza su producto cruzado para calcular el vector unitario "up".

tangent[{lon_, lat_}, {a_, b_}] := Module[{north, east, up},
  north = {-a Cos[lon] Sin[lat], -a Sin[lat] Sin[lon], b Cos[lat]};
  north = north / Norm[north];
  east = {-a Cos[lat] Sin[lon], a Cos[lat] Cos[lon], 0};
  east = east / Norm[east];
  up = Cross[east, north];
  {ellipsoidToCartesian[{lon, lat}, {a, b}], {east, north, up}}
  ]

Calcular las coordenadas cartesianas del satélite

La altitud y el acimut, junto con la distancia al observador, son esencialmente coordenadas esféricas locales. Convirtiéndolas en un desplazamiento vectorial se obtienen las coordenadas del satélite en relación con la ubicación del observador.

La altitud y el acimut son relativos a un sistema de coordenadas local basado en (este, norte, arriba). A continuación se supone que el acimut se da al este del norte.

satellite[{azimuth_, altitude_}, {lon_, lat_}, elevation_, {a_, b_}] := 
  Module[{origin, east, north, up},
   {origin, {east, north, up}} = tangent[{lon, lat}, {a, b}];
   elevation (east Sin[azimuth] Cos[altitude] + 
       north Cos[azimuth] Cos[altitude] + up Sin[altitude]) + origin
   ];

La solución

Una vez obtenidas las coordenadas cartesianas del satélite, ya hemos terminado, porque sólo hay que volver a convertirlas en coordenadas geodésicas.

position[{azimuth_, altitude_}, {lon_, lat_}, elevation_, {a_, b_}] := 
 Module[{xyz},
  xyz = satellite[{azimuth, altitude}, {lon, lat}, elevation, {a, b}];
  cartesianToEllipsoid[xyz, {a, b}]
  ]

Una visualización dinámica confirma que todo esto funciona como está previsto. Un segmento de línea trazado desde el origen (gris) hasta el satélite (azul) coincide efectivamente con las coordenadas calculadas del satélite (rojo).

Para la comodidad de aquellos que puedan tener acceso a Mathematica El código de la visualización es el siguiente.

Module[{world, azimuth, altitude, lon, lat, elevation, a = 1, s, p, o},
 Manipulate[
  world = 
   ParametricPlot3D[
    ellipsoidToCartesian[{lon, lat}, {1, b}], {lon, -\[Pi], \[Pi]}, {lat, -\[Pi]/2, \[Pi]/2}, 
    PlotStyle -> Opacity[0.5], MeshFunctions -> {#4 &, #5 &}, 
    MeshStyle -> Directive[Opacity[0.5], Gray], Boxed -> False];
  s = satellite[{azimuth, altitude}, {lon, lat}, elevation, {a, b}];
  p = ellipsoidToCartesian[cartesianToEllipsoid[s, {a, b}], {a, b}];
  o = ellipsoidToCartesian[{lon, lat}, {a, b}];
  Show[Graphics3D[{
     Thick, Line[{{0, 0, 0}, s}], Line[{o, s}], 
     PointSize[0.02], Gray, Point[{0, 0, 0}], Green, Point[o], Blue, Point[s], Red, Point[p]
     }, Boxed -> False], world],
  {{azimuth, 0, "Azimuth"}, 0, 2 \[Pi]},
  {{altitude, \[Pi]/2, "Altitude"}, 0, \[Pi]/2},
  {{lon, 0, "Longitude"}, -\[Pi], \[Pi]},
  {{lat, \[Pi]/3, "Latitude"}, -\[Pi]/2, \[Pi]/2},
  {{elevation, 1, "Distance"}, 0, 2},
  {{b, 1, "b/a"}, 0, 1}
  ]]

Answer 2

Yo trataría de traducir esto a coordenadas cartesianas: El observador está en $(R \cos \phi \cos \lambda, R\cos\phi R \sin \lambda, R \sin \phi)$ con $R=$ Radio de la Tierra, $\lambda=$ longitud, $\phi=$ latitud. A continuación, el acimut y la altitud (medidos con respecto a la tangente en el observador) determinan una recta y con la elevación se obtiene una cuadrática con dos soluciones (una está "detrás" del observador en la dirección equivocada), por lo que la posición $(x,y,z)$ del satélite en el espacio. Finalmente, $\frac z{H+R}$ es el seno de la latitud y $\frac y x$ es la tangente de la longitud.