Análisis y diseño de métodos de filtrado de conjunciones espaciales aplicados a poblaciones realistas en órbita baja terrestre
Abstract
Dado el cada vez mayor riesgo de colisiones entre cuerpos en órbita alrededor de la Tierra, se decide realizar una aproximación al problema del análisis del mismo y la elaboración de herramientas que lo traten. En ese sentido, se ha hecho un repaso del problema de los dos cuerpos, junto con algunos apuntes al respecto de la manera en la que se recoge la información sobre la situación de los satélites, los TLEs. A continuación, se elabora un propagador de órbitas basado en la ecuación del movimiento, al que incorporaremos los efectos de las perturbaciones vinculadas a la no esfericidad de la Tierra (J2 ) y la resistencia atmosférica. Además, se incluirán comprobaciones del funcionamiento de este propagador. Con el propagador completo, se van a comentar los desarrollos hechos para elaborar herramientas de filtrado de posibles conjunciones. Con estas en mano, también se harán algunas pruebas de comprobación para verificarlas. Finalmente, se va a realizar una prueba, usando todas las herramientas creadas sobre una población de objetos en órbita baja y analizando las posibles conjunciones y el riesgo de colisiones en la zona. Como conclusiones, la elaboración de las herramientas necesarias para cumplir los objetivos propuestos nos muestra las limitaciones propias de los modelos físicos empleados y la necesidad de ajustarlos al ámbito de aplicación que vaya a tener. Por otra parte, los resultados finales obtenidos de la prueba ejecutada nos muestran claramente como existen casos de conjunciones incluso en poblaciones relativamente reducidas, de lo que se intuye la necesidad de mejorar las herramientas de filtrado y detección para manipular una población cada vez mayor y más expuesta a posibles colisiones. Las futuras líneas de trabajo deberían enfocarse a mejorar las herramientas ya disponibles, mediante modelos más precisos y filtros más rápidos, además de continuar en el desarrollo de herramientas que, para conjunciones verificadas, analicen en profundidad la posibilidad de colisión, además de abordar el estudio de maniobras de evasión.
Full text
Proyecto Fin de Carrera Ingeniería de Telecomunicación Formato de Publicación de la Escuela Técnica Superior de Ingeniería Autor: F. Javier Payán Somet Tutor: Juan José Murillo Fuentes Dep. Teoría de la Señal y Comunicaciones Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2013 Trabajo de Fin de Grado Grado en Ingeniería Aeroespacial Análisis y diseño de métodos de filtrado de conjunciones espaciales aplicados a poblaciones realistas en órbita baja terrestre Autor: César Sánchez Puente Tutores: Rafael Vázquez Valenzuela Ana Sánchez Rivero Dpto. Ingeniería Aeroespacial Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2023
Trabajo de Fin de Grado Grado en Ingeniería Aeroespacial Análisis y diseño de métodos de filtrado de conjunciones espaciales aplicados a poblaciones realistas en órbita baja terrestre Autor: César Sánchez Puente Tutores: Rafael Vázquez Valenzuela Profesor Titular Ana Sánchez Rivero Dpto. Ingeniería Aeroespacial Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2023
Trabajo de Fin de Grado: Análisis y diseño de métodos de filtrado de conjunciones espaciales aplicados a poblaciones realistas en órbita baja terrestre Autor: César Sánchez Puente Tutores: Rafael Vázquez Valenzuela Ana Sánchez Rivero El tribunal nombrado para juzgar el trabajo arriba indicado, compuesto por los siguientes profesores: Presidente: Vocal/es: Secretario: acuerdan otorgarle la calificación de: El Secretario del Tribunal Fecha:
Agradecimientos Supongo que es algo tópico por mi parte, pero reconozco que está siendo una sección bastante "delicada" de redactar dentro del documento, igual que es algo tópico decir que no soy una persona muy dada a escribir este tipo de textos. También intentaría decir que es de las partes a las que más tiempo le he dedicado, pero no es cuestión de faltarle el respeto a mi ordenador (por primera vez en este trabajo usaré esta expresión: véanse las Secciones 5.3. y 5.4.). En cualquier caso, si reconozco que le he dado bastantes vueltas a cómo organizar en unas líneas unos agradecimientos que más que dirigidos a lo que ha supuesto este trabajo, irían dirigidos en su conjunto a estos 4 años en la Escuela. Las cosas hay que hacerlas (o al menos eso intento con mayor o menor éxito) con orden, y es por eso que, antes de nada, quiero agradecer profundamente a mis tutores en este proyecto, Rafael y Ana, por confiar en mí para hacer este trabajo y por su gran ayuda estos meses. Sin ellos hubiera sido imposible orientarlo y completarlo. También les quiero agradecer todas las lecciones que he podido aprender sobre esta disciplina, pero también sobre cómo hacer y redactar este tipo de documentos, buscar cierta información y también por los avisos sobre cuestiones mejorables del trabajo. Dicen que el Diablo está en los detalles, espero que con su ayuda lo haya espantado. Pero como dije al principio, también quería aprovechar esta página para agradecer por estos 4 años. Por un lado, quería empezar agradeciendo a los profesores, ya no solo por los conocimientos de incalculable valor que me han proporcionado, sino también por lecciones de vida que, de alguna forma, haya podido aprender gracias a ellos. Muchas gracias también a los que me habéis acompañado entre tantas clases, prácticas, exámenes y ratos en la cafetería. A algunos os llevo viendo desde el mismo día de la presentación del Grado, con otros he perdido algo el contacto con el tiempo y a otros os he conocido más bien al final. Con todos ha sido más fácil pasar este tiempo en la Escuela, 4 años más tarde le doy la razón al profesor que me avisó de que con gente al lado se hacía más sencilla esta experiencia. Por supuesto, también quería agradecer a mis amigos, los que habéis estado de siempre y los que llegasteis algo más tarde, me habéis hecho más llevadero este tiempo, y siempre ha sido un apoyo saber que, bien estuvierais aquí, bien a 400 km, bien pudiéramos hablar en persona o bien fuera por videollamada, estabais ahí si lo necesitaba. Y he querido dejar para el final a aquellas personas sin las que nada de esto hubiera sido posible. Si tuviera que agradecer a mi familia todas las cosas que les debo agradecer, como mínimo triplicaría la extensión del documento, así que me limitaré a este espacio. Gracias a mis padres y a mi hermano por estar siempre ahí, por toda la ayuda que me habéis dado, por aguantarme todo este tiempo, por acercarme al examen cuando la línea del autobús me ha dejado tirado y por tantas y tantas cosas que no caben en apenas estas líneas. César Sánchez Puente 2023 I
Resumen Dado el cada vez mayor riesgo de colisiones entre cuerpos en órbita alrededor de la Tierra, se decide realizar una aproximación al problema del análisis del mismo y la elaboración de herramientas que lo traten. En ese sentido, se ha hecho un repaso del problema de los dos cuerpos, junto con algunos apuntes al respecto de la manera en la que se recoge la información sobre la situación de los satélites, los TLEs. A continuación, se elabora un propagador de órbitas basado en la ecuación del movimiento, al que incorporaremos los efectos de las perturbaciones vinculadas a la no esfericidad de la Tierra ( J2 ) y la resistencia atmosférica. Además, se incluirán comprobaciones del funcionamiento de este propagador. Con el propagador completo, se van a comentar los desarrollos hechos para elaborar herramientas de filtrado de posibles conjunciones. Con estas en mano, también se harán algunas pruebas de comprobación para verificarlas. Finalmente, se va a realizar una prueba, usando todas las herramientas creadas sobre una población de objetos en órbita baja y analizando las posibles conjunciones y el riesgo de colisiones en la zona. Como conclusiones, la elaboración de las herramientas necesarias para cumplir los objetivos propuestos nos muestra las limitaciones propias de los modelos físicos empleados y la necesidad de ajustarlos al ámbito de aplicación que vaya a tener. Por otra parte, los resultados finales obtenidos de la prueba ejecutada nos muestran claramente como existen casos de conjunciones incluso en poblaciones relativamente reducidas, de lo que se intuye la necesidad de mejorar las herramientas de filtrado y detección para manipular una población cada vez mayor y más expuesta a posibles colisiones. Las futuras líneas de trabajo deberían enfocarse a mejorar las herramientas ya disponibles, mediante modelos más precisos y filtros más rápidos, además de continuar en el desarrollo de herramientas que, para conjunciones verificadas, analicen en profundidad la posibilidad de colisión, además de abordar el estudio de maniobras de evasión. III
XÍndice 3.3.1 Modelo del efecto de la resistencia atmosférica 34 3.3.2 Validación del propagador 35 4 Filtros de conjunciones 39 4.1 Estructura general del programa 39 4.2 Apogee-Perigee Filter 41 4.3 Orbit Path Filter 43 4.4 Time Filter 47 5 Puesta en funcionamiento: Análisis y resultados 51 5.1 Alcance del análisis y población seleccionada 51 5.2 Aplicación del filtro de apogeo-perigeo: Resultados 52 5.3 Aplicación del filtro Orbit Path Filter: Resultados 52 5.4 Aplicación del filtro Time Filter: Resultados 54 5.5 Propagación de un caso de posible conjunción 55 6 Conclusión y líneas de trabajo futuras 59 6.1 Conclusión 59 6.2 Líneas futuras de trabajo 60 Índice de Figuras 63 Índice de Tablas 65 Bibliografía 67
1 Introducción Múltiples causas que van desde el creciente números de actores interesados y con capacidad de enviar equipos al espacio hasta la acumulación de sistemas desplegados décadas atrás vuelven cada vez más necesaria la gestión y análisis del riesgo de colisiones de dichos cuerpos en órbita. Gracias a la observación y seguimiento se han podido propagar la trayectoria de dichos cuerpos y detectar posibles colisiones, sin embargo, esto puede plantear problemas desde el punto de vista computacional, haciéndose necesario optimizar los procesos usados. En este contexto se encuadra el desarrollo de este trabajo. 1.1 Motivación 1.1.1 El problema de la basura espacial El 4 de octubre de 1954, se lanzaba desde el cosmódromo de Baikonur (actual Kazajistán) el satélite Sputnik-1. Con su breve misión, de apenas 22 días, inauguró la ocupación del espacio extraterrestre, una ocupación que fue seguida por la Sputnik-2, la Explorer-1, el Project SCORE y tantos y tantos otros que, bien destinados a permanecer en órbita terrestre, bien dirigidos al espacio profundo, reclamaron para sí porciones del cielo. Este proceso, que ha aumentado de forma exponencial con los años (como podemos ver en la Figura 1.1.), ha llevado aparejado una consecuencia irremediable: la cada vez mayor posibilidad de que se produzcan coincidencias entre las trayectorias de objetos en órbita con el riesgo de colisión inherente. El año 2022 fue el que tuvo la mayor cantidad de lanzamientos en la historia de la humanidad (con 174 lanzamientos exitosos) batiendo al 2021 (con 132), que a su vez superó el record conseguido en 2020 (de 103) como podemos ver en la siguiente tabla: Tabla 1.1 Lanzamientos anuales exitosos desde el inicio de la década pasada hasta el año 2022 [1]. Año NºLanzamientos 2010 70 2012 71 2014 89 2016 82 2018 112 2020 103 2021 132 2022 174 1
2Capítulo 1. Introducción Este incremento se ha debido a las cada vez mayores posibilidades de negocio en el despliegue de sistemas de satélites en órbita, además de cuestiones tanto científicas en la investigación y observación de la Tierra y el Espacio como a cuestiones de seguridad y militares de diversos estados. Sin embargo, este aumento ha tenido como consecuencia lo siguiente: •Por un lado una cada vez mayor cantidad de satélites en órbita. • Por otro una cada vez mayor población de basura espacial, asociada a los lanzamientos de satélites ya mencionados, falta de regulación en el pasado respecto a esta cuestión y a la mala praxis de algunas entidades. Figura 1.1 Número de objetos lanzados al espacio anualmente, por paises y acumulado para todo el mundo. Observamos como se ha disparado el número desde 2018, estando el aumento liderado por EEUU [2]. En la actualidad diversas medidas se llevan a cabo para prevenir o minimizar la huella dejada por un objeto tras una misión. Así, por ejemplo, los satélites deben llevar una cantidad extra de combustible para realizar maniobras de desorbitado una vez su vida útil haya terminado, y también existe, en regiones del espacio de la órbita GEO, una órbita reservada para localizar objetos que ya no puedan dar servicio, la órbita cementerio, y así no estorbar o amenazar a otros que si estén en funcionamiento. Con estas medidas se busca reducir el riesgo de colisión, pero, aunque bajo, no es nulo, como muestra el caso de la colisión del Iridium-33 con el COSMOS-2251 en 2009. Es por ello que hay que mantener un cierto grado de observación para prevenir incidentes con graves consecuencias, pues por un lado se produce un incremento de basura espacial pero también se deteriora o, en el peor de los casos, se pierde el servicio que proporcionara el objeto u objetos implicados.
1.1 Motivación 3 Figura 1.2 Número de objetos en la órbita terrestre clasificados por su tipo. La basura espacial estaría representada por las curvas rosa y naranja [3]. 1.1.2 El interés de las partes implicadas La detección de posibles conjunciones entre las órbitas de satélites es una labor que implica múltiples agentes. Entre ellos, podemos encontrar a los propios operadores de dichos objetos (obviamente interesados en que estos no resulten accidentados) y organizaciones gubernamentales (interesadas en garantizar la seguridad de los objetos en órbita). Es por tanto un campo donde estos actores tienen motivos para buscar mejoras en los procesos que los hagan más eficientes y rápidos, así como preocupaciones por la actual forma de funcionar de estos. Estos intereses y preocupaciones se ven reflejados en la encuesta hecha por miembros de la compañía DEIMOS Space UK, con motivo de la 8ºConferencia Europea sobre la Basura Espacial en 2021, a operadores de satélites, proveedores de servicios de detección y control del entorno del satélite y desarrolladores de métodos de detección [4]. En sus respuestas a los encuestadores, los operadores afirmaron en su mayoría emplear métodos propios de detección aunque desearían contar con un software desarrollado o certificado por alguna entidad. Sobre la ejecución de maniobras ninguno de los encuestados declaró realizar más de 3 maniobras al año por satélite, e incluso el 60% afirmó que de media apenas realizaban 1. Sin embargo, un 60% también afirmó que en la actualidad detectan un incremento en el número de eventos en los que sus satélites se ven implicados. También resultó significativo que los encuestados no declararon temores al respecto del empleo de satélites cada vez más pequeños, pero si lo hicieron al aumento de la población de objetos en la órbita terrestre. Por su parte, los proveedores se mostraron a favor de mejorar los catálogos de objetos espaciales, así como el desarrollo de soportes técnicos para gestionar el tráfico espacial y en operación. También cuentan con sus propios métodos de análisis y, a diferencia de los operadores, no presentan un consenso en cuanto al futuro. Solo un 40% mostró preocupación por el aumento de la población,
4Capítulo 1. Introducción considerando el 60% restante que no iba a ser un problema o no contestando. No obstante, el 60% si mostró preocupación por las cada vez mayores constelaciones de satélites. Finalmente, entre desarrolladores actualmente la mayor parte de trabajos se dirigen a cálculo de probabilidades de colisión y recomendación de maniobras. Este grupo se encuentra dividido respecto al tipo de programas desarrollados, habiendo trabajado desde métodos de análisis numéricos a analíticos a propagadores de TLEs. De nuevo, y como los agentes anteriores, priorizan la mejora de catálogos y de los medios de soporte en operación para controlar el entorno del satélite. Como se puede ver en estos resultados existe, en mayor o menor medida, constancia de que en el futuro se van a dar fenómenos que pueden afectar a la operación y seguridad de los objetos en órbita, habiendo un amplio consenso en cuanto a la mejora de los medios y técnicas usadas para detectar las conjunciones y controlar los elementos en órbita entorno a la Tierra. 1.2 Estado del arte del análisis de conjunciones En esta sección se va a describir en qué consisten en la actualidad los métodos de análisis de conjunciones en órbita, como se realizan, con qué criterios y qué medidas se toman, en caso de ser necesario, para evitarlas. 1.2.1 Detección de colisiones El proceso de acción contra una posible colisión en órbita comienza con la detección de una conjunción entre las órbitas de cuerpos situados en el espacio. Esta detección se realiza usando los catálogos de objetos en órbita que conocemos. Con esa información, dicho a groso modo, comparamos las órbitas propagadas para cada satélite con las del resto. Como es evidente esto resulta muy costoso computacionalmente, de ahí que se lleve a cabo una criba de aquellos objetos que esté claro que no suponen una amenaza para el satélite analizado. Entre los criterios que, para un cierto satélite, nos permite descartar posibles riesgos se encuentran [4]: •Que el apogeo de la órbita de menor semieje sea menor que el perigeo de la de mayor. • Que, en la línea de intersección entre los planos orbitales de los dos cuerpos, sus distancias a la Tierra sean distintas. • Que, en caso de que las órbitas si tengan puntos de corte entre ellas, el desfase entre los cuerpos haga imposible que se crucen en dichos puntos. Con estos métodos conseguimos reducir la lista de cuerpos a valorar en nuestro análisis, y una vez hecho esto tendremos que propagar nuestro satélite estudiado, así como todos aquellos que hayamos determinado que pueden suponer una amenaza. No obstante, a la hora de propagar, y como se comentó en el apartado anterior, las imperfecciones de nuestros modelos dan lugar a unos ciertos errores e incertidumbres de nuestra medida, que valoraremos recurriendo a matrices de covarianzas para las distintas variables de la posición del satélite. Si como resultado del filtrado anterior llegamos a una o varias posibles conjunciones, entonces tenemos dos posibles medios para comprobar como de posible es que en dichas situaciones haya una colisión [4]: •Volumen de exclusión: En este caso, definimos en el entorno de nuestro satélite una región de seguridad. El criterio consistiría simplemente en comprobar si alguno de los objetos que potencialmente pudieran colisionar con este entra realmente en dicha región, de detectarse algún cuerpo que si entonces se procedería a realizar una maniobra de evasión. Este criterio, no obstante, resulta muy conservador, dando lugar a maniobras innecesarias.
1.2 Estado del arte del análisis de conjunciones 5 •Probabilidad de colisión o POC: En este caso definimos un cierto máximo de probabilidad aceptable, y solo actuaremos cuando se supere. El proceso consiste en definir una distribución de probabilidades de colisión para una cierta distancia ∆r del satélite que pueda colisionar de nuestro satélite y, sea sobre esa distancia una circunferencia de radio el llamado “Hard-Body Radius” (véase la Figura 1.3.), integramos la curva de probabilidad sobre toda esa región. El área obtenida será nuestra probabilidad [5]. Figura 1.3 Definición del Hard-Body Radius (HBR), que delimita el área de integración de la curva de probabilidad [5]. Las curvas de probabilidad reflejarán el grado de incertidumbre existente en la estimación de la posición de nuestro objeto, conforme menos segura sea esa predicción más distribuida estará la probabilidad a lo largo de las distancias (véase la Figura 1.4.) Figura 1.4 Representación de la curva de probabilidad de colisión para distintas distancias y el área Ac que usamos para integrar [6]. Aparte de este método, que requiere conocer las covarianzas, la probabilidad de colisión puede obtenerse por otras vías como pueden ser estimando dichas covarianzas, e incluso se pueden incluir factores que las escalen, para tener en cuenta como varían conforme se acerca el evento. También es posible la estimación mediante el método de Monte-Carlo. De esta manera nosotros podemos, para nuestro satélite estudiado y cada uno de los objetos que son susceptibles de colisionar con él, conocer el riesgo de que, en efecto, ocurra esa colisión. Aquellas conjunciones que si puedan ser un riesgo se denominan “High Interest Events” (HIE), y el umbral de riesgo a partir del cual una conjunción se denomina como tal es definido por cada operador de satélites en sus protocolos de seguridad. Por ejemplo, la ESA lo pone en 10´4 y la
6Capítulo 1. Introducción NASA los clasifica según su probabilidad en verdes (menor de 10´7 ), amarillo (entre 10´7 y 10´4 ) y rojo (mayor de 10´4 ). Una vez identificado un HIE, se procede a un análisis más profundo del fenómeno y, finalmente, si se siguen superando los límites aceptables, se procede a actuar [4]. 1.2.2 Maniobras de evasión Una vez evaluado en profundidad el HIE y verificado que existe un riesgo real de colisión, queda en manos del operador u operadores implicados definir las acciones a tomar para evitar el impacto. En este sentido existen múltiples variables que se deben considerar: •Que la maniobra sea lo más óptima posible y con el menor consumo de combustible. • Que, en la medida de lo posible, no altere o inhabilite el servicio que dicho satélite proporciona en tierra. •Que la maniobra no implique el riesgo de otra colisión. Tener en cuenta todos estos parámetros de forma simultánea supone un problema muy amplio, y es por ello que se simplifica, primero buscando la solución óptima en términos de combustible y, más adelante, incluyendo el resto de variables. De manera general, la maniobra óptima supone más riesgos que el resto, y por ello se opta en su lugar por una maniobra suficientemente buena. De cualquier forma, dado que estas maniobras suponen un gasto extra para la misión, cuando hay que realizarlas se suelen usar para realizar otras operaciones de mantenimiento, como puede ser el “station keeping". La maniobra más frecuente para ejecutar, incluso si el catálogo de opciones es más amplio, se suele reducir a pequeños impulsos que hagan cambiar los parámetros de la órbita lo suficiente como para evitar el impacto, seguido de una corrección de la órbita que nos devuelva a la posición original. Estos impulsos se ejecutan o bien tangenciales a la trayectoria seguida u oblicuas a esta, y simplifican los cálculos. De todas formas, la determinación de la dirección e instante adecuados para ejecutar la maniobra suponen un complejo problema [7]. Figura 1.5 Diagrama que muestra una posible colisión entre dos satélites y la maniobra de evasión realizada por el segundo de ellos [8].
1.3 Alcance y objetivos 7 1.3 Alcance y objetivos Con este trabajo se busca realizar un acercamiento al desarrollo de técnicas de análisis de conjunciones, así como a la comprobación de eventuales colisiones, para poder abordar los cada vez mayores problemas como los expuestos en las páginas anteriores. Con esta idea en mente lo primero que se ha llevado a cabo es la elaboración de un propagador de órbitas de satélites, basado en la ecuación del movimiento. Con este hecho, luego se ha mejorado, incluyendo los efectos perturbadores más importantes sobre la órbita de los satélites ubicados en la región de estudio: la órbita Baja Terrestre (LEO). Estos efectos son los vinculados a la no esfericidad de la Tierra (en nuestro caso a través del armónico J2 asociado al achatamiento terrestre) y a la resistencia que la atmósfera plantea al movimiento del satélite. Con el propagador completo, luego se han construido algoritmos de filtrado que permitan cribar las posibles parejas de satélites susceptibles de tener una conjunción. Se han elaborado versiones de los tres tipos de filtros comentados en la sección anterior, realizando pruebas para verificar su funcionamiento. Así, se busca un primer tanteo al problema de la detección de conjunciones, aunque el estudio en profundidad de posibles colisiones e incluso el estudio de maniobras de evasión ya quedan fuera del alcance de este trabajo. 1.4 Estructura Con el objetivo expuesto, este trabajo se ha organizado de la siguiente manera: •Capítulo 2: Antes de abordar en profundidad la elaboración de las herramientas necesarias para ejecutar la tarea señalada, se van a presentar algunas de las bases teóricas de la mecánica orbital, haciendo especial énfasis en la resolución del problema de los 2 cuerpos, las perturbaciones tanto consideradas como no consideradas en nuestro modelo y como le afectan, y el funcionamiento del Two-Line Element set (TLEs). •Capítulo 3: Con las bases sentadas, se mostrarán las claves del funcionamiento del propagador elaborado, incluyendo las bases del sistema de referencia usado, el modelado de las perturbaciones incorporadas y algunas muestras de su funcionamiento. •Capítulo 4: Se procederá a continuación a mostrar los filtros programados. Como en el Capítulo 3, se definirá la forma de funcionar de cada algoritmo, junto con muestras y ejemplos de la ejecución de los mismos. •Capítulo 5: Una vez desarrollado el funcionamiento de estas herramientas, se mostrará un caso práctico de uso. Se comentará la población objeto de estudio, la ejecución de cada filtro con el resultado obtenido y analizando cómo ha funcionado el filtro en su ejecución, y un caso de propagación para una de las parejas susceptible de tener una conjunción. •Capítulo 6: Para acabar, un último capítulo donde se expondrán las conclusiones de todo lo realizado con anterioridad y se mostrarán las líneas futuras de trabajo recomendadas.
2 Introducción a la mecánica orbital y otros conceptos Para este capítulo se expondrán una serie de conceptos e ideas básicas de la mecánica orbital y sobre el estudio de órbitas de satélites y posibles colisiones que pueden resultar de interés de cara a comprender el resto del siguiente trabajo. 2.1 Descripción del problema de los 2 cuerpos Base para todo el desarrollo posterior, este problema se trata de uno de los más relevantes de la mecánica clásica. Consiste en el estudio de los movimientos en un sistema formado por dos cuerpos con masa, cada uno de ellos sometido a la atracción gravitatoria del otro. Será Newton el cual, a finales del siglo XVII y demostrando las Leyes de Kepler, consiga definir unas ecuaciones que describieran de forma adecuada dicho movimiento. 2.1.1 La obtención de la ecuación del movimiento Para resolver dicho problema, se recurren a las siguientes hipótesis: 1. El sistema se encuentra completamente aislado del resto del Universo, de modo que no existe interacción entre ambos. 2. Todas las masas del problema se considerarán puntuales y localizadas en el centro de masas de cada cuerpo. Al respecto de esta segunda hipótesis, y a través de la comparación con una masa puntual, Newton mostró que la atracción generada por esta era equivalente a la de un cuerpo de geometría esférica y con su masa distribuida homogéneamente. Por lo tanto, la segunda hipótesis se puede sustituir por este teorema al tratar cuerpos celestes, de asumir que estos son esféricos [9]. Sean estas hipótesis, comenzamos a trabajar sobre la expresión de la fuerza gravitatoria obtenida por el propio Newton: F“GMm r2 r r(2.1) Donde G se corresponde con la constante de la gravitación universal, M y m las masas del problema y r el vector de la distancia entre ambos (y r su módulo). En la Ecuación 2.1 se puede particularizar para los dos cuerpos implicados en nuestro sistema como: 9
16 Capítulo 2. Introducción a la mecánica orbital y otros conceptos con v la velocidad del cuerpo y vatm la velocidad de la atmósfera en el espacio, que se puede modelar como: vatm “ ωCˆr(2.22) Finalmente, ρ se trata de la densidad atmosférica, que podemos aproximar mediante el siguiente modelo exponencial [15]: ρ“ρ0e˜RC´r H¸(2.23) con ρ0la densidad del aire al nivel del mar y Huna constante de valor 8500m •Presión de la radiación solar: Al impactar las partículas de la luz solar contra un cuerpo, generan una presión mecánica en este. La magnitud de esta presión dependerá de nuestra cercanía al Sol (pues define la intensidad de dicha luz) y el área sobre la que incida. Además, no es necesario que el satélite se encuentre frente al sol para sufrir esta presión, la luz reflejada por otros cuerpos como la Tierra o la Luna también genera una presión similar. Este fenómeno es la base de la propulsión de velas solares. La acción perturbadora para la radiación solar se describe como [9]: : rpert “´ F mv e@(2.24) Siendo F la fuerza ejercida y que para una placa plana se modela como: F“pAp1`εqcospφ@q(2.25) donde A es el área, ε el coeficiente de reflectividad (que va desde 1 para cuerpos totalmente reflectantes a -1 para totalmente transparentes), φ@ el ángulo de incidencia y p la presión que se calcula como: p“I c(2.26) siendo Iel flujo de radiación solar y cla velocidad de la luz. Así, estas son las principales causas de las perturbaciones, pero cada una de ellas tiene unos efectos distintos sobre nuestros elementos orbitales. Dependiendo de su continuidad en el tiempo se pueden dividir entre dos: •Seculares: Perturbaciones cuyos efectos en los elementos orbitales son permanentes. Estos, realmente, son generalmente los más importantes, pues son aquellos cuyos efectos permanecerán en el tiempo y crecerán de forma monótona, frente a los periódicos que irán oscilando con el tiempo •Periódicas: Perturbaciones cuyos efectos en los elementos orbitales se van repitiendo con un cierto período. Según la duración del mismo se pueden dividir a su vez en: – De largo periodo: Son aquellos cambios de los elementos orbitales que se suceden con un periodo de orden mayor al periodo orbital. – De corto periodo: Son aquellos cambios de los elementos orbitales que se suceden con un periodo de orden similar al periodo orbital. Todos estos fenómenos se pueden modelar aproximando la curva de la anomalía a una serie de Fourier [9]. Dentro de esta, el primer coeficiente será el asociado a los efectos seculares, y el resto a
2.1 Descripción del problema de los 2 cuerpos 17 los periódicos. Así, por ejemplo, si definiéramos nuestro elemento orbital e mediante la siguiente curva: eptq“0.3`0.00005t`0.00001 sinpω1tq`0.00001 cospω2tq(2.27) Con ω2ąω1 , entonces el primer término sería el valor inicial, el segundo el secular (el cambio monótono de la excentricidad con el tiempo), el tercero un cambio periódico de periodo largo, y el último uno de periodo corto. Figura 2.5 Tabla con el tipo de efectos (S para seculares y P para periódicas) para cada uno de los tipos de perturbaciones [9]. 2.1.5 Concepto de propagador y algunos ejemplos Por propagar entenderemos el cálculo de los elementos orbitales (y por tanto de las características del movimiento y de las coordenadas) de un determinado objeto a partir de esos mismos datos para un instante de tiempo conocido. De cara a propagar las órbitas se recurre al uso de propagadores, algoritmos que, dado el conjunto de elementos orbitales de partida, nos dan los datos actualizados para los instantes temporales que deseemos. Así, y asumiendo las hipótesis básicas del problema de los dos cuerpos, podemos recurrir al modelo más sencillo, el kepleriano, ya comentado anteriormente. En este modelo es importante mencionar que el set de elementos orbitales iniciales permanece constante para todo instante de tiempo, y que el recorrido de la órbita se puede resolver sin necesidad de recurrir a la resolución numérica de la Ecuación 2.14. mediante las leyes horarias. Con ellas podemos obtener la posición de un cuerpo en su órbita para un instante de tiempo determinado, aunque para nuestro caso se expondrá la ley horaria para órbitas elípticas. Con este modelo podemos resolver el problema en dos sentidos, tanto dado el tiempo transcurrido determinar la posición sobre la órbita como dada la posición final determinar el tiempo transcurrido. Para la resolución, no obstante, recurriremos a dos conceptos: •Velocidad orbital media: Se trata de una velocidad coincidente con la velocidad angular que tendría el objeto si la órbita fuera circular. Se define como: n“cµ a3(2.28)
18 Capítulo 2. Introducción a la mecánica orbital y otros conceptos •Anomalía media: Asociado a la velocidad anterior sería el ángulo recorrido en esa órbita circular en un cierto intervalo de tiempo. M“n∆t(2.29) Una vez descritos estos, el proceso para obtener el ángulo recorrido dentro de la órbita a partir de la posición inicial (de la que conocemos la anomalía verdadera) requiere de dar un pequeño rodeo. Esto se debe a que en las leyes horarias el tiempo ∆t viene referido al tiempo transcurrido desde que el satélite paso por el perigeo de la órbita, de esta manera el primer paso tendrá que ser determinar el tiempo que ha pasado desde que el satélite cruzó este punto hasta que llega a la posición inicial. En ese sentido, y a partir de la anomalía verdadera θ inicial, debemos obtener el valor de la anomalía excéntrica, E, de la ecuación: tan θ 2“c1`e 1´etan E 2(2.30) Una vez obtenida E calculamos la anomalía media de nuestra posición inicial mediante la ecuación de Kepler: M“E´esin E(2.31) Y una vez definida esta, y aplicando la definición de la anomalía media y la velocidad orbital media vistas antes en las Ecuaciones 2.17. y 2.16. respectivamente, podemos llegar al tiempo transcurrido entre que cruzamos el perigeo y llegamos al punto inicial conocido. A continuación, añadimos dicho tiempo al tiempo transcurrido entre la posición inicial y la posición final que deseamos conocer. Con estos datos calculamos M a través de la Ecuación 2.17., con M resolvemos la ecuación de Kepler y obtenemos E y, finalmente, con E entramos en la Ecuación 2.18. y determinamos θ . De esta manera podremos resolver de forma no numérica la órbita del satélite, para nuestro modelo sin perturbaciones de forma que sus elementos orbitales sean constantes. Según la forma en que se tengan en cuenta estas perturbaciones que afecten a nuestro sistema podemos encontrar múltiples propagadores distintos [16]. Un caso sería el propagador básico o kepleriano, que no tiene en cuenta ninguna perturbación, pero además también podemos encontrar otros, como el J2 medio que incluye una dependencia funcional del tiempo que tiene en cuenta el efecto del achatamiento de la Tierra en nuestra órbita (simple de usar y muy útil especialmente en órbita baja terrestre). En este modelo, los elementos orbitales se definen de la siguiente manera: a“a0(2.32) e“e0(2.33) i“i0(2.34) Ω“Ω0´3 2nR2 C p2J2cos ipt´t0q(2.35) ω“ω0`3 4nR2 C p2J2p5 cos2i´1qpt´t0q(2.36)
2.2 Two-Line Element set: Estructura y ejemplo 19 M“M0`ˆn`3 4nR2 C p2J2a1´e2p2´3 sin2iq˙pt´t0q(2.37) J2 es el armónico asociado al achatamiento de la Tierra y RC el radio terrestre. Nótese como en lugar de la anomalía verdadera se usa la media [9]. Como se puede apreciar en las ecuaciones, se incluyen términos que varían con t , los cuales complican el algoritmo, pero incluyen las variaciones de los elementos orbitales que en nuestro modelo básico se omitían, aumentando su precisión. No obstante, estos modelos de las perturbaciones siguen sin ser completamente exactos, dando lugar, con el paso del tiempo, a nuevos errores acumulados que nos obligarán a actualizar nuestros datos. 2.2 Two-Line Element set: Estructura y ejemplo En el siguiente apartado se van a abordar otra serie de conceptos que pueden resultar de utilidad para entender los desarrollos posteriores, como pueden ser en que consiste el TLEs o que son los propagadores y algunos ejemplos. 2.2.1 Diferencia entre elementos orbitales medios y osculantes Antes de hablar sobre como se estructuran los TLEs, es conveniente señalar la diferencia entre los elementos orbitales medios (los que aparecen en los mismos) y los llamados osculantes. Los elementos orbitales osculantes o instantáneos son aquellos que nos permiten para cada instante de tiempo obtener la posición real de nuestro cuerpo. En este caso se considerarían, implícitamente, los efectos de todas las perturbaciones comentadas anteriormente y, por tanto, la órbita reproducida será una órbita no elíptica (no kepleriana) y que dependerá del propagador usado para obtenerla. Debido a la acción de todas estas perturbaciones sobre la resolución de la órbita, si nosotros representáramos como cambian los elementos orbitales veríamos una curva irregular, con muchas variaciones (algunas periódicas y otras permanentes) y la órbita descrita sería diferente para cada revolución dada por el satélite. De esta manera, esta clase de elementos orbitales reproducen una posición específica del satélite en un instante concreto de tiempo y solo serán válidos para este y no para ningún otro momento, ni anterior ni posterior. También es de mención que en ocasiones para este tipo de elementos orbitales se incluyen otros rasgos que puedan ser necesarios para considerar esas perturbaciones, como los armónicos de la Tierra o el coeficiente balístico. Por su parte, los elementos orbitales medios resultan de una aproximación del movimiento del satélite en el cual se desprecien los efectos periódicos de corto periodo. Así, estos serán unos elementos orbitales que recogerán el cambio de los mismos por la acción de las distintas perturbaciones, pero solo tendrá en cuenta las variaciones seculares. Al representar la órbita, las posiciones no cambian de forma tan rápida entre una revolución y otra y, a diferencia de los osculantes, los medios servirán para representar la posición del satélite en un instante de tiempo específico, pero además también serán próximos a los de otros instantes de tiempo. Como contraparte de no ser tan diferentes los valores para cada revolución, las posiciones obtenidas por estos elementos no serán completamente correctas, y tendrán un error propio de haber despreciado los efectos periódicos. Como podemos ver en la Figura 2.7, la solución osculante (azul) presenta múltiples picos y cambios que se repiten con cierta periodicidad, pero la tendencia general se corresponde con la mostrada por la solución media (rojo).
20 Capítulo 2. Introducción a la mecánica orbital y otros conceptos Figura 2.6 Variaciones de los elementos orbitales debidos a la acción del J2 , en rojo los resultados para el modelo de J2medio y en azul la solución mediante integración numérica. 2.2.2 TLEs El “Two-Line Element set” consta de dos líneas que recogen la información básica sobre la identidad y los elementos orbitales de un objeto en órbita para una cierta época. Ambas líneas constan de 69 caracteres alfanuméricos, estando la primera destinada a recoger la identidad, los parámetros necesarios para modelos que incluyan perturbaciones y la época el objeto, y la segunda a los elementos orbitales [17]. Los 69 caracteres de la primera línea se distribuyen entre 14 campos, comenzando por un 1 que simplemente nos indica que estamos en la primera línea del formato. A continuación, les sigue el código identificador del objeto, su clasificación como objeto clasificado, secreto o cuya información es pública y los dos últimos dígitos del año de lanzamiento, el número de lanzamiento de ese año y la pieza a la que corresponde dentro de este. Posteriormente, se especifica la época con los datos de los dos últimos dígitos del año y el día y fracción del mismo al que corresponden. Finalmente se incluyen coeficientes de resistencia aerodinámica y derivadas, una serie de datos sobre el número de TLEs generados para ese cuerpo y un último dígito que nos avisa de posibles errores. Como ejemplo de una primera línea de un TLEs se expone el siguiente caso de un satélite GPS:
2.2 Two-Line Element set: Estructura y ejemplo 21 1 24876U 97035A 23027.22472234 .00000021 00000+0 00000+0 0 9992 Primero vemos el identificativo de la línea 1, el satélite está identificado con el código 24876, siendo pública su información como denota el “U” final de este campo (de “Unclassified”). El satélite fue lanzado en el 1997 y fue el 35ºde aquel año. Los dos últimos dígitos del año de la época son 23 y fue lanzado en el día 27 en la fracción 0.22472234 del mismo. La primera derivada es 0.00000021 y su segunda es nula, los términos asociados al coeficiente de resistencia también están en blanco. Para completar podemos ver que se han generado 999 TLEs para este satélite. A continuación, se aborda la segunda línea. En este caso comenzamos de nuevo con el indicador de la línea donde estamos, acto seguido, de nuevo, se nos muestra el código identificador del catálogo. Posteriormente comienza la sucesión de elementos orbitales (todos los ángulos en grados): inclinación, RAAN, excentricidad (se colocan directamente los decimales pues el 0. se da por hecho), argumento del perigeo, anomalía media, revoluciones por día y finalmente el número de dicha revolución en esa época. De nuevo, la línea se completa con un número que nos avisa de posibles errores. Como ejemplo se expone la segunda línea para el mismo satélite GPS anterior: 2 24876 55.5453 146.9403 0065760 52.8122 307.7679 2.00563955187142 Comenzamos con un indicador de que estamos en la fila 2 y luego el código identificador 24876 (obviamente el mismo que antes). Para este satélite la inclinación i es de 55.5453˝ , el RAAN Ω es de 146.9403˝ , la excentricidad e de 0.006576, el argumento de perigeo ω52.8122˝ , la anomalía media 307.7679˝ y realiza 2.00563955 revoluciones por día, el número de revoluciones en esa época es de 18714.
3 Propagador de órbitas Vistos ya una serie de conceptos básicos de la mecánica orbital, en el siguiente capítulo se va a abordar el desarrollo de la primera de las herramientas elaboradas para este trabajo: un propagador. Este se desarrollará para un ámbito de aplicación concreto, lo cual nos llevará a implementar unas ciertas perturbaciones en nuestro modelo. Además, se abordarán pruebas de validación de la herramienta para garantizar que funciona correctamente. 3.1 Consideraciones iniciales 3.1.1 Área de estudio del propagador Sean los satélites y objetos en órbita que vamos a tener que manipular, como parte de nuestro trabajo vamos a tener que estudiar su posición y velocidad en distintos instantes, siendo necesario para esto propagar la órbita a partir de un momento inicial dado. En ese sentido (y como se desarrollará en las próximas páginas) se han elaborado herramientas que, a través de una resolución numérica de la ecuación del movimiento, nos puedan aportar la información necesaria. Sin embargo, y como ya se mencionó en el Capítulo 2, existen fenómenos perturbadores que hacen el modelo basado en la ecuación del movimiento inexacto, siendo necesario incluir otros términos en el mismo que introduzcan estas perturbaciones. Es por ello que el modelo en que se base nuestro propagador se debe ajustar, atendiendo a las necesidades del tipo de órbita y el tipo de perturbación que experimentarán los cuerpos objeto de estudio, para poder mejorar la solución obtenida en última instancia. Idealmente nuestro modelo incluiría el efecto de todas las perturbaciones posibles, sin embargo esto aumentaría demasiado la complejidad del mismo, por lo que debemos acotar cuales de ellas conviene incluir. Como criterio para seleccionar dichas perturbaciones, vamos a tomar aquellas que sean de especial relevancia para los objetos contenidos en nuestra región de estudio: la órbita baja terrestre o LEO, y en ese sentido seleccionamos dos: • La perturbación asociada a lo no esfericidad de la Tierra, centrándonos concretamente en la más fuerte de todas (la vinculada al armónico J2 ), pues para los lapsos de tiempo que vamos a analizar es la que va a tener efectos más apreciables. •La resistencia que la atmósfera ejerce a los cuerpos en órbita junto a la Tierra. Como último punto recordar que los modelos usados para representar ambas perturbaciones tampoco son exactos, por lo que a largo plazo podrían dar lugar a la acumulación de un cierto error, pero para la precisión del trabajo y los periodos de tiempo que se van a manipular resultan suficiente. 23
24 Capítulo 3. Propagador de órbitas 3.1.2 Sistema de referencia A la hora de referir las coordenadas de posición y velocidad de un cuerpo, es necesario referirlas a un cierto sistema de referencia. En nuestro caso recurriremos al sistema de referencia geocéntrico ecuatorial, pero para llegar a este empezaremos por el sistema perifocal. En este último sistema, los ejes x e y definen un plano que contiene la órbita y donde el eje x apunta al perigeo, mientras el z se dirige hacia "arriba" en relación a la órbita. Por otra parte, de este es posible pasar al sistema geocéntrico en el cual el eje x se orienta hacia el primer punto de Aries, con el z en la dirección del polo norte del eje de rotación terrestre y el y completando el triedro a derechas. Ambos ejes se pueden apreciar en la Figura 3.1. Figura 3.1 Sistema de referencia geocéntrico (izquierda) y perifocal (derecha). A partir de los elementos orbitales, se pueden obtener las coordenadas de posición y velocidad en el sistema de referencia perifocal, y de este podemos pasar al geocéntrico mediante la siguiente matriz [9]: CR F“» – cos Ωcos ω´sin Ωsin ωcos isin Ωcos ω`cos Ωsin ωcos isin ωsin i ´cos Ωsin ω´sin Ωcos ωcos i´sin Ωsin ω`cos Ωcos ωcos icos ωsin i sin Ωsin i´cos Ωsin icos ifi fl(3.1) donde Ω es la ascensión recta del nodo ascendente, ω el argumento de perigeo e i la inclinación. 3.1.3 Conversión de vectores de posición y velocidad a elementos orbitales Una vez definido nuestro sistema de referencia se van a abordar los siguientes dos problemas: •Partiendo de la posición y velocidad del móvil obtener sus elementos orbitales. •Partiendo de sus elementos orbitales obtener su posición y velocidad. A continuación, se muestra el algoritmo para resolver el primero de ellos. primero determinaremos el momento angular específico ha partir de la posiciónry velocidadv: h“rˆv(3.2) Así como el vector de excentricidadea partir de los mismos datos:
3.1 Consideraciones iniciales 25 e“ rˆv´´µ r¯r µ(3.3) donde µ es la constante gravitacional para la Tierra y de cuyo módulo obtenemos la excentricidad e. Una vez obtenidos estos el semieje mayor ase puede calcular como: ε“v2 2´µ r(3.4) con εla energía específica de la órbita, y de modo que el semieje mayora termina siendo: a“´µ 2ε(3.5) La inclinación ise obtendrá de la ecuación: cos i“ h¨ kr h(3.6) Con kr el vector normal al plano ecuatorial. Para obtener finalmente tanto el RAAN como el argumento de perigeo, primero tendremos que sacar la línea de nodos que forman los planos orbital y ecuatorial: n“ krˆ kr | krˆ h|(3.7) Donde n es el vector que defina la dirección de la línea de nodos. Una vez con este vector el RAAN se obtiene como: cos Ω“nx(3.8) con nxla componente x del vector n. Y por último, el argumento de perigeo se obtiene de: cos ω“e¨n e(3.9) Un problema que se presenta a la hora de definir los ángulos del RAAN, argumento de perigeo y de anomalía verdadera es que el arcocoseno ofrece soluciones en el rango entre 0˝ y 180˝ , descartando por defecto otras posibles soluciones. Para corregir esto se aplican los siguientes criterios [9]: •θ θ θ : Si r¨vă0 entonces el ángulo de trayectoria se vuelve negativo. En consecuencia, nos encontramos viajando de apoapsis a periapsis, lo que supone que nuestra anomalía queda en el rango 180˝a360˝. •Ω Ω Ω : En este caso se usa el vector de la línea de nodos n que forma el plano orbital con el ecuatorial. Si su componente en el eje y ny es negativa nos encontramos a la “izquierda” del primer punto de Aries, por lo que el RAAN está entre 180˝y360˝. •ω ω ω : Para este último caso usamos como criterio la componente z del vector de excentricidad. Si es negativa, el perigeo está debajo del plano, lo que supone que el valor de ω queda entre 180˝y360˝. 3.1.4 Conversión de elementos orbitales a vectores de posición y velocidad Una vez vista la resolución de este primer problema, abordamos el problema en sentido opuesto, partiendo de los elementos orbitales obtener la posición y la velocidad. Para ello recurrimos primero a calcular el módulo del radio, empezando por obtener el parámetro de la cónica p:
32 Capítulo 3. Propagador de órbitas Figura 3.9 Variación de la excentricidad del objeto 2 durante un periodo de 1 día. Figura 3.10 Variación de la inclinación del objeto 2 durante un periodo de 1 día.
3.2 Propagador para la perturbación J233 Figura 3.11 Variación del RAAN del objeto 2 durante un periodo de 1 día. Figura 3.12 Variación del argumento de perigeo del objeto 2 durante un periodo de 1 día. De nuevo, nuestro propagador sigue la tendencia marcada por el J2 medio a lo largo del periodo analizado. Como detalle en estas comparativas del J2 medio con los elementos osculantes podemos ver como la media de los osculantes no acaba de coincidir con los medios. Mediante unas ecuaciones que eliminen los efectos periódicos podemos ser capaces de corregir estas ligeras diferencias de ambas soluciones.
34 Capítulo 3. Propagador de órbitas 3.3 Propagador para resistencia aerodinámica 3.3.1 Modelo del efecto de la resistencia atmosférica De manera análoga a lo expuesto en la sección anterior, para incorporar el efecto de la perturbación asociada a la resistencia atmosférica tendremos que incluir un término que la modele en la ecuación del movimiento. En cuanto a este término, ya expuesto en la Sección 2.1.4., se trata de [14]: D“1 2Bρv2 rel (3.28) Donde D será la resistencia experimentada por el objeto, ρ la densidad del aire y vrel se trata de la velocidad relativa del satélite respecto a la atmósfera y que se ha calculado según lo visto en la Sección 2.1.4. en las Ecuaciones 2.21. y 2.22. Por su parte, el modelo de densidad empleado se muestra en la Ecuación 2.23. Procede hacer una aclaración, en la expresión podemos identificar un término que contiene cuanto influye la atmósfera en el móvil, que se trata del coeficiente balístico B , no obstante, en los TLEs no se nos proporciona este, sino una variante, B˚ o BSTAR, que está asociado con el coeficiente balístico de la ecuación anterior, pero que, en caso de trabajar con TLEs, debemos convertir. Esta relación viene dada como [16][17]: B˚“ρ0B 2(3.29) Siendo esa densidad de referencia ρ0de valor: ρ0“0.157 kg m2¨RC (3.30) Es decir, kilos por metro cuadrado y radio terrestre. Para poder simplificar nuestro modelo, nosotros en lugar de B usaremos un coeficiente Cb , el cual definiremos simplemente como la inversa del balístico: Cb“1 B(3.31) Con la cual podremos sacar el coeficiente que finalmente introduciremos en nuestro modelo. De forma similar a como se hizo con la perturbación de J2 , este se incorporará en nuestro sistema de ecuaciones diferenciales como: : rpertAer “» — — — — — – ´1 2ρ¨Cb¨vrelx ´1 2ρ¨Cb¨vrely ´1 2ρ¨Cb¨vrelz fi ffi ffi ffi ffi ffi fl (3.32) Siendo vrelx,vrely,vrelz las componentes de la velocidad relativa en los ejes x y z . Como último detalle del funcionamiento del propagador, dado el modelo de densidad empleado se puede ver como conforme bajamos aumenta rápidamente la densidad y con ella la resistencia, haciendo poco fiable la solución de la ecuación del movimiento a baja altitud. Para corregir esto, se ha programado el propagador de manera que, en caso de bajar de una altitud de 122 km (6500 km del centro de la Tierra) se detendría la propagación al asumirse que ha reentrado en la atmósfera, evitando así las consecuencias de que se dispare la resistencia.
3.3 Propagador para resistencia aerodinámica 35 3.3.2 Validación del propagador Para poder comprobar lo adecuado del propagador desarrollado, se ha optado por hacer una comparativa entre los resultados obtenidos por el mismo y los resultados internos del departamento para unos ciertos satélites. El análisis se ha centrado en las variaciones del semieje mayor y de la excentricidad, pues son estos los parámetros afectados. Asimismo, se ha eliminado el efecto del armónico J2 a fin de poder apreciar mejor la tendencia. Finalmente, y como criterio, se va a comparar el tiempo que nuestro satélite permanece en órbita antes de reentrar con el tiempo obtenido en el departamento. Figura 3.13 Variación del semieje mayor del objeto NORAD 46950 en un período de 5 días. Figura 3.14 Variación de la excentricidad del objeto NORAD 46950 en un período de 5 días. En las Figuras 3.13. y 3.14. vemos el resultado obtenido para el propagador en el caso del objeto NORAD 46950. Sus elementos orbitales al inicio de la propagación son a“6630 km, e“0.0014 , i“97˝,Ω“92˝,ω“0˝,M“0˝y tiene Cb“5.67 ¨10´8m2/kg.
36 Capítulo 3. Propagador de órbitas En estimaciones internas del departamento el tiempo que se mantiene en órbita antes de bajar de 6500 kilómetros de distancia al centro de la Tierra se estima en unos 2.7 días, en el caso de nuestro propagador observamos 2.8 días. Figura 3.15 Variación del semieje mayor del objeto NORAD 46935 en un período de 3 días. Figura 3.16 Variación de la excentricidad del objeto NORAD 46935 en un período de 3 días. En las Figuras 3.15. y 3.16. vemos el resultado obtenido para el propagador en el caso del objeto NORAD 46935. Sus elementos orbitales al inicio de la propagación son a“6590 km, e“0.0015 , i“97.413˝,Ω“86.2025˝,ω“0˝,M“0˝y tiene Cb“3.504 ¨10´8m2/kg. En estimaciones internas del departamento el tiempo que se mantiene en órbita antes de bajar de 6500 kilómetros de distancia al centro de la Tierra se estima en unos 1.5 días, en el caso de nuestro propagador observamos 1.1 días.
3.3 Propagador para resistencia aerodinámica 37 Figura 3.17 Variación del semieje mayor del objeto NORAD 38270 en un período de 5 días. Figura 3.18 Variación de la excentricidad del objeto NORAD 38270 en un período de 5 días. En las Figuras 3.17. y 3.18. vemos el resultado obtenido para el propagador en el caso del objeto NORAD 38270. Sus elementos orbitales al inicio de la propagación son a“6640 km, e“0.0029 , i“97.33˝,Ω“235.4565˝,ω“100˝,M“´100˝y tiene Cb“9.23 ¨10´8m2/kg. En estimaciones internas del departamento el tiempo que se mantiene en órbita antes de bajar de 6500 kilómetros de distancia al centro de la Tierra se estima en unos 2.06 días, en el caso de nuestro propagador observamos 2.2 días.
4 Filtros de conjunciones Como ya se mencionó en las últimas secciones del Capítulo 2 de este trabajo, la labor de detectar las posibles colisiones entre dos objetos en órbita requeriría de estudiar las órbitas para cada pareja posible de cuerpos en el espacio. Obviamente, esto supone un gasto de recursos gigantesco desde un punto de vista computacional [4], algo aún más sangrante dado que en algunas de esas parejas, por las características propias de sus órbitas, dichas colisiones son imposibles sin necesidad de propagarlas para comprobarlo. Es por esto que se diseñan algoritmos que sean capaces de detectar que parejas si son susceptibles de generar una colisión, descartando buena parte de los pares de objetos en órbita por el camino, y reduciendo sensiblemente la cantidad de trabajo necesario a realizar. Estos algoritmos son los filtros. Así, y ya comentado el funcionamiento del propagador empleado para recrear las órbitas de los cuerpos que se estudien, a continuación se van a abordar los algoritmos empleados para discriminar posibles colisiones. 4.1 Estructura general del programa Antes de entrar a explicar el funcionamiento de dichos filtros y su validación, conviene comentar el programa en el que se van a incorporar. Este programa consta de dos bloques, en uno primero barrimos una base de datos de satélites dada (en nuestro caso el catálogo del sitio web SpaceTrack actualizado al 2022) obteniendo para cada uno tres datos: •Su código identificativo NORAD. •El tipo de órbita que traza (LEO, GEO, Molniya, ...). •Sus radios de apogeo y perigeo (necesarios para iniciar el filtrado). De cara a este filtrado resulta fundamental hacer este barrido de la base de datos y es para ello para que se usa una adaptación de un procesador de TLEs. En esta adaptación, lo que hacemos es entrar en la base de datos de Excel con un año que es en el que trabajamos (esta base de datos incluye el código NORAD de cada objeto, su fecha de lanzamiento, de caída si ya ha reentrado, . . .) aportada por SpaceTrack para, posteriormente, entrar en unos documentos de texto en formato txt los cuales incluyen TLEs para todos los objetos de ese año. El programa selecciona cada uno de los objetos de la base de datos, verifica que, por las fechas de lanzamiento y reentrada, se encuentra en órbita para el año analizado y, de estarlo, extrae sus TLEs y de ellos deduce sus elementos orbitales más su periodo juliano. 39
40 Capítulo 4. Filtros de conjunciones Figura 4.1 Diagrama conceptual con el funcionamiento del programa. En la Figura 4.1. podemos ver de forma esquemática como funciona este código. Como anotación mencionar que en una matriz se almacenan el identificador de cada satélite, el radio de apogeo Ra y el de perigeo Rp . Esta información se recoge para aquellos satélites en órbita que se encuentren activos y cumplan que sus valores de Ra y Rp son no nulos. Por último, las funciones que se recorren en la columna de la derecha se tratan de los filtros, que se recorren empezando por el de apogeo-perigeo, luego el Orbit Path Filter y finalmente el Time Filter. Una vez hecho esto, comenzamos el filtrado, donde pasaremos por tres algoritmos: • Un primer filtro que descarte parejas de objetos cuya colisión es imposible por las posiciones de su apogeo y perigeo. • Un segundo que nos descarte parejas donde los puntos de sus órbitas en la intersección de sus planos orbitales no coincidan. • Un tercero que, finalmente, nos descarte aquellos donde nunca pasan a la vez por las intersecciones de las órbitas. Para trabajar con los filtros, no obstante, será necesario volver a sacar los TLEs, es por ello que se ha hecho otra adaptación al código del programa clasificador de satélites, solo que en este caso únicamente nos extrae los elementos orbitales con los cuales trabajaremos los filtros.
4.2 Apogee-Perigee Filter 41 4.2 Apogee-Perigee Filter Como ya se ha comentado, un primer criterio de filtrado se basaría en hacer una comparativa entre los perigeos y apogeos de las órbitas de los distintos objetos abordados. La idea consistiría en que sean dos cuerpos, si el apogeo del primero queda más próximo a la Tierra que el perigeo del segundo entonces no puede existir ninguna posición donde ambas órbitas coincidan. Figura 4.2 Representación del funcionamienteo del filtro de apogeo-perigeo. A la izquierda tenemos una conjunción imposible, y a la derecha una posible. Como podemos apreciar en la Figura 4.2., en el caso de la izquierda tanto el radio de perigeo como el de apogeo para la órbita azul son menores que el perigeo para la órbita verde, luego la colisión es imposible. No ocurre así en el caso de la derecha, donde el perigeo de la órbita verde es menor que el de la azul y viceversa para los apogeos La primera idea que se nos podría ocurrir para hacer este filtrado, y encontrar una primera serie de parejas cuyas órbitas sí podrían cortarse según este criterio, sería comparar cada elemento con todos los demás y así para toda nuestra población. Sin embargo, por lo mismo que se comentó en la sección anterior, esto resultaría muy costoso computacionalmente. Es por ello que se opta por un algoritmo que funcionaría de la siguiente manera: • Primero entramos con una matriz de dos columnas, en la primera tenemos los radios de perigeo y apogeo de los distintos objetos, en la segunda los códigos identificativos del satélite al que le pertenece esa posición. •Luego ordenamos las filas de esa matriz según los radios. • Finalmente, vamos recorriendo la matriz de la siguiente forma: Cada vez que pasamos por un radio de perigeo (es decir, cada vez que vemos un código NORAD por primera vez) activamos un puntero. Por cada radio de perigeo distinto por el que pasamos creamos una pareja entre los NORAD de los punteros activos.
48 Capítulo 4. Filtros de conjunciones derecha, donde los dos cuerpos llegarán en instantes de tiempo similares Por tanto, y para completar nuestro filtrado, nos quedará un pequeño trabajo de propagación de la órbita, de manera que, dada la fecha juliana de cada uno de los satélites, podamos ver en qué instantes pasan por este punto. Si en ninguno de esos instantes están los dos a la vez en él, podemos descartar la colisión. Figura 4.7 Diagrama conceptual con el funcionamiento del filtro Time Filter. De nuevo, y como en el caso del filtro anterior, es conveniente resaltar que este filtro tampoco se ha diseñado para incluir el caso de órbitas coplanarias. Estas requerirían un trato específico, el cual no se ha considerado, para poder identificar posibles conjunciones al ser distinto el modo en el que intersecan las órbitas. Una vez explicado cómo funciona nuestro filtro, queda probarlo y comprobar que actúa correctamente. En ese sentido, se ha decidido continuar con el ejemplo ya ilustrado en la Sección 4.3. basado en dos satélites que se definían de forma que pasaran por el mismo punto en un mismo instante de tiempo. Ambos van a seguir dados por la misma serie de elementos orbitales expuestos en la Tabla 4.3., solo que incluiremos ahora como dato de entrada la fecha juliana (que se define
4.4 Time Filter 49 como el número de días desde el mediodía del 1 de enero del 4713 a.C.) de los elementos orbitales, omitida en el Orbit Path Filter por no ser necesaria. Para garantizar que se dé la coincidencia se ha dado el mismo valor de esta fecha para la serie de elementos orbitales que definen el punto de intersección, este valor se ha seleccionado de forma arbitraria, pudiendo haber valido cualquier otro, y se muestra a continuación: Tabla 4.9 Fechas Julianas característicos de los objetos de prueba empleados, identificados como 1 y 2. Objeto JD 1 100 2 100 Así, garantizamos que en la primera pasada por el punto de intersección planteado para nuestra prueba se produzca la coincidencia. Efectivamente, si lanzamos la simulación el filtro detecta que nuestros objetos de prueba cruzan por el punto de coordenadas medidas según el sistema de referencia geocéntrico ecuatorial: r“» – 6700 0 0fi flkm (4.6) en el mismo instante definido como JD “100 . Efectivamente la respuesta del programa es un vector con los identificadores de nuestra pareja de objetos, demostrando que ha sido capaz de captar la posible colisión. Con esta primera comprobación se ha podido verificar que el filtro es capaz de detectar objetos que coinciden en un instante, sin embargo también nos interesa comprobar si en futuros pasos por los puntos de intersección entre las órbitas se producirían potenciales coincidencias. Es por ello que realizamos la siguiente comprobación: • Comprobamos los instantes de paso por el punto de intersección r para ambos objetos y seleccionamos un instante de paso por este para cada uno. En nuestro caso tomamos el instante del cuarto paso del objeto 1 y del tercero del objeto 2, dados por las fechas julianas siguientes: JD1“100.0964354 (4.7) JD2“100.2131020 (4.8) Por tanto, el satélite 1 pasa por él en su cuarto barrido antes de que el satélite 2 haga su tercer paso. • Una vez hecho esto, calculamos la diferencia de tiempo entre ambos pasos por dicho punto y se lo añadimos al satélite 1: JD2´JD1“0.1166666 (4.9) • Finalmente, con la fecha juliana del primer objeto recalculada volvemos a lanzar la prueba, obteniendo el resultado esperado que será el mismo vector de detección de una conjunción entre ambos objetos, solo que ahora lo hace en un instante de tiempo posterior.
5 Puesta en funcionamiento: Análisis y resultados Terminada la presentación y explicación sobre el funcionamiento de las herramientas elaboradas con objeto de este trabajo, se van a poner a prueba con la ejecución de un análisis de conjunciones sobre una población seleccionada a tal efecto. En las próximas páginas se irán detallando tanto la muestra tomada, como los resultados obtenidos como una serie de conclusiones que se pueden extraer sobre el funcionamiento de cada uno de los filtros. 5.1 Alcance del análisis y población seleccionada El primer paso antes de iniciar las tareas de filtrado y propagación, que se desarrollarán en las secciones posteriores, ha consistido en la selección de una muestra de satélites con la cual trabajar. La base para seleccionar dicha población ha sido el ya mencionado anteriormente (al inicio de la Sección 4.1.) catálogo que el sitio web SpaceTrack pone a nuestra disposición. Idealmente se hubiera buscado desarrollar el análisis a la base de datos entera. Sin embargo, limitaciones técnicas respecto de la capacidad de los medios empleados para procesar han obligado a adoptar un enfoque más "modesto". En ese sentido, se ha hecho un primer "filtrado" de la base de datos, excluyendo la mayor parte de la misma para centrarnos en una porción del espacio que por sus características fuera representativa de los problemas futuros de sobrepoblación y aumento del riesgo de colisiones entre objetos. Con esta filosofía en mente, se ha buscado acotar el análisis a una órbita LEO, concretamente a una franja de la misma, por ser una región muy expuesta a la "sobrepoblación" del espacio y proporcionar ya de por sí una variedad importante de satélites, junto a la posibilidad de poder ver el efecto de las perturbaciones asociadas al armónico J2 y la resistencia atmosférica. Se han valorado dos candidatas como población de estudio: •Población 1: Objetos cuyas órbitas cumplen que su altitud de perigeo es superior a 600 km y su altitud de apogeo es inferior a 750 km. •Población 2: Versión reducida de la anterior, solo que en este caso el apogeo se encontraría a una altitud inferior de 700 km. Por cuestiones de tiempo se ha optado por este segundo caso que, tras barrer la base de datos, nos ha dado una población de 801 satélites. Esta selección presenta las siguientes limitaciones: • Se excluyen satélites que rebasen o se desplacen en los márgenes de los límites impuestos pero que podrían ser susceptibles de cruzar órbitas de otros cuerpos si incluidos en nuestra región. 51
52 Capítulo 5. Puesta en funcionamiento: Análisis y resultados • Se excluyen satélites que presenten elevadas excentricidades y tenga partes de sus órbitas en la región estudiada, aun estando su perigeo y/o apogeo fuera de la misma. Si bien esto va a dar lugar a que el resultado final sea incompleto, la muestra tomada sigue siendo válida para poner en funcionamiento todas las herramientas elaboradas, además de poder tener una idea de cómo es el estado de la región objeto de análisis. 5.2 Aplicación del filtro de apogeo-perigeo: Resultados El primer paso realizado ha consistido en la aplicación del filtro de apogeo-perigeo ya comentado en el capítulo anterior (Sección 4.2.). Una vez extraída la población de la base de datos según los criterios comentados anteriormente, se han obtenido los valores de sus altitudes de apogeo y de perigeo con las cuales se ha construido una matriz válida como entrada para el programa del filtro usado. Sea ese conjunto de 801 satélites con sus perigeos y apogeos, el total de parejas únicas posibles se puede calcular como: NPare jas “n¨pn´1q 2(5.1) Siendo n el número total de elementos en nuestro conjunto, para nuestro caso 801. Aplicando la Ecuación 5.1. resultan un total de 320400 parejas posibles, sin embargo, y tras pasar el conjunto por el filtro, vemos como de todas esas parejas las que pasan el corte y por lo tanto son realmente susceptibles de tener una conjunción son 151386, un 47.25 % del conjunto original. También resulta interesante ver los tiempos de ejecución del filtro para distintos tamaños de población. Como era de esperar, y al igual que el número de parejas posibles crece de forma exponencial con el número de satélites, también lo hace el tiempo de ejecución. A continuación, se muestra el tiempo necesario para procesar distintas cantidades de objetos: Tabla 5.1 Tiempos de ejecución para distintos tamaños de las poblaciones. NºObjetos Tiempo 10 0.000304 segundos 50 0.004808 segundos 100 0.174448 segundos 300 29.203408 segundos 500 3.967 minutos 801 26.445 minutos Como último apunte, con estos tiempos delante podemos ver como tampoco es uniforme el tiempo de procesado por pareja, aumentando junto al tamaño de la muestra usada. Esto se debe a que el tamaño de la población implica un aumento del número de candidatos con los que el filtro tendrá que cotejar posibles parejas, dando lugar a mayores tiempos de procesado por unidad. Una vez el filtro termina de ejecutarse, obtendremos como resultado la matriz de parejas susceptibles de tener una conjunción entre sus órbitas, identificadas estas por los NORAD de sus miembros. Es con esta información con la que ingresaremos al siguiente filtro. 5.3 Aplicación del filtro Orbit Path Filter: Resultados Una vez terminamos de ejecutar el filtro anterior, hemos obtenido una matriz con los NORAD de las parejas susceptibles de tener conjunciones entre sus órbitas, pero en muchas de estas realmente
5.3 Aplicación del filtro Orbit Path Filter: Resultados 53 no las llegará a haber, pues, como se comentó en la Sección 4.3., en muchos casos las posiciones de dichas órbitas sobre la línea de nodos relativa a ambas no coinciden. Es aquí donde entra en juego el Orbit Path Filter. Sin embargo, exigir una coincidencia absoluta de dichas posiciones resulta muy riguroso, principalmente por dos razones: • Porque los modelos usados para definir la órbita no son completamente exactos, al no modelar completamente el efecto de las perturbaciones sobre las órbitas. • Porque los métodos para ir localizando y trazando la órbita del satélite cuentan con un cierto error, existiendo una diferencia entre la posición real del satélite y la que podemos medir. Es por ello que nuestro filtro no se programa de manera que filtre todas las parejas que no tengan una coincidencia total en dichos puntos, sino que se ha programado para que permita una cierta tolerancia considerada suficiente como para que siga habiendo riesgo de conjunción. A la hora de definir esa tolerancia, se ha usado como referencia los órdenes de magnitud en los errores en las estimaciones y medidas. Así, el criterio usado ha sido el siguiente: ∆ră“5km (5.2) De esta forma, si sobre la línea de nodos relativa a ambas órbitas el módulo de la distancia entre los puntos por los cuales pasan uno y otro satélite es menor o igual a 5 kilómetros, se asumirá que existe riesgo de conjunción. Con este criterio, se ha procedido al barrido de las 151386 parejas que superaron el filtro anterior. Aplicando el criterio comentado anteriormente, el conjunto de pares susceptibles de tener una conjunción se reduce significativamente hasta 57058 parejas, lo que supone un 37.69% del total original. Como podemos ver, el filtrado ha sido más riguroso que en el caso anterior, aunque evidentemente se podría haber hecho más laxo o más duro modificando el criterio de filtrado convenientemente. Una vez comentado el resultado obtenido de la ejecución del filtro, se van a realizar unos apuntes sobre el rendimiento y el funcionamiento del mismo. En este caso se observa cómo, a diferencia del filtro anterior, el incremento de tiempo destinado a procesar una cierta población no aumenta de forma exponencial, pues en este caso el número de operaciones a realizar por el programa viene a ser el mismo para cada pareja, aumentando de forma lineal, como podemos ver en la siguiente tabla: Tabla 5.2 Tiempos de ejecución para distintos tamaños de las poblaciones. NºParejas Tiempo 10 46.56 segundos 50 3.46 minutos 100 6.81 minutos 300 20.24 minutos 500 33.66 minutos A través de esta tabla, podemos además concluir que el tiempo que nuestro filtro dedica a cada pareja es de unos 4 segundos. Dado que el número original de parejas obtenido del filtro de apogeo perigeo fue 151386, resulta rápido obtener las horas que han sido necesarias para poder completar el filtrado: 168.2 horas, o lo que es lo mismo, un poco más de una semana de tiempo de procesado. Como último detalle a comentar, si calculamos el tiempo dedicado a cada pareja para cada uno de los casos mostrados en la Tabla 5.2., vemos como disminuye conforme aumenta la cantidad hasta converger en los 4 segundos ya mencionados. Esto se debe a que: • Por un lado, cuantas más veces se repite una cierta operación, menos tiempo tarda el ordenador en completarla.
54 Capítulo 5. Puesta en funcionamiento: Análisis y resultados • Por otro, al principio del programa y antes de iniciar el barrido de las parejas, se lleva a cabo una lectura de todos los TLEs asociados al año analizado. Esta operación requiere un cierto tiempo, que distorsiona el resultado en poblaciones más reducidas pero que se vuelve despreciable cuando esta alcanza una cierta magnitud. Finalmente, una vez se complete la ejecución del filtro, este nos devolverá tres matrices: •Una matriz con las parejas que han pasado el filtro. • Dos matrices más con las anomalías verdaderas de cada órbita asociadas a sus puntos sobre la línea de nodos relativa (una matriz para cada uno de los dos nodos). Además, devolverá un vector que, para cada una de las parejas que pasen el filtro, nos dirá en cuál de las dos mitades de la línea de nodos relativa ocurre la conjunción o si ocurre en ambas. Con estas 3 matrices y ese vector es con lo que entraremos al tercer y último filtro, el Time Filter. 5.4 Aplicación del filtro Time Filter: Resultados Ya pasados los dos filtros anteriores hemos llegado a una lista de parejas que verifican tener órbitas con puntos suficientemente próximos entre sí como para considerar posible una conjunción, sin embargo, esto no implica que necesariamente deba de haberla (como se comentó en la Sección 4.4.). Si los satélites no llegan a coincidir en ningún instante sobre estos puntos considerados de gran proximidad, no se podrá producir la conjunción ni una eventual colisión. Para poder verificar si se produce esta coincidencia temporal recurriremos al Time Filter. Sin embargo, y por motivos similares a los comentados anteriormente respecto de las tolerancias y márgenes de error en los medios para medir y predecir la trayectoria del satélite, no exigiremos en este filtro una coincidencia total en el tiempo, si no que se permitirá una cierta tolerancia. En este caso, no obstante, en lugar de usar como referencia para la tolerancia solo el error de la medición, se va a combinar este con las velocidades típicas de los satélites en órbita. De esta manera el criterio resultante ha sido: ∆tă“1segundo (5.3) De modo que, si la diferencia entre los instantes de paso de los cuerpos estudiados por alguna de las intersecciones con la línea de nodos relativa baja en algún momento de 1 segundo, se considerará que si existe conjunción para dicha pareja. Una vez definido el criterio según el cual filtrar, se ha llevado a cabo finalmente el barrido del conjunto de parejas obtenido tras la ejecución del segundo filtro. De las 57058 parejas iniciales, el resultado obtenido ha sido de 1783 parejas que han cumplido nuestro segundo requisito, es decir, que cruzan la posición de la conjunción con un segundo o menos de diferencia. Esto supone "apenas" un 3.12% del grupo original, una reducción mucho mayor que la vista en los dos filtros anteriores y que da cuenta de lo exigente de este en comparación. Ya vistos los resultados obtenidos con este filtro, y al igual que se ha hecho en los casos anteriores, se va a proceder a comentar como es la ejecución del mismo. Dada la estructura tan similar que posee el código del Time Filter con respecto al Orbit Path Filter, no es sorprendente que la relación de aumento del tiempo de ejecución del filtro con el número de parejas sea igualmente lineal. Esto se puede comprobar con los siguientes tiempos obtenidos tras pasar el filtro a muestras de distintos tamaños:
5.5 Propagación de un caso de posible conjunción 55 Tabla 5.3 Tiempos de ejecución para distintos números de parejas. NºParejas Tiempo 10 46.51 segundos 50 3.46 minutos 100 6.83 minutos 300 20.54 minutos 500 33.72 minutos Con estos resultados, resulta sencillo deducir que el programa tarda unos 4.1 segundos en procesar cada pareja. Resulta, además, notorio como no solo la relación es igualmente lineal como en el caso del Orbit Time Filter, sino que los tiempos de ejecución son muy similares (compárese la tabla anterior con la Tabla 5.2.). Esto apunta a que los algoritmos destinados a verificar el criterio realmente no suponen la mayor parte del tiempo de ejecución, sino que la mayor parte de los 4.1 segundos que se tarda en procesar cada pareja en este filtro (y de los 4 segundos en el caso del Orbit Path Filter) se dedican a ejecutar la función que nos extrae los elementos orbitales necesarios para los cálculos posteriores (véase la función común a ambos filtros donde se extraen los elementos orbitales en las Figuras 4.5. y 4.7.) ya mencionada al final de la Sección 4.1. Por último, una vez obtenidos los segundos que tardamos en cada pareja, y sabiendo el total de parejas, podemos determinar que, para que el filtro pueda procesar la totalidad de parejas que han superado el Orbit Path Filter, el tiempo total necesario será de 63.4 horas o, lo que es lo mismo, 2.64 días. Como resultado de este filtro, obtendremos aquellas parejas que, finalmente, han demostrado ser susceptibles de tener una conjunción entre sus órbitas ateniéndonos a los criterios que hemos definido. Ya después de esto, habríamos completado la tarea de detectar posibles conjunciones, con una lista de casos de riesgo, tocando abordar caso a caso cual será el riesgo real de una colisión y, de haberlo, diseñar una maniobra de evasión. 5.5 Propagación de un caso de posible conjunción Tras la ejecución de los filtros mostrada en las secciones anteriores, se ha llegado finalmente a un conjunto de pares de satélites que, desde el punto de vista de los requisitos que hemos impuesto, son susceptibles de presentar una conjunción. Es por ello que, para cerrar este capítulo, se ha decidido hacer una comprobación más exhaustiva de una posible conjnunción para uno de estos pares de satélites, ayudándonos para ello del propagador comentado en el Capítulo 3. Esta propagación se hará en dos casos, primero para una comprobación inicial en un escenario sin perturbaciones y luego con los efectos de la resistencia atmosférica y de la perturbación asociada al armónico J2 incluidas, para hacer un estudio más preciso y dar cuenta del error que se estaría cometiendo entre ambas situaciones para un intervalo de tiempo breve. Con esto en mente, se ha entrado en el vector de posibles conjunciones obtenido del Time Filter, y se ha extraído la pareja dada por los códigos NORAD: Tabla 5.4 Identificadores NORAD de la pareja analizada. 40055 43262 Correspondiéndose el primero con el satélite canadiense NLS 7.1/CANX 4 y el segundo con el chino GAOFEN 1-04. Para esta pareja se van a analizar sus movimientos entorno al 31 de diciembre del 2021, siendo sus elementos orbitales para el final de este día los siguientes:
56 Capítulo 5. Puesta en funcionamiento: Análisis y resultados Tabla 5.5 Elementos orbitales característicos de los satélites identificados por los NORAD 40055 y 43262. Satélite a[km] e ir˝sωr˝sΩr˝sθr˝s 40055 7029.04 0.0017 98.33 124.20 152.64 -124.20 43262 7024.70 0.00097 97.40 84.34 80.95 -84.34 Como último apunte, estos TLEs están vinculados a unos instantes específicos. El TLE asociado al satélite 40055 se corresponde con el 31 de diciembre de 2021 a las 21:08:01, mientras que el del satélite 43262 se corresponde con el 31 de diciembre de 2021 a las 21:11:13. Analizando las dos órbitas, podemos identificar dos posibles puntos donde ocurra la conjunción, ateniéndonos a lo comentado en la Sección 4.3. centrada en el Orbit Path Filter, que se corresponden con sus pasos por la línea de nodos relativa. Comprobando los instantes en los que los barrerían, encontramos los siguientes instantes de interés en los que cruzarían por dicha recta, en la mitad asociada al sentido opuesto al vector de dirección de la línea de nodos relativa Nod , con la suficientemente poca diferencia de tiempo: Tabla 5.6 Primeros instantes de paso de interés de los objetos estudiados. Satélite Instante de paso 40055 2021/12/31 21:33:56.74 43262 2021/12/31 21:33:56.55 que se corresponden con una diferencia de tiempo ∆t de 0.2 segundos, cumpliendo el requisito del Time Filter. Aplicando el propagador para el objeto 40055 desde la fecha con la que se corresponde cada TLE hasta el instante en el que pasa por la intersección, obtenemos que cruza la siguiente posición al no aplicar perturbaciones: r“» – 1066.45 587.13 6912.29fi flkm (5.4) y con estas aplicadas cruza por: r“» – 1067.62 587.17 6905.40fi flkm (5.5) Por su parte, propagando el satélite 43262 en su periodo correspondiente obtenemos que, sin perturbaciones, pasa por la intersección en: r“» – 1066.33 587.08 6911.53fi flkm (5.6) y con ellas: r“» – 1066.34 585.47 6905.63fi flkm (5.7)
5.5 Propagación de un caso de posible conjunción 57 Si para estas posiciones ahora calculamos la diferencia ∆r , obtendremos finalmente como de cerca llegan a pasar: Tabla 5.7 Distancias entre los puntos de paso de los satélites analizados por la intersección. ∆r[km] Sin perturbación 0.77 Con perturbación 2.14 Como se puede comprobar en la Tabla 5.7., en el caso sin perturbaciones pasan muy cerca (a menos de 1 km de distancia), lo cual teniendo en cuenta las tolerancias de medición de los TLEs hace aún más peligrosa la situación. No obstante, al no haber computado las perturbaciones cometemos un error y dado que se trata de un potencial caso de riesgo se hace recomendable hacer un estudio más preciso. Así, al incluirlas vemos como la distancia sigue siendo muy reducida (aunque mayor que antes), por lo tanto tendríamos una posible colisión entre ambos objetos en esta posición que requeriría entrar en el uso de herramientas como el Volumen de Exclusión o la Probabilidad de Colisión (comentadas en la Sección 1.2.1.) para acabar de determinar si el riesgo es lo suficientemente alto como para lanzar una maniobra de evasión. Como se mencionó antes, además de este punto donde ya hemos verificado que existe una posible conjunción, habría además una segunda posición de la línea de nodos relativa en la cual, por la escasa diferencia de tiempo entre los tránsitos de los satélites, podría haber riesgo, y que se corresponde con el tránsito por la línea de nodos relativa por la mitad en la que apunta el vector Nod . En este caso los instantes de tiempo se corresponden con: Tabla 5.8 Segundos instantes de paso de interés de los objetos estudiados. Satélite Instante de paso 40055 2021/12/31 22:22:46.14 43262 2021/12/31 22:22:46.20 En este caso la diferencia de tiempo ∆t es aún menor, de 0.06 segundos, de nuevo se cumple el requisito usado en el Time Filter. Siguiendo el procedimiento usado para analizar el posible punto de conjunción anterior, volvemos a calcular los puntos de paso del satélite 40055 por la intersección, que sin perturbación será: r“» – ´1068.60 ´588.87 ´6932.70fi flkm (5.8) y con perturbación: r“» – ´1095.88 ´576.05 ´6923.65fi flkm (5.9) Como vemos, en este caso la diferencia es significativa. Repitiendo el proceso para el satélite 43262, obtenemos sin perturbaciones: r“» – ´1068.39 ´588.21 ´6924.89fi flkm (5.10)
64 Índice de Figuras 3.14 Variación de la excentricidad del objeto NORAD 46950 en un período de 5 días 35 3.15 Variación del semieje mayor del objeto NORAD 46935 en un período de 3 días 36 3.16 Variación de la excentricidad del objeto NORAD 46935 en un período de 3 días 36 3.17 Variación del semieje mayor del objeto NORAD 38270 en un período de 5 días 37 3.18 Variación de la excentricidad del objeto NORAD 38270 en un período de 5 días 37 4.1 Diagrama conceptual con el funcionamiento del programa 40 4.2 Representación del funcionamienteo del filtro de apogeo-perigeo. A la izquierda tenemos una conjunción imposible, y a la derecha una posible 41 4.3 Diagrama conceptual con el funcionamiento del filtro de apogeo-perigeo programado 42 4.4 Representación de la intersección entre dos órbitas con la correspondiente línea de nodos relativa a ambas. A la izquierda tenemos una conjunción imposible, y a la derecha una posible 44 4.5 Diagrama conceptual del funcionamiento del filtro Orbit Path Filter 45 4.6 Representación de dos satélites en situación de potencial colisión. A la izquierda tenemos una conjunción imposible, y a la derecha una posible 47 4.7 Diagrama conceptual con el funcionamiento del filtro Time Filter 48
Índice de Tablas 1.1 Lanzamientos anuales exitosos desde el inicio de la década pasada hasta el año 2022 [1] 1 4.1 Objetos de prueba junto con sus radios de perigeo y apogeo para la prueba realizada. Estos se miden en unidades de distancia genéricas, que se han denominado UD 43 4.2 Parejas de objetos susceptibles de tener conjunciones en sus órbitas ateniéndonos al filtro de apogeo-perigeo 43 4.3 Elementos orbitales característicos de los objetos de prueba empleados, identificados como 1 y 2 46 4.4 Distancias al origen del sistema de referencia para los satélites de prueba 1 (denotado por r1 ) y 2 (denotado por r2 ) en su tránsito por la mitad de la línea de nodos relativa en el sentido del vector Nod 46 4.5 Distancias al origen del sistema de referencia para los satélites de prueba 1 (denotado por r1 ) y 2 (denotado por r2 ) en su tránsito por la mitad de la línea de nodos relativa en el sentido opuesto del vector Nod 46 4.6 Elementos orbitales característicos del objeto de prueba identificado como 3 46 4.7 Distancias al origen del sistema de referencia para los satélites de prueba 1 (denotado por r1 ) y 3 (denotado por r3 ) en su tránsito por la mitad de la línea de nodos relativa en el sentido del vector Nod 47 4.8 Distancias al origen del sistema de referencia para los satélites de prueba 1 (denotado por r1 ) y 3 (denotado por r3 ) en su tránsito por la mitad de la línea de nodos relativa en el sentido opuesto del vector Nod 47 4.9 Fechas Julianas característicos de los objetos de prueba empleados, identificados como 1 y 2 49 5.1 Tiempos de ejecución para distintos tamaños de las poblaciones 52 5.2 Tiempos de ejecución para distintos tamaños de las poblaciones 53 5.3 Tiempos de ejecución para distintos números de parejas 55 5.4 Identificadores NORAD de la pareja analizada 55 5.5 Elementos orbitales característicos de los satélites identificados por los NORAD 40055 y 43262 56 5.6 Primeros instantes de paso de interés de los objetos estudiados 56 5.7 Distancias entre los puntos de paso de los satélites analizados por la intersección 57 5.8 Segundos instantes de paso de interés de los objetos estudiados 57 5.9 Distancias entre los puntos de paso de los satélites analizados por la intersección 58 65
Bibliografía [1] Statista. Number of orbital space launches worldwide from 1957 to 2022. Statista, 2022. Extraido el 22/6/2023 de https://www.statista.com/statistics/1343344/orbital-space-launches- global/. [2] Edouard Mathieu and Max Roser. Space exploration and satellites. Our World in Data, 2022. Extraído el 23/6/2023 de https://ourworldindata.org/space-exploration-satellites. [3] NASA ARES division. Legend : 3d/od evolutionary model. NASA ORBITAL DEBRIS PROGRAM OFFICE, 2022. [4] Emma Kerr and N Sanchez. State of the art and future needs in conjunction analysis methods, processes and software. In Proceedings of 8th European Conference on Space Debris (virtual edition, SDC8), pages 20–23, 2021. [5] Alinda K Mashiku and Matthew D Hejduk. Recommended methods for setting mission conjunction analysis hard body radii. In 2019 AAS/AIAA Astrodynamics Specialist Conference, number GSFC-E-DAA-TN71115-1, 2019. [6] Matthew D Hejduk, DE Snow, and LK Newman. Satellite conjunction assessment risk analysis for “dilution region” events: issues and operational approaches. In Space Traffic Management Conference, 2019. [7] Javier Hernando Ayuso. Collision Avoidance Maneuver Optimization: A Fast and Accurate Semi-analytical Approach. PhD thesis, Aeronauticos, 2014. [8] A Lapusneanu. Use of thrusters for collision avoidance maneuver. StackExchange, 2022. extraído el 26/2/2023 de https://space.stackexchange.com/questions/33870/use-of-thrusters- for-collision-avoidance-maneuver. [9] R Vázquez Valenzuela. Apuntes de clase. extraídos de http://aero.us.es/astro/desc.html, 2022. [10] Elbaz I Abouelmagd, Daniele Mortari, and Hadia H Selim. Analytical study of periodic solutions on perturbed equatorial two-body problem. International Journal of Bifurcation and Chaos, 25(14):1540040, 2015. [11] Yuchul Shin, Seyoung Yoon, Yongmyung Seo, Ho Jin, and Jongho Seon. Radiation effect for a cubesat in slow transition from the earth to the moon. Advances in Space Research, 55(7):1792–1798, 2015. 67
68 Bibliografía [12] E Sinem Ince, Franz Barthelmes, Sven Reißland, Kirsten Elger, Christoph Förste, Frank Flechtner, and Harald Schuh. Icgem–15 years of successful collection and distribution of global gravitational models, associated services, and future plans. Earth system science data, 11(2):647–674, 2019. [13] David A Vallado. Fundamentals of astrodynamics and applications, volume 12. Springer Science & Business Media, 2001. [14] Delores Knipp and Yihua Zheng. Space weather in the thermosphere: Satellite drag. Space Weather Training for Mission Operators and Engineers, University of Colorado Boulder, https://ccmc.gsfc.nasa.gov/RoRWWW{SWREDI{2015{SatDragYZheng060415.pd f ,2014. [15] R Fitzpatrick. Effect of atmospheric drag on artificial satellite orbits. University of Texas at Austin, https://farside.ph.utexas.edu/teaching/celestial/Celestial/node94.html, 2016. [16] Felix R Hoots. Spacetrack report no. 3, models for propagation of norad element sets. http://www. itc. nl/-bakker/orbit. html, 1980. [17] TS Kelso. “frequently asked questions: Two-line element set format. See http://CelesTrak. com/columns/v04n03/, 2004.