1. На главную
  2. Электронные версии
  3. Исследование Земли из Космоса
  4. № 2, 2022

Исследование Земли из Космоса, 2022, № 2, стр. 88-94

Двухканальный алгоритм определения аэрозольной оптической толщины над морской поверхностью по данным прибора МСУ-МР спутника “Метеор-М” № 2

М. О. Кучма a, *, Е. И. Холодов a, Ю. А. Амельченко a

a Дальневосточный центр Федерального государственного бюджетного учреждения “Научно-исследовательский центр космической гидрометеорологии "Планета”
Хабаровск, Россия

* E-mail: m.kuchma@dvrcpod.ru

Поступила в редакцию 13.09.2021

Полный текст (PDF)

Аннотация

В настоящей работе рассматривается задача определения аэрозольной оптической толщины над морской поверхностью по данным коротковолновых каналов многозонального сканирующего устройства малого разрешения космического аппарата “Метеор-М” № 2. Представленный алгоритм основан на использовании оптической модели, в качестве которой выступает шестимерная справочная таблица, содержащая информацию о коэффициентах спектральной яркости каналов спутникового прибора для различных атмосферных условий и геометрии наблюдений. Проведена валидация результатов, показавшая корреляцию 80% с эталонными данными сети измерений AERONET.

Ключевые слова: дистанционное зондирование Земли, аэрозоль, аэрозольная оптическая толщина, МСУ-МР, справочная таблица, AERONET

ВВЕДЕНИЕ

Аэрозоли играют важную роль в формировании регионального и глобального климата за счет взаимодействия с солнечным излучением и облачностью. Сложность исследований аэрозоля обусловлена разнообразием его дисперсного и химического состава, процессов возникновения, рассеяния и поглощения. Для лучшего понимания свойств аэрозоля активно разрабатываются системы регулярных долгосрочных измерений аэрозольной оптической толщины (далее – АОТ) во всех континентальных регионах планеты. АОТ является важной составляющей при оценке климатических изменений, например, изменений температуры подстилающей поверхности (Асарин и др., 2014), а также приходящей суммарной солнечной радиации (Stocker at al., 2013).

Важную роль в увеличении количества и качества регулярных измерений АОТ играют разработанные методы его восстановления по спутниковым наблюдениям, которые необходимо развивать, в том числе, и для новых космических аппаратов (далее – КА). В настоящее время российская группировка полярно-орбитальных КА состоит из двух аппаратов – “Метеор-М” № 2 и “Метеор-М” № 2-2, а до 2024 года предполагается наличие на орбите уже шести КА этой серии.

АОТ, получаемая по спутниковым данным КА серии “Метеор-М”, необходима, в том числе, для проведения атмосферной коррекции коротковолновых каналов многозонального сканирующего устройства малого разрешения (МСУ-МР) (Kuchma, Bloshchinskiy, 2020), установленного на этих КА.

Проблема в настоящее время заключается в отсутствии подхода к определению параметра АОТ по данным МСУ-МР. В рамках настоящей работы впервые представлен метод определения параметра АОТ над морской поверхностью по данным МСУ-МР российского КА “Метеор-М” № 2.

ОБЗОР ИССЛЕДОВАНИЙ

Наиболее распространенными методами расчета АОТ над морской поверхностью на сегодняшний день являются алгоритмы для спутниковых приборов VIIRS и MODIS, установленных на КА Suomi NPP/NOAA-20 и TERRA/AQUA соответственно (Ray, 2014; Levy et al., 2009). Для расчета АОТ над морской поверхностью по данным VIIRS используются данные каналов с центральными длинами волн 0.67, 0.74, 0.86, 1.24, 1.62 и 2.25 мкм, по данным MODIS – 0.47, 0.55, 0.66, 0.86, 1.2, 1.6 и 2.12 мкм. Ввиду отсутствия у прибора МСУ-МР некоторых из перечисленных каналов, эти алгоритмы невозможно использовать для решения данной задачи.

Для реализации методики восстановления АОТ применимо использование алгоритма (Zhao, 2017), который используется в оперативной практике National Centers for Environmental Information для построения суточных композитных карт АОТ над морской поверхностью по данным прибора AVHRR, установленного на КА серии NOAA и MetOp. Функции спектральной чувствительности коротковолновых каналов прибора МСУ-МР и прибора AVHRR имеют сходство (рис. 1), что также подтверждает возможность применения этой методики определения АОТ применительно к данным прибора МСУ-МР.

Рис. 1.

Функции спектральной чувствительности коротковолновых каналов МСУ-МР и AVHRR.

Суть данного алгоритма заключается в использовании справочных таблиц (LUT), представляющих собой смоделированные измерения коэффициентов спектральной яркости (КСЯ) в каналах спутникового прибора с центральными длинами волн 0.63, 0.83 и 1.61 мкм для различных значений АОТ морского аэрозоля с учетом атмосферных условий, геометрии наблюдения и освещения. Смоделированные значения КСЯ из LUT для конкретной схемы наблюдений сопоставляются с измерениями в каналах спутникового прибора. С помощью метода минимизации среднеквадратической невязки между рассчитанными и измеренными КСЯ осуществляется конкретный выбор значений АОТ.

СОЗДАНИЕ СПРАВОЧНОЙ ТАБЛИЦЫ

Так как функции спектральной чувствительности коротковолновых каналов прибора МСУ-МР широкие (см. рис 1), то они перекрывают часть электромагнитного спектра, в котором молекулы озона и водяного пара имеют сильные полосы поглощения (Тимофеев, Васильев, 2007) (рис. 2). В результате чего происходит ослабление регистрируемого в каналах спутникового прибора электромагнитного излучения.

Рис. 2.

Спектры солнечного излучения, приходящего на верхнюю границу атмосферы и на поверхность Земли.

Из этого следует, что для учета общего содержания этих параметров необходимо использовать шестимерную LUT таблицу, в которой присутствуют различные значения водяного пара и озона, задаваемые в модели переноса излучения 6S (Vermote et al, 1997).

С использованием быстрой радиационной модели на основе программного кода DISORT (Buras et al., 2011), входящего в состав библиотеки libradtra (Emde et al., 2016), был произведен расчет КСЯ в коротковолновых каналах прибора МСУ-МР для морской поверхности. Расчет проводился для различной геометрии наблюдения при зенитном угле Солнца 40 градусов. Согласно полученным результатам (рис. 3), при совпадении угла наблюдения спутникового прибора с зенитным углом Солнца, наблюдается возрастание значения КСЯ, самое большое из которых в первом канале – почти на 100% больше, чем средние значения по всем углам. Поэтому в данной работе при восстановлении АОТ использовались только два канала МСУ-МР – каналы № 2 и № 3.

Рис. 3.

Расчет КСЯ для морской поверхности для каналов прибора МСУ-МР: 1 – канал № 1; 2 – канал № 2; 3 – канал № 3.

Для моделирования значений КСЯ, которые должны быть зарегистрированы спутниковым прибором в каналах МСУ-МР, с помощью 6S задавался ряд параметров: функция спектральной чувствительности канала с шагом 0.025 мкм, КСЯ чистой водной поверхности, тип аэрозоля “Maritime”, предопределенный в 6S, углы наблюдения и освещения, содержание озона и водяного пара, а также значения АОТ, которая в модели 6S определена только на длине волны 550 нм. С использованием 6S было установлено, что КСЯ чистой водной поверхности составляет 0.00379 и 0.0001 для каналов № 2 и № 3 прибора МСУ-МР соответственно. Для всех этих параметров в определенном диапазоне составлялись дискретные наборы значений (табл. 1), для каждой возможной комбинации которых запускалось моделирование. Результатом моделирования 6S является значение КСЯ на верхней границе атмосферы.

Таблица 1.  

Параметры при создании LUT-таблицы

Параметр Количество Значения
Зенитный угол Солнца (град.) 9 0, 25, 40, 50, 60, 70, 75, 80, 85
Зенитный угол спутника (град.) 8 0, 10, 20, 30, 40, 50, 60, 65
Относительный азимутальный угол (град.) 4 0, 60, 120, 180
Озон (DU) 6 220, 270, 300, 350, 400, 600
Водяной пар (cm) 8 0.0, 0.25, 0.5, 0.75, 1.0, 1.5, 2.5, 10.0
АОТ 11 0.0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.5, 0.7, 1.3, 5.0

АЛГОРИТМ ВЫЧИСЛЕНИЯ АОТ

При реализации алгоритма получения АОT данные об углах Солнца и наблюдения рассчитываются по имеющейся информации о положении КА “Метеор-М” № 2 в момент съемки. Данные о содержании водяного пара и озона извлекаются из прогностической модели GFS (Global Forecast System). Для совмещения времени съемки и данных модели применяется интерполяция по времени с использованием оптического потока (Симоненко и др., 2017). АОТ рассчитывается только для безоблачных пикселей, определенных по маске облачности (Андреев и др., 2019). В зимний период, помимо маски облачности, для отделения чистой воды ото льда используется маска ледяного покрова (Кучма и др., 2021).

Ввиду того, что коэффициент зеркального отражения солнечного света от поверхности океана (ρgl) влияет на точность восстановления АОТ (чем меньше угол отражения, тем больше регистрируемое излучение), необходимо производить его коррекцию для каждого из коротковолновых каналов прибора МСУ-МР, прошедших предварительную интеркалибровку (Филей и др., 2016). Стоит отметить, что ρgl является одним из составляющих коэффициента отражения поверхности океана ρos (λ) (Koepke, 1984):

$\begin{gathered} {{\rho }_{{os}}}\left( {{{\theta }_{s}},{{\theta }_{v}},\varphi ,\lambda } \right) = \,\,~{{\rho }_{{wc}}}\left( \lambda \right) + \left( {1 - W} \right) \times \\ \times \,\,{{\rho }_{{gl}}}\left( {{{\theta }_{s}},{{\theta }_{v}},\varphi ,\lambda } \right) + \left( {1 - {{\rho }_{{wc}}}\left( \lambda \right)} \right){{\rho }_{{sw}}}\left( {{{\theta }_{s}},{{\theta }_{v}},\varphi ,\lambda } \right), \\ \end{gathered} $
где ${{\theta }_{s}}$ – зенитный угол Солнца, ${{\theta }_{v}}$ – зенитный угол спутника, φ – относительный азимутальный угол (φ = φs – φv, где φs – азимутальный угол Солнца, φv – азимутальный угол спутника), λ – длина волны, ${{\rho }_{{wc}}}\left( \lambda \right)$ – коэффициент отражения гребня волны на водной поверхности, ${{\rho }_{{gl}}}\left( \lambda \right)$ – коэффициент зеркального отражения солнечного света от поверхности океана, ${{\rho }_{{sw}}}\left( \lambda \right)$ – коэффициент отражения морской воды, W – относительная площадь, покрытая гребнем волны.

Коэффициент ρgl рассчитывается по формуле:

(1)
${{\rho }_{{gl}}}\left( {{{\theta }_{s}},{{\theta }_{v}},\varphi } \right) = {{\left( {\pi P(Z_{x}^{'},~Z_{y}^{'})R\left( {n,{{\theta }_{s}},{{\theta }_{v}},\varphi ,\lambda } \right)} \right)} \mathord{\left/ {\vphantom {{\left( {\pi P(Z_{x}^{'},~Z_{y}^{'})R\left( {n,{{\theta }_{s}},{{\theta }_{v}},\varphi ,\lambda } \right)} \right)} {\left( {4\cos \left( {{{\theta }_{s}}} \right)\cos \left( {{{\theta }_{v}}} \right){{{\cos }}^{4}}\left( \beta \right)} \right)}}} \right. \kern-0em} {\left( {4\cos \left( {{{\theta }_{s}}} \right)\cos \left( {{{\theta }_{v}}} \right){{{\cos }}^{4}}\left( \beta \right)} \right)}},$
где $Z_{x}^{'}$ и $Z_{y}^{'}$ – компоненты наклона поверхности (Cox, Munk, 1954), $P(Z_{x}^{'},~Z_{y}^{'})$ – распределение наклона поверхности, выраженное рядом Грама-Шарлье (Cox, Munk, 1954), n – комплексный показатель преломления морской воды, $R\left( {n,{{\theta }_{s}},{{\theta }_{v}},\varphi ,\lambda } \right)$ – коэффициент отражения Френеля (Born, Wolf, 1975), β – угол наклона поверхности.

Значение ρgl вычисляется с помощью выражения (1) для каждого пикселя каналов № 2 и № 3, и это значение вычитается из КСЯ канала. Также стоит отметить, что алгоритм выполняет расчет параметра АОТ над поверхностью воды, исключая области с искаженными значениями на спутниковых изображениях, вызванных солнечными бликами. Для этого рассчитывается угол солнечного блика, и учитываются значения только больше 40° (Zhao, 2017).

Алгоритм для каждого пикселя находит ближайшие значения всех параметров из LUT и считывает значения смоделированных АОТ и значений КСЯ (которые должны быть зарегистрированы на спутниковом приборе), соответствующих им. После этого производится линейная интерполяция (3) значений КСЯ относительно каждого из параметров:

(2)
$\begin{gathered} {\text{Alb}} = {\text{Al}}{{{\text{b}}}_{1}} + ({\text{Al}}{{{\text{b}}}_{2}} - {\text{Al}}{{{\text{b}}}_{1}}) \times \\ \times \,\,\left( {{{({\text{Va}}{{{\text{l}}}_{{curr}}} - {\text{Va}}{{{\text{l}}}_{1}})} \mathord{\left/ {\vphantom {{({\text{Va}}{{{\text{l}}}_{{curr}}} - {\text{Va}}{{{\text{l}}}_{1}})} {({\text{Va}}{{{\text{l}}}_{2}} - {\text{Va}}{{{\text{l}}}_{1}})}}} \right. \kern-0em} {({\text{Va}}{{{\text{l}}}_{2}} - {\text{Va}}{{{\text{l}}}_{1}})}}} \right), \\ \end{gathered} $
где Alb – результирующее значение АОТ, Val1 и Val2 – граничные значения параметра, Alb1 и Alb2 – значения АОТ для Val1 и Val2, Valcurr – текущее значение параметра в пикселе.

Для каждого пикселя после вычисления (2) получается массив с 11 значениями КСЯ, каждое из которых соответствует значению АОТ (табл. 1). Этот массив интерполируется относительно массива значений АОТ с шагом 0.001, то есть для каждого значения АОТ от 0 до 5 с шагом 0.001 рассчитываются значения КСЯ каждого канала с помощью (3), при этом в качестве Val используются значения АОТ из LUT. После этого, с использованием реальных значений КСЯ каждого из каналов, выполняется поиск ближайшего расстояния до результирующего значения АОТ:

$distance = \sqrt {{{{\left( {{\text{КСЯ}}_{{curr}}^{1} - {\text{КСЯ}}_{{{\text{LUT}}}}^{1}} \right)}}^{2}} + {{{\left( {{\text{КСЯ}}_{{curr}}^{2} - {\text{КСЯ}}_{{{\text{LUT}}}}^{2}} \right)}}^{2}}} ,~$
где ${\text{КСЯ}}_{{curr}}^{1},~\,{\text{КСЯ}}_{{curr}}^{2}$ – измеренные значения КСЯ каналов № 2 и № 3 прибора МСУ-МР соответственно, ${\text{КСЯ}}_{{{\text{LUT}}}}^{1},{\text{\;КСЯ}}_{{{\text{LUT}}}}^{2}~$ – смоделированные значения КСЯ каналов № 2 и № 3 прибора МСУ-МР соответственно.

Кроме наименьшего расстояния между значениями КСЯ двух каналов значение distance представляет собой индекс элемента массива, по которому и выбирается результирующее значение АОТ.

ВАЛИДАЦИЯ

Валидация полученных значений АОТ проводилась по данным сети измерений AERONET, а именно по 14 островным станциям и данным корабельных измерений. Информация отбиралась в период с 2014–2017 гг. с усреднением по времени в пределах 1-часового интервала от времени пролета КА “Метеор-М” № 2, при этом значения АОТ по данным МСУ-МР также были пространственно усреднены в квадрате 25 × 25 км. Таким образом, для валидации было отобрано 380 наблюдений уровня 2.0. Результаты сравнения показаны на рис. 4.

Рис. 4.

Результат сравнения восстановленных значений АОТ по данным МСУ-МР и сети AERONET.

Для оценки точности восстановления АОТ над морской поверхностью были рассчитаны среднеквадратичная (RMSE), средняя абсолютная (MAE) ошибки и коэффициент корреляции Пирсона (R). RMSE составил 0.137, MAE – 0.09, а R – 80.26%. Результаты валидации показали, что разработанный метод определения АОТ по данным МСУ-МР не уступает по точности зарубежным алгоритмам по данным приборов VIIRS (Ray, 2014) и AVHRR (Zhao, 2017). Таким образом, по результатам проведенной валидации можно считать, что разработанный алгоритм может быть применен для расчета параметра АОТ на длине волны 550 нм по данным прибора МСУ-МР, поскольку он отвечает предъявляемым требованиям к качеству продукта.

ЗАКЛЮЧЕНИЕ

В ходе работы был реализован алгоритм расчета параметра АОТ над морской поверхности по данным прибора МСУ-МР. В качестве исходного был использован метод восстановления АОТ по данным прибора AVHRR, в основе которого лежит оптическая модель аэрозоля, опирающаяся на LUT таблицы, размерность которых была увеличена с четырехмерной до шестимерной. В частности, в LUT таблицу были добавлены общие содержания в столбе атмосферы водяного пара и озона. Анализ результатов проведенной валидации по данным AERONET показал, что разработанный алгоритм не уступает по точности зарубежным аналогам по восстановлению параметра АОТ и может применяться для оперативного расчета параметра.

Список литературы

  1. Андреев А.И., Шамилова Ю.А., Холодов Е.И. Применение сверточной нейронной сети для детектирования облачности по данным прибора МСУ-МР спутника “Метеор-М” № 2 // Метеорология и гидрология. 2019. № 7. С. 44–53.

  2. Асарин А.Е., Балонишникова Ж.А., Битков Л.М., Булыгина О.Н., Бугров Л.Ю., Виноградова В.В., Гаврилова С.Ю., Ганюшкин Д.А., Гинзбург А.И., Георгиевский М.В., Глазовский А.Ф., Говоркова В.А., Голованов О.Ф., Голод М.П., Гребенец В.И., Гройсман П.Я., Груза Г.В., Губарьков А.А., Данченков М.А., Дацюк Т.А., Дмитриева Е.В., Дроздов Д.С., Елистратов В.В., Жиркевич А.Н., 3амолодчиков Д.Г., Иванов В.В., Иванова Т.Б., Ильин Б.М., Инсаров Г.Э., Карклин В.П., Кароль И.Л., Катцов В.М., Кириллов С.А., Киселев А.А., Клепиков А.В., Ковалевский Д.В., Кононова Н.К., Комаров А.С., Коршунова Н.Н., Краев Г.Н., Курганова И.Н., Лебедев С.А., Лейбман М.О., Лурье П.М., Лучин В.А., Максимов А.А., Малеев В.В., Махоткина Е.Л., Мезенцева Л.И., Мелешко В.П., Минин А.А., Мирвис В.М., Муравьев А.Я., Мухортова Л.В., Надежина Е.Д., Неушкин А.И., Носенко Г.А., Осипов Э.Ю., Осипова О.П., Павлова Т.В., Панов В.Д., Пересыпкин В.И., Платова Т.В., Поздняков Л.А., Пономарев В.И., Попов И.О., Пчелкин А.В., Радионов В.Ф., Разуваев В.Н., Ранькова Э.Я., Решетников А.И., Романенков В.А., Рочева Э.В., Рябченко В.А., Самохина О.Ф., Санина А.Т., Семенов В.А., Семенова И.В., Сирин А.А., Соколов О.В., Спорышев П.В., Степанов А.Л., Страшная А.И., Стрелецкий Д.А., Тимохов Л.А., Титкина С.Н., Трусенкова О.О., Устинова Е.И., Ушакова Р.Н., Хабаров Н.В., Хен Г.В., Хомутов А.В., Хромова Т.Е., Черенкова Е.А., Шахгеданова М.В., Швиденко А.З., Шеремет Н.А., Шерстюков А.Б., Шилин М.Б., Шиятов С.Г., Школьник И.М., Шкорба С.П., Щепащенко Д.Г. Второй оценочный доклад Росгидромета об изменениях климата и их последствиях на территории Российской Федерации. Федеральная служба по гидрометеорологии и мониторингу окружающей среды. М.: ААНИИ, 2014. 1009 с.

  3. Кучма М.О., Лотарева З.Н., Слесаренко Л.А. Определение маски ледяного покрова дальневосточных морей по данным прибора МСУ-МР спутника “Метеор-М” № 2 // Исслед Земли из космоса. 2021. № 2. С. 31–41. https://doi.org/10.31857/S0205961421020032

  4. Симоненко Е.В., Чудин А.О., Давиденко А.Н. Дифференциальный метод вычисления векторов перемещения облачности // Метеорология и гидрология. 2017. № 3. С. 29–39.

  5. Тимофеев Ю.М., Васильев А.В. Основы теоретической атмосферной оптики: Учебно-методическое пособие. СПб.: Физический ф-т СПбГУ, 2007. 152 с.

  6. Филей А.А., Рублев А.Н., Зайцев А.А. Радиометрическая интеркалибровка коротковолновых каналов многоканального спутникового устройства КА “Метеор-М” № 2 по радиометру AVHRR КА “Metop-A” // Соврем. пробл. дист. зондир. Земли из космоса. 2016. Т. 13. № 6. С. 251–263. https://doi.org/10.21046/2070-7401-2016-13-6-251-263

  7. Born M., Wolf E. Principles of Optics - fifth edition. Oxford: Pergamon Press, 1975. P. 852. https://doi.org/10.1016/0030-3992(75)90061-4

  8. Buras R., Dowling T., Emde C. New secondary-scattering correction in DISORT with increased efficiency for forward scattering // J. Quant. Spectrosc. Radiat. Transfer. 2011. V. 112. № 12. P. 2028–2034. https://doi.org/10.1016/j.jqsrt.2011.03.019

  9. Cox C., Munk W. Statistics of the sea surface derived from sun glitter // J. Marine Research. 1954. V. 13. P. 198–227.

  10. Emde C., Buras-Schnell R., Kylling A., Mayer B., Gasteiger J., Hamann U., Kylling J., Richter B., Pause C., Dowling T., Bugliaro L. The libradtran software package for radiative transfer calculations (version 2.0.1) // Geoscientific Model Development. 2016. V. 9. № 5. P. 1647–1672. https://doi.org/10.5194/gmd-9-1647-2016

  11. Koepke P. Effective Reflectance of Oceanic Whitecaps // Applied Optics. 1984. V. 23. № 11. P. 1816–1824. https://doi.org/10.1364/ao.23.001816

  12. Kuchma, M.O., Bloshchinskiy, V.D. Algorithm for the Atmospheric Correction of Shortwave Channels of the MSU-MR Radiometer of the Meteor-M No. 2 Satellite // Izv. Atmos. Ocean. Phys. 2020. V. 56. P. 909–915. https://doi.org/10.1134/S0001433820090145

  13. Ray G. Joint Polar Satellite System (JPSS) VIIRS Aerosol Optical Thickness (AOT) and Particle Size Parameter Algorithm Theoretical Basis Document (ATBD) // GSFC JPSS CMO. 2014. Rev. B. P. 64.

  14. Levy R., Remer L., Tanre D., Mattoo S., Kaufman Y. Algorithm for Remote Sensing of Tropospheric Aerosol Over Dark Targets from MODIS: Collections 005 and 051: Product ID: MOD 04 / MYD 04 // NASA. 2009. Rev. 2. P. 87.

  15. Stocker T.F., Qin D., Plattner G.-K., Tignor M., Allen S.K., Boschung J., Nauels A., Xia Y., Bex V., Midgley P.M. Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change / Eds. Midgley P.M. Cambridge: Cambridge University Press, 2013. P. 1535.

  16. Vermote E. F., Tanre D., Deuze J.L., Herman M., Morcrette J.J. Second Simulation of the Satellite Signal in the Solar Spectrum, 6S: An overview // IEEE Trans. Geosci. Remote Sens. 1997. V. 35. № 3. P. 675–686. https://doi.org/10.1109/36.581987

  17. Zhao X. NOAA Climate Data Record (CDR) of AVHRR Daily and Monthly Aerosol Optical Thickness (AOT) over Global Oceans, Version 3.0 // NOAA National Centers for Environmental Information. 2017. P. 29. https://doi.org/10.7289/V5BZ642P

Дополнительные материалы отсутствуют.

Инструменты

следующая статья выпуска предыдущая статья выпуска содержание выпуска

Исследование Земли из Космоса

Архивы выпусков Информация о журнале Отправить рукопись в журнал