Мониторинг ледовых покровов с помощью ГНСС-рефлектометрии

14
АВТОМЕТРИЯ. 2025. Т. 61, № 5
УДК 528.88; 629.7.066.1
МОНИТОРИНГ СОСТОЯНИЯ ЛЕДОВЫХ ПОКРОВОВ ПРЕСНЫХ
И СОЛЁНЫХ ВОДОЁМОВ С ПОМОЩЬЮ МЕТОДОВ
ГНСС-РЕФЛЕКТОМЕТРИИ
c М. Ю. Реушев1, 3 , Д. С. Макаров1 , Д. В. Харламов1 ,
В. Г. Подопригора2, 3 , К. В. Краснощеков1
1
Федеральный исследовательский центр СО РАН,
660036, г. Красноярск, Академгородок, 50
2
Институт физики им. Л. В. Киренского СО РАН,
660036, г. Красноярск, Академгородок, 50/38
3
Сибирский федеральный университет,
660041, г. Красноярск, просп. Свободный, 79
E-mail: reuqem@mail.ru
Проведены рефлектометрические измерения состояния ледовых поверхностей водоёмов
Красноярского края и республики Хакасия с помощью радиосигналов от навигационных
спутников, излучающих в диапазоне L1. В ходе обработки полученных рефлектограмм применена математическая обработка амплитудно-временны́х зависимостей с помощью методики быстрого преобразования Фурье, позволяющая экспериментально определять толщины ледовых покровов с различными макрофизическими характеристиками.
Ключевые слова: метод ГНСС-рефлектометрии, ледовые поверхности водоёмов, методика быстрого преобразования Фурье, сигналы навигационных спутников, диапазон L1.
DOI: 10.15372/AUT20250502
EDN: LPBXJB
Введение. В настоящее время спутниковые технологии находят всё большие применения в дистанционном зондировании поверхности Земли [1, 2]. Особый интерес представляют направления исследований, использующих сигналы диапазона L1 от навигационных спутников (НС), входящих в состав глобальной навигационной спутниковой системы
(ГНСС). Использование сигналов диапазона L1 от навигационных спутников позволяет
получать информацию о рассеянии, отражении и поглощении электромагнитных волн в
атмосфере и земных покровах. Изменения сигналов, принимаемые навигационными приёмниками, обусловлены различными электрофизическими характеристиками исследуемых
сред. Наиболее известный инструмент изучения таких изменений — это методы многолучевой ГНСС-рефлектометрии. С помощью таких методов можно восстанавливать значимые параметры исследуемых сред, таких, к примеру, как толщина и плотность слоёв,
земных покровов, коэффициенты ослабления и рассеяния сигналов, влажность, масштаб
неровностей поверхности и т. д. [3, 4].
В основе рефлектометрических измерений с помощью сигналов от навигационных
спутников лежит явление многолучевой интерференции. Навигационные сигналы диапазона L1 представляют собой непрерывно-фазомодулированные сигналы типа BPSK (двоичная фазовая манипуляция: один бит на символ) и занимают частотную полосу порядка
∆f = 1 МГц в области центральной частоты 1575 МГц [5]. Длина когерентности радиосигнала определяется выражением lког = 2c/∆f . Полоса частот сигналов в диапазоне L1
составляет величину порядка 1 МГц, вследствие чего сигналы можно считать пространственно когерентными с радиусом когерентности в пределах 300 м. В пределах радиуса когерентности можно наблюдать интерференцию отражённых сигналов от различных
участков слоистой структуры ледяного покрова.
М. Ю. Реушев, Д. С. Макаров, Д. В. Харламов и др.
15
В том случае, когда среднеквадратичное отклонение высот неровностей отражающей
поверхности меньше длины волны сигнала от НС, отражённое излучение испытывает малое диффузное рассеяние и, следовательно, незначительную флуктуацию фазы отражённого сигнала. Результат интерференции регистрируется приёмником в виде амплитудновременно́й зависимости (АВЗ) A(t), имеющей периодический характер.
Явление многолучевой интерференции проявляется в том, что на фоне гармонического
колебания основной частоты f1 , которая обусловлена сложением прямого сигнала от НС
и сигнала, отражённого от поверхности льда, появляется дополнительные частоты f2 ,
обусловленные проникновением части сигналов в слои снега и льда [6].
Несмотря на значительное количество публикаций исследований в области
ГНСС-рефлектометрии, остаются актуальными вопросы, связанные с методиками обработки результатов рефлектометрических измерений, пригодными для практического применения в качестве инструментальной базы мониторинга состояния ледовых покрытий.
Известно, что при решении ряда типовых задач обработки сигналов (таких как фильтрация, определение гармонического состава и спектральной плотности и др.) наряду с
данными об изменении сигналов во временно́й области АВЗ (A(t)) используются также
их частотные характеристики (A(f )). Эти две формы представления данных дают дополняющую друг друга информацию об исследуемых сигналах [7]. Математической операцией перехода от временно́й области к частотной является дискретное преобразование Фурье (ДПФ). В настоящее время ДПФ широко применяется при обработке результатов исследований в различных областях науки и техники. С целью развития методов
ГНСС-рефлектометрии были проведены рефлектометрические измерения ледовых поверхностей пресных и солёных водоёмов Красноярского края и Республики Хакасия. При обработке АВЗ применялась методика быстрого преобразования Фурье (БФП). Итоговые
результаты измерений и их обработка были использованы для определения связи между толщиной льда d, высотой антенны H, показателями преломления льда и воды n1 ,
n2 , длиной волны навигационного сигнала λ, скоростью угла возвышения навигационного
спутника над горизонтом ∆θ/∆t, частотами f2 и f1 .
Исследования и результаты. В период с 09.01.2025 г. по 30.03.2025 г. была проведена серия рефлектометрических измерений состояния ледовых поверхностей на озёрах Белё,
Шира (Республика Хакасия, Россия), а также на ледовой переправе в районе села Новосёлово (Красноярский край, Россия). Регистрация навигационных сигналов производилась
приёмником НСРП-04 с помощью антенны с правой круговой поляризацией. Обновление
данных об амплитуде сигналов происходила с частотой 1 Гц. Временно́й интервал сеанса
записи составлял около 90 мин. В основных экспериментах фазовый центр антенны находился на высоте 2,4 м от поверхности ледового покрова водоёмов. Параллельно во время
каждого сеанса производились измерения толщин и состояний структуры ледового покрова ручным способом. В качестве дополнительной информации использовались данные о
минерализации (TDS) водоёмов [8, 9]. В таблице приведены сводные данные об условиях
проведения исследований.
В качестве иллюстрации проведённых исследований на рис. 1, 2 даны примеры АВЗ
и БФП, полученные на всех исследуемых водоёмах 29.01.2025 г. и 30.01.2025 г.
В описании рисунков для НС приняты следующие обозначения: R — навигационная
система ГЛОНАСС; G — система GPS; α — азимутальный угол для НС относительно
точки фазового центра приёмной антенны; θ — угол возвышения НС над горизонтом;
∆θ/∆t — скорость изменения угла возвышения над горизонтом от времени.
Как можно заметить, из приведённых выше примеров рефлектометрических измерений разница частот ∆f = f2 − f1 отличается при разной толщине льда d. Общее решение
определения толщины ледового покрова является достаточно сложной задачей. В полном
16
АВТОМЕТРИЯ. 2025. Т. 61, № 5
Сводные данные условий рефлектометрических сеансов измерений
Дата проведения/
толщина льда, см
Минерализация
TDS, г/л
П/п
Место проведения
1
Ледовая переправа,
с. Новосёлово
09.01.2025
38
29.01.2025
54,5
19.02.2025
68
29.03.2025
76,8
0,08
2
О. Белё
09.01.2025
59
29.01.2025
74
19.02.2025
77,5
29.03.2025
89,5
8,8
3
О. Шира
10.01.2025
63
30.01.2025
76,5
29.02.2025
83
30.03.2025
95
11
a
A(t), äÁ
b
A(t), äÁ
4
f1 = 0,003153 Ãö
30
3
20
2
a, ãðàä
o, ãðàä
C/N, äÁ
10
1
Df = 0,003153 Ãö
f2 = 0,006307 Ãö
0
2000
3000
4000
5000
6000 t, c
0
0,005
0,010
0,015 f, Ãö
Рис. 1. Результаты рефлектометрических исследований (ледовая переправа, с. Новосёлово): a — АВЗ, b — БФП (R17, толщина льда d = 54,5 cм,
∆θ/∆t = −0,00764t + 46,555 град/с, R2 = 0,99991)
варианте решения должны быть учтены множественные параметры и условия проведения рефлектометрических измерений, такие как поляризация падающего и отражённого
излучения, макрофизическое состояние ледовой поверхности, диаграмма направленности
приёмной антенны, концентрация минеральных солей воды, температура воздуха, скорость ветра и т. д. Тем не менее для практического применения допустима упрощённая
версия решения в приближении геометрического подхода.
Рассмотрим условия задачи.
Навигационный спутник движется по орбите, излучая сигнал с длиной волны λ в диапазоне L1. Угол возвышения спутника над горизонтом изменяется с постоянной скоростью:
∆θ/∆t = const. Фазовый центр навигационного приёмника расположен на высоте H над
поверхностью льда толщиной d. Лёд имеет показатель преломления n1 , вода n2 .
После применения процедуры быстрого БФП к A(t) в спектре появляются характеристические частоты: f1 — соответствующая интерференции прямого сигнала и сигнала,
отражённого от поверхности льда; f2 — соответствующая интерференции прямого сигнала и сигнала, прошедшего через лёд и отразившегося от границы лёд — вода. Требуется
найти связь между толщиной льда d и другими параметрами: H, n1 , n2 , λ, ∆θ/∆t, f2 − f1 .
На рис. 3 приведена схема распространения лучей, используемая при решении задачи.
Здесь RF — радиус первой зоны Френеля как наиболее эффективный участок отражения лучей для получения интерференционной картины высокой контрастности [10].
17
М. Ю. Реушев, Д. С. Макаров, Д. В. Харламов и др.
a
70
f1 = 0,002648 Ãö
60
3
50
40
Df = 0,002648 Ãö
2
30
20
a, ãðàä
10
o, ãðàä
C/N, äÁ
0
b
A(f ), ãÖ
4
1000 2000 3000 4000 5000 6000 7000 8000 t, c
1
f2 = 0,005296 Ãö
0
0,005
0,010
0,015 f, Ãö
Рис. 2. Результаты рефлектометрических измерений (о. Шира): a — АВЗ,
b — БФП (G8, толщина льда d = 76,5 см, ∆θ/∆t = −0,00667t + 48,64 град/с,
R2 = 0,9995)
a
b
Do/Dt
HC
r0
r1
o r1
r2
o
H
d
o
o
0
n1
d
n2
RF
Рис. 3. Схема распространения лучей: a — общая геометрия лучей, b — геометрия распространения лучей в толщине льда
В первую очередь определим разность хода различных путей распространения сигнала.
1. Прямой сигнал от НС до приёмника. Обозначим длину пути как r0 .
2. Сигнал, отражённый от поверхности льда. Длину пути обозначим как r1 . Если
угол возвышения спутника над горизонтом θ, а высота фазового центра приёмника H, то
дополнительный путь, который проходит этот сигнал по сравнению с прямым сигналом,
составляет величину
∆r1 = r1 − r0 = 2H sin θ.
(1)
3. Сигнал, прошедший через лёд, отражается от границы лёд — вода и возвращается
через лёд к приёмнику. В связи с тем, что показатели преломления льда n1 ≈ 1,76 [11]
и n2 ≈ 9,2 [12] на частоте излучения навигационного сигнала 1,5 ГГц (L1) существенно
отличаются, упростим решение задачи, считая отражение Френеля от границы раздела
сред зеркальным. При падении луча на ледовую поверхность под углом θ он преломляется.
Угол падения связан с углом преломления θ0 по закону Снеллиуса: sin θ/ sin θ0 = n1 или
sin θ0 = sin θ/n1 .
18
АВТОМЕТРИЯ. 2025. Т. 61, № 5
Дополнительный путь в слое льда толщиной d равен
∆r2 = r2 − r1 = 2H sin θ + rd ,
(2)
где rd — оптический путь через лёд и обратно, который равен физической длине пути,
умноженной на показатель преломления n1 .
Физическая длина пути через лёд (в одну сторону) составляет величину d/ cos θ0 . Полный оптический путь через лёд и обратно определяется выражением dr =p2n1 d/ cos θ0 .
1 − (sin θ0 )2 =
Используем
закон Снеллиуса
для преобразования cos θ0 : cos θ0q=
q
q
1 − (sin θ)2 /n21 =
n21 − (sin θ)2 /n1 . Откуда dr = 2(n1 d/( n21 − (sin θ)2 /n1 ) =
q
= 2n21 d/ n21 − (sin θ2 ). Окончательное выражение для ∆r2 :
=
q
∆r2 = 2H sin θ + 2n21 d/ n21 − (sin θ)2 .
(3)
Известно, что частота интерференционных полос связана с изменением разности оптического хода во времени. При этом разность фаз между сигналами пропорциональна
разности хода

∆ϕ1 = 2π ∆r1 /λ;

q
(4)
∆ϕ2 = 2π ∆r2 /λ = 4πH sin θ/λ + 4πdn21 /(λ n21 − (sin θ)2 ). 
При изменении угла θ эти фазы будут меняться, создавая периодическую интерференцию. Частоты интерференции пропорциональны скорости изменения фазы
f1 ∝
d (∆ϕ1 )
d (∆ϕ1 ) dθ
cos θ dθ
dθ
=
= 4πH
= 2H(cos θ/λ) 2π ,
dt
dθ
dt
λ dt
dt
f2 ∝
d (∆ϕ2 )
d (∆ϕ1 )
=
+
dt
dt
q
d (4πn21 )/(λ n21 − (sin θ)2 )
dt
.
(5)
(6)
После вычисления производной от второго слагаемого для f2 и подстановки в уравнение (5) получим выражение
f2 = f1 +
2n21 d sin θ · cos θ dθ/dt
.
λ(n21 − (sin θ)2 )3/2
(7)
Из выражений (7) и (5) следует выражение для разности частот
2dn21 sin θ · cos θ dθ/dt
∆f = f2 − f1 =
.
λ(n21 − (sin θ)2 )3/2
(8)
Выражая d из (8), получим уравнение
d=
(f2 − f1 )λ(n21 − (sin θ)2 )3/2
.
2n21 sin θ · cos θ dθ/dt
(9)
Формула (9) определяет связь между толщиной льда d и параметрами: n1 , λ, ∆θ/∆t,
f2 − f1 . Выражение (9) не содержит значение H и n2 , так как является упрощённым
19
М. Ю. Реушев, Д. С. Макаров, Д. В. Харламов и др.
и ограниченным в рамках геометрического подхода при решении поставленной задачи.
Для точного решения необходимо учитывать множество других факторов и моделировать
разности хода и фазовые сдвиги с учётом кривизны земной поверхности, минерализации
водоёма, физического состояния ледовой поверхности и др.
Особенности применения формулы (9):
1. Показатель преломления воды n2 не входит в формулу (9), но влияет на коэффициент отражения от границы лёд — вода и, следовательно, на амплитуду сигнала.
2. Высота приёмника H не влияет на разность частот, но определяет абсолютные
значения f2 и f1 .
3. При практических измерениях необходимо учитывать угол возвышения θ в момент
измерения, что требует синхронизации с данными о положении спутника.
4. Для повышения точности измерений следует проводить калибровку системы на
участках с известной толщиной льда. Следовательно, формулу (9) можно записать в виде
d=K
f2 − f1
.
dθ/dt
(10)
В этом случае при создании оборудования для контроля толщины ледовых покровов
коэффициент K можно определить в ходе практических полевых испытаний как величину,
равную отношению
K=d
dθ/dt
,
f2 − f1
(11)
где значения dθ/dt и f2 − f1 измеряются при известных толщинах льда.
С целью определения общих закономерностей соотношений между толщиной льда d
и величинами ∆θ/∆t, f2 − f1 , а также влияния на них минерального состава воды TDS
была проведена процедура обобщения всех результатов рефлектометрических измерений
в виде последовательной зависимости толщины льда d от отношения k = (f2 − f1 )/(dθ/dt)
по возрастанию. Результаты такой процедуры приведены на рис. 4.
Как следует из рис. 4, зависимость толщины ледового покрова d от отношения коэффициента K аппроксимируется линейной функцией
d = −0, 07421k + 2, 59 м.
(12)
d, ì
1,0
d, ì
0,9
Ëèíèÿ èíòåðïîëÿöèè
0,8
0,7
0,6
0,5
0,4
0,3
23
24
25
26
27
28
29
30 k, Ãö . ñ/ðàä
Рис. 4. Обобщённые результаты обработки рефлектометрических измерений
при d = −0,0742k + 2,59 м; R2 = 0,8894
20
АВТОМЕТРИЯ. 2025. Т. 61, № 5
При аппроксимации зависимости (12) прямой линией (линейная интерполяция) среднеквадратичное отклонение не превышало 15 %.
Значительная погрешность не позволила выявить влияние минерального состава TDS
на рефлектометрические измерения. Такое влияние может проявляться, когда оптические
длины лучей, проходящих через лёд, существенно отличаются вследствие разницы в показателях преломления. Для оценки значений диэлектрической проницаемости льда были
выполнены расчёты с использованием результатов исследований, приведённых в работе [11]. Расчёты позволяют утверждать о том, что значения TDS заметно влияют на
затухание излучения, но слабо влияют на оптическую длину и соответственно на фазовые
сдвиги интерферирующих лучей (4).
Преобразуем уравнение (12) к следующему виду:
d∗ = d − 2, 59 = Kk,
(13)
где K, рад · м — коэффициент пропорциональности в уравнении (11).
Таким образом, для оценки толщины льда требуется коэффициент K умножить на
отношение k = (f2 − f1 )/(dθ/dt) и из полученного значения вычесть 2,59 м. Необходимые
условия применения эмпирической формулы (13) следующие:
1. Фазовый центр приёмной антенны должен быть расположен на уровне 2,4 м над
поверхностью льда;
2. Среднеквадратичное отклонение неровностей поверхности льда не должно превышать λ — длины волны излучения от НС.
Заключение. В ходе рефлектометрических измерений состояния ледовых поверхностей на озёрах Белё, Шира (Республика Хакасия, Россия), а также на ледовой переправе в районе села Новосёлово (Красноярский край, Россия) были измерены амплитудновременны́е зависимости A(t) сигналов от НС.
Результаты измерений были подвергнуты обработке с помощью математической процедуры БФП (A(t) → A(f )).
Для объяснения зависимостей была предложена модель многолучевой интерференции
в приближении геометрического подхода.
Результаты работы могут быть использованы в качестве методики для мониторинга
состояния ледовых поверхностей пресных и солёных водоёмов, в частности определения
толщины льда.
В связи с этим предложенная методика допускает погрешности определения толщины
льда порядка 15 %, требуется её совершенствование путём дополнительных экспериментальных и теоретических исследований.
СПИСОК ЛИТЕРАТУРЫ
1. Jin S., Camps A., Jia Y. et al. Remote sensing and its applications using GNSS reflected
signals: Advances and prospects // Satellite Navigation. 2024. 5. 19. DOI: 10.1186/s43020-02400139-4.
2. Козик В. И., Нежевенко Е. С. Сокращение размерности системы признаков при классификации гиперспектральных данных дистанционного зонирования Земли с помощью нейронных сетей // Автометрия. 2022. 58, № 1. С. 5–12. DOI: 10.15372/AUT20220101.
3. Прикладные аспекты исследования геосфер с использованием спутниковых технологий /В. Б. Кашкин, Т. В. Рублёва, К. В. Симонов и др. Красноярск: СФУ, 2023. 256 с.
4. Jin S., Cardellach E., Xie F. GNSS Remote Sensing. New York, London: Springer Dordrecht
Heidelberg, 2014. 286 p.
М. Ю. Реушев, Д. С. Макаров, Д. В. Харламов и др.
21
5. Перов А. И. Основы построения спутниковых радионавигационных систем. М.: Радиотехника, 2012. 236 с.
6. Reushev M. Yu., Podoprigora V. G., Makarov D. S. et al. Interaction of satellite navigation
signals with forest and snow – ice layers // Russ. Phys. Journ. 2023. 66, N 6. P. 706–711.
DOI: 10.1007/s11182-023-02996-z.
7. Айфичер Э. C., Джервис Б. У. Цифровая обработка сигналов: практический подход: Пер.
с англ. М.: Издательский дом «Вильямс», 2004. 992 с.
8. Zadereev E. S., Drobotov A. V., Tolomeev A. P. et al. The effect of salinity and nutrient
load on the ecosystems of selected lakes in the south of Siberia // Journ. Sib. Fed. Univ. Biol.
2021. 14, N 2. P. 133–153. DOI: 10.17516/1997-1389-0343.
9. Озёрский А. Ю. Геохимические особенности водной миграции в Красноярском водохранилище // Сб. ст.: «Вопросы географии Сибири». Томск: ТГУ, 2006. Вып. 26. С. 117–126.
10. Чёрный Ф. Б. Распространение волн. М.: Сов. радио, 1972. 463 с.
11. Vant M. R., Ramseier R. O., Makios V. The complex-dielectric constant of sea
ise at frequencies in the range 0.1-40 GHz // Journ. Appl. Phys. 1978. 49. P. 126–128.
DOI: 10.1063/1.325018.
12. Kiein L. A., Swift C. T. An improved model for the dielectric constant of sea water at microwave
frequencies // IEEE Trans. Antennas Propagation. 1977. AP-25, N 1. P. 104–111.
Поступила в редакцию 28.07.2025
После доработки 15.08.2025
Принята к публикации 15.09.2025