Космическая наука с Python
Космическая наука с Python — бурные времена кометы
В 13-й части учебной серии будут рассмотрены элементы орбиты 67P/Чурюмов–Герасименко. Это простая орбитальная механика. Правильно?
Предисловие
Это 13-я часть моей серии руководств по Python Космические науки с Python. Все коды, которые показаны здесь, загружены на GitHub. Наслаждайтесь!
Пользователи предлагали добавлять небольшие задачи в конце каждого блокнота Jupyter Notebook для решения научных задач или вопросов по программированию. Начиная с этой сессии, я добавлю эти задания и предоставлю решение через некоторое время на GitHub. Удачи!
Введение
На последнем занятии (Космическая наука с Python — комета в 3D) мы начали рассматривать очень конкретную комету: 67P/Чурюмова–Герасименко. Комета, которую посетила миссия ЕКА Rosetta/Philae. С 2014 по 2016 год орбитальный аппарат исследовал этот ледяной мир и развернул свой посадочный модуль Philae в ноябре 2014 года. Через несколько часов после сближения с кометой человечество собрало больше данных и информации о кометах, чем все предыдущие миссии и наблюдения вместе взятые. В 2004 г. была опубликована книга под названием Comets 2, в которой собрана вся необходимая информация об этих испаряющихся мирах [1]. Теперь у нас есть еще несколько научных открытий, которые частично находятся в свободном доступе в Интернете: Специальные выпуски по астрономии и астрофизике.
Наш последний сеанс дал нам подробный обзор трехмерной структуры кометы, полученной с помощью одной из систем камер Rosetta. Сегодня и на следующем занятии мы хотим узнать больше о его траектории; его движение вокруг Солнца. Как комета вращается вокруг Солнца? Существуют ли различные информационные ресурсы? И какие выводы мы можем из этого извлечь?
Нечего ожидать, верно?
Учитывая простые задачи с двумя телами, орбитальную механику довольно просто понять. Две массы вращаются вокруг общего центра масс. На предыдущем занятии (Космическая наука с Python: кометы — гости издалека) мы создали базу данных SQLite, содержащую орбитальные элементы всех известных комет. Данные получены из Центра малых планет (ЦМП). Эти данные могут быть использованы для расчета траектории кометы. Но так ли это просто?
Во-первых, мы импортируем наши необходимые модули. На этот раз нам понадобится Инструментарий SPICE от NASA и оболочка Python spiceypy. Также используются другие модули, такие как numpy, pandas и matplotlib.
Мы устанавливаем соединение с нашей базой данных в строке 2. Из-за небольшого размера базы данных, около 180 КБ, база данных также хранится в моем репозитории GitHub, поэтому вам не нужно беспокоиться о создании этой базы данных с нуля.
Какие данные нам нужны для описания траектории 67P? Очевидно, нам нужны элементы орбиты, такие как перигелий (линия 5), большая полуось (линия 6), афелий и эксцентриситет (строка 7), а также аргумент перигелия и наклонения (строка 8). Названия комет начинаются с аббревиатуры номера и типа (C, P, X и т. д.), за которыми следует фактическое название. Чтобы получить данные о комете, указав только номер, мы используем команду SQL NAME LIKE 67P%. Оператор where ищет имена, строки которых похожи на 67P с присоединенной любой (под) строкой (%). Данные сохраняются в кадре данных pandas с помощью команды read_sql.
Давайте посмотрим на данные из базы данных. Мы просто печатаем результаты с помощью команды iloc внутри f-строки.
Результат показан ниже. Мы успешно извлекли информацию о комете 67P. Перигелий находится на расстоянии 1,2 а.е., то есть в пределах диапазона, в котором нагрев Солнца заставляет кометы испаряться [2]. Афелий находится на расстоянии 5,7 а.е.; а небольшой наклон около 4 градусов указывает на то, что комета пересекает плоскость орбиты Юпитера, что делает ее кометой семейства Юпитера (JFC).
NAME 67P/Churyumov-Gerasimenko PERIHELION_AU 1.20213 SEMI_MAJOR_AXIS_AU 3.45513 APHELION_AU 5.70813 ECCENTRICITY 0.652073 ARG_OF_PERIH_DEG 21.8178 INCLINATION_DEG 3.8868 Name: 0, dtype: object
С этими данными мы могли легко визуализировать и проанализировать орбиту кометы. Но рассмотрим статью о параметре Тиссерана (Космическая наука с Пифоном — рандеву с Юпитером): На орбиту малых тел влияют планетарные гравитационные возмущения (как Юпитер). Возмущения изменяют элементы орбиты, однако база данных предоставляет только статические результаты в определенную эпоху!
Нужно ли нам создавать сложные и ресурсоемкие численные модели? Это не обязательно. Пересмотрите уже используемые ядра SPICE. Возможно, мы найдем ядро, которое уже содержит некоторые более сложные соображения. Давайте посмотрим.
В общедоступном репозитории ядра есть подпапка для миссии Rosetta. Там мы находим подпапки kernels, а затем spk. Мы подробнее рассмотрим старые данные в подпапке former_versions. В aareadme.txt указано:
Name Comments
--------------------------------- ---------------------------------
67P_CHURY_GERAS_2004_2016.BSP Contains ephemeris data for the
Comet Churyumov-Gerasimenko/67P.
This file has been originated from
JPL Horizons On-Line Ephemeris
System and IT IS NOT an official
Rosetta Project SPK. It is
provided as support to study
cases until an official SPK is
released. Spans from January 2004
to January 2016.
Хотя это не официальный проект spk, мы используем это ядро из-за охвата времени между 2004 и 2016 годами. Мы добавляем путь к ядру (с некоторыми другими необходимыми ядрами) в kernel_meta.txt. и загрузите файл через furnsh.
Для элементов орбиты, которые были загружены ранее, нам понадобится гравитационная постоянная, умноженная на массу Солнца. Мы извлекаем эту информацию в строке 5 и 6 с помощью команды bodvcd.
Чего мы хотим достичь? Данные MPC предоставляют элементы орбиты, которые не меняются со временем. С другой стороны, ядро SPICE может содержать более сложный набор данных о траектории кометы. Итак, давайте проанализируем временной диапазон в 10 лет; с 2004 по 2014 год. Во-первых, нам нужно создать массив, содержащий шаги даты и времени между обоими годами. Строки 2 и 3 создают панды Timestamp(s) на 1 января каждого года. Строка 6 создает numpy массив с 1000 временными шагами, а строка 9 преобразует его обратно в объект даты и времени с помощью pandas to_datetime .
Теперь мы вычисляем орбитальные элементы 67P для каждого временного шага на основе соответствующего ядра SPICE. В строке 2 создается пустой фрейм данных, в котором будут храниться все результаты. Строки времени UTC хранятся в строке 5 и преобразуются в эфемеридное время (ET) в строках с 8 по 9 с помощью функции SPICE utc2et. Затем вычисляется вектор состояния 67P, видимый со стороны Солнца (строки 12–16). Для этого вычисления используется функция spkgeo, а ECLIPJ2000 устанавливается в качестве системы отсчета. Идентификационный код НАИФ для 67P — 1000012 [*]. Наконец, строки с 19 по 23 вычисляют соответствующие орбитальные элементы векторов состояния с помощью функции oscltx.
Орбитальные элементы хранятся в виде массивов в столбце STATE_VEC_ORB_ELEM. Теперь мы извлекаем необходимые орбитальные элементы из массива и назначаем их отдельным столбцам фрейма данных. Порядок элементов орбит в массивах объясняется в документации oscltx:
elts are equivalent conic elements describing the orbit
of the body around its primary. The elements are,
in order:
RP Perifocal distance.
ECC Eccentricity.
INC Inclination.
LNODE Longitude of the ascending node.
ARGP Argument of periapsis.
M0 Mean anomaly at epoch.
T0 Epoch.
MU Gravitational parameter.
NU True anomaly at epoch.
A Semi-major axis. A is set to zero if
it is not computable.
TAU Orbital period. Applicable only for
elliptical orbits. Set to zero otherwise.
Пространственная информация дается в км, а угловая информация дается в радианах. Мы конвертируем значения в астрономические единицы (используя функцию convrt в строках с 3 по 7 и с 26 по 30) и градусы (используя numpy degrees в строках 14–15, 18–19 и 22–23) соответственно. В строках с 33 по 36 мы дополнительно вычисляем афелий, используя большую полуось и эксцентриситет.
Теперь можно отобразить результаты ядра SPICE и статический набор данных из набора данных MPC. Мы используем matplotlib для визуализации и интерпретации результатов. Мы устанавливаем темную тему в строке 4 и увеличиваем размер шрифта для лучшей читаемости в строке 7. Цикл for теперь перебирает 3 разных параметра: перигелий, эксцентриситет и аргумент перигелия (строки 11–50). Данные SPICE отображаются в строках с 22 по 24, а данные MPC отображаются в виде горизонтальной опорной линии в строках кода с 28 по 29 с использованием функции matplotlib hlines. Процедура построения графика требует статического значения y и необходимого диапазона оси x (xmin и xmax). Задаются сетка (строка 35) и метки осей (строки 38 и 39). В строках с 43 по 47 настраиваются параметры легенды (прозрачность строк удалена для лучшей читаемости), и рисунок, наконец, сохраняется (строка 50).
Давайте посмотрим на получившиеся графики. Первые два рисунка показывают перигелий и эксцентриситет относительно даты и времени, указанных в формате UTC. Сплошная линия представляет значения из ядра SPICE, а горизонтальная пунктирная линия — справочные данные из набора данных MPC. Как видите, оба результата различаются. Кроме того, результаты SPICE со временем меняются. Например, между 2004 и 2008 годами перигелий изменился примерно с 1,29 а.е. до менее чем 1,25 а.е., то есть около 6 миллионов километров! Наименьшая разница между данными MPC и результатом ядра SPICE составляет около 0,04 а.е. для данных перигелия и 0,0125 для эксцентриситета. Учитывая, что мы не учли ошибку из данных MPC (ядро SPICE не предоставило ошибку для значений), разница не является существенно большой.
На последнем рисунке показано отношение перигелия к дате и времени в формате UTC. Вы можете видеть, что положение перигелия кометы повернулось примерно на 2 градуса за 2 года. Это "много"? Наша родная планета Земля также подвержена влиянию гравитационных возмущений. Однако из-за большего расстояния до других планет и собственной массы эффекты довольно малы. Аргумент перигелия Земли поворачивается на 360° за более чем 100 000 лет [3]!
Вывод
Как можно доверять источникам данных? Есть ли способ получить настоящие данные? На этот вопрос нельзя ответить. Это всегда зависит от вашего варианта использования и требуемой степени детализации. Для анализа большой популяции всех комет (Космическая наука с Python — все ли мы наблюдали?) данных MPC достаточно. Даже для краткосрочных прогнозов и расчетов траекторий статические результаты подходят. Однако исследования по проектированию космических аппаратов и детальный научный анализ конкретных комет требуют более сложного и скорректированного взгляда. Численный анализ, учитывающий гравитационные возмущения и выделение газа из кометы, является наиболее подробным. Предоставление ядер SPICE в общедоступном репозитории всегда является первым хорошим подходом к обеспечению хорошего качества данных.
Но можно ли объяснить возмущения количественно? Без слишком сложных вычислений? Да! И это будет сделано в следующий раз, когда мы будем использовать SPICE для получения собственного решения истории 67P.
Томас
Задания
По совету некоторых пользователей я буду добавлять небольшие задачи в конце каждого урока. Таким образом, вы сможете применить полученные знания. Не стесняйтесь обращаться ко мне, если у вас возникнут какие-либо проблемы.
- Имеет ли параметр Тиссерана w.r.t. Юпитер меняется со временем?
- Визуализируйте расстояние между Юпитером и 67P между начальным и конечным временем. Используйте для этой цели spkgps и подумайте о параметрах targ и obs. Преобразуйте результаты в AU.
Сноски
Как узнать идентификаторы NAIF, не проверяя веб-сайт NAIF ID? Каковы временные окна или охваты ядер SPICE? На эти вопросы можно ответить с помощью функций SPICE, которые будут показаны в следующем руководстве, где я представлю дополнительную статью об извлечении информации о ядре SPICE.
использованная литература
[1] Мишель С. Фестоу (редактор), Х. Уве Келлер (редактор), Гарольд А. Уивер-младший (редактор). (2004). Кометы II. Издательство Аризонского университета. стр. 780. ISBN-10: 0816524505. ISBN-13: 978–0816524501
[2] Комби, М. Р.; Харрис, WM; Смит, WH (2004). Газовая динамика и кинетика кометной комы: теория и наблюдения. В: Кометы II [1]. п. 523–552
[3] van den Heuvel, EPJ (1966). О прецессии как причине плейстоценовых вариаций температуры воды Атлантического океана. Международный геофизический журнал. 11: 323–336. Бибкод: 1966GeoJ…11..323V. doi: 10.1111/j.1365–246X.1966.tb03086.x
