Обрезка Landsat

Nule
Завсегдатай
Сообщения: 276
Зарегистрирован: 26 дек 2011, 12:26
Репутация: 3

Обрезка Landsat

Сообщение Nule » 03 авг 2012, 09:39

Здравствуйте!
Задача - пакетно обрезать снимки Landsat, с целью убрать артефакты с боков.
Т.е. надо наложить маску такого вида:
Изображение
Причем маска должна быть "универсальной", т.е. не зависеть от координат обрабатываемой сцены.
Можно ли это реализовать в GRASS?
Спасибо!

Александр Мурый
Гуру
Сообщения: 5173
Зарегистрирован: 26 сен 2009, 16:26
Статьи: 3
Проекты: 5/1
Репутация: 793
Ваше звание: званий не имею
Откуда: Москва

Re: Обрезка Landsat

Сообщение Александр Мурый » 03 авг 2012, 11:26

Может, и можно. Нужны примеры снимков с "артефактами" (устно не совсем понятно).
Редактор материалов, модератор форума

Nule
Завсегдатай
Сообщения: 276
Зарегистрирован: 26 дек 2011, 12:26
Репутация: 3

Re: Обрезка Landsat

Сообщение Nule » 03 авг 2012, 11:32

amuriy писал(а):Может, и можно. Нужны примеры снимков с "артефактами" (устно не совсем понятно).
Примером может служить абсолютно любой снимок landsat.
Изображение
Изображение

ericsson
Гуру
Сообщения: 3321
Зарегистрирован: 27 июл 2009, 19:26
Репутация: 746
Ваше звание: Вредитель полей

Re: Обрезка Landsat

Сообщение ericsson » 03 авг 2012, 11:52

Nule писал(а):Здравствуйте!
Причем маска должна быть "универсальной", т.е. не зависеть от координат обрабатываемой сцены.
А кто сказал, что маска может быть универсальной? Вот browse-файлы из архива GLCF, одного спутника и одного класса обработки. Похоже, но не совпадает:
ftp://ftp.glcf.umd.edu/glcf/Landsat/WRS ... browse.jpg
ftp://ftp.glcf.umd.edu/glcf/Landsat/WRS ... browse.jpg
ftp://ftp.glcf.umd.edu/glcf/Landsat/WRS ... browse.jpg

Я бы при решении задачи обрезки мусора сделал бы следующее:
Для каждой сцены извлек бы незначащие точки из каждого канала в композите (можно даже из уменьшенной в два или четыре раза копии), инвертировал, нашел бы пересечение получившихся масок значащих данных, инвертировал лезультат, (увеличил бы во столько раз во сколько уменьшал в начале), добавил бы некоторый threshold и по получившемуся - обрезал бы композит. GRASS не владею, но раз уж это можно сделать в ImageMagick (коим пользуюсь для обычной графики в т.ч. в подобных задачах), то и в GRASS, наверняка, можно.

Nule
Завсегдатай
Сообщения: 276
Зарегистрирован: 26 дек 2011, 12:26
Репутация: 3

Re: Обрезка Landsat

Сообщение Nule » 03 авг 2012, 12:03

В ERDAS IMAGINE есть универсальная маска в инструментах Mosaic Pro (маску можно нарисовать вручную) и Mosaic Express (можно задать, сколько процентов срезать). Но это инструменты для создания мозаик и использовать их для обрезания отдельного снимка неудобно.
Я надеялся, что может есть аналогичные, но более гибкие инструменты в других пакетах.

ericsson
Гуру
Сообщения: 3321
Зарегистрирован: 27 июл 2009, 19:26
Репутация: 746
Ваше звание: Вредитель полей

Re: Обрезка Landsat

Сообщение ericsson » 03 авг 2012, 12:51

Ну так а что может быть более гибким, чем пара-тройка булевых операций над растрами?

Александр Мурый
Гуру
Сообщения: 5173
Зарегистрирован: 26 сен 2009, 16:26
Статьи: 3
Проекты: 5/1
Репутация: 793
Ваше звание: званий не имею
Откуда: Москва

Re: Обрезка Landsat

Сообщение Александр Мурый » 03 авг 2012, 14:03

ericsson писал(а): Я бы при решении задачи обрезки мусора сделал бы следующее:
Для каждой сцены извлек бы незначащие точки из каждого канала в композите (можно даже из уменьшенной в два или четыре раза копии), инвертировал, нашел бы пересечение получившихся масок значащих данных, инвертировал лезультат, (увеличил бы во столько раз во сколько уменьшал в начале), добавил бы некоторый threshold и по получившемуся - обрезал бы композит.
Здорово было бы увидеть этот метод в действии с помощью IM. Хотя бы на примере приведённого выше увеличенного куска снимка. А самое интересное -- совмещение модулей GRASS + утилит IM.
Редактор материалов, модератор форума

ericsson
Гуру
Сообщения: 3321
Зарегистрирован: 27 июл 2009, 19:26
Репутация: 746
Ваше звание: Вредитель полей

Re: Обрезка Landsat

Сообщение ericsson » 03 авг 2012, 16:35

Постараюсь сделать в ближайшие дни.

ericsson
Гуру
Сообщения: 3321
Зарегистрирован: 27 июл 2009, 19:26
Репутация: 746
Ваше звание: Вредитель полей

Re: Обрезка Landsat

Сообщение ericsson » 05 авг 2012, 15:30

Ну собственно, как обещал...
Требуемые средства - listgeo, geotifcp (из Libgeotiff), convert (свежий ImageMagick).

1) Сделать дамп метаданных любого канала в отдельный текстовик через listgeo

2) Для каждого канала, который предполагается использовать в композите, выполняем команду такого вида:

Код: Выделить всё

convert L71169028_02820060620_B10.TIF -alpha set -channel RGBA -fuzz 5% -fill none -floodfill +0+0 black L71169028_02820060620_B10.withalpha.TIF
Параметр fuzz отвечает за то чтобы захватить не совсем черные края.
На выходе получаем тот же файл с убитыми метаданными (привязкой) но двухканальный - восьмибитный greyscale, плюс маска там, где был черный ±5%.
Примечание: если используются "дурные" данные с вырезанными полосами и прочим безобразием (которые с поврежденного спутника, если я не путаюсь в причине их возникновения), можно использовать другой код:

Код: Выделить всё

convert L71169028_02820060620_B10.TIF  -transparent black  L71169028_02820060620_B10.withalpha.TIF
Он пропишет прозрачность черным точкам не только по краям растра, но и изолированным незначащим точкам среди значащих.

3) Отделяем маску командой

Код: Выделить всё

convert.exe L71169028_02820060620_B10.withalpha.TIF -alpha extract L71169028_02820060620_B10.alpha.TIF
И получаем однобитные маски в виде отдельных файлов.

4) Попарно (1+2)+3 складываем маски командами типа

Код: Выделить всё

convert L71169028_02820060620_B10.alpha.TIF L71169028_02820060620_B20.alpha.TIF -compose Darken -composite L71169028_02820060620_B10_B20.alpha.TIF

convert L71169028_02820060620_B10_B20.alpha.TIF L71169028_02820060620_B30.alpha.TIF -compose Darken -composite L71169028_02820060620_B10_B20_B30.alpha.TIF
Наложение растров Darken эквивалентно логическому пересечению, если незначащие части - черные, значащие - белые.

5) Убиваем получившиеся файлы *_withalfa.TIF - они больше не нужны.

6) Через geotifcp вернуть файлу L71169028_02820060620_B10_B20_B30.alpha.TIF привязку из дампа, который был сделан на первом шаге. Далее эту маску можно использовать где угодно вместе с исходными растрами.

Александр Мурый
Гуру
Сообщения: 5173
Зарегистрирован: 26 сен 2009, 16:26
Статьи: 3
Проекты: 5/1
Репутация: 793
Ваше звание: званий не имею
Откуда: Москва

Re: Обрезка Landsat

Сообщение Александр Мурый » 05 авг 2012, 16:05

Спасибо большое, ericsson. Будем пробовать :)
Редактор материалов, модератор форума

Александр Мурый
Гуру
Сообщения: 5173
Зарегистрирован: 26 сен 2009, 16:26
Статьи: 3
Проекты: 5/1
Репутация: 793
Ваше звание: званий не имею
Откуда: Москва

Re: Обрезка Landsat

Сообщение Александр Мурый » 06 авг 2012, 08:31

Создание подобной маски для всех каналов собственно в GRASS (взято с отсюда):

Код: Выделить всё


BASE=L71074092_09220040924

g.region rast=$BASE.1
r.series in=`g.mlist pat="$BASE.[0-8]*" sep=,` -n out=$BASE.thresh method=threshold thresh=1
r.mapcalc "$BASE.mask = if(isnull($BASE.thresh))"
g.remove $BASE.thresh

Проверено -- работает, причём довольно шустро. :)
На заметку: надо следить за именами растров (КО: они могут быть не как в примере).

P.S. А вот как вообще обрезать "бахрому" на автомате?
Редактор материалов, модератор форума

ericsson
Гуру
Сообщения: 3321
Зарегистрирован: 27 июл 2009, 19:26
Репутация: 746
Ваше звание: Вредитель полей

Re: Обрезка Landsat

Сообщение ericsson » 06 авг 2012, 08:42

Ну еще можно вооружиться полигонами сетки path/row (где-то в недрах архива GLCF такие шейпы были), извлекать соответствующие номера из имен файлов или сопутствующих метаданных и обрезать по этим полигонам.

Nule
Завсегдатай
Сообщения: 276
Зарегистрирован: 26 дек 2011, 12:26
Репутация: 3

Re: Обрезка Landsat

Сообщение Nule » 06 авг 2012, 09:06

Странно, казалось бы для такой массовой вещи как ландсат должна быть куча разнообразнейших тулзовин.
Что-нибудь типа: находим самый верхний пиксель, отличный от noData, потом нижний, правый и левый, по ним строим векторный прямоугольник, уменьшаем его на 10%, обрезаем по шейпу.

Аватара пользователя
rhot
Гуру
Сообщения: 1727
Зарегистрирован: 25 янв 2011, 17:50
Статьи: 1
Репутация: 194
Ваше звание: доктор
Откуда: Архангельск

Re: Обрезка Landsat

Сообщение rhot » 06 авг 2012, 09:24

amuriy писал(а):Создание подобной маски для всех каналов собственно в GRASS (взято с отсюда):

Код: Выделить всё


BASE=L71074092_09220040924

g.region rast=$BASE.1
r.series in=`g.mlist pat="$BASE.[0-8]*" sep=,` -n out=$BASE.thresh method=threshold thresh=1
r.mapcalc "$BASE.mask = if(isnull($BASE.thresh))"
g.remove $BASE.thresh

Проверено -- работает, причём довольно шустро. :)
На заметку: надо следить за именами растров (КО: они могут быть не как в примере).

P.S. А вот как вообще обрезать "бахрому" на автомате?
Объясните, пожалуйста, что это за метод threshold?
Nule писал(а):Странно, казалось бы для такой массовой вещи как ландсат должна быть куча разнообразнейших тулзовин.
А зачем? Таких задач бывает так же много, как пользователей ландсат. ИМХО, программа не должна быть перегружена всякими примочками и тулзовинами. GRASS в этом смысле очень позитивно смотрится; исходники - чуть больше 20 Мб.
Хочешь, качай и устанавливай доп. модули, а хочешь, сам пиши на питоне - свобода :) 8)
___________(¯`·.¸(¯`·.¸ Scientia potentia est _/ {SILVA}:::{FOSS}:::{GIS} \_ Знание сила ¸.·´¯)¸.·´¯)___________

Александр Мурый
Гуру
Сообщения: 5173
Зарегистрирован: 26 сен 2009, 16:26
Статьи: 3
Проекты: 5/1
Репутация: 793
Ваше звание: званий не имею
Откуда: Москва

Re: Обрезка Landsat

Сообщение Александр Мурый » 06 авг 2012, 11:25

rhot писал(а): Объясните, пожалуйста, что это за метод threshold?
Честно говоря, для меня это тоже не до конца прояснилось. Кроме мануала, пояснения от автора модуля (Glynn Clements) есть, например, здесь.
Редактор материалов, модератор форума

Ответить

Вернуться в «GRASS»

Кто сейчас на конференции

Сейчас этот форум просматривают: нет зарегистрированных пользователей и 1 гость