Астрономический журнал, 2021, T. 98, № 11, стр. 935-943

1. ВВЕДЕНИЕ

Задача сличения и синхронизации с высокой точностью шкал времени, разнесенных на большие расстояния, была и остается актуальной по настоящее время. В масштабах Земли эта задача решается с помощью работающей спутниковой навигационной системы ГЛОНАСС, с помощью методов РСДБ. На космических расстояниях далее орбиты Луны предлагается использовать сличение и синхронизацию шкал времени с использованием гигантских импульсов (ГИ) пульсара. Этот метод был предложен в работе [1].

Гигантские импульсы пульсаров не подходят для хронометрирования и построения шкал времени, поскольку их появление внутри пульсарного импульса носит случайный характер. Но гигантские импульсы имеют очень малую длительность, как правило, определяемую шириной частотной полосы регистратора, их амплитуда и интервал повторения имеют уникальный характер в заданном промежутке времени, они имеют очень высокую спектральную плотность излучения, в импульсе на 4–6 порядков превышающую среднее значение. Такие свойства позволяют использовать для их регистрации антенны с относительно небольшой эффективной площадью, которые могут быть развернуты на других телах Солнечной системы.

Освоение Луны является первоочередной задачей для космических агентств многих стран как по причине самостоятельной ценности Луны, например, в плане добычи полезных ископаемых, так и по причине отработки космических технологий для полетов к более далеким телам Солнечной системы. В связи с этим сличение шкал времени представляет собой одну из важнейших частей радиотехнического обеспечения межпланетной навигации, связи, распространения сигналов точного времени. Потенциальная точность сличения временных шкал по ГИ оценивается величиной <1 нс (см. ниже), которая сравнима с точностью сличения посредством глобальных навигационных спутниковых систем (ГНСС) на Земле.

Основной целью данной работы является демонстрация возможностей технологии синхронизации шкал времени по ГИ на примере сличения лунной и земной шкал времени. Для этого в статье разработана модель задержки между станциями на Луне и Земле, изложен алгоритм сличения шкал по гигантским импульсам на основе их кросс-корреляции, приведена оценка точности сличения шкал.

2. СВОЙСТВА ГИГАНТСКИХ ИМПУЛЬСОВ И ОСОБЕННОСТИ ИХ НАБЛЮДЕНИЙ

Выделим отдельно пульсары, излучающие гигантские радиоимпульсы, пиковая плотность потока которых может превышать пиковую плотность потока средних импульсов от нескольких десятков до сотен тысяч раз. Основными свойствами ГИ являются [2–4]:

• аномально высокая пиковая плотность потока радиоизлучения (для пульсара в Крабовидной туманности B0531+21 на частотах 600–1400 МГц – до нескольких миллионов янских);

• короткая длительность (от нескольких наносекунд до нескольких микросекунд);

• степенной закон распределения по пиковой плотности потока/энергии в импульсе;

• высокая направленность излучения, которая следует из узкой долготной локализации ГИ;

• высокая степень поляризации ГИ.

Общее число известных пульсаров с ГИ составляет 16 [5], включая 11 нормальных и 5 миллисекундных пульсаров, в том числе “классические” представители: пульсар в Крабовидной туманности B0531+21 (P = 33 мс) и одиночный миллисекундный пульсар B1937+21 (P = 1.56 мс), а также ряд пульсаров, наблюдавшихся только на низких частотах и обладающих пиковой плотностью потока до нескольких сотен Янских.

Индивидуальные импульсы пульсара PSR 0531+21 в Крабовидной туманности полностью состоят из коротких (обычных и гигантских) импульсов. Распределение их амплитуд соответствует степенному закону. Частота появления особо мощных ГИ этого пульсара была оценена в результате его длительных наблюдений на радиотелескопе РТ-64 в Калязине на частотах 600 МГц и 4.85 ГГц. На частоте 600 МГц в течение часа регистрировалось 6–8 ГИ с пиковой плотностью потока >300 кЯн, на частоте 4.85 ГГц регистрировалось 2–3 ГИ с потоком >10 кЯн [6]. Наименьшая длительность ГИ была зарегистрирована на радиотелескопе в Аресибо на частоте 9.25 ГГц в полосе 2.5 ГГц от пульсара 0531+21 в Крабовидной туманности [2]. Ширина самой мощной детали микроструктуры импульса составила 0.4 нс, пиковая плотность потока – 2.2 МЯн. Отметим, что наблюдавшаяся деталь микроструктуры не разрешена и ширина ее ограничена скоростью выборки (обратной к полосе частот) приемной системы. То есть в данном случае при корреляции ГИ, принятых в разных пунктах, мы, по сути, используем метод РСДБ, лишенный своего основного недостатка – неопределенности периода. Кроме этого, нет необходимости передавать в центр обработки с борта КА весь массив данных наблюдений на длительном интервале, достаточно передать данные на интервале появления ГИ, что на порядки снижает объем передаваемых данных. В случае наблюдений ГИ пульсаров в диапазоне 400 МГц с полосой регистрации 50 МГц точность сличения шкал времени по ГИ пульсаров будет лучше 10 нс. Если принять очень консервативную оценку отношения сигнал/шум ~10, то погрешность синхронизации будет порядка 1 нс. В случае реализации на борту КА додетекторного метода компенсации дисперсионного запаздывания импульса в межзвездной среде объем передаваемых с борта КА данных в сеансе связи может быть ограничен размером 20 Мб.

3. МОДЕЛЬ ЗАДЕРЖКИ МЕЖДУ НАБЛЮДАТЕЛЯМИ НА ЛУНЕ И ЗЕМЛЕ

Привязка двух шкал времени основана на измерении разности времени прихода импульсов короткой продолжительности относительно этих шкал. Как уже было сказано, в роли импульсов короткой продолжительности ранее в работе [1] было предложено использовать гигантские импульсы от пульсаров. Мощность этих импульсов такова, что позволяет регистрировать их на антеннах относительно небольшой площади.

Для корректной привязки шкал необходимо разработать модель задержки, которая позволит точно рассчитывать время распространения между двумя приемными антеннами, расположенными на Луне и Земле. Расчет задержки производится в несколько этапов (см. рис. 1):

1) Рассчитывается вектор центра масс системы Земля–Луна относительно барицентра Солнечной системы $\overrightarrow {BG} (t)$ на необходимый момент времени t. Для этого используются современные планетные эфемериды EMP2XXX (ИПА РАН) или же DE4XX (JPL). Вызов процедуры расчета координат барицентра Земля–Луна производится с параметром NPL = 3.

2) Рассчитывается положение Луны относительно геоцентра $\overrightarrow {TL} (t)$. Это производится также с помощью эфемерид с параметром NPL = 11.

3) При вычисленных компонентах векторов $\overrightarrow {BG} (t)$ и $\overrightarrow {TL} (t)$ рассчитываются вектора $\overrightarrow {GL} (t)$ и $\overrightarrow {GT} (t)$ по следующим формулам:

где ${{m}_{L}}$, ${{m}_{T}}$ – массы Луны и Земли соответственно.

4) Вектор положения наблюдателя на Луне $\overrightarrow {LR} (t)$ определяется через повороты на углы либрации, которые, в свою очередь, также вычисляются из эфемерид с параметром NPL = 13.

где x, y, z – координаты радиотелескопа на поверхности Луны в селеноцентрической системе координат, определяемой относительно главных осей инерции Луны, ${\mathbf{\Psi }}(t)$, ${\mathbf{\Theta }}(t)$, ${\mathbf{\Phi }}(t)$ – углы либрации, ${{{\mathbf{R}}}_{x}}$, ${{{\mathbf{R}}}_{z}}~$ – матрицы поворота вокруг осей x и z соответственно.

5) Вектор положения наблюдателя на Земле $\overrightarrow {TO} (t)$ вычисляется через ряд последовательных поворотов вектора наблюдателя в системе ITRF, учитывающих суточное вращение Земли, прецессию, нутацию и движение полюса:

где $u$, ${v}$, $w$ – координаты наблюдателя на Земле в системе ITRF, P, N – матрицы прецессии и нутации, ${{{\mathbf{R}}}_{x}}$, ${{{\mathbf{R}}}_{y}}$, ${{{\mathbf{R}}}_{z}}$ – матрицы поворота вокруг осей x, y, z, ${{x}_{p}}$, ${{y}_{p}}$ – координаты полюса Земли, θ – звездное время.

6) Барицентрические вектора наблюдателей $\overrightarrow {BO} $ и $\overrightarrow {BR} $ вычисляются как соответствующие суммы

7) Вычисляется время распространения сигнала вдоль направления на пульсар (задержка Рёмера):

где $\vec {r}$ – барицентрический радиус-вектор наблюдателя ($\overrightarrow {BO} $ или $\overrightarrow {BR} $), $\vec {k}$ – барицентрический единичный вектор в направлении на пульсар, c – скорость света.

8) Вычисляется задержка сигнала, вызванная конечностью расстояния до пульсара по формуле:

где $\vec {r}$ – барицентрический радиус-вектор наблюдателя ($\overrightarrow {BO} $ или $\overrightarrow {BR} $), R – расстояние до пульсара, $\vec {k}~$ – барицентрический единичный вектор в направлении на пульсар, c – скорость света.

9) Вычисляется релятивистская поправка распространения сигнала в гравитационном поле тел Солнечной системы по следующей формуле:

где $\vec {k}~$ – барицентрический единичный вектор в направлении на пульсар, ${{\vec {r}}_{p}}$ – вектор “наблюдатель – гравитирующее тело” (в нашем случае это большие планеты Солнечной системы), ${{R}_{{gp}}}$ – гравитационный радиус p-й планеты, который вычисляется по формуле ${{R}_{{gp}}} = \frac{{2G{{M}_{p}}}}{{{{c}^{3}}}}$. График релятивистской поправки показан на рис. 2.

10) Разность моментов прихода импульса к наблюдателям на Земле и Луне $\Delta {{t}_{c}}$ вычисляется по окончательной формуле

где моменты t₁ и t₂ выражены в шкале времени TB. График модели задержки показан на рис. 3.

На этом алгоритм вычисления разности времен прихода пульсарного импульса к наблюдателю на Земле и Луне завершается.

Необходимо отметить несколько особенностей при использовании этого алгоритма:

1) В нем отсутствует поправка за распространение сигнала в межзвездной среде. Это связано с тем, что при вычислении разности моментов прихода импульсов между Землей и Луной эти задержки эффективно взаимоуничтожаются, так как путь распространения сигнала при вычислении этих задержек с достаточной точностью можно считать одинаковым. Единственной ионизованной средой между Землей и Луной является земная ионосфера, которая вносит дополнительную задержку на наблюдательных частотах 400–1400 МГц порядка 2–50 нс, что с большим запасом укладывается в точность метода пульсарного тайминга. Если в будущем точность измерений возрастет, то не составит труда вычислить ионосферную поправку через задержку в зените и учет наклонного пути распространения сигнала между наблюдателями на Луне и Земле.

2) Необходимо учитывать разный ход шкал времени на Земле и Луне. При этом можно придерживаться того же подхода, что и в шкалах TT и TB: средний темп сделать одинаковым, а различие будет только в периодических членах. График взаимного хода шкалы TT на Земле и шкалы на Луне приведен на рис. 4.

3) Для учета светового промежутка (времени распространения сигнала между Луной и Землей) используется метод итераций. Процесс сходится очень быстро. На практике достаточно одной итерации.

4. МОДЕЛИРОВАНИЕ ПРИВЯЗКИ ДВУХ ШКАЛ ВРЕМЕНИ ПО ГИГАНТСКИМ ИМПУЛЬСАМ ПУЛЬСАРА PSR 0531+21

Для проведения моделирования были использованы наблюдения пульсара в Крабовидной туманности 02.01.2018 и 02.02.2018 г. Наблюдения проводились на антенне БСА ФИАН на частоте 111 МГц в полосе 2.5 МГц с временным разрешением 2.4576 мс. Записывались 6-секундные файлы, показанные на рис. 5.

В работе [7] была обоснована наиболее оптимальная частота для наблюдений пульсаров 400–500 МГц. Для перевода наблюдений с частоты 111 на 500 МГц используем следующие соображения. Так как на частоте 111 МГц импульсы от пульсара PSR B0531+21 показывают заметное экспоненциальное рассеяние, а рассеяние в среде распространения изменяется с частотой как ${{f}^{{ - 4}}}$, то для масштабирования этого эффекта достаточно изменить шкалу времени в (500/111)⁴ = 412 раз. Таким образом, новый интервал выборки составит величину примерно 6 мкс. Это не изменит форму импульса, но изменит его характерную ширину. В дальнейшем удобно будет оперировать именно дискретами, а приводить к единицам времени в конце расчетов.

Для моделирования наблюдений того же импульса второй антенной на Луне исходные импульсы сдвигались на величину, определяемую моделью задержки, и на них накладывался гауссовский белый шум с амплитудой эквивалентной результирующему отношению сигнал/шум = 3 и 10. На величину С/Ш влияет много факторов. Здесь мы будем считать, что в первую очередь она зависит от эффективной площади антенны. В случае фазированной решетки, сходной с БСА, для получения вышеупомянутых отношений С/Ш антенна на Луне должна состоять из 1.5 или 5 тыс. диполей.

При определении величины расхождения шкал использовалась процедура кросс-корреляции импульсов с последующим сглаживанием шаблоном экспоненциальной формы (рис. 7) и вписыванием в максимум кросс-корреляционной функции полинома 4-й степени. Именно эта процедура позволяет определить положение максимума на субдискретном уровне. Кросс-корреляционные функции показаны на рис. 8, вписывание полинома показано на рис. 9.

Результаты измерения привязки шкал времени сведены в таблице 1. Основной вывод: метод работает даже при относительно слабом сигнале с Луны (С/Ш = 3), так как основной сигнал, записанный на Земле, уверенно “вытягивает” максимум кросс-корреляционной функции. Добавим, что формальная ошибка измерения сдвига не сильно зависит от отношения сигнал/шум, так как вершина кросс-корреляционной функции по-прежнему остается гладкой и хорошо описывается полиномом, чего нельзя сказать про фактический сдвиг, который ожидаемо отличается сильнее от модельного при ухудшении отношения сигнал/шум.

Теперь рассмотрим влияние ширины частотного канала на точность измерения положения импульса и, соответственно, измерения привязки двух шкал. Приведем формулу уширения импульса с мерой дисперсии DM в полосе $\Delta f$ [Гц] на частоте f [Гц]:

При используемых параметрах f = 500 МГц, $B = 2.5$ МГц, число каналов = 512 получаем $\Delta f = 4873$ Гц. Подставляя эти значения в формулу (1), получим величину уширения 18 мкс, что равно 3 дискрета.

Далее моделирование проводилось для ширины канала 53.6 кГц. С точки зрения математики использование более широкого частотного канала, приводящего к уширению импульса, эквивалентно свертке суммарного импульса $x(\tau )$ с П-образным ядром, ширина которого равна уширению внутри канала

В дискретном случае ширина канала 53.6 кГц эквивалентна прямоугольному окну шириной 11 отсчетов. Свертка с окном такой ширины была применена к исходному импульсу.

Также при расчетах использовалось два значения отношения сигнал/шум = 3 и 10. Результаты этой части моделирования показаны в нижних четырех строчках таблицы 1. Из них хорошо видно, что уширение импульса приводит к ухудшению точности определения положения максимума кросс-корреляционной функции, и отличие от фактического сдвига шкал может достигать нескольких отсчетов. Ситуация ожидаемо усугубляется, если один из импульсов имеет низкое значение С/Ш.

Еще раз специально отметим, что формальная ошибка вписывания полинома слабо зависит от С/Ш и уширения в канале. Более того, для широких импульсов максимум кросс-корреляционной функции становится более гладким, и это приводит даже к уменьшению формальной ошибки вписывания. Таким образом, для реалистичных оценок точности необходимо использовать величину фактического отклонения максимума от наперед заданного при моделировании положения.

5. ЗАКЛЮЧЕНИЕ

Модель задержки построена на основе разности координат лунного и земного пунктов в барицентрической системе отсчета с учетом эффектов распространения сигнала в гравитационном поле Солнечной системы и кривизны фронта. Погрешность модели оценивается величиной ~30 нс. Путем численного интегрирования посчитан взаимный ход лунной и земной шкал времени. Показано, что амплитуда вариаций после вычитания векового члена составляет величину ~150 мкс. Размер вариаций таков, что их необходимо учитывать при синхронизации шкал.

Результаты моделирования точности привязки двух шкал укладываются в ожидаемые. Частота выборки импульсов должна соответствовать уширению импульсов в канале. При ухудшении С/Ш одного из импульсов метод продолжает работать, так как другой импульс с высоким отношением С/Ш “вытягивает” корреляционную функцию, у которой положение максимума определяется с субдискретной точностью. В едини-цах времени при дискрете 6 мкс и С/Ш > 10 точность синхронизации шкал достигает величины ~0.5 мкс.

Ситуация существенно ухудшается, если используется аппаратура с широкими частотными каналами, в которых уширение импульса равно нескольким дискретам. В таком случае точность сличения падает более чем на порядок до величины ~20 мкс.

При использовании высоких частот (несколько гигагерц) и широких полос (несколько сотен мегагерц) потенциальная точность синхронизации шкал с помощью предложенного метода достигает субнаносекундного уровня на интервалах более суток, что позволяет сличать шкалы времени на телах Солнечной системы на самом современном уровне точности. Из-за нерегулярного и относительно редкого следования ГИ предложенный метод, по нашему мнению, не подходит для потребителей, которым нужны координаты и (или) метки шкалы в реальном времени. Им лучше использовать чисто радиотехнические методы.