CC BY
149
НАУЧНОЕ ИЗДАНИЕ МЕТУ ИМ. Н. Э. БАУМАНА
НАУКА и ОБРАЗОВАНИЕ
Эл № ФС77 - 48211. Государственная регистрация №0421200025. КБМ 1994-0408
электронный научно-технический журнал
Контроль состояния сельскохозяйственных полей на основе прогнозирования динамики индекса NDVI по данным космической мультиспектральной и гиперспектральной съёмки
# 07, июль 2013 DOI: 10.7463/0713.0577991
Майорова В. И., Банников А. М., Гришко Д. А., Жаренов И. С., Леонов В. В., Топорков А. Г., Харлан А. А.
УДК 528.88, 528.714.1, 528.852.4
Россия, МГТУ им. Н.Э. Баумана [email protected] alexm.fighter@ gmail. com [email protected] [email protected] [email protected] toporkov@, student.bmstu. ru
a.harlan@hotmail. com
Введение
Космические средства наблюдения представляют собой быстро развивающееся техническое средство эффективного контроля и управления природными ресурсами. Данные, поступающие с космических аппаратов (КА) дистанционного зондирования Земли (ДЗЗ), характеризуются следующими свойствами:
— объективностью, выражающейся в представлении данных независимо от воли заинтересованных лиц;
— масштабностью, определяемой техническими возможностями камеры, существенно большей по сравнению с потенциалом аэрофотосъёмки;
— оперативностью, зависящей от состава орбитальной группировки КА ДЗЗ и составляющей срок от одного дня до недели.
Среди всего многообразия хозяйственной деятельности, осуществляемой в Российской Федерации, сельское хозяйство (c/х) является приоритетным направлением [1]. Огромная площадь и требования к повышению качества обработки земель определяют прямую необходимость внедрения космического сегмента в процесс контроля состояния посевов.
В настоящей работе представлены результаты проведённых исследований, направленных на установление возможности применения космической мультиспектральной и гиперспектральной съёмки для прогноза состояния с/х полей России и дана оценка точности полученных результатов моделирования.
Объекты исследования
Для отработки моделей прогнозирования были выбраны четыре экспериментальных района (рис. 1): юг Московской области (р-н с. Мочилы), центр Московской области (р-н п. Михнево), Новосибирская область (р-н Бердска, с. Улыбино), Ростовская область (р-н г. Морозовска).
Рисунок 1. Районы прогнозирования
Все эти 4 района принципиально отличаются по условиям вегетации: как по рельефу и климату, так и по уровню антропогенного воздействия в виде уровня выброса вредных веществ в атмосферу, что может повлиять как на динамику вегетации с/х культур, так и на качество атмосферной коррекции.
Экологическая обстановка в районе с. Мочилы оценивается, как благоприятная, а степень загрязнённости почвы и атмосферы - низкая по сравнению с территориями, расположенными ближе к Москве [2]. Экологическая ситуация в районе п. Михнево характеризуется большей напряженностью по сравнению с экологической ситуацией в районе с. Мочилы, степень загрязненности и изменения природной среды здесь выше. К югу от Новосибирска экологическая ситуация неплохая. Уровень вредных выбросов здесь низкий, а уровень модернизации промышленности позволяет говорить о низкой степени ее влияния на экологию региона [3]. Выбранная территория в Ростовской области подвержена среднему антропогенному воздействию, экологическое состояние оценивается как напряжённое, что ниже среднего по области [4].
Ростовская область характеризуется ярко выраженным степным ландшафтом, представляющим собой многочисленные поля при отсутствии развитых лесных массивов. Таким образом, с/х посевы расположены на открытой местности с континентальным климатом, что обуславливает высокую чувствительность посевов к солнечной активности, так как устойчивость концентрации хлорофилла в листе невелика. Кроме того, открытость территории и соответствующий почвенный состав [5] приводят к быстрому поглощению осадков почвой.
Выбранный район Новосибирской области характеризуется ярко выраженным лесостепным ландшафтом, представляющим собой наличие с/х полей на распаханных территориях лесных массивов. Таким образом, с/х посевы расположены на локально укрытой от действия ветров местности, что обуславливает переменную чувствительность посевов к солнечной активности, устойчивость концентрации хлорофилла в листе хорошая. Кроме того, укрытость территории и соответствующий почвенный состав [6] (подзолистые, серые лесные и чернозёмные) приводит к достаточному поглощению осадков почвой.
Московская область (ближнее Подмосковье) характеризуется высокой индустриальной насыщенностью. Таким образом, с/х посевы расположены на локальноукрытой местности рядом с промышленными и инфраструктурными объектами, что обуславливает переменную чувствительность посевов к солнечной активности, устойчивость количества хлорофилла в листе хорошая.
Методика исследования
При выполнении исследования была предложена модель прогнозирования вегетационной активности, основанная на использовании относительного разностного индекса вегетации (ЫБУГ). Общий подход к моделированию динамики ЫБУГ состоит в установлении функциональной зависимости между метеоданными и индексом. После установления однозначной связи эти зависимости предполагается использовать для предсказания значений ЫБУГ на некотором интервале времени. При этом в качестве аргументов выступают начальное состояние системы, прогнозные значения метеовеличин и характеристик, определяющих специально вводимые параметры влияния. С целью апробации разработанного методического и алгоритмического обеспечения был проведён спутниковый мониторинг опытного участка - с. Мочилы, который выбирался из соображений удалённости от крупных промышленных объектов, с тем, чтобы свести влияние промышленных выбросов на результаты моделирования к минимуму. Полученные предварительные результаты описаны в [7].
Снимки
LandSat5, LandSat7
Атмосферная
коррекция
Индекс NDVI
\ f
Прогнозное значение
LST, NUS
VI, NDVI '
і 4
Снимки Terra
Рисунок 2. Различные подходы к формированию начальных данных для моделирования
На данном этапе исследований были сформулированы три подхода к формированию начальных данных для моделирования (рис. 2). Изначально использовались данные мультиспектральных сенсоров системы LandSat, а климатические данные вводились на основании информации с наземных метеостанций (схема М-Н). Наземные метеостанции имеют некоторую географическую дискретность расположения и, таким образом, не совсем объективно отображают реальное климатическое состояние интересующего района. Кроме того, эти станции предоставляют информацию о состоянии приземного воздуха, которое характеризуется быстрой изменчивостью и не отражает в полной мере суммарный эффект воздействия климата на процесс вегетации, хотя нельзя не отметить, что значения параметров воздуха являются критичными для c/х культур в случае погодных аномалий. В качестве альтернативы было предложено использовать данные с гиперспектрального сенсора в виде температуры подстилающей поверхности (ошибка определения около 250 м) и индекса влагосодержания, устойчивого к атмосферному влиянию. Полный переход на гиперспектральный источник информации (индекс NDVI также рассчитывался по гиперспектральным снимкам) не принёс положительных результатов (схема Г-Г), однако комбинация мультиспектральных данных по NDVI и гиперспектральных данных по климатическим параметрам (схема М-Г) позволила заметно улучшить качество прогнозирования.
В исследовании использовались гиперспектральные снимки с КА Aqua и Terra. Снимки скачивались с сайта NASA [8] за даты, совпадающие с датами снимков с КА LandSat, по которым ранее рассчитывались значения индекса вегетации NDVI. Район съёмки, в котором лежат с/х поля: между 38° и 40° в.д. и между 54° и 55° с.ш.
Стандартные средства программы ScanEx Image Processor® [9] дают возможность автоматически получать маску LST - значения температуры подстилающей поверхности.
Включаемая географическая сетка позволяет безошибочно определить район наблюдений и увеличить изображение до требуемого функционального уровня. Ориентируясь на цветовую раскраску, можно выделить области примерно равных температур, а исходя из относительной малости района наблюдения и его хозяйственной однородности, можно считать колебания температуры поверхности незначительными. Для приведения их к среднему значению предлагается определять квазиминимальную и квазимаксимальную температуру (у каждого участка поверхности, соответствующего пикселю, своё значение температуры и просмотреть каждый не представляется возможным и не имеет смысла), а затем, пользуясь небольшими различиями между ними, вычислять значение температуры подстилающей поверхности для района как среднее арифметическое квазиминимального и квазимаксимального значений. Вышеописанный алгоритм приведёт к получению температуры подстилающей поверхности именно в интересующей точке на поверхности Земли. Физический смысл получаемого параметра говорит о том, что при отсутствии всякого рода термальных аномалий он должен находиться в некоторой взаимосвязи (рис. 3б) с температурой приземного воздуха, данные по которой поступают с наземных метеостанций.
Наземные метеорологические станции предоставляют данные о влажности воздуха в относительных единицах (иногда для удобства переводимых в проценты). Среди готовых масок для прибора ЫОБШ не предусмотрено определение относительной влажности воздуха, вместе с тем, для подстилающей поверхности и для растений существуют специальные индексы влажности и влагосодержания, которые по своему физическому смыслу также должны коррелировать с фактом выпадения и количеством осадков.
0,8
0,6
0,4
0,2
0
-0,2Z
-0,4
-0,6
-0,8
-1
40
35
30
25
20
15
10
5
О
0,8
0,7
0,6
0,5
0,4
0,3
0,2
ОД
О
NDVI а)
/
_ _ --- ^ _
Дата
7.4 17 .5 6 У/ 6 26 .6 16 .7 5 8 21 .8 14 .9 24
\
\
\
27.4
17.5
6.6
26.6
16.7
5.8
Поле 1.1 ■Поле 1.2 ■Поле 1.3 ■Поле 1.4 Поле 1.5 ■ Поле 1.6 Поле 1.7 ■Поле 1.8 •Ср1
Т, °С б)
/Н
н Г">
— — л V ' N
Дата
ТВ_2011
ТП_Ак_11
*ТП_Тер_11
25.8
14.9
NDWI е)
л
/7^
/ /
Дата
ВЛ_11 NDWI_A_11 NDWI Т 11
27.4
17.5
6.6
26.6
16.7
5.8
25.8
14.9
Рисунок 3. Предварительные данные по с. Мочилы: а) изменение во времени индекса NDVI;
б) изменение во времени температуры приземного воздуха и температуры подстилающей поверхности по данным с КА Terra и Aqua; в) изменение во времени влажности приземного воздуха и индекса NDWI, рассчитанного по данным с КА Terra и Aqua
Одним из таких доступных для расчёта по данным MODIS индексов, отражающих влагосодержание в растениях и в почве, является Нормализованный разностный водный
индекс (Normalized Difference Water Index, NDWI), введенный впервые в 1996 г. Гао (Gao) и определяемый по аналогии с NDVI как [10]:
где
NIR - участок ближнего инфракрасного диапазона с длинами волн в интервале 0.841 - 0.876 нм; SWIR - участок диапазона с длинами волн в интервале 1.628-1.652 нм.
Функциональность формулы объясняется следующими соображениями: вместо использования красного диапазона, интенсивность отражения в котором определяется наличием хлорофилла, применяется коротковолновый ближний инфракрасный диапазон (SWIR), в котором происходит высокое поглощение света водой. В общем случае возможно применение более широкого интервала 1500-1750 нм. Использование же ближнего инфракрасного диапазона (NIR) как и в случае с NDVI связано с тем, что вода не поглощает этот участок электромагнитного спектра, тем самым индекс является устойчивым к атмосферному воздействию, что выгодно отличает его от NDVI [11]. Отметим, что применительно к лесным массивам известно, что для индекса NDWI характерно более устойчивое снижение значений после достижения некоторых критических значений уровня антропогенной нагрузки, что может служить его характеристикой, как более чувствительного, по сравнению с NDVI, индикатора экологического состояния лесов. Из сравнения динамики обоих спектральных индексов можно также предположить, что реакция лесов на возрастание антропогенной нагрузки, проявляющаяся в снижении концентрации хлорофилла происходит с некоторой задержкой после начала процессов обезвоживания древесной растительности [12].
Из рис. 3в видно, что наиболее точное совпадение с динамикой данных, полученных с наземных метеостанций, в случае с. Мочилы показывает индекс NDWI, рассчитанный по снимкам КА Terra. График индекса повторяет все узловые точки графика относительной влажности воздуха. Индекс, полученный с КА Aqua совпадает с наземными данными лишь частично. Это может быть связано с тем, что в летние месяцы влага, содержащаяся в растениях и в почве, в течение дня испаряется в условиях достаточной температуры, но к моменту пролёта КА Terra процесс ещё не успевает набрать полную силу [13].
Математическая модель прогнозирования
Учитывая то обстоятельство, что эмпирически найденная кривая изменения индекса NDVI во времени носит нелинейный характер, для осуществления прогнозирования необходимо использовать аппарат дифференциального исчисления. В случае прогнозирования на короткий срок участок кривой принято в самой простейшей схеме заменять на отрезок касательной, считая функцию локально-линейной на рассматриваемом интервале прогнозирования, для чего изначально и подбирались аналитические зависимости с использованием программного приложения Formulize.exe [7]. Однако подобный подход является неприемлемым для данного случая, так как интервалы наблюдений объекта с использованием мультиспектральных сенсоров системы LandSat составляют от недели до двух, что делает шаг по времени сравнимым со всем промежутком наблюдений. В этом случае выходом является либо использование имеющих более высокий порядок слагаемых ряда разложения, либо добавление в модель возмущающих факторов с весовыми коэффициентами. Вместе с тем, влияние возмущающих факторов (в данном исследовании -влажность и температура), может быть заложено как постоянной величиной (весовой
коэффициент), так и специальной зависимостью (весовая функция). В данном исследовании предполагается, что изменения температуры и влажности пропорциональны влиянию этих факторов на NDVI с точностью до некоторых констант.
Примем справедливой зависимость следующего вида:
где
*
t - момент времени прогноза,
^ - начальный момент времени,
Щ^о) - априорно известное значение NDVI в момент времени to,
*
N(t) - искомое значение NDVI к моменту времени прогноза,
: Е Г;,, г' - текущее время,
Т^) - прогнозное значение среднесуточной температуры,
H(t) - прогнозное значение среднесуточной влажности,
^(Г(£), £), Рн (Н (£), с)- некоторые функции для определения зависимости N1)VI от среднесуточной температуры и среднесуточной влажности:
дЫ(Т,Н,€)
1уСГ(£),£) = ^(Я(£)Д)
дТ
дЩХМХ)
дН
(3)
(4)
Методика и алгоритм определения функций ^(Г(£), £) и Рн (Н (£), с)- следующие. На выбранном интервале времени £ Е [£<>£"] получают данные о N1)VI и метеоданных. Предпочтительной является информация на небольшом интервале времени [£1, £2] (не более двух недель), с большим количеством (п) измерений параметра по назначению (не менее 4 измерений). Измерению I Е [1, гг] соответствует значение М, в момент времени I,.
В силу того, что исходных данных по ограниченное количество, для определения зависимостей Р} (], £), где / Г, Н, в первом приближении делается допущение о линейном
характере ^ (/, £) в диапазоне метеоданных (текущие значения Г, и Н1 определяются на момент 13:00 местного времени):
^(г(0,0 = ст + ахдт, г е (//(О, О = сн + анАН, г е
(5)
(6)
Зависимость функций (?"(£), £) и Гп (Н(£), £)только от одного соответствующего метеопараметра объясняется использованием допущения о суперпозиции. Если в выражении
*
(2) для ) стандартного ряда разложения изначально прогнозируемый параметр по назначению зависит формально только от времени прогноза, то с учётом (5) и (6) N есть уже функция от Т, Н и t, где Т, Н - неявные функции времени.
Для четырёх измерений N (¡=4), то есть трёх интервалов времени, составляется система линейных уравнений относительно Су + С^, и ССн, которая решается одним из известных методов:
Описанные выше алгоритмы были реализованы в программе «ПУСК» (Прогнозирование урожайности сельскохозяйственных культур), предназначенной для вычисления прогнозируемых значений индекса NDVI и соответствующих ему в зависимости от типа исследуемой сельскохозяйственной культуры предполагаемых значений урожайности.
Необходимость модификации математической модели прогнозирования
В процессе выполнения работы в математическую модель прогнозирования были внесены изменения структурного и функционального уровня, корректировались как концепция её составления, так и типы исходных данных, форма их представления.
На начальном этапе предполагалось, что функциональная зависимость индекса NDVI от времени будет получена по имеющимся измерениям с космических снимков с помощью обработки в программных приложениях, позволяющих осуществлять подбор аналитической зависимости по экспериментальным данным. Были проведены предварительные исследования и показана теоретическая возможность подбора таких зависимостей с приемлемым качеством аппроксимации, причём с использованием стандартных полиномов, со степенью не выше шестой [7]. Построенные кривые частично подтвердили предположение о возможности классификации с/х культур по геометрическим характеристикам кривой NDVI (по углам наклона восходящей и нисходящей ветвей графика и по величине максимума) [14]. Вместе с тем, основной идеей выполнения исследований являлось прогнозирование состояния локальной экосистемы, и в этом случае использование аппроксимирующих кривых затруднительно, так как их поведение в межузловых точках и при экстраполяции зачастую вовсе не соответствует участку с опорными данными [15], следовательно, необходим переход на численные соотношения между достоверными точками и использование этих соотношений в прогнозировании.
В качестве первого варианта модели предлагалась схема Эйлера, в которой производная по времени представлялась в виде частных производных NDVI по температуре и влажности в предположении линейной зависимости каждой из частных производных от соответствующего параметра, зависимость частных производных от времени принималась неявной. Результаты моделирования показали, что при прогнозировании на недельный срок наблюдается сильное расхождение между истинными и рассчитанными данными (рис. 4). Приближенный к реальности результат показало прогнозирование на срок в 3 дня (рис. 5), вместе с тем, было отмечено, что имеют место вылеты значений индекса NDVI за допустимый диапазон значений, а также локальные колебания относительно истинного значения с большой амплитудой (рис. 6). По итогам спутникового мониторинга наилучший результат показало моделирование, проведённое по схеме М-Г для Ростовской области с интервалом прогноза в 3 дня (рис. 7).
0.200
0.400
Рисунок 4. Истинная и прогнозная кривые значений ЫВУ1 при прогнозировании на недельный срок с использованием изначальной математической модели по схеме М-Н
------ист
Прогн
0.11
Рисунок 5. Истинная и прогнозная кривые значений ЫВУ1 при прогнозировании на срок в 3 дня с использованием изначальной математической модели по схеме М-Н (Московская
область)
-----ист
-----ПРОГН
Дата
.9
Рисунок 6. Истинная и прогнозная кривые значений ЫВУ1 при прогнозировании на срок в 3 дня с использованием изначальной математической модели по схеме М-Н (Ростовская
область)
-----ист
Прогн
Дата
.9
Рисунок 7. Истинная и прогнозная кривые значений NDVI при прогнозировании на срок в 3 дня с использованием изначальной математической модели по схеме М-Г (Ростовская
область)
Атмосферная коррекция мультиспектральных изображений с КА LandSat
Следующим шагом с целью улучшения достоверности начальных данных стало применение к мультиспектральным снимкам атмосферной коррекции (АК), которая позволила получить более адекватную картину состояния рассматриваемых c/х полей.
Данные дистанционного зондирования с КА серии LandSat, предоставляемые конечным пользователям, проходят определенную обработку, включающую радиометрическую коррекцию и масштабирование [16] полученных значений яркости пикселов на шкалу возможных значений яркости элемента изображения. Эти данные зависят от радиометрического разрешения матрицы (количества уровней яркости) и представляют собой величины, пропорциональные количеству приходящей радиации (т.н. DN - digital numbers).
Таким образом, для того, чтобы провести АК, необходимо предварительно подготовить снимок. Самым важным этапом подготовки снимка к АК является процедура конвертации исходных значений яркости DN снимка (которые изменяются от 0 до 255 для LandSat 5 TM) в реальные значения приходящего излучения. Значения DN безразмерны и пропорциональны количеству попадающего на сенсор излучения, которое измеряется в Вт/ (м стерадианмкм).
Перед запуском КА определяется соотношение между DN и измеряемым потоком энергии. Этот процесс носит название калибровки сенсора. Перевод данных из DN в реальные значения приходящего излучения осуществляется с помощью специальных формул [17], главная из которых (8) представлена ниже:
L, = L'1 - О" (Qcal - Qamin ) + Lm|ni, (8)
cal max s^cal min
где
Lx - количество приходящего излучения;
Qcalmax - максимальное калибровочное значение DN (=255);
Qcaimin - минимальное калибровочное значение DN (=0 или 1);
LmaxX - количество приходящего излучения, которое после масштабирования становится Qcaimw¿; LminX - количество приходящего излучения, которое после масштабирования становится Qcalmin; Qcal - калибровочное значение DN.
Значения соответствующих величин, которые нужно подставить в (8) и другие важные параметры съёмки обычно распространяются с самими данными LandSat и находятся в так называемом файле метаданных. Пример части файла метаданных LandSat 5 TM представлен ниже
GROUP = MIN_MAX_RADIANCE LMAX_B AND 1 = 193.000 LMIN_BAND1 = -1.520 LMAX_B AND2 = 365.000 LMIN BAND2 = -2.840
END_GROUP = MIN_MAX_RADIANCE GROUP = MIN_MAX_PIXEL_VALUE QCALMAX_BAND1 = 255.0 QCALMIN_BAND1 = 1.0 QCALMAX_BAND2 = 255.0 QCALMIN_BAND2 = 1.0
END GROUP = MIN MAX PIXEL VALUE
Зная параметры и используя формулу (8) можно пересчитать полученные данные с КА LandSat. Важно то, что расчёт проводится для каждого канала отдельно. При пересчете
следует учитывать, что для данных TM (LandSat 5) Qca¡тт = 0 Также нужно обратить внимание, что в зависимости от даты съёмки, коэффициенты разные (см. Таблицу 1) [18].
Таблица 1
Значения калибровочных коэффициентов для проведения АК данных КА LandSat
Дата съёмки 01.03.1984 - 04.05.2003 с 05.05.2003
Канал Lmin Lmax Lmin Lmax
1 -1,52 152,10 -1,52 193,0
2 -2,84 296,81 -2,84 365,0
3 -1,17 204,30 -1,17 264,0
4 -1,51 206,20 -1,51 221,0
5 -0,37 27,19 -0,37 30,2
6 1,2378 15,303 1,2378 15,303
7 -0,15 14,38 -0,15 16,5
Для АК снимков с КА LandSat 5/TM, LandSat 7/ETM+ применялся программный модуль FLAASH программного комплекса ENVI [19], при работе с которым использовались значения представленных выше калибровочных коэффициентов и констант. В ходе обработки снимков выяснилось, что в целом динамика кривой NDVI не меняется, зато абсолютные значения индекса увеличиваются на квазипостоянную величину, повышая устойчивость индекса к влиянию климатических изменений (рис. 8). Это происходит вследствие устранения замутнённости снимка и повышения его контрастности. Эффект влияния процессов рассеяния и поглощения излучения атмосферой состоит в том, что на изображении появляется некий аналог атмосферной дымки, приводящий к уменьшению контраста изображений. При прохождении атмосферы некоторое количество излучения рассеивается и преломляется, то есть уменьшается количество излучения, попадающего на датчик оптико-электронной аппаратуры, вследствие чего сцена на снимке является слегка размытой.
На рис. 9 представлены результаты моделирования динамики NDVI с применением АК для того же с/х поля, результаты моделирования для которого были показаны на рис. 7.
1.0
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0.0
-0.1
NDVI
Дата
17.5 6.6 26.6 16.7 5.8 25.8 14.9
----Среднее значение №)У1 (С применением атмосферной коррекции)
----Среднее значение Ж)У1 (Без атмосферной коррекции)
Рисунок 8. График зависимости значений NDVI от времени до АК и после АК
Рисунок 9. Реальная и прогнозная кривые значений NDVI при прогнозировании на срок в 3 дня с использованием изначальной математической модели по схеме М-Г (Ростовская
область) с применением АК
Таким образом, используя атмосферно откорректированные снимки, можно получить более точную классификацию объектов исследования, а применение их в задаче прогнозирования динамики вегетационного индекса позволяет существенно повысить качество прогноза, что видно из рис. 9.
Радиометрическая коррекция
В результате проведённых исследований была разработана комплексная методика радиометрической коррекции космических снимков по средствам наблюдения, включающая математическую модель процесса искажения изображений, моделирование искажающей функции оптической системы и алгоритмы восстановления снимков.
Процесс искажения изображения может быть математически представлен в следующем виде [20]:
где к(х, у) - искажающая функция, д (х, у) - искажённое изображение на выходе из
оптической системы, f(xíy) идеальное (неискажённое) изображение, г](х,у) - функция
аддитивного шума. Приведённое уравнение модели искажения (9) в пространственной области можно записать в эквивалентном представлении в частотной области:
В случае работы с растровыми изображениями функции д(х,у) и f(x,y) являются
матрицами размерностью М*К, представляющими соответственно искажённое и идеальное изображение в цифровом виде. к(х,у) будет являться искажающей функцией оптической
системы, представленной в виде квадратной матрицы, размерностью Б*Б, где должны
выполняться условия К<М и К<1М. Размерность Ых, у) определяет размер окрестности
пикселей идеального изображения, участвующих в формировании значения одного пикселя искажённого изображения. При этом степень вклада каждого пикселя идеального изображения, принадлежащего этой окрестности, в формирование значения пикселя искажённого изображения целиком определяется видом искажающей функции. В качестве примера на рис.10 приведены графики искажающих функций размытия (а) и смаза (б) размерностью 30x30.
Рисунок 10. Искажающие функции размытия (а) и смаза (б) размерностью 30х30
Причина появления смаза заключается в смещении изображения во время накапливания заряда на светочувствительной матрице. Поэтому для того, чтобы определить размерность искажающей функции смаза оптической системы, можно, например, вычислить количество пикселей на матрице, на которое происходит «сдвиг» изображения за время интервала дискретизации (периода) оптической системы. Зная период орбиты космического аппарата, временной интервал дискретизации (период) оптической системы, пространственное разрешение датчика ДЗЗ, можно определить размерность Б искажающей функции как:
(11)
1рт
где й - радиус Земли, р - пространственное разрешение датчика ДЗЗ, Т - период орбиты
спутника ДЗЗ, At - временной интервал дискретизации (время выдержки) датчика ДЗЗ.
Правая часть равенства в формуле (11) округляется до ближайшего целого натурального значения, поскольку О £ М.
Фильтрация Винера
Данный метод восстановления изображений соединяет в себе учёт свойств искажающей функции и статистических свойств шума. Фильтрация Винера основана на рассмотрении изображения и шума как случайных процессов. При этом задача алгоритма состоит в том, чтобы найти такую оценку / для идеального (неискажённого) изображения/
чтобы среднеквадратическое отклонение этих величин было минимальным. Среднеквадратическое отклонение а задаётся следующей формулой [20]:
М
{(/-Я2}.
(12)
где М{х}- функция математического ожидания аргумента х. Предполагается, что выполнены следующие условия:
- шум и неискажённое изображение не коррелированны между собой;
- либо шум, либо неискажённое изображение имеют нулевое среднее значение;
- оценка линейно зависит от искажённого изображения.
При выполнении перечисленных условий минимум среднеквадратического отклонения достигается для функции, которая задаётся в частотной области следующим выражением [21]:
(13)
где V) = |^(и, V)2 |- энергетический спектр шума, Sf(UI I?) = |/?(ц, V)2 |-
энергетический спектр неискажённого изображения.
Регуляризация Тихонова
Решение задачи минимизации сглаживающего функционала со связью (регуляризация Тихонова) в частотной области представляется следующим выражением [20, 21]:
Н (и, V)
у) =
|Я(гх,17)|2 +у|Р(и,17)|2у
(7 (и, г),
где
у - параметр регуляризации, обычно принимающий значения, близкие к 10"; Р(п,\) - преобразование Фурье для матрицы регуляризирующих параметров
(14)
(15)
Тестирование алгоритмов
Для проведения численного эксперимента по восстановлению размытого и искажённого изображения с помощью параметрической фильтрации Винера и регуляризации Тихонова был взят фрагмент снимка ДЗЗ, смоделировано его размытие и зашумление с использованием искажающей функции с гауссовым шумом с дисперсией 10"5, а затем
выполнена процедура восстановления изображения с помощью алгоритмов, реализованных в специальной программе [22] (рис. 11).
Рисунок 11. Пример размытия изображения и восстановления его двумя методами: а) исходное изображение; б) зашумлённое изображение; в) Результат восстановления изображения методом фильтрации Винера; г) восстановление изображения методом
регуляризации Тихонова
Для осуществления радиометрической коррекции космических снимков необходимо обладать информацией об искажающей функции датчика ДЗЗ, для построения которой требуется знание параметров конкретного датчика (пространственное разрешение и временной интервал дискретизации), а также параметров орбиты КА. Наиболее приближённый к оригиналу результат восстановления размытого зашумлённого изображения (дисперсия шума 10"5) показали алгоритмы восстановления изображений методом минимизации сглаживающего функционала (регуляризация Тихонова). Качество
работы алгоритмов восстановления зависит от значения дисперсии и математического ожидания самого изображения и шума, а также от соотношения сигнал/шум.
Изменения структуры модели
Несмотря на внесённые изменения в процесс подготовки начальных данных, не была в полной мере достигнута цель - получение адекватных данных прогноза на интервале 3-7 дней. Для участков Ростовской области, где результаты моделирования показали наилучшие результаты, было выдвинуто предположение, что это во многом объясняется степным ландшафтом местности [23], при котором высоко влияние климата на вегетационные процессы. Однако результаты по остальным участкам из Московской и Новосибирской областей свидетельствовали о необходимости доработки модели. В связи с этим были проведены следующие изменения:
1) введено ограничение на допустимые значения прогнозных величин индекса ЫБУГ.
ЩБЩ < 1;
2) введён контроль скорости изменения значений ЫВУ1 с целью получения адекватной реакции системы на изменение климатических параметров в формате ед./день:
^ ^—- = 0.04 ед/день;
(Введение регулирования также возможно не только по абсолютному значению спрогнозированной величины, но и по скорости её изменения, то есть по полной производной индекса по времени. Это значение легко получается с учётом геометрического смысла производной. Критический угол, при котором необходимо осуществлять переключение модели, или соответствующая ему величина относительного рассогласования предыдущего и спрогнозированного значения определяется экспериментально. Важно отметить, что ось времени также должна быть безразмерной, приведение к безразмерным величинам можно осуществить ко всему интервалу наблюдений);
3) добавлена компонента, описывающая изменения аэрозольной оптической толщины (АОТ) атмосферы:
Рисунок 12. Представление величины АОТ на сайте NASA в виде палитры значений от 0 до 1
АОТ является безразмерной величиной и определяет количество света, которое рассеивается или поглощается атмосферными частицами во время прохождения света через атмосферу. Goddard Space Flight Center NASA [24] предоставляет ежесуточные данные MOD04 L2 об АОТ с пространственным разрешением 1000 м (рис. 12). Этот тип данных относится как ко второму уровню обработки гиперспектральных изображений, так и, будучи дополненным, к третьему уровню [24]. АОТ представляет собой параметр, чувствительный к оптической прозрачности атмосферы и является индикатором содержания в воздухе различных веществ природного и техногенного происхождения. В атмосфере постоянно присутствуют различные природные частицы, в том числе водяной пар, концентрация которого (АОТ, соответственно) может изменяться в силу климатических вариаций. Вместе с тем, постоянная составляющая аэрозоля в атмосфере и колебания её величины с учётом климата относительно малы по сравнению с антропогенным влиянием, выражающимся, в частности, в выбросах продуктов горения и в распылении химических веществ при обработке полей. Отметим, что и некоторые природные факторы (например, недавнее крупное извержение вулкана в Исландии) могут оказывать аналогичное воздействие на АОТ, но эти случаи оговариваются отдельно;
4) предложено повышение порядка системы за счёт времени или за счёт климатических параметров:
Вариант 1. Линейная модель с добавлением АОТ. Основание прогноза - 5 временных точек, система линейных алгебраических уравнений (СЛАУ) 4-го порядка:
Вариант 2. Нелинейная по климатическим составляющим модель. Основание прогноза - 7 временных точек, СЛАУ 6 порядка:
Вариант 3. Нелинейная по времени модель. Основание прогноза - 8 временных точек, СЛАУ 7 порядка:
:: л : :_1: ■ _ г. ;.г(16г.)
Основные результаты
В ходе моделирования было установлено, что повышение порядка по времени в варианте №3 (16в) приводит к быстрому расхождению прогнозных данных с реальными. По факту повышение порядка этим методом приводит к фиксированному усилению (так как постоянен шаг по времени) изменения температуры и влажности.
Важным моментом является существенное повышение точности шестидневного прогноза (рис. 13, 14), полученное после всех изменений как для линейной (вар. 1), так и для нелинейной (вар. 2) моделей. В отличие от графика, представленного на рис. 9, прогнозные зависимости на рис. 13 характеризуются большей достоверностью. Нельзя не отметить, что моделирование и оценка расхождений с реальной кривой проводились для так называемой кривой тенденции, отражающей отфильтрованную от каждодневных колебаний динамику изменения значений индекса. При определении опорных точек графика ЫВУ1 по имевшимся в наличии мультиспектральным снимкам не было получено неадекватных данных (в отличие
от каждодневных данных в случае использования гиперспектрального сенсора КА Тегга), что позволило исключить ту часть ошибки прогнозирования, которая определяется некорректностью задания начальных условий. Хорошие результаты показало повышение порядка климатических компонент, вместе с тем заметно повышение колебательности прогнозной кривой (рис. 14). Дать оценку достоверности этой колебательности не представляется возможным, так как используемые мультиспектральные снимки имели период обновления от 7-ми до 21 дня.
Рисунок 13. Результаты моделирования по варианту №1 на 3 дня и на 6 дней: а) Михнево б) Мочилы в) Новосибирск г) Ростов
Рисунок 14. Результаты моделирования по варианту №2 на 3 дня и на 6 дней: а) Михнево
б) Мочилы в) Новосибирск г) Ростов
Оценка точности прогнозирования
В ходе выполнения исследований была детально проработана методика оценки точности прогноза. В частности, был выполнен ряд численных экспериментов, в которых исследовалась зависимость погрешности прогнозирования от задаваемых исходных данных.
В первой серии экспериментов дата прогноза оставалась неизменной, а основание прогноза (интервал исходных дат) постепенно сдвигалось назад во времени. Затем (во второй серии экспериментов) основание прогноза оставалось неизменным, а дата прогноза сдвигалась вперед. И, наконец, в третьей серии экспериментов последняя из исходных дат (опорная) и дата прогноза оставались неизменными, а все остальные даты сдвигались назад во времени, расширяя период накопления статистики. Был введен параметр т/;, характеризующий разреженность интервала прогнозирования с целью уменьшения зависимость величины разреженности от типа модели:
где ДО - интервал прогноза, количество дней между опорной и начальной датой,
N - количество исходных дат, равное 5 для варианта 1 и равное 7 для варианта 2.
Численные эксперименты были проведены с использованием серий №1 и №2. Далее в зависимости от интервала прогнозирования значения погрешностей для каждой из серий
были сведены к одному графику. Как видно из рис. 15, все значения погрешностей расположены ниже некоторой кривой максимальной погрешности:
Рисунок 15. Кривая максимальной погрешности, полученная по сериям №1 (А) и №2 (Б{)
После осреднения значений погрешности по всем экспериментам для каждого значения ширины интервала прогнозирования и аппроксимации полученной зависимости экспонентой для каждой серии были получены кривые ожидаемой погрешности (рис. 16):
Рисунок 16. Кривые ожидаемой погрешности для линейной (вариант №1) и нелинейной
моделей (вариант №2) прогнозирования
Эти законы были записаны в виде эмпирических зависимостей:
где
Е];:_ч;- максимально возможная нормированная погрешность прогнозирования;
Г?:-: - параметр модели, равный 1 для линейной модели, и 0,87 для нелинейной модели (вариант 2);
(.3 - климатический коэффициент. /? = 0.085 ... 0.1 : меньшие значения соответствуют плавным изменениям климата от опорной даты к дате прогноза, большие - резким локальным изменениям климата;
I - величина интервала прогнозирования в днях;
к ■ {■!' — - ) - слагаемое, учитывающее разреженность основания прогноза.
Закон изменения средней (ожидаемой) погрешности выглядит аналогично (16):
Выводы
1. Разработаны схемы оптимального сочетания мультиспектральной и гиперспектральной космической съёмки с различным пространственным и временным разрешением в составе методики мониторинга с/х угодий;
2. Реализован алгоритм атмосферной коррекции мультиспектральных космических снимков, позволивший повысить устойчивость и достоверность получаемых значений индекса ЫБУ1 во времени;
3. Проведена серия численных экспериментов по радиометрической коррекции изображений и определён наиболее подходящий для данных задач метод восстановления изображений;
4. Проанализированы и отобраны основные факторы, влияющие на процесс вегетации;
5. Получены достоверные данные прогнозирования состояния с/х угодий с 6-дневным интервалом прогноза;
Работа выполнена при финансовой поддержке Министерства образования и науки Российской Федерации в рамках мероприятия 1.9 Федеральной целевой программы «Исследования и разработки по приоритетным направлениям развития научно-технического комплекса России на 2007-2013 годы», государственный контракт от 18 августа 2011 года № 11.519.11.5002, шифр «2011-1.9-519-004».
Список литературы
1. Приоритетные национальные проекты. Совет при Президенте России по реализации приоритетных национальных проектов: Официальный сайт. Режим доступа: http://www.rost.ru/news/2007/10/091614 11170^Ш1 (дата обращения 15.04.2013).
2. Экологические карты Московской области. Официальный сайт посёлка Загорянка. Режим доступа: http://zagoryansky.com/environmenta1-maps (дата обращения 15.06.2012).
3. Экологическая карта Новосибирской области. Режим доступа: http://www.mapnso.ru/ (дата обращения 15.03.2013).
4. Экологический мониторинг Ростовской области. Донской экологический центр.
Режим доступа: http://ektor.ru/pages/mon1.asp?Idr=2&id=56 (дата обращения 12.09.2012).
5. Нормативы и методика применения побочной продукции сельскохозяйственных культур для обеспечения бездефицитного баланса органического вещества в почвах на землях сельскохозяйственного назначения. Донской НИИСХ. Режим доступа: www.don-agro.ru/fi1es/nauka/gosznr/So1oma.doc (дата обращения 23.11.2012).
6. Природа Новосибирской области. Режим доступа: http://bober.ru/research/novosibirsk.htm (дата обращения 15.01.2013).
7. Майорова В.И., Гришко Д.А., Муравьёв В.В., Топорков А.Г. Использование космических средств наблюдения для мониторинга локальных экосистем // Вестник МГТУ им. Н.Э. Баумана. Сер. Машиностроение. 2012. Спец. вып. «Крупногабаритные трансформируемые космические конструкции и материалы для перспективных ракетнокосмических систем». С. 182-192.
8. Центр оперативной обработки данных ДЗЗ с прибора ЫОБШ. NASA. Режим доступа: http://rapidfire.sci.gsfc.nasa.gov/rea1time/ (дата обращения 10.02.2013).
9. Программное обеспечение. Информационно-технологический центр СканЭкс. Режим доступа: http://www.scanex.ru/ru/software/index.html (дата обращения 20.03.2012).
10. Remote Sensing and GIS in Agriculture. Vegetation Indices. Science Education through Earth Observation for High Schools (SEOS). Режим доступа: http://www.seos-project.eu/modules/agriculture/agriculture-c01-s03.html (дата обращения 25.06.2012).
11. Gao B.C. NDWI - A normalized difference water index for remote sensing of vegetation liquid water from space // Remote sensing of environment. 1996. Vol. 58, no. 3. P. 257-266.
12. Разработка методики региональной экологической оценки состояния лесов по данным спутниковых наблюдений. ГИС-анализ экологического состояния лесов Московской области. Режим доступа: http://www.wonder-lands.net/places-284-1.html (дата обращения
15.07.2012).
13. Прогноз погоды. Режим доступа: http://www.rp5.ru (дата обращения 06.05.2012).
14. Пугачева И.Ю., Шевырногов А.П. Изучение динамики NDVI посевов сельскохозяйственных культур на территории Красноярского края и республики Хакасия. Режим доступа: http://d33.infospace.ru/d33 conf/2008 pdf/2/46.pdf (дата обращения 20.03.13).
15. Аппроксимация функций. Факультет физики РГПУ им. А.И. Герцена. Режим доступа: http://physics.herzen.spb.ru/library/01/01/nm labs/approximation.htm (дата обращения
11.12.2012).
16. Конвертация данных TM, ETM+ в показатели излучения на сенсоре. GIS-Lab («ГИС Лаборатория»): Географические информационные системы и дистанционное зондирование. Режим доступа: http://gis-lab.info/qa/dn2radiance.html (дата обращения 26.09.2012).
17. Vermote E.F., Vermeulen A. Atmospheric correction algorithm: spectral reflectances (M0D09). Version 4.0. NASA, April 1999. Режим доступа:
http://modis.gsfc.nasa.gov/data/atbd/atbd mod08.pdf (дата обращения 21.04.2012).
18. Chander G., Markham B. Revised Landsat-5 TM Radiometric Calibration Procedures and Postcalibration Dynamic Ranges // IEEE Transactions on geoscience and remote sensing. 2003.
Vol. 41, no. 11. P. 2674-2677. DOI: 10.1109/TGRS.2003.818464
19. Руководство по эксплуатации модуля ENVI-FLAASH. Лаборатория ГИТ и ДЗ ИГМ СО РАН. Режим доступа: http://www.nrcgit.ru/metod/nrcgit flaash.pdf (дата обращения 26.07.2012).
20. Гонсалес Р., Вудс Р. Цифровая обработка изображений. М.: Техносфера, 2005. 1072 с.
21. Шовенгердт Р.А. Дистанционное зондирование. Модели и методы обработки изображений : пер. с англ. М.: Техносфера, 2010. 560 с.
22. Майорова В.И., Банников А.М., Зайцев К.И. Математическое моделирование процесса радиометрической коррекции снимков ДЗЗ и сравнительная характеристика алгоритмов восстановления изображений // Вестник МГТУ им. Н.Э. Баумана. Сер. Машиностроение. 2013. Спец. вып. С. 110-121.
23. Смагина Т.А., Кутилин В.С. Ландшафты Ростовской области. Режим доступа: http://stepnoy-sledopyt.narod.ru/geologia/landshafts.htm (дата обращения: 24.01.2013).
24. Ежесуточные данные MOD04 L2 Центра оперативной обработки данных ДЗЗ с прибора MODIS. NASA. Режим доступа: http://ladsweb.nascom.nasa.gov/data/search.html (дата обращения 10.01.2013).
SCIENTIFIC PERIODICAL OF THE BAUMAN MSTU
SCIENCE and EDUCATION
EL № FS77 - 48211. №0421200025. ISSN 1994-0408
electronic scientific and technical journal
Monitoring condition of agricultural fields based on prediction of NDVI with the use of multi-spectral and hyper-spectral data from space imagery
# 07, July 2013 DOI: 10.7463/0713.0577991
Maiorova V.I., Bannikov A.M., Grishko D.A., Jarenov I.S., Leonov V.V., Toporkov A.G., Harlan A. A.
Bauman Moscow State Technical University, 105005, Moscow, Russian Federation
[email protected] alexm.fighter@ gmail. com [email protected] [email protected] [email protected] toporkov@ student.bmstu. ru
a.harlan@hotmail. com
Problems of condition monitoring of agricultural fields based on prediction of NDVI dynamics considering climate and atmospheric influence were studied. Different algorithms for preparing initial data were proposed. The research was conducted for the regions of Russia which are notable for different relief, climate and environment. Initial data for simulation were formed by combining the information from the multi-spectral and hyper-spectral space imagery. The conducted mathematical simulation showed a principal opportunity to use Earth observing systems for appropriate short-time predicting of vegetation processes.
Publications with keywords: monitoring, remote Earth sensing, vegetation dynamics predicting, local ecosystems, agriculture, atmospheric correction, hyperspectral and multispectral imagery
Publications with words: monitoring, remote Earth sensing, vegetation dynamics predicting, local ecosystems, agriculture, atmospheric correction, hyperspectral and multispectral imagery
References
1. Prioritetnye natsional'nyeproekty [The priority national projects]. Sovet pri Prezidente Rossii po realizatsii prioritetnykh natsional'nykh proektov [Council under RF President on implementation of priority national projects]. Official website. Available at: http://www.rost.ru/news/2007/10/091614 11170.shtml , accessed 15.04.2013.
2. Ekologicheskie karty Moskovskoy oblasti. [Ecological maps of the Moscow region]. Official website of settlement Zagoryanka. Available at: http://zagoryansky.com/environmental-maps , accessed 15.06.2012.
3. Ekologicheskaya karta Novosibirskoy oblasti [Ecological map of Novosibirsk region.]. Available at: http://www.mapnso.ru/ , accessed 15.03.2013.
4. Ekologicheskiy monitoringRostovskoy oblasti. [Ecological monitoring of the Rostov region]. Don Environmental Centre. Available at: http://ektor.ru/pages/mon1.asp?Idr=2&id=56 , accessed 12.09.2012.
5. Normativy i metodika primeneniya pobochnoy produktsii sel'skokhozyaystvennykh kul'tur dlya obespecheniya bezdefitsitnogo balansa organicheskogo veshchestva v pochvakh na zemlyakh sel'skokhozyaystvennogo naznacheniya [Norms and methods of of application of by-products of agricultural crops to ensure the sufficient balance organic substance in soils on agricultural lands]. Don Scientific Research Institute of Agriculture. Available at: www.don-agro.ru/files/nauka/gosznr/Soloma.doc , accessed 23.11.2012.
6. PrirodaNovosibirskoy oblasti [Nature of the Novosibirsk region]. Available at: http://bober.ru/research/novosibirsk.htm , accessed 15.01.2013.
7. Mayorova V.I., Grishko D.A., Murav'ev V.V., Toporkov A.G. Ispol'zovanie kosmicheskikh sredstv nablyudeniya dlya monitoringa lokal'nykh ekosistem [Using of space-based observation aids for monitoring of local ecosystems]. VestnikMGTUim. N.E. Baumana. Ser. Mashinostroenie [Bulletin of the Bauman MSTU. Ser. Mechanical Engineering], 2012. Spets. vyp.
“Krupnogabaritnye transformiruemye kosmicheskie konstruktsii i materialy dlya perspektivnykh raketno-kosmicheskikh system ” [Spec. iss. “Large-sized transformable space designs and materials for advanced space rocket systems”], pp. 182-192.
8. Center-line processing of remote sensing data with MODIS instrument. NASA. Available at: http://rapidfire.sci.gsfc.nasa.gov/realtime/ , accessed 10.02.2013.
9. Programmnoe obespechenie [Software]. Informatsionno-tekhnologicheskiy tsentr SkanEks [Information Technology Center ScanEx]. Available at: http://www.scanex.ru/ru/software/index.html , accessed 20.03.2012.
10. Remote Sensing and GIS in Agriculture. Vegetation Indices. Science Education through Earth Observation for High Schools (SEOS). Available at: http://www.seos-project.eu/modules/agriculture/agriculture-c01-s03.html , accessed 25.06.2012.
11. Gao B.C. NDWI - A normalized difference water index for remote sensing of vegetation liquid water from space. Remote sensing of environment, 1996, vol. 58, no. 3, pp. 257-266.
12. Razrabotka metodiki regional'noy ekologicheskoy otsenki sostoyaniya lesovpo dannym sputnikovykh nablyudeniy. GIS-analiz ekologicheskogo sostoyaniya lesov Moskovskoy oblasti [Development of methodology of regional environmental assessment of forest from data of satellite observations. GIS analysis of the ecological condition of forests of Moscow region]. Available at: http://www.wonder-lands.net/places-284-1.html , accessed 15.07.2012.
13. Prognozpogody [Weather forecast]. Available at: http://www.rp5.ru , accessed 06.05.2012.
14. Pugacheva I.Yu., Shevyrnogov A.P. Izuchenie dinamiki NDVIposevov sel'skokhozyaystvennykh kul'tur na territorii Krasnoyarskogo kraya i respubliki Khakasiya [Study of the dynamics of NDVI sowing of agricultural crops on the territory of Krasnoyarsk territory and
the Republic of Khakassia]. Available at: http://d33.infospace.ru/d33 conf/2008 pdf/2/46.pdf , accessed 20.03.13.
15. Approksimatsiyafunktsiy [Approximation of functions]. Fakul'tet fiziki RGPU im. A.I. Gertsena [Department of Physics, Herzen State Pedagogical University.]. Available at: http://physics.herzen.spb.ru/library/01/01/nm labs/approximation.htm , accessed 11.12.2012.
16. Konvertatsiya dannykh TM, ETM+ vpokazateli izlucheniya na sensore [Data Conversion TM, ETM + into indicators of radiation at the sensor]. GIS-Lab : Geograficheskie informatsionnye sistemy i distantsionnoe zondirovanie [GIS-Lab : Geographic information systems and remote sensing]. Available at: http://gis-lab.info/qa/dn2radiance.html , accessed 26.09.2012.
17. Vermote E.F., Vermeulen A. Atmospheric correction algorithm: spectral reflectances (MOD09). Version 4.0. NASA, April 1999. Available at: http://modis.gsfc.nasa.gov/data/atbd/atbd mod08.pdf , accessed 21.04.2012.
18. Chander G., Markham B. Revised Landsat-5 TM Radiometric Calibration Procedures and Postcalibration Dynamic Ranges. IEEE Transactions on geoscience and remote sensing, 2003, vol. 41, no. 11, pp. 2674-2677. DOI: 10.1109/TGRS.2003.818464
19. Rukovodstvopo ekspluatatsii modulya ENVI-FLAASH [Operation manual of module ENVI-FLAASH]. Laboratoriya GIT i DZ IGM SO RAN [Laboratory of GIS and remote sensing. Sobolev Institute of Geology and Mineralogy. Siberian Branch of Russian Academy of Sciences]. Available at: http://www.nrcgit.ru/metod/nrcgit flaash.pdf , accessed 26.07.2012.
20. Gonzales R.C., Woods R.E. Digital Image Processing. 2nd ed. Prentice Hall, 2001. (Russ. ed.: Gonzales R., Vuds R. Tsifrovaia obrabotka izobrazhenii. Moscow, Tekhnosfera, 2006. 1072 p.).
rd
21. Showengerdt R.A. Remote Sensing: Models and Methods for Image Processing. 3 ed. Academic Press, 2007. 515 p. (Russ. ed.: Shovengerdt R.A. Distantsionnoe zondirovanie. Modeli i metody obrabotki izobrazheniy. Moscow, Tekhnosfera, 2010. 560 p.).
22. Mayorova V.I., Bannikov A.M., Zaytsev K.I. Matematicheskoe modelirovanie protsessa radiometricheskoy korrektsii snimkov DZZ i sravnitel'naya kharakteristika algoritmov vosstanovleniya izobrazheniy [Mathematical modeling of process of radiometric correction of remote sensing images and comparative characteristics of algorithms of image reconstruction]. VestnikMGTU im. N.E. Baumana. Ser. Mashinostroenie [Bulletin of the Bauman MSTU. Ser. Mechanical Engineering], 2013, Spec. iss., pp. 110-121.
23. Smagina T.A., Kutilin V.S. Landshafty Rostovskoy oblasti [Landscapes of the Rostov region]. Available at: http://stepnoy-sledopyt.narod.ru/geologia/landshafts.htm , accessed 24.01.2013.
24. The daily data MOD04L2 Center operational processing of remote sensing data with MODIS instrument. NASA. Available at: http://ladsweb.nascom.nasa. gov/data/search.html , accessed 10.01.2013.
