CC BY
16
Моделирование задачи низкоорбитальной спутниковой УФ-томографии ионосферы
И. А. Нестеровa, А.М. Падохинb, Е.С. Андреева, С. А. Калашникова
Московский государственный университет имени М. В. Ломоносова, физический факультет, кафедра физики атмосферы. Россия, 119991, Москва, Ленинские горы, д. 1, стр. 2.
E-mail: a nia2002@yandex.ru, b padokhin@physics.msu.ru
Статья поступила 02.12.2015, подписана в печать 14.01.2016.
Представлены результаты моделирования прямой и обратной задач низкоорбитальной спутниковой ультрафиолетовой (УФ) томографии скорости объемной эмиссии OI 135.6 нм в ионосфере. Для решения прямой задачи использовалась геометрия орбиты спутников DMSP, блок 5D3, несущих среди прочей научной аппаратуры УФ-спектрометры SSUSI и SSULI, реальные параметры работы этих приборов (скорость и интервал углов сканирования), а также набор модельных распределений скорости объемной эмиссии, содержащих неоднородности различных масштабов. Полученные в результате решения прямой задачи данные об интенсивности УФ-излучения на длине волны 135.6 нм использовались в качестве входных данных для восстановления исходных модельных распределений скорости объемной эмиссии. При этом для решения полученной системы линейных уравнений (СЛУ) применялись хорошо зарекомендовавшие себя в задачах низкоорбитальной радиотомографии ионосферы итерационные алгоритмы ART (algebraic reconstruction technique) и SIRT (simultaneous iterative reconstructive technique). Показано, что для успешного восстановления исходных модельных распределений необходимо учитывать условие неотрицательности решения, использовать весовые функции для уменьшения решения в областях, где оно априори мало, а также использовать межитерационное сглаживание для устранения влияния погрешности аппроксимации, причем параметр сглаживания должен уменьшаться в ходе итераций. При выполнении этих условий вычислительные затраты для реализации алгоритмов решения СЛУ ART и SIRT оказываются сопоставимыми, а погрешность реконструкции составляет около 6%. Проанализировано влияние случайной и систематической погрешности данных на результаты реконструкций и показано, что, зная уровень погрешности исходных данных, на этапе моделирования можно заранее подобрать параметры алгоритмов реконструкции так, чтобы эффективно подавить влияние на решение задачи шумов с относительными величинами порядка 2-3%.
Ключевые слова: ионосфера, томография, УФ-излучение, моделирование.
УДК: 528.81. PACS: 94.20.Tt.
Введение
Возможность использовать собственное свечение ночного неба в ультрафиолетовом и видимом диапазонах длин волн для исследования распределения электронов и ионов О + в Р-области ионосферы была описана в конце 1960-х — начале 1970-х гг. Подобные измерения были проведены в видимом диапазоне с Земли [1, 2] и в ультрафиолетовом диапазоне из космоса [3, 4]. Эти эксперименты выявили два возможных механизма, вызывающих собственное свечение ионосферы в ультрафиолетовом диапазоне: излучающую рекомбинацию ионов О + и электронов и нейтрализацию ионов О + окружающими их ионами О_ [5], причем свечение на длине волны 91.1 нм порождается только реакцией излучающей рекомбинации [6]. Более поздние измерения на средних широтах [7] подтвердили, что излучающая рекомбинация является основным механизмом для излучения и на длине волны 135.6 нм.
Использование собственных ночных эмиссий кислорода для дистанционного зондирования электронной концентрации Р -слоя ионосферы описано в работах [8, 9]. Определение профилей электронной
концентрации по данным УФ-спектрометрии было успешно продемонстрировано в работе [9], в которой использовались данные излучения на длинах волн 91.1 и 135.6 нм со спектрографа LORAAS (low-resolution airglow and aurora spectrograph) на спутнике ARGOS.
В последнее десятилетие для восстановления двух- и трехмерных распределений электронной концентрации в F -слое ионосферы по данным спутниковой УФ-спектрометрии стали активно применятся томографические методы [10-13]. Общей чертой данных работ является использование УФ-спек-трометров следующего поколения типа GUVI (global ultraviolet imager) на спутнике TIMED, SSUSI (spécial sensor ultraviolet spectrographic imagers) и SSULI (special sensor ultraviolet limb imagers) на спутниках DMSP 5D3 и попытка восстановления непосредственно распределения электронной концентрации, что требует привлечения дополнительной информации о фотохимии ионосферы. Кроме того, в большинстве подобных работ используется приближение чепменовского слоя для описания высотного профиля электронной концентрации
и простые алгоритмы алгебраической реконструкции (ART, MART...) для решения СЛУ задачи. Все это накладывает ограничения на качество получаемых реконструкций. Поэтому в последние годы все чаще решается задача томографического восстановления скорости объемной эмиссии на выбранной длине волны в ионосфере, не требующая дополнительной информации о физических механизмах, вызывающих данное свечение. Данной проблематике посвящена и настоящая работа, в которой проведено моделирование и оптимизация алгоритмов задачи низкоорбитальной спутниковой УФ-томографии скорости объемной эмиссии OI 135.6 нм в ионосфере для реальной геометрии орбит спутников DMSP 5D3 и режимов работы установленных на них УФ-спек-трометров SSUSI и SSULI.
1. Описание инструментов и параметров орбит спутников
Американские метеорологические спутники последнего поколения DMSP 5D3 со спектрометрами SSUSI и SSULI на борту начали запускать с 2003 г. К настоящему времени запущено 4 спутника, последний — в апреле 2014 г. Все спутники работают на практически круговых (высотой ~ 850 км) солнечно-синхронных полярных орбитах с высоким наклонением (~ 99°).
Первый прибор, параметры которого используются в моделировании, — УФ-спектрометр SSULI, разработанный в Naval Research Laboratory, позволяющий измерять излучение лимба Земли в УФ-диа-пазоне 80-170 нм с разрешением 1.5 нм. Сканирование производится вдоль плоскости пролета спутника с интервалом 90 с в угловом диапазоне от 10 до 27° ниже траектории движения спутника.
Вторым источником данных является SSUSI — УФ-спектрометр пространственного сканирования, разработанный в Johns Hopkins University Applied Physics Laboratory (JHU/APL). Сканирование производится поперек траектории полета спутника в дальнем УФ-диапазоне и позволяет получать измерения для пяти заданных интервалов длин волн: линия Н + 121.6 нм, линии О + 130.4 и 135.6 нм и два интервала N2 Лаймана-Бриджа-Хопфилда 140-150 и 165-180 нм. Результатом каждого сканирования являются 2 профиля светимости: в направлении лимба (угловой диапазон — 72.8... —63.2° от вертикали) и от горизонта до горизонта. Пространственное разрешение при сканировании лимба составляет 0.4°, а при сканировании в направлении Земли — 0.8°. При моделировании в работе были использованы значения светимости, полученные строго вертикально (под спутником). Эти данные доступны раз в 22 с.
2. Постановка задачи
Рассмотрим постановку задачи для 2D-томо-графии по УФ-данным. Реконструкции подлежит
высотный разрез скорости объемной эмиссии 01 135.6 нм в ионосфере вдоль траектории спутника. Исходя из параметров орбиты спутников ЭМБР 5Э3 и режимов работы спектрометров ББиБ! и ББии геометрия лучей, вдоль которых осуществляется сканирование области реконструкции в координатах широта-высота за время полета спутника 33 мин, представлена на рис. 1 (вверху).
В качестве модельного распределения скорости объемной эмиссии возьмем распределение, представленное на рис. 1 (внизу), качественно отражающее высотный профиль, широтный ход, наличие экваториальной аномалии (ЭА), провалов, а также добавленные в него дополнительные локальные структуры: узкие слои, располагающиеся как ниже, так и выше максимума слоя (здесь и далее по тексту для всех графиков высотно-широтных распределений скорости объемной эмиссии используется единая цветовая шкала, нормированная на максимальное значение модельного распределения).
Интенсивность собственного свечения ионосферы ¡к на длине волны 135.6 нм вдоль направления сканирования 1к определяется как
е(ф(г), h(r))exp —
— p(r' ) dl^j
p(r')dl') dl = 4nIk, (1)
где е(ф(г), Н(т)) — скорость объемной эмиссии излучения 01 135.6 нм, р(г') — коэффициент поглощения УФ-излучения на длине волны 135.6 нм в ионосфере.
Отметим, что скорость объемной эмиссии на длине волны 135.6 нм связана с распределением электронной концентрации [е], концентраций ионов [0+] и атомов кислорода [0] [2, 14]
k1k2^1356[O][e][O+] k2[O+]+ k3[O]
+ û!1356[e][O+],
где коэффициенты реакций к1, к2 и к3 приблизительно равны 1.3 • 10~15 см3/с, 1.0 • 10~7 см3/с и 1.4 • 10~10 см3/с соответственно; вклад реакции нейтрализации в возмущенное состояние (05Б), соответствующее излучению на длине волны 135.6 нм, Азб6 = 0.54 [2]; коэффициент для реакции излучающей рекомбинации а1356 зависит от температуры и при 1160 К приблизительно равен 7.5-10~13 см3/с. Принимая во внимание характерные значения концентраций электронов, ионов и атомов кислорода на ионосферных высотах, скорость объемной эмиссии на длине волны 135.6 нм пропорциональна квадрату электронной концентрации. Таким образом, можно рассматривать томографическую задачу (1) относительно неизвестного распределения электронной концентрации, как делается в ряде работ [10-13]. Однако сам коэффициент пропорциональности между скоростью объемной эмиссии и квадратом электронной концентрации может как зависеть от высоты, так и меняться в соответствии с лоальным временем и сезоном [14]. В связи
20 0 Широта
Модель
Рис. 1. Геометрия лучей в координатах широта-высота (вверху) и комплексная модель для реконструкции (внизу), содержащая ионосферные структуры различных типов
с этим разумным представляется томографическое восстановление распределения скорости объемной эмиссии е(ф(г), Н(т)), не требующее дополнительной информации о физических механизмах, вызывающих данное свечение, с последующим восстановлением распределения электронной концентрации по полученному распределению е(ф(г), Н(г)).
Для получения модельных данных интегралы (1) исходно рассчитываются численно вдоль всех направлений сканирования (рис. 1, вверху) для заданной модели е = ето^(ф, Н) (рис. 1, внизу). Интегрирование производится методом прямоугольников вдоль известных лучей с разбиением каждого луча на Ыт равных интервалов в декартовой системе координат. Полученные данные у0 = 4п1т°А являются входными для реконструкции распределения скорости объемной эмиссии е(ф(г), Н(г)).
Будем рассматривать искомую реконструируемую функцию как разложение по некоторому конечному базису функций {В} с неизвестными числовыми коэффициентами {х}:
(2)
е(ф, = х^В^(ф, Н).
ч
Будем использовать систему базисных функций, обеспечивающих билинейную интерполяцию значений реконструируемой функции в узлах сетки, размерности Ыф х ЫН по широте и высоте соответственно. Тогда (1) можно представить в виде
е(ф(г), Н(г)) ехр -
р(г') сС/^ (И-,
«Яф.
1<ЫН
П<Мт
Вч Ш), 4^)) ехр-екп) Д/0 = уи,
(3)
п=0
гП = * +
где вП « ) р( ги0) Д/ + £ П^ р( Л) ; + (Во\ - (п + 0.5)/Ыт, п = 0,..., Ят-1; Д/0 = = 1*2 - *0|/Ыт, а Я\, *2 — начальная и конечная точки луча /0.
Сопоставляя каждому узлу сетки уникальный индекс (', 1) ^ т, получим СЛУ относительно коэффициентов разложения {хт}
£А0тХ т = у0 ,
(4)
или, в векторной записи, с учетом погрешности аппроксимации и исходных данных (:
Ах = у + С.
(5)
Не имея строгих моделей, позволяющих учесть эти погрешности, явно модифицировав форму СЛУ, будем ставить задачи на минимизацию невязки (4), обращая, однако, существенное внимание на выбор алгоритмов минимизации. При этом для сравнения результатов будем следить за успешностью реконструкции модельного распределения, качественно и количественно контролируя отклонение результата реконструкции от модели. Для отработки алгоритмов на первом этапе можно пренебречь вкладом поглощения в (3) и учесть его на последнем этапе, используя оптимально подобранные на первом этапе параметры алгоритмов.
3. Используемые алгоритмы решения СЛУ
Для решения СЛУ задачи (4) в данной работе применялись итерационные алгоритмы ART и SIRT, хорошо зарекомендовавшие себя в задачах радиотомографии ионосферы [15, 16] и способные эффективно (с точки зрения объема вычислений) решать СЛУ с разреженными матрицами, не имеющими специальной структуры (диагональной, ленточной).
Идея алгоритма ART заключается в том, чтобы последовательно перебирать уравнения системы и для каждого из них обращать в ноль невязку уравнения путем добавления к вектору решения добавочного вектора, пропорционального строке матрицы системы, задающей данное уравнение:
yk - (ak, xn)
xn+1 = xn + ak-
= xn + pkak--
(au, ak)
(7)
Выбрав рк = ррт, получим для алгоритма БШТ наиболее простой и прозрачный вариант
xn+1 = xn +
An
NAH2
AT (y - Axn
(8)
результат реконструкции с точки зрения работы алгоритмов ART и SIRT. Использование малых значений весов позволяет подавить «артефакты» решения в тех областях, где решение заведомо должно быть мало. В то же время в тех областях, где могут содержаться важные реконструируемые структуры, веса должны быть близки к единице, чтобы не вносить искажений в реконструкцию. Пример подобной весовой функции приведен на рис. 2.
( ) (6) (ak, ak)
(из чего непосредственно следует, что (ak, xn+1) = = yk). Если последовательность итераций для алгоритма ART сходится, то она сходится к нормальному (относительно начального приближения) решению задачи на квазирешение системы (4).
Теперь рассмотрим в той же ситуации применение итерационного алгоритма SIRT, который имеет усредняющий характер. Добавки к решению, полученные для каждого из уравнений в методе ART, в данном случае вначале суммируются с некоторыми весами pk и лишь затем прибавляются к итоговому решению:
yk - (ak, xn)
соответствующим естественному итерационному алгоритму поиска нормального решения для квазирешения СЛУ min NAx - yN2.
Скорость сходимости SIRT на несколько порядков меньше, чем ART (сходимость является теоретически гарантированной для 0< An <2). С целью получения более быстрой сходимости итераций параметр tn = jA^ можно выбирать из соображений минимизации невязки СЛУ на каждой итерации.
Важной особенностью итерационных алгоритмов является возможность учета условия неотрицательности решения (скорости объемной эмиссии) на каждом шаге итераций. Для этого будем занулять все отрицательные элементы x после каждой итерации. Еще одна удобная возможность при использовании итерационных алгоритмов — введение в постановку задачи весовых коэффициентов, позволяющих скорректировать высотный профиль скорости объемной эмиссии. Тогда от набора переменных x можно перейти к набору переменных 1 таких, что xm = 1m • wm, где веса wm могут модифицировать
0.4 0.6 0.8
тЩ
Рис. 2. Квадрат весовой функции, w2(h)
Таким образом, в интервале высот от 150 до 450 км никаких искажений не вносится. Ниже этого интервала решение подавляется при подходе к уровню земли (что физически совершенно оправдано), а выше — плавно уменьшается с приближением к высоте пролета спутника (что качественно соответствует ходу профиля электронной концентрации в ионосфере).
4. Результаты моделирования
Приведем результаты восстановления модельного распределения скорости объемной эмиссии, показанной на рис. 1 (внизу), по набору данных сканирования вдоль направлений, показанных на рис. 1 (вверху). На данном этапе пренебрежем поглощением. Размерность сетки х N для задания базисных функций выберем равной 60 х 60, что, с одной стороны, достаточно для разрешения всех неоднородностей, включенных в модель, и, с другой стороны, плотность сетки не превышает плотности заполнения реконструируемой области лучами (между лучами каждого пучка нет «пустых» ячеек сетки, не пересекаемых хотя бы одним из них). Начальное приближение выберем нулевым. Погрешность реконструкции относительно исходного модельного распределения будем рассматривать в нормах С (\\х\\с = шах,- \х,\), ¿1 (\\х\\и = ¿\ \х,\)
и ¿2 (\\х\к = ).
Результаты реконструкций представлены на рис. 3. Сходимость решений по невязке СЛУ изображена на рис. 4. Отклонения решения от модели
ART (условие неотрицательности, с весами)
0 -20 Широта
SIRT (условие неотрицательности, с весами)
0 -20 Широта
Рис. 3. ART (вверху) и SIRT (внизу) с условием неотрицательности и весовой функцией, представленной
на рис. 2
80 120 Итерация
4000 6000 8000 10000 Итерация
Рис. 4. Изменение невязки СЛУ в ходе итераций для методов ART и SIRT с условием неотрицательности
и весовой функцией, представленной на рис. 3
для ART и
составляют
||е - етобУс/\\emodУс = 14.3%,
IIе - етобУ^/ ||етоб\\ь{ = 13.2%,
\\е - етоб\\ь2/||етоб\\^2 = 12.7%
\\е - етоб\\с/||етоб\\с = 16.°%, \\е - етоб У^/||етоб \= 13.3%, \\е - етоб\\ь2/||етоб\\^2 = 12.4%
для SIRT. Видно, что хорошо восстанавливается общий уровень решения, положение структур (слои, положение горбов ЭА, положение провала в северном полушарии). Отметим, что провал в южном полушарии не восстанавливается из-за отсут-
ствия в данной широтной области вертикальных лучей сканирования (модельные данные спектрометра SSUSI) (см. рис. 1, вверху), что говорит о важности этих лучей для разрешающей способности метода.
4.1. Межитерационное сглаживание решения
Поскольку, при неплохом количественном соответствии, основным качественным недостатком реконструкции на рис. 3 является наличие «артефактов», связанных с грубостью аппроксимации задачи («шахматная» структура, через которую проявляется сетка), следующий естественный ход по улучшению реконструкции — межитерационное сглаживание решения. Между последовательными итерациями алгоритмов ART и SIRT будем применять к промежуточному решению сглаживающий фильтр,
производящий усреднение решения по ближайшей окрестности каждого узла сетки:
¿+1 j+i
xi,j = pii-i,ii-ixii,ii, ¿i=i-i ji=j-i
где
p ii- i,ji-j
ц
¿i- i,ji-j
(9) (i0)
^¿2 = —1 ^/2=/-1 Р,2-,,12-1 а семейство симметричных по индексам и по координатам сглаживающих фильтров с параметром сглаживания р зададим как
= p'i(p)pj (p), p'i(p) =
i, i' = 0,
(ii)
[p, V\ = 1,
что автоматически обеспечит нормировку (9).
ART (условие неотрицательности, с весами, сглаживание с изменяющимся параметром)
Отметим, что параметр сглаживания p разумно уменьшать от итерации к итерации, чтобы вклад сглаживания в решение не превышал приращения к решению в алгоритмах ART и SIRT на данной итерации.
На рис. 5 показаны результаты реконструкции модельного распределения с использованием межитерационного сглаживания с экспоненциально убывающим параметром p . Сходимость решений по невязке СЛУ изображена на рис. 6. Визуально реконструкции для ART и SIRT в области сравнения практически неразличимы с моделью. Более того, конечная невязка СЛУ для SIRT уменьшилась до того же порядка величины, который дает ART. Количественные оценки при этом дают следующие значения.
О
Широта
SIRT (условие неотрицательности, с весами, сглаживание с изменяющимся параметром)
20 О Широта
Рис. 5. ART (вверху) и SIRT (внизу) с межитерационным сглаживанием решения с переменным параметром
сглаживания
SIRT
80 120 Итерация
200
400 600 Итерация
800 1000
Рис. 6. Изменение невязки СЛУ в ходе итераций для методов ART (слева) и SIRT (справа) с межитерационным сглаживанием решения с переменным параметром сглаживания
Для SIRT:
Для ART:
\\£ - Emodllc/||£mod||C = 6.78%, № - £mod\ii/\emod\ii = 6.71%, \\£ - £mod\L2/\£mod\L2 = 5.92%.
\\£ - £mod\c/\\emod\c = 7.94%, \\s - £mod\ii/\\£mod\k = 6.64%,
\\£ - EmodH^/\\£mod\L2 = 5.85%.
Таким образом, для получения хорошего результата реконструкции в задаче УФ-томографии (при применении стандартных для томографии ионосферы алгоритмов) необходимо использовать:
1) условие неотрицательности решения;
2) весовые функции для уменьшения решения в областях, где электронная концентрация априори мала;
3) межитерационное сглаживание для устранения влияния погрешности аппроксимации (при аппроксимации первого порядка), причем параметр сглаживания должен уменьшаться в ходе итераций.
При выполнении этих условий вычислительные затраты для реализации алгоритмов решения СЛУ ART и SIRT оказываются сопоставимыми, а погрешность реконструкции составляет около 6%.
4.2. Влияние шумов данных на точность реконструкции
Предыдущие результаты получены без специального добавления шумов в исходные данные задачи.
0.000
0.020
Рис. 7. Ошибка реконструкции методом ART модельного распределения в зависимости от уровня шума в данных
После выбора оптимальных алгоритмов реконструкции исследуем влияние на их работу дополнительных шумов, вносимых в данные. Рассмотрим задачу на примере алгоритма ART. Зафиксируем размерность сетки 120 х 120 и число итераций 2000 и будем вносить в исходные данные шум с уровнем а, меняющимся в диапазоне 0-2.0%, анализируя при этом погрешность решения в нормах C, Lj, L2. Результаты представлены на рис. 7. Можно сделать вывод, что случайные погрешности данных слабо влияют на результаты работы тестируемых алгоритмов при относительной величине погрешности, не превышающей 0.5%. При увеличении погрешности до 2% на реконструкциях появляется сильный шум (рис. 8, вверху), накладывающийся поверх исходных структур модели, хотя сами эти структуры остаются при этом ясно различимыми.
ART, CT=2.0%
20 0 -20 Широта
-60 -80 ART, о = 2.0%
Рис. 8. Реконструкция модельного распределения по зашумленным с а = 2.0% данным методом ART (вверху)
и методом ART с подавлением шумов (внизу)
Идеи сглаживания могут быть применены также и для подавления влияния шумов в данных. Для этого необходимо модифицировать сглаживающие алгоритмы таким образом, чтобы согласовать степень сглаживания решения с погрешностью исходных данных. Пусть в процессе итераций n = 1,...,N^r параметр сглаживания pn меняется от P1 до P2 в логарифмическом масштабе:
Pn = exp(ln(Pi ) + (ln(P2) - ln(Pi ))(n - 1) j(Nter - 1)).
(12)
Кроме того, поскольку сглаживание шумов необходимо производить на пространственных масштабах больших, чем расстояние между соседними лучами, сглаживающий фильтр к решению необходимо применять многократно, а также производить сглаживание не только между итерациями, но и после осуществления последней итерации.
Тогда зафиксировав P1, число итераций и число последовательных применений сглаживающего фильтра, варьируя P2, можно определить значение этого параметра, обеспечивающее минимальные ошибки реконструкции. Результаты применения данного подхода для P1 = 0.5, P2 = 0.03 представлены на рис. 8 (внизу).
Этот пример показывает, что, зная уровень погрешности исходных данных, на этапе моделирования можно заранее подобрать параметры алгоритмов реконструкции так, чтобы эффективно подавить влияние на решение задачи случайного шума с относительными величинами порядка 2-3%.
4.3. Учет поглощения
Основным источником поглощения собственного свечения ионосферы на длине волны 135.6 нм является поглощение молекулярным кислородом (полоса Шумана-Рунге). Таким образом, в первом приближении можно считать, что коэффициент поглощения р зависит только от высоты, за счет высотной зависимости концентрации O2 и зависимости сечения поглощения от температуры. Для моделирования примем, что форма высотного профиля коэффициента поглощения р соответствует представленной на рис. 9. Максимальное значение коэффициента поглощения pmax будет параметром, задающим величину поглощения.
Рассмотрим влияние pmax на качество реконструкции методом ART с межитерационным сглаживанием решения (рис. 10). Использовались следующие параметры: сетка 120 х 120, P1 = 1,, P2 = 3 • 10-6, 2000 итераций. График зависимости погрешностей решения в нормах C, L1, L2 от максимального значения коэффициента поглощения pmax показан на рис. 11.
Можно видеть, что искажения решения, вызванные влиянием поглощения:
а) незначительны при pmax < 5 • 10-7 м -1;
б) становятся заметными при pmax = 5 • 10-7 м -1;
1000 -1
800 -
2 и 600 -
400 -
200 -
0 -
PW/Pwbjl
Рис. 9. Модельное задание профиля коэффициента поглощения
в) становятся велики при pmax = 2 • 10-6 м-1 (ЭА сильно несимметрична).
Влияние поглощения на результат томографической реконструкции становится заметным тогда, когда для квазигоризонтальных лучей вклад интеграла по второй, восходящей части луча становится мал по сравнению с вкладом интеграла по первой, нисходящей части луча. Помимо потери части информации это приводит к тому, что геометрия задачи становится несимметричной, что порождает, в частности, асимметрию при реконструкции ЭА: горб ЭА, наклоненный в ту же сторону, что и нисходящие части лучей, восстанавливается лучше, чем горб, наклоненный в ту же сторону, что и восходящие части лучей.
Отметим, что при pmax = 2 • 10-6 м-1 уже для первой половины луча происходит ослабление излучения в 25 раз, а при pmax = 1 • 10-6 м-1 — в 5 раз. Понятно, что такое сильное поглощение ставит под вопрос саму идею томографической реконструкции и получение для этих значений поглощения хотя бы качественно адекватных результатов демонстрирует большие возможности разработанных методов. Слабым поглощение можно считать лишь для pmax = 1 • 10-7 м-1, а в этом случае точность реконструкции остается практически такой же, как и при отсутствии поглощения.
Заключение
Проведенное численное моделирование задачи УФ-томографии ионосферы для реальной геометрии орбит спутников DMSP 5D3, и режимов работы установленных на этих аппаратах спектрометров SSUSI и SSULI показало, что при использовании хорошо известных и зарекомендовавших себя в задачах радиотомографии ионосферы итерационных алгоритмов типа ART и SIRT возможно успешно восстанавливать параметры распределения скорости объемной эмиссии OI 135.6 нм в ионосфере. При этом постановка томографической задачи относительно скорости объемной эмиссии позволяет избежать дополнительных трудностей, связанных с учетом фотохимических реакций, вызывающих
ART
800
600
пГ
8 о 400
Я
PQ 200
0
20 0 Широта
Рис. 10. ART для задачи c поглощением: pmax = 1 • 10~7 (а), 5 • 10~7 (б), 1 • 10~6 (в), 2 • 10~6 (г)
собственное свечение ионосферы, что необходимо в случае постановки задачи относительно электронной концентрации. Вместе с тем результаты моделирования показали, что алгоритмы ART и SIRT требуют модификации для использования в задаче УФ-томографии ионосферы: необходимо учитывать условие неотрицательности решения, использовать весовые функции для уменьшения решения в областях, где оно априори мало, а также использовать межитерационное сглаживание для устранения влияния погрешности аппроксимации, причем параметр
сглаживания должен уменьшаться в ходе итерации. Кроме того, анализ влияния шумов на результаты реконструкции показал, что, зная уровень шумов исходных данных, на этапе моделирования можно заранее подобрать параметры алгоритмов реконструкции так, чтобы эффективно подавить влияние на решение задачи случаиного шума с относительными величинами порядка нескольких процентов. Разработанные алгоритмы также показали свою эффективность в случае учета поглощения УФ-излучения в ионосфере.
РмЯКС ' М
Рис. 11. Зависимости погрешности решения (в нормах C, Li, L2) от максимального значения коэффициента поглощения pmax
Работа выполнена при финансовой поддержке РФФИ (грант 15-35-21065-мол_а_вед).
Список литературы
1. Tinsley B.A., Christensen A.B., Bittencourt J.A. et al. // J. Geophys. Res. 1973. 78. P. 1174.
2. Tinsley B.A., Bittencourt J.A. //J. Geophys. Res. 1975. 80. P. 2333.
3. Hicks G.T., Chubb T.A. // J. Geophys. Res. 1970. 75. P. 6233,
4. Barth C.A., Schaffner S. // J. Geophys. Res. 1970. 75. P. 4299,
5. Knudsen W.C. // J. Geophys. Res. 1970. 75. P. 3862.
6. Tinsley B.A. // Ann. Geophys. 1972. 28. P. 155.
7. Brune W.H., Feldman P.D., Anderson W.G. et al. // Geophys. Res. Lett. 1978. 5. P. 383.
8. Meier R.R. // Space Sci. Rev. 1991. 59.
9. Dymond K.F., Thonnard S.E., McCoy R.P., Thomas R.J. // Radio Sci. 1997. 32. P. 1985.
10. Dymond K.F., Thomas R.J. // Radio Sci. 2001. 36. P. 1241.
11. DeMajistre R., Paxton L.J., Morrison D. et al. // J. Geophys. Res. 2004. 109. A05305.
12. Comberiate J.M., Kamalabadi F., Paxton L.J. // Radio Science. 2007. 42, N 2. RS2011.
13. McMahon E.M., Comberiate J.M., Kelly M.A., Paxton L.J. // 1st AIAA Atmospheric and Space Environments Conference, 22-25 June 2009, San Antonio, Texas, 2009.
14. Rajesh P.K., Liu J.Y., Hsu M.L. et al. // J. Geophys. Res. 2011. 116. A02313.
15. Kunitsyn, V., Tereshchenko, E. Ionospheric Tomography. B.; Heidelberg: Springer-Verlag, 2003.
16. Nesterov I.A. Kunitsyn V.E. // Adv. Space Res. 2011. 47, N 10. P. 1789.
Modeling the problem of low-orbital satellite UV-tomography of the ionosphere I. A. Nesterova, A.M. Padokhinb, E.S. Andreeva, S. A. Kalashnikova
Department of Physics of Atmosphere, Faculty of Physics, Lomonosov Moscow State University, Moscow 119991, Russia.
E-mail: a nia2002@yandex.ru, b padokhin@physics.msu.ru.
The results of modeling the direct and inverse problems of low-orbital satellite ultraviolet (UV) tomography of the ionospheric 135.6 OI volume emission rate are presented. The direct problem was solved with the orbital geometry of DMSP block 5D3 satellites with SSUSI and SSULI UV spectrographs among the other payloads, the real operating parameters of these instruments (the scan rate and the interval of scan angles), and the set of the model distributions of the volume emission rate that contain irregularities on various scales. The solution of the direct problem yields the radiation intensities in the 135.6 nm line, which is used as the input data for reconstructing the initial (prototype) model distributions of the volume emission rates. The obtained system of linear equations (SLE) was solved using the Algebraic Reconstruction Technique (ART) and Simultaneous Iterative Reconstructive Technique (SIRT) algorithms, which are highly efficient in problems of the low-orbit radio tomography of the ionosphere. It is shown that the initial model distribution can be successively reconstructed if one takes the non-negativity condition of the solution into account, uses weighting functions to decrease the solution in the regions where it is known to be a priori small, and applies inter-iteration smoothing to eliminate the effects of the approximation errors. Here, the averaging parameters should decrease in the course of the iterations. With these constraints fulfilled, the computational costs of the ART- and SIRT-based solutions are similar, while the reconstruction error is approximately 6%. The influence of random errors and bias in the data on the results of the reconstruction is explored. It is shown that with a given error level of the initial data the parameters of the reconstruction algorithms can be adjusted in such a way as to efficiently suppress the influence of the noise with a relative amplitude of 2-3% on the solution.
Keywords: ionosphere, tomography, UV radiation, modeling.
PACS: 94.20.Tt.
Received 02 December 2015.
English version: Moscow University Physics Bulletin. 2016. 71, No. 3. Pp. 329-338.
Сведения об авторах
1. Нестеров Иван Анатольевич — канд.физ.-мат. наук, доцент; тел.: (495) 939-15-41, e-mail: nia2002@yandex.ru.
2. Падохин Артем Михайлович — канд. физ.-мат. наук, ст. науч. сотрудник; тел.: (495) 939-20-89, e-mail: padokhin@physics.msu.ru.
3. Андреева Елена Станиславовна — канд. физ.-мат. наук, доцент; тел.: (495) 939-20-89, e-mail: es_andreeva@mail.ru.ru.
4. Калашникова Светлана Андреевна — аспирант; тел.: (495) 939-20-89, e-mail: kalashnikovasv.a@gmail.com.
