Релятивистские поправки в GPS / Хабр

В этой статье разберем, как релятивистские поправки влияют на точность определения координат по GPS‑спутникам. Для наших тестов будем использовать симулятор сигнала спутников gps‑sdr‑sim и утилиту gnss‑sdr для геолокации.

Устанавливаем утилиту gnss‑sdr:

sudo apt install gnss-sdr

Скачиваем и собираем gps‑sdr‑sim:

git clone https://github.com/osqzss/gps-sdr-sim.git
cd gps-sdr-sim
gcc gpssim.c -lm -O3 -o gps-sdr-sim

В исходном коде файла gpssim.h — никаких релятивистских поправок нет. Теперь посмотрим файл gpssim.c, здесь находим поправку (строка 434):

       relativistic = -4.442807633E-10*eph.ecc*eph.sqrta*sek;

далее эта поправка добавляется к часам спутника:

       clk[0] = eph.af0 + tk*(eph.af1 + tk*eph.af2) + relativistic - eph.tgd;

Чтобы собрать «нерелятивистский» симулятор, достаточно обнулить эту поправку и пересобрать gps‑sdr‑sim:

 sed-i "434s|.*|\trelativistic =0;// -4.442807633E-10*eph.ecc*eph.sqrta*sek;|"
gcc gpssim.c -lm -O3 -o gps-sdr-sim-no-rel

Теперь можно записать файлы с моделью сигналов, которые бы принял GPS‑ приемник. Координаты приемника задаются через параметр ‑l. Для запуска симулятора нужен файл эфимерид (он содержит параметры орбиты и поправки часов спутников). Файл эфимерид можно скачать с официального сайта НАСА, или использовать файл brdc0010.22n из архива gps‑sdr‑sim:

./gps-sdr-sim -e brdc2580.26n -l 56.46958800001232,84.94677100122212 -o gpssim_true.iq -s 2500000 -d 1800 -b 8
./gps-sdr-sim-no-rel -e brdc2580.26n -l 56.46958800001232,84.94677100122212 -o gpssim_norel.iq -s 2500000 -d 1800 -b 8

Для работы с gnss‑sdr будем использовать конфигурационные файлы. Создаем файл «gps_analyze_true.conf» и копируем туда настройки:

[GNSS-SDR]

;######### GLOBAL OPTIONS ##################
GNSS-SDR.internal_fs_sps=2500000


;######### SIGNAL_SOURCE CONFIG ############
SignalSource.implementation=File_Signal_Source
# Укажите путь к вашему .bin файлу
SignalSource.filename=gpssim_true.iq
# Тип данных ishort означает знаковые 8-битные целые (int8)
SignalSource.item_type=ibyte
# Укажите реальную частоту дискретизации вашего файла 
SignalSource.sampling_frequency=2500000
SignalSource.samples=0

;######### SIGNAL_CONDITIONER CONFIG ############
SignalConditioner.implementation=Signal_Conditioner
# Преобразование чередующихся int8 (I и Q) в комплексный формат float
DataTypeAdapter.implementation=Ibyte_To_Complex
InputFilter.implementation=Pass_Through
InputFilter.item_type=gr_complex

;######### CHANNELS GLOBAL CONFIG ############
Channels_1C.count=8
Channels.in_acquisition=8
Channel.signal=1C

;######### ACQUISITION GLOBAL CONFIG ############
Acquisition_1C.implementation=GPS_L1_CA_PCPS_Acquisition
Acquisition_1C.item_type=gr_complex
Acquisition_1C.pfa=0.01
# Диапазон поиска Доплера увеличен до +-10000 Гц на случай сильного сдвига частоты
Acquisition_1C.doppler_max=10000
Acquisition_1C.doppler_step=250
Acquisition_1C.blocking=true

;######### TRACKING GLOBAL CONFIG ############
Tracking_1C.implementation=GPS_L1_CA_DLL_PLL_Tracking
Tracking_1C.item_type=gr_complex
# Полоса PLL (фаза) — 40 Гц, полоса DLL (код) — 4 Гц для стабильного удержания
Tracking_1C.pll_bw_hz=40.0
Tracking_1C.dll_bw_hz=4.0
Tracking_1C.dump=false
# Путь, куда будут сохраняться файлы (папка должна существовать)
;######### TELEMETRY DECODER GPS CONFIG ############
TelemetryDecoder_1C.implementation=GPS_L1_CA_Telemetry_Decoder

;######### OBSERVABLES CONFIG ############
Observables.implementation=Hybrid_Observables
;######### PVT CONFIG ############
PVT.enable_monitor=true
PVT.monitor_client_addresses=127.0.0.1
PVT.monitor_udp_port=1111
PVT.implementation=RTKLIB_PVT
PVT.positioning_mode=Single
PVT.output_rate_ms=100
PVT.display_rate_ms=500
PVT.iono_model=Broadcast
PVT.trop_model=Saastamoinen
PVT.flag_rtcm_server=true
PVT.flag_rtcm_tty_port=false
PVT.rtcm_dump_devname=/dev/pts/1
PVT.rtcm_tcp_port=2101
PVT.rtcm_MT1019_rate_ms=5000
PVT.rtcm_MT1077_rate_ms=1000

Monitor.enable_monitor=true
Monitor.decimation_factor=1
Monitor.client_addresses=127.0.0.1
Monitor.udp_port=1112

GNSS-SDR.telecommand_enabled=true
GNSS-SDR.telecommand_tcp_port=3333

Создаем файл «gps_analyze_norel.conf», копируем те же настройки, изменив имя файла на gpssim_norel.iq(10 строка).

Запускаем утилиту gnss‑sdr:

gnss-sdr --config_file='gps_analyze_true.conf'
gnss-sdr --config_file='gps_analyze_norel.conf'

Координаты записываются в форматах.kml,.gpx и.geojson с частотой обновления 10 Гц (каждые 0,1 с). На рисунке представлены результаты сравнения двух симуляторов за период 30 минут.

Слева локация по данным gps-sdr-sim, справа - по данным «нерелятивистского» симулятора. — Слева локация по данным gps‑sdr‑sim, справа — по данным «нерелятивистского» симулятора.

Как видно из сравнения, отсутствие релятивистских поправок в симуляторе приводит к ошибке геолокации: облако точек смещается на несколько метров относительно эталонной позиции.

Теперь попробуем разобраться с физическим смыслом релятивистской поправки. В работе Нила Ашби «Relativity in the Global Positioning System» (2003) приводится полное выражение для относительного сдвига частоты:

$\begin{equation}\frac{\Delta f}{f} = -\frac{1}{2}\frac{v^2}{c^2} - \frac{GM_E}{rc^2} - \frac{\Phi_0}{c^2}\label{eq:ashby_full}\end{equation}$

— скорость спутника в инерциальной системе отсчета (ECI), м/с;
— расстояние от центра Земли до спутника, м;
$GM_E = 3.986004418 \times 10^{14}$ м/с — гравитационный параметр Земли;
$\Phi_0 \approx -6.264 \times 10^7$ Дж/кг — гравитационный потенциал на геоиде;
м/с — скорость света.

Исключая расстояние через уравнение живой силы $v^2 = GM_E\left(\frac{2}{r} - \frac{1}{a}\right)$ — большая полуось орбиты, преобразуем исходное выражение:

$\left( \frac{\Delta f}{f} \right) = \left( \frac{\Delta f}{f} \right) _{const}+\left( \frac{\Delta f}{f} \right)_{fluct}$

Постоянная часть сдвига:

$\left( \frac{\Delta f}{f} \right)_{const} = -\frac{3GM_E}{2ac^2} -\frac{\Phi_0}{c^2} \approx 4.46*10^{-10}$

Переменная:

$\left( \frac{\Delta f}{f} \right)_{fluct} = \frac{GM_E}{ac^2} - \frac{v^2}{c^2}= \frac{v_s^2-v^2}{c^2}$ $v_s= \sqrt{\frac{GM_E}{a}}$ — средняя скорость спутника.

Таким образом, $\left( \frac{\Delta f}{f} \right)_{fluct}$ можно рассматривать как эффект Доплера второго порядка

Получить поправку relativistic = $-\frac{2\sqrt{GM_E \cdot a}}{c^{2}} \cdot e \cdot \sin (E)$ , где

— эксцентриситет орбиты;
— эксцентрическая аномалия,угол, который определяет текущее положение спутника на орбите.

можно путем интегрирования $\left( \frac{\Delta f}{f} \right)_{fluct}$ с момента времени ${t_0}$ , когда спутник находится в

апогее орбиты. Как видно из определения, поправка relativistic периодична (период соответствует обороту спутника вокруг орбиты) и имеет нулевое среднее.

Постоянная часть сдвига $\left( \frac{\Delta f}{f} \right)_{const}$ в симуляторе gps‑sdr‑sim не учитывается, так как

считается, что она уже скомпенсирована настройкой часов спутника. Именно $\left( \frac{\Delta f}{f} \right)_{const}$

принято «обвинять» в накоплении ошибок по 10 км в день.