На правах рукописи

Алексеева Ольга Михайловна

Интерполяционная модель спектральной яркости объектов для задач

имитационного моделирования излучения земной поверхности при

наблюдении из космоса

Специальность:25.00.34 - Аэрокосмические исследования Земли,

фотограмметрия

Автореферат на соискание ученой степени

кандидата технических наук

Москва - 2013 2

Работа выполнена в Московском государственном университете геодезии и картографии на кафедре аэрокосмических съемок

Научный руководитель: член-корреспондент РАН, доктор технических наук профессор, кафедра аэрокосмических съемок, Московский государственный университет геодезии и картографии Савиных Виктор Петрович

Официальные оппоненты: Ильин Юрий Александрович доктор технических наук, профессор, кафедра физики, Московский государственный университет геодезии и картографии, профессор Старостенко Алексей Михайлович кандидат технических наук, ОАО «Межгосударственная акционерная компания «Вымпел», старший научный сотрудник

Ведущая организация: Федеральное государственное бюджетное учреждение науки Институт космических исследований РАН

Защита состоится «» _ 2013 г. в часов на заседании диссертационного совета Д 212.143.01 при Московском государственном университете геодезии и картографии в зале заседаний Учного совета по адресу: 105064, г. Москва, Гороховский пер. д.4.

С диссертацией можно ознакомиться в библиотеке Московского государственного университета геодезии и картографии.

Автореферат разослан «» _ 2013 г.

1

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

Актуальность темы.При создании космических систем дистанционного зондирования Земли (ДЗЗ), предназначенных для решения широкого спектра задач в области мониторинга окружающей среды, исследований природных ресурсов (ИПР) и контроля ситуаций, связанных с катастрофами, требуется достоверная информация о полях сигналов, которые регистрируются средствами ДЗЗ.

Адекватное описание полей таких сигналов представляет собой сложную научную задачу, так как их параметры зависит от большого количества внутренних и внешних факторов, природа которых не всегда ясна.

Для решения этой задачи часто применяется многомодельный подход, при котором для описания отдельных аспектов формирования полей сигналов применяются так называемые «частные модели», а описание поля сигналов в целом, формируется как результат взаимодействия этих моделей. Примерами частных моделей являются, например, для случая наблюдения из космоса, модели излучения подстилающей поверхности Земли, модели оптико-электронной аппаратуры, модели атмосферы и т.д.

При решении задач проектирования и испытаний космических систем методами математического моделирования и статистических испытаний модели полей сигналов представляются в виде имитаторов, воспроизводящих пространственное распределение случайных амплитуд в соответствии со свойствами наблюдаемой ситуации.

Для решения задач ДЗЗ в оптическом диапазоне спектра наиболее важным аспектом является зависимость регистрируемых сигналов от оптических свойств объектов наблюдения и атмосферы, которые в свою очередь зависят от индивидуальных свойств и физического состояния объектов и изменяются в зависимости от условий наблюдения и освещения. Так, спектры излучения в области зеркального отражения для некоторых объектов (облачность, водная поверхность) существенно отличаются от спектров тех же объектов при «стандартных» дневных условиях; при больших углах линии визирования существенно увеличивается вклад в интенсивность излучения объекта многократного рассеянного атмосферой солнечного излучения и т.д. Соответственно, адекватный учет спектральных свойств объектов и их изменчивости является важной задачей при моделировании полей сигналов.

В то же время применение ресурсоемких методов расчета яркости в имитационных моделях полей излучения нецелесообразно. Это обусловлено методологией применения таких моделей: на стадии проектирования систем– для оценки сигналов от объектов при характерных сочетаниях условий, на стадии создания и испытаний – для формирования массивов изображений, покрывающих весь диапазон условий функционирования и позволяющих получать статистические обоснованные оценки качества решения задачи ДЗЗ. В связи с этим организация корректного и быстрого расчета спектральных характеристик объектов для широкого диапазона изменения условий наблюдения и освещения является актуальной задачей при построении моделей полей излучения для проектирования и тестирования аппаратуры систем ДЗЗ и систем в целом.

Степень разработанности темы В существующих к настоящему времени моделях полей излучения расчет спектральной яркости базируется либо на теоретических решениях уравнения переноса излучения, в частности методом Монте-Карло, требующих для достижения необходимой точности расчета больших вычислительных ресурсов (в качестве примера таких моделей можно привести модели, разрабатываемые в ОАО МАК «Вымпел», ИОА СО РАН, AirForceResearchLaboratory и др., либо на эмпирических зависимостях, которые в связи с определенными трудностями при получении экспериментальных данных ограничены по условиям применения, и в таком случае при имитации наблюдаемой обстановки прибегают к «частным» моделям и как следствие возникают проблемы с последующей их «стыковкой». К таким моделям можно отнести модель, разработанную в ГОИ им. С.И.Вавилова и основанную на полуэмперической формуле пересчета экспериментально полученных спектров при конкретных условиях яркости на другие условия.

Однако применение эмпирических методик имеет ограничения, связанные с ростом погрешности вычислений при экстраполяции экспериментальных данных на широкий диапазон условий функционирования космических средств.

Таким образом, для решения задач имитации полей излучения в целях проектирования и тестирования систем ДЗЗ, на данный момент существует необходимость в модели расчета яркости, которая позволила бы производить вычисления в широком диапазоне условий наблюдения и освещения в рамках одной расчетной схемы и в тоже время требовала минимальные затраты вычислительных ресурсов.

Цель работы и задачи исследований.

Целью работы является разработка модели расчета спектральной яркости объектов для имитации подстилающей поверхности при проектировании и тестировании систем ДЗЗ, позволяющая использовать ее при моделировании излучения объектов при глобальном обзоре, при этом затрачивая минимальное количество вычислительных ресурсов и сохраняя адекватность получаемых результатов, что имеет существенное значение для отрасли аэрокосмических исследований Земли.Эту цель можно достигнуть путем разработки интерполяционной модели, которая выражает в явном (как правило, в линейном) виде связь между переменными модели, определяющими значения спектральной яркости объектов.

Для этого необходимо решить следующие задачи:

1. Разработать структуру и принципы работы модели расчета спектральной яркости объектов, с учетом ее направленности на применение для проектирования систем ДЗЗ на основе анализа существующих моделей расчета спектроэнергетических характеристик объектов.

2. Определить параметры разрабатываемой модели, удовлетворяющие требованиям к моделям расчета яркости на основе анализа принципов формирования яркости различных объектов на входе бортовой аппаратуры КА в зависимости от условий наблюдения и освещения и на его основе.

3. Провести исследование разработанной модели для оценки ее эффективности, быстродействия и возможности использования при проектировании систем ДЗЗ.

Научная новизна результатов диссертационной работы заключается в том, что разработанная модель, в отличие от существующих моделей позволяет:

получать спектральную яркость объектов в широком диапазоне условий наблюдения и освещения (день, ночь, терминатор, пригоризонтная область, область зеркального отражения) в рамках единой схемы расчета, что делает возможным ее использование в задачах имитации наблюдаемой обстановки при глобальном обзоре и исключает проблемы совмещения «частных» моделей; при этом простота модели позволяет доработку и включение новых данных с минимальными изменениями в программном коде имитационной модели;

производить одновременно быстрый и корректный расчет спектральной яркости объектов и тем самым моделировать излучение объектов в условиях «реального времени» или близких к нему.

Теоретическая и практическая значимость работы.

Результаты работы использованы в НИР «Листва-2020» в части разработки модели расчета спектральной яркости объектов и в ОКР «КИУС» в части расчета отраженного и диффузного излучения объектов. Результаты работы могут использоваться при решении различных частных задач, связанных с априорной оценкой уровня регистрируемых сигналов и имитацией полей яркости, в частности, при разработке принципов построения бортовой аппаратуры, выборе ее характеристик, создании бортовых и наземных программных и технических средств обработки информации, а также в учебном процессе в дисциплинах, связанных с физическими принципамидистанционного зондирования Земли.

Методология и методы исследования, применяемые в настоящей работе, основаны на использовании теории переноса излучения в неоднородной атмосфере, физики природной среды, метеорологии, методов дешифрирования космических снимков, методов математического моделирования, методов классификации объектов, теории вероятностей и математической статистики. Методология исследования заключалась в проведении анализа существующих моделей, на основе которого разрабатывалась структура модели, принципы работы модели и выбирались способы расчета спектральной яркости объектов. На основе результатов анализа физических принципов формирования яркости различных объектов на входе бортовой аппаратуры КА в зависимости от условий наблюдения и освещения были выбраны параметры модели. Проведены исследования характеристик разрабатываемой модели.

Положения, выносимые на защиту 1. Структура и принципы работы модели расчета спектральной яркости объектов, позволяющие ее использование при имитации наблюдаемой обстановки системой ДЗЗ глобального обзора.

2. Методика выбора количества значений параметров информационных массивов модели для различных условий наблюдения и освещения, основанная на особенностях физических принципов формирования яркости объектов на входе БА и обеспечивающая быстрый и корректный расчета спектральной яркости объектов при имитации излучения от объектов в условиях «реального времени» или близких к нему.

3. Упрощенная методика расчета спектральной яркости объектов в области зеркального отражения, необходимая для оценки максимального уровня входного сигнала.

4. Методика расчета спектральной яркости при наблюдении сложной ситуации с учетом затенения, позволяющая при моделировании подстилающей поверхности имитировать "смешанные" пиксели.

Степень достоверности и апробация результатов работы.

Достоверность результатов работы подтверждается корректным применением методов расчета и математического моделирования зависимостей энергетических характеристик объектов, обоснованным выбором исходных данных, сравнением полученных параметров с результатами публикаций и экспериментальными данными.

Изложенные в диссертации результаты представлены в трех статьях, опубликованных в научно-техническом журнале «Геодезия и аэрофотосъемка», в четырех научно-технических отчетах и документах.

Апробация результатов диссертационных исследований проводилась на четырех научнотехнических конференциях, в том числе: на научно-технической конференции молодых специалистов ОАО «Корпорация «Комета» в 2013 г. и на трех конференциях студентов, аспирантов и молодых ученых МИИГАиК в 2010, 2012 и 2013 гг.

Структура и объем диссертацииДиссертация состоит из введения, 4 глав, заключения и списка литературы. Общий объем работы составляет 151страниц.Диссертация включает 72 рисунка и 13 таблиц. Библиография содержит 103 наименования.

2 СОДЕРЖАНИЕ РАБОТЫ

Во введении обоснованы актуальность поставленной научной задачи, показана степень разработанности темы решаемой научной задачи, сформулирована цель и задачи работы, показана научная новизна и практическая значимость полученных результатов.

Первая главапосвящена анализуособенностей моделирования излучения подстилающей поверхности в задачах создания систем ДЗЗ.

При моделировании подстилающей поверхности значения яркости объектов должны быть заданы в абсолютных величинах для решения задач оценки амплитуд реального сигнала, поступающего на вход БА. На примерах имитационных моделей, разработанных в ОАО «Корпорация «Комета» и ОАО МАК «Вымпел» показано, что для получения абсолютных величинможно использовать различные модели расчета спектральной яркости объектов.

Анализ применимости моделей расчета абсолютных значений яркостипоказал, что существующиемоделине позволяютимитироватьполе излучения в широком диапазоне условий наблюдения и освещения в рамках одной расчетной схемы, и в тоже время затрачивать минимальное количество вычислительных ресурсов, что является основнымтребованием при имитационном моделировании.

Так какподстилающая поверхность Земли, наблюдаемая с КА, характеризуется большим разнообразием, то без предварительной классификации объектов создание модели не представляется возможным.На основе таких факторов как назначение разрабатываемой аппаратуры, ее параметры и наличие исходных данных, а также на основе принципов агломеративного иерархического алгоритма классификации автором были выбраны 16 классов однотипных объектов, которые считаются «различимыми» в настоящей модели: пустыня, почва, зеленная растительность, водная поверхность, горы и облачность(Ci, Cc, Ac, As, St, Sc, Ns, Cu, Cb).

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

Геометрические условий, которые встречаются при наблюдении Земли из космосапоказаны на рисунке 1.

Рисунок 1 – Пример существующих условий наблюдения и освещения:

А) 1 - стандартные условия, 2- условия перспективной съемки, 3- условия терминатора, 6 – условия зеркального отражения,5- ночные условия;

Диапазоны изменения значений геометрических условий приведены в таблице 1.

Таблица 1- Диапазоны изменения геометрических условий наблюдения Условия наблюдения и Условия перспективной съемки Пригоризонтный ореол Зенитный угол наблюдения и зенитный угол Солнца отсчитываются от местной вертикали в наблюдаемой точке, азимутальный угол Солнца – от плоскости наблюдения до плоскости солнечного вертикала по часовой стрелке.

Проблемными в плане получения энергетических характеристик объектов по экспериментальным данным и по теоретическим расчетам являются условия зеркального отражения, пригоризонтного ореола Земли, перспективной съемки и терминатора из-за таких особенностей как: сложности полной интерпретации картины физических явлений из-за суперпозиции отраженной, рассеянной и тепловой составляющих яркости и существенные изменения в излучении объектов, вносимые большими оптическими толщами атмосферы.

Анализоптико-географических условий показал, что пропускание атмосферы и как следствие значения яркости объектов, может существенно зависеть от таких параметров как: тип аэрозоля, широта и сезон; особенно их влияние сказывается с увеличением оптической толщи атмосферы, что вызывает необходимостьих учета в задачахмоделирования наблюдаемой обстановки, особенно при глобальном обзоре.

Третья глава посвящена описанию разработанной модели расчета спектральной яркости объектов и методик расчета яркости входящих в нее.

Основная идея предлагаемой модели расчета яркости объектов состоит в замене сложных вычислительных процедур линейной интерполяцией между фиксированными значениями в информационных массивах (ИМ), представленных в виде таблиц. Информационные массивы могут формироваться при этом как с использованием эмпирическихданных, так и теоретических расчетов, покрывающих широкий диапазон условий. При выполнении простых требований эта модель позволит оптимизировать затраты вычислительных ресурсов и при этом обеспечить адекватность получаемых результатов.

Так как модель предполагается применять в составе имитационно-моделирующего комплекса полей излучения, то на выходе формируется в зависимости от типа аппаратуры наблюдения спектр яркости или/и интегральная яркость наблюдаемого объекта.

В ИМ хранятся спектры яркости, соответствующие условиям наблюдения, сгруппированные по типам объектов в соответствии с их классификацией. Выбор количества дискретных значений яркости в ИМ и соответственно шага по аргументам интерполяции определяется тем, что погрешность линейной интерполяции не должна превышать естественных вариаций яркости.

Функциональная схема предлагаемой модели приведена на рисунке 2.

На вход модели поступают исходные данные, полный перечень которых приведен таблице 2.

Анализ входных параметров начинается с определения типа объекта – простой(в элементе расчета присутствует излучение от одного объекта) или сложный(в элементе расчета присутствует излучение от двух или более объектов), если простой то переходят к анализу условий наблюдения, если сложный то определяется состав простых объектов и их соотношение. Здесь же определяются объекты, для которых характерно зеркальное отражение.

Во время анализа условий наблюдения проверяется взаимное расположение наблюдателя и источника освещения на соответствие условиям зеркального отражения ("блика").

Таблица 2 - Входные параметры модели * для удобства входных данных можно ввести Аsun, который пересчитывается в угол рассеяния cos cos Zobs cos Zsun + sin Zobs sin Zsun cos Asun ** в зависимости от типа облачности Если заданные условия не соответствуют условиям зеркального отражения или для объектов оно не характерно, то запрос поступает в ИМ1, при этом несмотря на различие в физических механизмахформирования яркости в зависимости от условий, в модели они объединены в единый ИМ, что позволяет исключить проблемыиспользования различных схем расчета в зависимости от условий наблюдения и освещения.

Порядок считывания информации сИМ по атрибутам начинается от самого быстро изменяющегося параметра к медленнее изменяющимся параметрам и имеет следующий порядок:

, Zsun,Zobs,, h.

Спектральная яркость объектов рассчитывается путем трилинейной интерполяции по зенитным углам Солнца, наблюдателя и углу рассеяния. Для облачности дополнительно производится линейная интерполяция по высоте облаков h. Высота наблюдателя не учитывается, так как исследования показали, что при наблюдении объектов из космоса вклад атмосферы выше 20-30 км пренебрежимо мал.

Если условия соответствуют условиям зеркального отражения, то яркость находится по методике расчета яркости в области зеркального отражения. Выделение условий блика отдельно от основной методики расчета яркости (см. рисунок 2) связано вариациями размера и интенсивности бликовой зоны в зависимости от состояния подстилающей поверхности, что приводит к необходимости ввода дополнительного параметра (например, скорость ветра для водной поверхности), который не предусмотрен в структуре ИМ1.

Если задан простой объект, то полученные значения спектральной яркости поступают на выход модели вместе с интегральной яркостью по всему заданному диапазону; если задан сложный объект, то рассчитанные спектральные яркости простых объектов в составе сложного поступают в блок расчета сложных объектов. Параллельно с этим вычисляются параметры затенения одним объектом другого, например облачность и подстилающая поверхность.

Неудобство предлагаемой схемы связано с различной генерализацией данных в зависимости от рассматриваемого диапазона - в окнах прозрачности спектры объектов более разнообразны, чем для полос поглощения, так как величина погрешности расчета или измерения в полосе поглощения превышают разность значений спектральной яркости объектов в пределах однородногокласса, то спектры сливаются, «укрупняя» классы. Это затрудняет получение спектра яркости каждого наблюдаемого объекта в относительно широком диапазоне, включающем и окна прозрачности, и крылья полос поглощения. Предлагаемая методика решает эту задачу путем ввода специальных правил перехода от одних спектральных поддиапазонов к другим для каждого класса объектов путем выбора репрезентативных спектров из узких классов для более широкого.

Так как одним из главных требований является быстродействие модели, была проведена оценка длительности расчета яркости по предложенной методике.Расчеты проводились на ПЭВМ со следующими характеристиками Intel®Pentium ® 2/4 Ghz, Windows 2000 в диапазоне длин волн 0,4-5,0 мкм.Длительность расчета по предлагаемой модели составила 4,21•10-9 с. Для сравнения также была измерена длительность расчета по методике, в основе которой лежит модель ГОИ им. С.И. Вавилова равна 1,516•10-4 с. Коэффициент ускорения предлагаемой модели относительно используемой в настоящее время в имитационной модели составил 36009,5. Естественно, что при имитации полноценного кадра время расчета существенно увеличивается, напримерпри моделировании одного кадра размером 8192х8192 методом линейной интерполяции время не превышает одной секунды, в то время как расчет методом, разработанным в ГОИ им. С.И. Вавилова составляет 2,8 ч. Расчет яркости путем решения УПИ методом Монте-Карло требует существенных временных ресурсов - расчет занимает от нескольких минут до суток в зависимости от различных параметров.

Для моделирования излучения объектов в области зеркального отражения, что позволит правильно выбрать уровень ограничения по уровню сигналов в системе обработки и динамический диапазон работы БА, предложена методика расчета яркости объектов в области зеркального отражения – «блика».

Расчет спектральной яркости в области «блика» производится только для объектов, для которых характерно зеркальное или зеркально-диффузное отражение: водная поверхность, облачность, снег.

Принцип расчета с учетом характера распределение яркости в области зеркального отражениясхематично показан на рисунке 3.

Рисунок 3 – Принцип расчета яркости в области зеркального отражения Для расчета яркости в области блика формируются дополнительный ИМ 2, в котором хранятся спектры яркости объектов непосредственно в точке зеркального отражения (ТЗО).

Фиксированные значения ИМ2 принимаются равными фиксированным значениям ИМ1.

Методом линейной интерполяции по ИМ2 находится яркость в ТЗО, которая определяется зенитным углом Солнца, по ИМ1 также методом линейной интерполяции находится яркость вне зоны блика с учетом максимально возможного размера бликовой области, которое отсчитывается от ТЗО. Под максимальным размером блика понимается максимально возможное расстояние от центра блика до его края, которое определяется максимально возможными наклонами отражающих площадок (угловой размер равен 2 ). Размер блика является фиксированным для каждого типа объекта и ограниченным, так как при удалении от ТЗО бликовая составляющая уменьшается настолько, что ее вклад становится незначительным и укладывается в естественные вариации яркости объекта.

Необходимо найти яркость для заданной произвольной точки М в области зеркального отражения по зенитному углу наблюдения Zobs (см. рисунок 3).

Для описания изменения яркости в области зеркального отражения выбрана аппроксимирующая функция:

где Z – зенитный угол линии визирования в искомой точке М;

параметры функции, определяющие характеристики бликовой области:

Z0 - зенитный угол линии визирования, при котором наблюдается точка зеркального отражения, B0 – спектральная яркость при Z0, Zn – зенитный угол линии визирования, при котором наблюдается конечная точка бликовой области, Bn - спектральная яркость при Zn;

m,, x0 - параметры функции, определяющие характер кривой.

B ( Z, Z 0, B0, Z n, Bn, m,, x0 ) для различных типов объектов и состояний их поверхности. Для определения параметров m,, x0 необходимо минимизировать невязку между измеренными или рассчитанными значениями яркости Bi (Zi ) и значениями B ( Z i, Z 0, B0, Z n, Bn, m,, x0 ), полученными по аппроксимирующей функции (1).

Для нахождения min L(m,, x0 ) была использована функция minimize, встроенная в паm,, x кет Mathcad, и определены параметры m,, x0. Так как в функции minimize используется метод сопряженных градиентов, существует вероятность пропуска экстремума функции, поэтому полученные значения использовались как приближенные и уточнялись для лучшего совпадения.

Для расчета яркости Bi ( Z i ) водной поверхности использовались экспериментальные данные, а для облачности и снежной поверхности использовались методики расчета яркости, в основе которых лежит расчет зеркальной составляющей по CoxC., MunkW, а диффузной по методике ГОИ им.И.С.Вавилова. Пример результатов аппроксимации показан на рисунке5.

Рисунок 5 – Аппроксимирующие кривые для водной поверхности при скорости ветра 10 км/ч Таким образом, получены значения параметров m,, x0 для водной поверхности(скорость ветра 2, 10 и 15 км/ч), снега и облачности, значения приведены в таблице 3.

Таблица 3– Параметры аппроксимации *- в связи с отсутствием доступных адекватных исходных данных, использовалось приближенные распределение яркости.

В разработанной модели предусмотрен расчет спектральной яркости сложных объектов, чтопозволяет получать значения энергетических характеристик для разнообразных наблюдаемых ситуаций, например, подстилающая поверхность с различными баллами облачности. Расчет яркости сложных объектов происходитпо следующей очевидной формуле:

где B 1 и B 2 - спектральная яркость первого и второго объектов соответственно; S1и S2– доли относительной площади первого и второго объектов в элементе расчета, задающиеся пользователем от 0 до 1.

При моделировании излучения объектов с различными уровнями высот (например, облачность и подстилающая поверхность) требуется поправка за затенение одного объекта другим.Для расчета спектральной яркостизатененного объекта использованы коэффициенты затененияk.

Значения kдля облачности получены по данным исследований В.И.Чапурского, а такжерассчитаны по гипрспектральным снимкам для зданий и древесной растительности,пример значений kприведен в таблице 4.

Таблица 4 – Значения коэффициентов затенения Процент площади затенения относительно всей наблюдаемой поверхности зависит от Zsun, Zobs, Аsun, а так же в случае облачности высоты нижней границы облачного поля и мощности облачности.

Учитывая, что в модели рассматриваются усредненные по пространству спектры яркости объектов, для определения площади тени относительно всей моделируемой площади,облачность можно представить как параллелепипед в основании которого – квадрат. Зная балл облачности, который задает пользователь, можно определить площадь, занимаемую облачностью и ее относительный размер.Тогда общая площадь тени от облака высотой h, длиной а и на высоте Н будет равна:

Необходимо определить долю площади тени, которую видит наблюдатель при данных условиях, так как часть тени может перекрываться самим облаком. Если зафиксировать наблюдателя и поворачивать Солнце по азимуту, то площадь перекрытия тени будет отличаться, от случая, когда Солнце зафиксировано, а наблюдатель движется. Для учета этого пользователь может во входных данных указать зафиксированный объект. Например, при движении Солнца, площадь перекрытия будет равна:

где l1 a H (tgZobs tgZsun cos Asun ) ;

Очевидно, что площадь затенения, которую видит наблюдатель, определяется как:

Для объектов, находящихся на земной поверхности, то есть для зданий и древесной растительности расчет затенения производится следующим образом.Общая площадь тени от зданий и лесов высотой h, длинойа равна: