Генератор сигналов GPS
В примере рассматривается работа блока gpsWaveformGenerator, а именно генерируемые с его помощью сигналы различных типов (P-, L1C-, L2C- и L5-коды).
Введение
GPS — это космическая радионавигационная инфраструктура, состоящая из группировки спутников, передающих навигационные сигналы, и сети наземных станций и станций управления спутниками для мониторинга и контроля этих спутников.
Спутники GPS используют три разные частоты для гражданских применений: L1 на 1575,42 МГц, L2 на 1227,60 МГц и L5 на 1176,45 МГц.
GPS использует пять различных гражданских сигналов: код грубого захвата (C/A) и точный код (P) в унаследованном диапазоне L1, модернизированный гражданский код L1 (L1C) в диапазоне L1, гражданский средний код (CM) и гражданский длинный код (CL) в диапазоне L2C, а также синфазный код (I5-код) и квадратурный код (Q5-код) в диапазоне L5.
P-код
P-код — дальномерный код, использует метод модуляции BPSK (Binary Phase Shift Keying) для передачи на частотах L1 и L2 с чиповой скоростью 10,23 Мчип/с. Более высокая чиповая скорость P-кода способствует повышению точности по сравнению с C/A-кодом.
P-код является общедоступным, однако его шифрование с помощью Y-кода для формирования P(Y)-кода ограничивает его использование военными целями.
Ниже представлен скрипт для генерации и отображения спектрограммы Р-кода.
SampleRate = 4 * 10.23e6 задает частоту дискретизации 40,92 МГц, то есть 4 отсчета на один чип P‑кода. Это удовлетворяет теореме Котельникова, поскольку ширина главного лепестка спектра P‑кода составляет 20,46 МГц (от –10,23 до +10,23 МГц), и позволяет увидеть несколько боковых лепестков.
navdata — вектор из двух случайных битов (0 или 1), имитирующий навигационное сообщение GPS (в реальности скорость передачи данных 50 бит/с, P‑код модулируется этими битами через BPSK).
p = welch_pgram(waveform_vec, 1024, 0;
fs=gpswaveobj.SampleRate,
window = hann,
onesided=false)
Метод Уэлча с окном Ханна дает сглаженную оценку, усредняя периодограммы отдельных сегментов.
Подключение вспомогательных файлов и библиотек
Для визуализации результатов работы блока необходимо подключить сторонние библиотеки.
let
installed_packages = collect(x.name for (_, x) in Pkg.dependencies() if x.is_direct_dep)
list_packages = ["DSP"]
for pack in list_packages
pack in installed_packages || Pkg.add(pack)
end
end
import DSP
using EngeeDSP, EngeeSatellites
# Создание и настройка генератора
gpswaveobj = gpsWaveformGenerator()
gpswaveobj.EnablePCode = true
gpswaveobj.SampleRate = 4 * 10.23e6
# Генерация сигнала
navdata = rand([0,1], 2, 1)
waveform = gpswaveobj(navdata)
# Преобразуем матрицу сигнала в одномерный вектор
waveform_vec = vec(waveform)
# nfft=1024 — размер БПФ (чем больше, тем выше разрешение)
# Задаем длину окна n=1024 и перекрытие noverlap=0 (0%)
p = DSP.welch_pgram(
waveform_vec, 1024, 0;
fs=gpswaveobj.SampleRate,
window = DSP.hann,
onesided=false
)
# Извлекаем частоты и спектральную плотность мощности
freqs = DSP.freq(p) # вектор частот в Гц
psd = DSP.power(p) # вектор PSD в Вт/Гц
# Центрируем спектр
freqs_shifted = EngeeDSP.Functions.fftshift(freqs)
psd_shifted = EngeeDSP.Functions.fftshift(psd)
# Построение спектра в dBW/Гц
plot(freqs_shifted/1e6, 10*log10.(psd_shifted .+ eps()),
xlabel="Частота (МГц)",
ylabel="Мощность (dBW/Гц)",
title="Спектр GPS P-кода",
legend=false)
Сигнал L1C
L1C — Этот сигнал открыт для гражданского использования и передается на частоте L1. Он использует схему модуляции с двоичным смещением несущей (BOC)
Глобальное отличие от предыдущего скрипта - наличие двух спутников и соответственно, двух сигналов. В основе лежат те же основные принципы, что и в генерации P-кода.
prn = [4, 70] # идентификаторы спутников
numsat = length(prn) # 2 спутника
numbits = 100 # длина сообщения в битах
navdata = rand([0, 1], numbits, numsat)
fs = 25e6 # частота дискретизации 25 МГц
gpswaveobj = gpsWaveformGenerator(; SignalType="l1c", PRNID=prn, SampleRate=fs)
waveform = gpswaveobj(navdata) # Это матрица размером N×2, где N - количесство отсчетов
# Настройка параметров для Welch (сглаживание спектра)
# Если сигнал короче 2048 точек, используем его полную длину
seg_len = min(2048, size(waveform, 1))
# Построение графика
# Создаем пустую область с подписями осей и легендой
plot_obj = plot(; xlabel="Частота (МГц)", ylabel="Мощность (dBW/Гц)",
title="Спектр GPS L1C (PRN 4 и 70)", legend=:topleft)
# Цикл по двум каналам (спутникам)
for i in 1:numsat
# Берем только i-й столбец (канал)
channel = waveform[:, i]
# Считаем спектр методом Уэлча
p_i = DSP.welch_pgram(channel, seg_len, 0; fs=fs, onesided=false)
# Извлекаем и центрируем частоты
freqs_i = DSP.freq(p_i)
psd_i = DSP.power(p_i)
freqs_shifted = EngeeDSP.Functions.fftshift(freqs_i)
psd_shifted = EngeeDSP.Functions.fftshift(psd_i)
# Добавляем график
plot!(plot_obj, freqs_shifted / 1e6, 10 * log10.(psd_shifted .+ eps()),
label="PRN $(prn[i])")
end
# Отображаем график
display(plot_obj)
Сигнал L2C
SignalType="l2c" – основой является гражданский сигнал L2C. Он передаётся на частоте 1227,60 МГц и содержит два кода, мультиплексированных по времени:
-
CM (Civil Moderate) – код длиной 10 230 чипов (период 20 мс), модулирован данными CNAV (25 бит/с).
-
CL (Civil Long) – код длиной 767 250 чипов (период 1,5 с), без данных (пилотная компонента).
Оба кода имеют чиповую скорость 511,5 кчип/с, но мультиплексируются таким образом, что суммарная чиповая скорость составляет 1,023 Мчип/с (как у C/A).
waveform = gpswaveobj(lnavdata, cnavdata)
- В отличие от предыдущего примера L1C, здесь данные передаются кортежем из двух матриц, потому что L2C использует два разных навигационных потока: LNAV и CNAV. Генератор использует их для модуляции соответствующих компонентов.
Число спутников увеличено до четырёх, но в остальном принцип работы остаётся тем же.
# Задаем параметры
prn = [7, 11, 20, 28]
numsat = length(prn)
numbits = 10
# Генерация случайных данных LNAV и CNAV (матрицы 10x4)
lnavdata = rand([0, 1], numbits, numsat)
cnavdata = rand([0, 1], numbits, numsat)
# Создаем генератор с настройками
gpswaveobj = gpsWaveformGenerator(;
SignalType="l2c",
PRNID=prn,
SampleRate=15e6,
EnablePCode=true
)
# Генерация сигнала
# ВНИМАНИЕ: Передаем данные в виде кортежа (аналог cell-массива MATLAB)
waveform = gpswaveobj(lnavdata, cnavdata)
# Проверка размера
println("Размер сгенерированного сигнала: ", size(waveform))
# Ожидаемый вывод: (3000000, 4)
# Настройка БПФ для сглаживания
seg_len = min(2048, size(waveform, 1))
# Построение графика
plot_obj = plot(; xlabel="Частота (МГц)", ylabel="Мощность (dBW/Гц)",
title="Спектр GPS L2C (PRN 7, 11, 20, 28)", legend=:topleft)
for i in 1:numsat
channel = waveform[:, i] # Берем канал i-го спутника
# Вычисляем спектр
p_i = DSP.welch_pgram(channel, seg_len, 0; fs=15e6, onesided=false, nfft=seg_len)
freqs_i = DSP.freq(p_i)
psd_i = DSP.power(p_i)
freqs_shifted = EngeeDSP.Functions.fftshift(freqs_i)
psd_shifted = EngeeDSP.Functions.fftshift(psd_i)
# Рисуем график для этого спутника
plot!(plot_obj, freqs_shifted / 1e6, 10 * log10.(psd_shifted .+ eps()),
label="PRN $(prn[i])")
end
display(plot_obj) # Отображаем график
Сигнал L5
Сигнал GPS L5 — это один из современных гражданских сигналов системы GPS, передаваемый на частоте 1176,45 МГц (диапазон L5). Он был введён в рамках программы модернизации GPS и предназначен в первую очередь для обеспечения высокой точности и целостности, необходимых для авиационных и других ответственных приложений (Safety-of-Life).
Основное тело скрипта остаётся схожим с предыдущими, однако к сигналу добавляются искажения, вызванные доплеровским смещением, возникающим из-за движения спутника и приемника. Введение разных сдвигов (+3 кГц и –4,6 кГц) разносит их спектры, что приводит к появлению двух близких максимумов внутри общей огибающей – интересный и наглядный эффект.
Внесение искажений: доплеровский сдвиг частоты
# Внесение искажений
# Частотный сдвиг
t = (0:size(waveform, 1) - 1) / fs # Вектор времени
freq_offsets = [3e3, -4.6e3] # Сдвиги для двух каналов
# Создаем комплексную экспоненту и умножаем на сигнал
fwave = waveform .* exp.(im * 2 * pi * t .* freq_offsets')
# Внесение временной задержки
delay_val = 5.23e-4
delay_samples = round(Int, delay_val * fs) # 13075 отсчетов
# Физическая задержка: добавляем нули в начало и отбрасываем конец
dwave = vcat(zeros(ComplexF64, delay_samples, numsat), fwave[1:end-delay_samples, :])
# Суммирование двух каналов
bbwave = sum(dwave, dims=2) # Размер (10_000_000 x 1)
bbwave_vec = vec(bbwave) # Превращаем в одномерный вектор
# Генерация сигнала GPS L5
prn = [185, 189]
numsat = length(prn)
numbits = 40
msg = rand([0, 1], numbits, numsat)
fs = 25e6
gpswaveobj = gpsWaveformGenerator(;
SignalType="l5",
PRNID=prn,
InitialTime=42123,
SampleRate=fs
)
# Генерируем сигнал. Получаем матрицу (10_000_000 x 2)
waveform = gpswaveobj(msg)
println("Размер исходного сигнала: ", size(waveform))
# Внесение искажений
# Частотный сдвиг
t = (0:size(waveform, 1) - 1) / fs # Вектор времени
freq_offsets = [3e3, -4.6e3] # Сдвиги для двух каналов
# Создаем комплексную экспоненту и умножаем на сигнал
fwave = waveform .* exp.(im * 2 * pi * t .* freq_offsets')
# Внесение временной задержки
delay_val = 5.23e-4
delay_samples = round(Int, delay_val * fs) # 13075 отсчетов
# Физическая задержка: добавляем нули в начало и отбрасываем конец
dwave = vcat(zeros(ComplexF64, delay_samples, numsat), fwave[1:end-delay_samples, :])
# Суммирование двух каналов
bbwave = sum(dwave, dims=2) # Размер (10_000_000 x 1)
bbwave_vec = vec(bbwave) # Превращаем в одномерный вектор
println("Размер результирующего сигнала (BB): ", size(bbwave_vec))
# Визуализация спектра
seg_len = min(2048, length(bbwave_vec))
# Вычисляем спектр
p = DSP.welch_pgram(bbwave_vec, seg_len, 0; fs=fs, onesided=false, nfft=seg_len)
freqs = DSP.freq(p)
psd = DSP.power(p)
# Центрируем спектр
freqs_shifted = EngeeDSP.Functions.fftshift(freqs)
psd_shifted = EngeeDSP.Functions.fftshift(psd)
# Рисуем график
plot(freqs_shifted / 1e6, 10 * log10.(psd_shifted .+ eps()),
xlabel="Частота (МГц)",
ylabel="Мощность (dBW/Гц)",
title="Спектр GPS L5",
legend=false)
Вывод
Представленные примеры позволяют оценить работу блока gpsWaveformGenerator, его функционал и влияние настроек на результирующий GPS-сигнал на выходе.