Результат интеллектуальной деятельности: СПОСОБ РАЗВЕРТЫВАНИЯ ФАЗЫ ПРИ ИНТЕРФЕРОМЕТРИЧЕСКОЙ ОБРАБОТКЕ ИНФОРМАЦИИ ОТ КОСМИЧЕСКИХ СИСТЕМ РАДИОЛОКАЦИОННОГО НАБЛЮДЕНИЯ ЗЕМЛИ

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

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

Отсчеты интерферограммы (которые далее будем называть интерферометрической фазой или просто фазой), как правило, содержат шумы. Степень зашумленности сигнала на интерферограмме описывается матрицей значений когерентности. Единичная когерентность означает полное отсутствие шумов, а нулевая - максимальную зашумленность сигнала на интерферограмме.

Диапазон истинных значений интерферометрической фазы определяется диапазоном высот рельефа и может выходить за пределы периода волны излучения, т.е. может составлять более чем 2π радиан. В то же время на интерферограмме представлены лишь свернутые значения фазы, лежащие в пределах от -π до π радиан. Для получения матрицы истинной фазы к каждому свернутому значению фазы на интерферограмме необходимо прибавить некоторое целое количество периодов. Этот процесс определен в литературе как двумерное развертывание фазы. Развернутая фаза путем простого преобразования может быть впоследствии пересчитана в высоту точек земной поверхности.

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

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

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

Авторами указанного выше способа (патент US 5424743 А) предложена его модификация, где минимизируется не квадратичная (критерий наименьших квадратов), а любая p-норма отклонения градиента результирующей поверхности от свернутого градиента фазы (см. Ghiglia D.C., Romero L.A. Minimum L^p-norm two-dimensional phase unwrapping // Journal of the Optical Society of America A. 1996. Vol. 13. Issue 10. P. 1999-2013). Минимизация p-нормы сводится к многократному развертыванию фазы по критерию наименьших квадратов с изменяющимися весами. Правило изменения весов определяется величиной р. При р<2 паразитный градиент распространяется по меньшей части интерферограммы, и точность развертывания фазы повышается. Тем не менее, недостатки, отмеченные для способа по патенту US 5424743 А, сохраняются и в данном способе, только проявляются в меньшей степени. В предельном случае, когда р=0, области распространения паразитного градиента вырождаются в линии, однако их положение может не совпадать с истинным положением линий разрыва фазы, которое определяется не минимальной длиной линий, а конкретным сюжетом. По этой причине фаза на значительной части интерферограммы может оказаться неправильно развернутой.

Другим подходом к двумерному развертыванию фазы является способ Голдштейна-Зебкера-Вернера (см. Goldstein R.M., Zebker Н.А., Werner C.L. Satellite radar interferometry: Two-dimensional phase unwrapping // Radio Science. 1988. Vol. 23. No. 4. P. 713-720). Он основан на соединении сингулярных точек так называемыми линиями отсечения. Считается, что проведенные линии отсечения и есть линии разрыва фазы. Свернутый градиент фазы в пикселях, лежащих на этих линиях, игнорируется при интегрировании поля градиента. Соединение сингулярных точек линиями отсечения выполняется путем квазиминимизации суммарной длины линий. Эта операция выполняется правильно лишь тогда, когда линии разрыва фазы коротки и расположены далеко друг от друга, что характерно для гладкого рельефа земной поверхности и редких импульсных шумов на интерферограмме. Основным недостатком способа является то, что в случае сложного рельефа и сильного фазового шума, когда линии разрыва фазы имеют большую протяженность и расположены близко друг к другу, возникают ошибки развертывания фазы, распространяющиеся на значительную часть интерферограммы. Еще одним недостатком является то, что в наиболее сложных случаях проведенные линии отсечения могут изолировать друг от друга целые области на интерферограмме, а способ Голдштейна-Зебкера-Вернера не предлагает никаких решений, позволяющих получить единую фазовую поверхность по всей интерферограмме, а не отдельно по каждой изолированной области.

Еще одним подходом к развертыванию фазы является способ Constantini М. A novel phase unwrapping method based on network programming // IEEE Transactions on Geoscience and Remote Sensing. 1998. Vol. 36. No. 3. P. 813-821 (далее способ Константини). В нем также выявляются линии разрыва фазы, но градиент фазы в пикселях, лежащих на этих линиях, не игнорируется, а корректируется таким образом, чтобы сингулярные точки на концах линий взаимно устранились. Коррекция градиента осуществляется путем прибавления или вычитания к нему некоторого целого неотрицательного количества периодов. При этом линии проводятся так, чтобы взвешенная сумма числа периодов, прибавляемых или вычитаемых из градиента фазы во всех пикселях, была минимальна. Для этого строится граф, представляющий собой транспортную сеть, узлы которой располагаются в решетке, смещенной на половину пикселя по строке и столбцу относительно решетки пикселей интерферограммы. Соседние узлы сети попарно соединяются вертикальными и горизонтальными дугами. Узлы, совпадающие по положению с сингулярными точками, объявляются источниками или стоками в зависимости от знака ротора свернутого градиента фазы в сингулярной точке. Каждому источнику сопоставляется сток, в который от источника по дугам сети проводится единица потока. Стоимость потока определяется суммарной длиной дуг сети, по которым он проводится. Источники и стоки разбиваются на пары таким образом, чтобы суммарная стоимость всех потоков была минимальной. Пути проведения потоков определяют линии разрыва фазы, а величина потока в дугах транспортной сети определяет, на сколько периодов нужно корректировать градиент в соответствующих пикселях интерферограммы. В результате восстанавливается поле градиента развернутой фазы, обладающее свойством потенциальности, из которого может быть получена фазовая поверхность путем простого интегрирования. При этом появляется возможность выявления линий разрыва фазы с учетом дополнительной информации (например, о направлении и величине наклона фазовой поверхности, когерентности и т.д.) путем управления длиной дуг транспортной сети, за счет чего результаты способа Константини оказываются в среднем точнее результатов способа Голдштейна-Зебкера-Вернера. Кроме того, в отличие от способа Голдштейна-Зебкера-Вернера, способ Константини гарантирует получение единой фазовой поверхности даже в тех случаях, когда пути проведения потоков изолируют друг от друга отдельные области на интерферограмме, поскольку градиент фазы вдоль этих путей восстанавливается, а не игнорируется. Тем не менее, недостатком способа Константини является то, что строгая минимизация суммарной стоимости всех потоков не гарантирует правильного выявления линий разрыва фазы, особенно в тех случаях, когда длина этих линий велика, что характерно для широких крутых склонов на земной поверхности (чем шире крутой склон, тем длиннее линия разрыва фазы). Поэтому, если на снятом участке земной поверхности имеется много широких крутых склонов, пути проведения потоков не соответствуют истинному положению линий разрыва фазы и возникают ошибки развертывания, занимающие значительную часть интерферограммы.

Тем не менее, авторы заявки на изобретение считают способ Константини наиболее близким к их предложениям.

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

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

Неправильное выявление длинных линий разрыва фазы приводит к распространению ошибок развертывания фазы на значительную часть интерферограммы, поэтому очень важно правильно выявить такие линии. Истинное положение длинных линий разрыва фазы зависит от конкретного сюжета и может не соответствовать минимуму того или иного критерия (например, суммарной длины линий или суммарной стоимости потоков в транспортной сети). Следовательно, в общем случае невозможно правильно выявить длинные линии разрыва фазы путем решения различных оптимизационных задач, на которых основаны указанные выше существующие способы (US 5424743 А, способ на основе минимизации p-нормы отклонения градиента, способы Голдштейна-Зебкера-Вернера и Константини). Поэтому, в отличие от способа Константини и других указанных выше способов, предлагается более надежно определять положение длинных линий разрыва фазы путем анализа сюжета интерферограммы. Это одно из наиболее принципиальных положений предлагаемого авторами способа. Под анализом сюжета авторы подразумевают выделение участков на интерферограмме, в которых наблюдаются наиболее резкие изменения величины и направления наклона фазовой поверхности. На выделенных участках как раз и располагаются длинные линии разрыва фазы.

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

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

Если выделенные и проигнорированные участки не разделяют интерферограмму на несколько изолированных частей, то для получения единой поверхности развернутой фазы достаточно выполнить интерполяцию фазы на проигнорированных участках. Однако в наиболее сложных случаях интерферограмма может распасться на несколько областей, изолированных друг от друга проигнорированными участками. Фаза в каждой изолированной области оказывается развернутой с точностью до постоянной, и для получения единой фазовой поверхности требуется определить, как соотносятся между собой средние уровни фазы в изолированных областях. В крупных изолированных областях средний уровень развернутой фазы предлагается восстанавливать с привлечением существующих низкодетальных цифровых моделей рельефа (SRTM, ASTER GDEM и др.) в качестве опорной информации. Под крупными областями понимаются такие, на которые приходится несколько отсчетов низкодетальной опорной цифровой модели рельефа. Восстановление среднего уровня фазы в таких областях выполняется очень надежно и не приводит к появлению глобальных ошибок развертывания фазы. Таким образом, привлечение опорной низкодетальной информации о рельефе для восстановления среднего уровня фазы в крупных изолированных областях является еще одним конструктивным предложением.

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

Таким образом, конструктивность предлагаемого способа сводится к следующим положениям:

- предварительное выделение участков, содержащих длинные линии разрыва фазы, на основе анализа сюжета интерферограммы и игнорирование градиента фазы на этих участках;

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

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

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

1) вводят ограничение на максимальную стоимость потоков, проводимых в транспортной сети, ассоциированной с интерферограммой, и восстанавливают истинный градиент фазы на основе проведенных потоков только в окрестности коротких линий разрыва фазы;

2) оценивают наклон фазовой поверхности по отфильтрованному свернутому градиенту фазы, выявляют пиксели на интерферограмме, в окрестности которых наиболее сильно изменяется модуль или направление наклона фазовой поверхности, и назначают данным пикселям нулевые веса; остальным пикселям веса назначают в соответствии с когерентностью;

3) выполняют восстановление поверхности развернутой фазы из частично восставленного ее градиента по критерию наименьших квадратов с учетом назначенных весов;

4) формируют разностную интерферограмму путем попиксельного вычитания свернутой фазы из восстановленной поверхности развернутой фазы и приведения результата вычитания к диапазону от -π до π радиан;

5) выделяют пиксели, в которых наблюдается локальный максимум градиента фазы в строке или столбце разностной интерферограммы (который и является паразитным градиентом), и обнуляют веса в них, если локальный максимум превышает среднее значение паразитного градиента;

6) повторяют действия 3-5 с уточненными весами, до тех пор пока паразитный градиент не станет пренебрежимо мал (обычно требуется от одной до десяти таких итераций в зависимости от сложности рельефа);

7) для каждой крупной области, изолированной от остальной интерферограммы пикселями нулевого веса, приводят средний уровень развернутой фазы максимально близко к среднему уровню фазовой картины опорного рельефа с ограничением, что развернутая фаза отличается от свернутой на целое количество периодов;

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

9) выполняют интерполяцию фазы в пикселях нулевого веса по значениям развернутой фазы в ближайших пикселях ненулевого веса.

Предложенный способ апробирован на натурной информации от космических радиолокаторов Radarsat-2 (Канада), COSMO-SkyMed (Италия), TerraSAR-X (Германия). На фиг. 1 и фиг. 2 в качестве примера приводятся пересчитанные в высоту результаты развертывания фазы на интерферограмме, полученной по паре изображений от радиолокатора на космическом аппарате Radarsat-2. Более высокие участки земной поверхности показаны более светлым оттенком. На фиг. 1 представлены результаты способа Константини, а на фиг. 2 - предлагаемого способа. Как видно, способ Константини привел к образованию крупных областей, в которых высота завышена или занижена на 27 м из-за неправильного развертывания фазы; эти области обведены белой линией на фиг. 1. В то же время предложенный способ развернул фазу без глобальных ошибок. Отсчеты высоты, полученные по способу Константини и по способу, предлагаемому авторами, сравнивались с цифровой моделью рельефа SRTM1 версии 3. В результате установлено, что среднеквадратическое отклонение оценок высот для способа Константини составило 9,1 м, а для предложенного способа - 4,7 м.