1
Geoinformatyka
Ćwiczenie 6
Podstawy przetwarzania zdjęć satelitarnych
Katarzyna Ostapowicz, Jacek Kozak
kostapowicz@gis.geo.uj.edu.pl
;
jkozak@gis.geo.uj.edu.pl
Zakład Systemów Informacji Geograficznej, Kartografii i Teledetekcji
Instytut Geografii i Gospodarki Przestrzennej UJ
Kraków 2011
2
Cel
Celem ćwiczenia jest zapoznanie się z możliwościami wykorzystania zdjęć satelitarnych jako źródła
informacji o środowisku przyrodniczym. Przedstawione zostaną zasady wstępnej interpretacji zdjęć
satelitarnych m.in. zmian w czasie, oraz przeprowadzania klasyfikacji wzorcowej (nadzorowanej)
zdjęcia satelitarnego.
Słowa klucze
Jaskrawość (DN), rozdzielczość: przestrzenna, czasowa, radiometryczna i spektralna, zakresy
spektralne, histogram, kompozycja barwna, sygnatury, analiza zmian (change detection), pola
treningowe, wzorce/sygnatury, klasyfikacja wzorcowa/nadzorowana
Po tych zajęciach powinniście potrafić:
rozróżniać cztery typy rozdzielczości,
rozróżniać obiekty występujące na zdjęciu,
przeprowadzić wstępną analizę zmian na zdjęciu satelitarnym,
wyznaczyć pola treningowe,
utworzyć sygnatury spektralne dla pól treningowych,
ocenić jakość sygnatur,
przeprowadzić procedurę klasyfikacji
zweryfikować wynik klasyfikacji.
Dane
Część 1:
•
ortofotomapa satelitarna SPOT 5 (fragment sceny z 24 sierpnia 2004, kanał PAN) zapisana w
formacie *.img, źródło: Zakład GIS, Kartografii i Teledetekcji, IGiGP UJ
•
ortofotomapa satelitarna EO - 1 ALI (fragment sceny 188/25 z 15 czerwca 2006, dziewięć
kanałów 1 – 9) zapisana w formacie *.img, źródło:
http://glovis.usgs.gov/
•
ortofotomapa satelitarna Landsat ETM+ (fragment sceny 187/26 z 20 sierpnia 2000 roku,
siedem
kanałów
1
-
7)
zapisana
w
formacie
*.img,
źródło:
http://edcsns17.cr.usgs.gov/EarthExplorer/
Część 2:
Dwie otofotomapy satelitarne Landsat (TM – z roku 1990, oraz ETM+ z roku 2000). Pochodzą one z
serwera udostępniającego bezpłatne zdjęcia satelitarne Landsat: Global Land Cover Facility
(Uniwersytet Maryland):
http://www.landcover.org
.
Część 3:
Dane podstawowe:
•
ortofotomapa satelitarna Landsat TM (fragment sceny 188/25 z 19 lipca 2006 roku, trzy
kanały 2, 3, 4) zapisana w formacie *.img, źródło:
http://edcsns17.cr.usgs.gov/EarthExplorer/
3
Dane dodatkowe: baza danych CORINE Land Cover dla roku 2000
1
, ortofotomapa z okolic Krakowa,
mapa topograficzna w skali 1:25 000 okolic Krakowa
Oprogramowanie
ArcGIS 9.3 z rozszerzeniem Spatial Analyst
Zestaw narzędzi ArcToolbox: Multivariate, wybrane elementy z pasków narzędzi Tools, Editor oraz
Effects oraz polecenia z menu kontekstowego warstw w Table of Contents aplikacji ArcMap.
Wstęp
Zdjęcia satelitarne mogą służyć m.in. jako źródło informacji o pokryciu terenu (land cover) czyli
fizycznym stanie powierzchni ziemi, który odnosi się do jej właściwości biofizycznych (rozmieszczenia
roślinności, pustyń, wód) lub właściwości związanych z działalnością człowieka (np. zabudowa,
infrastruktura techniczna).
Tok postępowania
Przegraj katalog cwiczenie_6 z danymi do dzisiejszego ćwiczenia z dysku studenci w katalogu
Geoinformatyka_zaoczne na dysk lokalny E:\
Uruchom program ArcGIS>ArcMap
CZĘŚĆ 1: PODSTWOWE WŁAŚCIWOŚCI ZDJĘĆ SATELITARNYCH
Rozdzielczość zdjęć satelitarnych
Rozróżniane są cztery rozdzielczości charakteryzujące systemy teledetekcyjne: rozdzielczość
przestrzenna, rozdzielczość spektralna, rozdzielczość radiometryczna i rozdzielczość czasowa.
W tej części ćwiczenia przeprowadzisz interpretację zdjęć satelitarnych o różnych rozdzielczościach.
Rozdzielczość przestrzenna jest miarą najmniejszego obiektu, który jest/może być rozróżniany przez
sensor lub obszarem na powierzchni ziemi reprezentowanym przez piksel. Rozdzielczość
przestrzenna zależy od chwilowego pola widzenia (instantaneous field of view – IFOV) sensora.
Wyświetl projekt teledetekcja1.mxd (katalog dane_cz1).
Na podstawie zdjęć o różnej rozdzielczości przestrzennej - Landsat ETM+ (30 m) i SPOT 5 (5 m) –
przeprowadź wizualną analizę zapisu różnego typu obiektów na zdjęciach. Jakiego typu obiekty
potrafisz wyróżnić na zdjęciach?
Rozdzielczość radiometryczna - liczba rozróżnianych poziomów promieniowania (liczba bitów, na
które podzielono zarejestrowaną energię w danym zakresie spektralnym). Rozdzielczość spektralna
określa przedział rejestrowanego promieniowania elektromagnetycznego (przedział długości fal
widma elektromagnetycznego), w których dany sensor może rejestrować energię.
Wywołaj opcję Properties, zakładkę Source i porównaj zapis zdjęć satelitarnych Landsat ETM+ i EO –
1 ALI, jakie kanały spektralne i ile poziomów rozróżnianych jest w tych zdjęciach satelitarnych? Czy
prawdą jest, że Landsat zapisywany jest jako obraz 8-bitowy (pojedynczy kanał: 256 poziomów - 2
8
) a
EO – 1 ALI jako obraz 16-bitowy (pojedynczy kanał: 65536 poziomów – 2
16
). Jakie wartości DN
zapisane są w poszczególnych kanałach spektralnych obydwu zdjęć? W tabeli 1 znajdziesz
zestawienie kanałów spektralnych rejestrowanych przez analizowane sensory.
1
W bazie danych CORINE Land Cover (CLC) znajdują się mapy wektorowe pokrycia terenu i użytkowania ziemi, wykonane na podstawie
interpretacji wizualnej zdjęć satelitarnych. Wyróżniane, w procesie tworzenia bazy danych powierzchnie o jednorodnym pokryciu terenu i
użytkowaniu ziemi przypisywano do jednej z 44 klas (załącznik 1) zgrupowanych w trójstopniowym systemie hierarchicznym. Utworzone
mapy odpowiadają dokładnością i szczegółowością mapie w skali 1:100 000 (minimalna jednostka kartowania wynosi 25 ha) i zostały do tej
pory wykonane dla trzech momentów czasowych, przełomu lat osiemdziesiątych i dziewięćdziesiątych XX wieku (CLC1990), przełomu XX i
XXI wieku (na podstawie zdjęć satelitarnych Landsat i SPOT) oraz roku 2006 (na podstawie zdjęć satelitarnych SPOT i IRS).
4
Tabela 1. Kanały spektralne rejestrowane przez sensory EO-1 ALI i Landsat ETM+.
Kolorystyka zdjęcia i schemat RGB
Zdjęcia satelitarne wykonane w wielu różnych zakresach spektralnych pokazuje się najczęściej tak, że
wybrane zakresy spektralne przyporządkowuje się do trzech kolorów wyświetlania monitora:
czerwonego (R), zielonego (G), i niebieskiego (B)
W przypadku naszego zdjęcia Landsat ETM+, przyporządkowanie wygląda tak:
•
rejestracja w zakresie podczerwonym (TM4) -> R
•
rejestracja w zakresie czerwonym (TM3) -> G
•
rejestracja w zakresie zielonym (TM2) -> B
ponieważ akurat roślinność jest ‘najjaśniejsza’ w podczerwieni, stąd dominującym kolorem jest
czerwony (bo ten jest przyporządkowany jako kolor wyświetlania do podczerwieni.
Przypisanie kolorów można łatwo zmienić dla warstwy (zdjęcia): Properties > Symbology. Zmień
ustawienia zdjęcia EO – 1 ALI tak aby kompozycja barwna odpowiadała tej ze zdjęcia satelitarnego
Landsat ETM+.
Oceń wizualnie na ile zmiana kompozycji barwnej zwiększyła/zmniejszyła rozróżnialność obiektów na
zdjęciu.
Porównanie zawartości informacyjnej różnych kanałów spektralnych i histogram
Promieniowanie elektromagnetyczne odbija się od obiektów terenowych w różny sposób, w
zależności od długości fali i fizycznych właściwości tych obiektów (Ryc. 1). Na wykładzie
przedstawiono zapis typowego obrazu satelitarnego – przypomnij sobie jak on wygląda.
Ryc. 1. Typowe charakterystyki spektralne wody (linia niebieska), roślinności zielonej (linia czerwona) i
odkrytego gruntu (linia brązowa) (źródło: Longley 2006).
5
Najprościej odczytać rozkład wartości (m.in. cechy statystyczne: wartości DN minimalne,
maksymalne, amplitudę, średnią, odchylenie standardowe lub wariancję) zakresu spektralnego
zdjęcia (a także dowolnej mapy rastrowej) za pomocą histogramu.
Histogram to wykres przedstawiający częstość występowania określonych wartości macierzy
rastrowej. Na osi poziomej mamy wartości (albo ich przedziały), na osi pionowej – liczbę wystąpień o
danej wartości, albo w danym przedziale.
Sprawdź jakie wartości pojawiają się dla kanałów spektralnych 2, 3, 4 zdjęcia satelitarnego Landsat
EMT+ oraz kanałów 3, 4, 5 zdjęcia satelitarnego EO - 1 ALI dla analizowanego fragmentu Karpat (Ryc.
2).
! Większość poszukiwanych wartości odczytasz dla każdego zakresu zdjęcia wywołując opcję
Properties, zakładkę Symbology i opcję Histogram (rycina poniżej).
Ryc. 2. Przykładowy histogram dla obrazu Landsat.
Zwróć uwagę, że maksimum histogramu to wartość 255 (Landsat) lub 65535 (możliwa dla ALI) dla
wszystkich trzech zakresów (dlaczego ?); jednakże analiza histogramów wskazuje, że na obrazie nie
ma zbyt wielu pikseli o wartościach np. > 150 dla zdjęcia Landsat, a większość pikseli skupia się w
dość wąskim zakresie wartości. Na podstawie histogramów spróbuj wskazać, który z zakresów
cechuje największe zróżnicowanie często występujących wartości (tj. takich, które mają co najmniej
500 reprezentantów [pikseli] – wartość ta jest dobrana zupełnie arbitralnie). Przesuwając kursor
poziomo przez obszar wykresu, obserwuj zmiany wartości ‘Input’ oraz ‘Count in’ – pokazują one, ile
piskeli (Count In) o danej wartości (Input) jest na obrazie. Obliczając różnicę między tak wyznaczoną
wartością maksymalną i minimalną, określ zakres wartości faktycznie występujących na zdjęciu
Który z kanałów ma największą pojemność informacyjną (niesie najwięcej informacji) ? Dlaczego ?
W dalszej części ćwiczenia będziesz korzystać tylko ze zdjęć Landsat.
Profile spektralne pikseli
Na podstawie wizualnej oceny zdjęć o dużej rozdzielczości potrafisz rozpoznać przynajmniej kilka
podstawowych typów pokrycia terenu: zabudowa, użytki rolne, lasy, woda.
Wykorzystaj narzędzie Identify (pasek narzędzi Tools) odczytaj wartości z trzech dowolnych pikseli
reprezentujących wymienione powyżej klas pokrycia terenu, odczytane wartości zapisz w Excelu i
utwórz wykresy przedstawiające profile (ryciny poniżej).
6
! Zakres TM2 przypisany jest do warstwy ‘Blue’, zakres TM3 – do warstwy ‘Green’, zakres TM4 – do
warstwy ‘Red’. W tabeli Identify pojawiają się tylko wartości z kanałów wybranych przy ustawianiu
kompozycji barwnej.
Zwróć uwagę na to że:
•
roślinność (lasy) cechują najwyższe wartości jaskrawości w zakresie TM4 (bliska podczerwień), a
niskie w TM2 i TM3
•
woda ma bardzo niskie wartości w tym zakresie
•
woda i lasy mają podobne wartości jaskrawości w zakresie TM2 i TM3
•
zabudowa i użytki rolne to klasy o dużym zróżnicowaniu wartości we wszystkich zakresach,
zabudowa często daje wartości najwyższe
Wynika to z fizycznych własności konkretnych materiałów, i tego, jak reagują na promieniowanie o
określonej długości fali.
CZĘŚĆ 2: ANALIZA ZMIAN W CZASIE
Do tej pory analizowałeś zdjęcia dla wybranych momentów czasowych nie porównując ich między
sobą. Teraz przejdziesz do prostej analizy zmian w czasie (analiza ilościowa).
Rozdzielczość czasowa (powtarzalność zdjęć) - określa częstotliwość pozyskiwania danych, to
częstotliwość z jaką dany satelita obrazuje ten sam fragment powierzchni kuli ziemskiej. Z jakiego
powodu jest to ważna informacja?
Automatyczna analiza zmian
Wyświetl projekt teledetekcja2.mxd (katalog dane_cz2)
7
W tej części ćwiczenia przeanalizujesz zmiany powierzchni leśnej na obszarze położonym na granicy
Paragwaju, Brazylii oraz Argentyny (nieco na zachód od znanych wodospadów Iguazu). Obszar ten
podlegał szybkiemu wylesianiu pod koniec XX wieku.
Użyj znanych Ci opcji Effects w celu wizualnego porównania obu zdjęć. Czy potrafisz zlokalizować
miejsca największych wylesień?
Więcej informacji o tym obszarze, w kontekście wylesiania można znaleźć tu:
http://na.unep.net/OnePlanetManyPeople/
(obszar opisany jest w rozdziale 3.4 „Forests”)
Analizę rozpocznij od ustalenia jakie są typowe wartości dla obszarów leśnych, użytków rolnych oraz
wód na danym obszarze. W interpretacji kieruj się kolorystyką zapamiętaną z poprzedniego ćwiczenia
– tak jak poprzednio, masz do dyspozycji zdjęcie Landsat, i tak jak poprzednio – kanały TM2, TM3 i
TM4, i podobnie jak wcześniej są one przypisane do warstw R, G, B. Ograniczymy się w tym ćwiczeniu
wyłącznie do zakresu TM3 – a więc dla obu zdjęć, zwróć uwagę na wartości w polu ‘Green’.
Znajomość tych wartości pozwoli Ci w późniejszym etapie analizy ustalić miejsca, w których nastąpiła
zmiana.
? Zastanów się czy różnią się od siebie warstwy tm.img i tm,img - layer2
Zwróć uwagę, że wartości różnią się między zdjęciami – wynika to z różnych charakterystyk
technicznych sensora ETM+ i TM.
Landsat TM
Forma pokrycia
Jaskrawości pikseli (od-do) [TM3]
Woda
13-15
Użytki rolne
20-45
Lasy
14-17
8
Landsat ETM+
Forma pokrycia
Jaskrawości pikseli (od-do) [TM3]
Woda
30-40
Użytki rolne
30-70
Lasy
22-27
Zastanów się jakie (mniej więcej) wartości uzyskasz, odejmując od wartości zdjęcia (ta informacja
pozwoli Ci zlokalizować miejsca, w których nastąpiły zmiany) z roku 2000 wartości zdjęcia z roku
1990, w zakresie TM3, w tych miejscach, gdzie nastąpiło wylesienie (patrz tabele poniżej) ? UR –
oznacza użytki rolne; W – wylesienie, BZ – bez zmian. Jakie wartości uzyskasz tam, gdzie zarówno w
2000, jak i 1990 roku był las ? A jakie tam, gdzie w obu momentach były użytki rolne?
ETM+
TM
różnica
UR
UR
UR
LAS
UR
LAS LAS LAS
BZ
W
W
BZ
UR
UR
UR
LAS
UR
LAS LAS LAS
BZ
W
W
BZ
UR
UR
UR
LAS
UR
LAS LAS LAS
BZ
W
W
BZ
UR
UR
LAS LAS
UR
LAS LAS LAS
BZ
W
BZ
BZ
Spróbuj wypełnić puste pola poniżej.
ETM+]
TM
różnica
UR
UR
UR
LAS
UR
LAS LAS LAS
UR
UR
UR
LAS
UR
LAS LAS LAS
UR
UR
UR
LAS
UR
LAS LAS LAS
UR
UR
LAS LAS
UR
LAS LAS LAS
Co do wylesień, wg tabeli wartości mogą mieć zakres od 13 (30-17) do 56 (70-14)
Operację tą powtórzymy w komputerze.
Upewnij się, że w menu głównym Tools > Extensions zaznaczona jest opcja ‘Spatial Analyst’
Uruchom Arc Toolbox.
Uruchom w ArcToolbox opcję Spatial Analyst Tools > Map Algebra > Raster Calculator. Uruchamia się
okno analizy danych, w którym można wpisywać różnego typu komendy analityczne (więcej na ten
temat – kolejne ćwiczenia).
Z pola Layers and variables wybierz z listy zakres TM3 obu zdjęć (w obu wypadkach jest to Layer_2).
Wskaż poprawną lokalizację i nazwę dla nowo tworzonej mapy (w katalogu roboczym, nazwa np.
‘roznica’)
9
Przeciągnij nazwy warstw do górnego okna, pomiędzy nimi wstaw znak odejmowania (-); tak aby
powstało działanie typu:
etm.img(layer_2) – tm.img(layer_2)
Okno powinno wyglądać tak:
Wybierz OK, i poczekaj, aż po wykonaniu działania dodany zostanie nowy obraz do projektu.
Zbadaj zakres wartości nowo powstałego obrazu. Czy wartości większe od 20 występują tam gdzie
wylesienia ?
W kolejnym kroku wykonamy mapę wylesień. W tym celu w oknie Raster Calculatora wpiszemy
działanie:
roznica > 20
Wykonaj to analogicznie jak poprzednio, wybierając obraz ‘roznica’, deklarując obraz wynikowy, a
następnie wpisując działanie.
10
Efektem będzie nowy obraz (np. wylesienia), który jako zero potraktuje wszystkie wartości obrazu
roznica <= 20, a jako 1 – wartości większe. Zgodnie z naszą analizą, wylesienia powinny mieć wartości
>20.
Przeanalizuj wynik. Czy wszystko wygląda tak, jak powinno ? Jakie mankamenty widzisz na mapie
wylesień (poniżej – potencjalne wylesienia zaznaczone są na żółto, na zdjęciu Landsat ETM) ? Jak
można je wyeliminować ?
Uwaga:
działanie można napisać raz:
(etm.img(layer_2) – tm.img(layer_2)) >20
można też eksperymentować z progami: 20, 25, 30 itd. – próg określa wartość, od której uznajemy że
na danym obszarze nastąpiła zmiana – w naszym przypadku wylesienie.
CZĘŚĆ 3: KLASYFIKACJA
W kolejnej części ćwiczenia przejdziesz do klasyfikacji zdjęć satelitarnych (analiza jakościowa)
wynikiem, której będzie mapa pokrycia terenu przedstawiająca cztery klasy pokrycia terenu:
obszarów zabudowanych, obszarów rolniczych, lasów i wód.
Wyświetl projekt teledetekcja3.mxd (katalog dane_cz3). Spróbuj zlokalizować na zdjęciu wymienione
powyżej cztery klasy pokrycia terenu.
11
Wykorzystując dodatkowa informacje w postaci mapy pokrycia terenu utworzonej z danych
pochodzących z bazy danych CORINE Land Cover oceń wizualnie na ile trafna była Twoja analiza
zdjęcia. (Przydatne narzędzie Swipe layer w pasku narzędzi Effects)
Mając do dyspozycji tylko jeden zakres zdjęcia, można dla potrzeb wyróżnienia typów pokrycia
terenu wykorzystać reklasyfikację (ArcToolbox>Spatial Analyst Tools). Jednak jaskrawości różnych
typów pokrycia terenu i relacje między nimi w różnych zakresach zdjęcia nie są identyczne. Oznacza
to, iż w celu precyzyjnego rozróżnienia pokrycia terenu na zdjęciu dobrze jest wykorzystać informacje
z kilku zakresów. Jeśli w jednym zakresie lasy i wody mają podobne jaskrawości, to w innym
jaskrawości mogą być różne. Ta dodatkowa informacja pozwala na dokładniejsze, niż w przypadku
jednego zakresu, wyróżnienie typów pokrycia terenu.
Etapy klasyfikacji
Wyróżnienie typów pokrycia terenu na zdjęciu satelitarnym, wykorzystujące informacje z wielu
zakresów spektralnych, jest możliwe poprzez jego klasyfikację. Wyróżnia się dwa rodzaje klasyfikacji:
wzorcową (supervised) i bezwzorcową (unsupervised). Ta druga w programie ArcGIS reprezentowana
jest przez opcje Iso cluster nie będzie omawiana na zajęciach
Klasyfikacja wzorcowa (nadzorowana) pozwala na wyróżnienie klas pokrycia terenu na zdjęciu
satelitarnym na podstawie przygotowanych wzorców. W jej trakcie określane jest podobieństwo
pikseli klasyfikowanego zdjęcia do posiadanych wzorców i przydzielenie ich na tej podstawie do
konkretnej klasy. Na proces klasyfikacji wzorcowej składają się następujące etapy (Jezioro, Kozak
2004):
wyznaczenie pól treningowych,
utworzenie sygnatur spektralnych dla pól treningowych,
ocena jakości sygnatur,
wybór i przeprowadzenie procedury klasyfikacyjnej,
weryfikacja wyniku klasyfikacji.
Wzorce (sygnatury) spektralne opisują sposób, w jaki obiekty danej klasy pokrycia terenu odbijają
promieniowanie w poszczególnych kanałach zdjęcia. Każdy wzorzec opisany jest przez statystykę
odpowiedzi spektralnych (np. średnie wartości odbicia, odchylenie standardowe itp.). Wzorce
wyliczane są przez program na podstawie pól treningowych, wyznaczanych np. poprzez digitalizację
(lub z istniejących bardzo dokładnych materiałów kartograficznych) na ekranie lub mierzonych w
terenie, w miejscach o znanym pokryciu terenu. Podczas procesu klasyfikacji każdy piksel zdjęcia
porównywany jest z wzorcami a następnie przypisywany według zdefiniowanych w metodzie
algorytmu (reguł decyzyjnych) klasyfikacji do odpowiedniej np. tak jak w naszym przypadku klasy
pokrycia terenu. Przygotowanie reprezentatywnych dla danego obszaru, wzorców spektralnych nie
jest proste i bardzo czasochłonne.
! Dobór odpowiednich pól treningowych determinuje dokładność klasyfikacji.
Podstawowa zasada, o której należy pamiętać tworząc (digitalizując) pola treningowe to taka, że nie
powinny one być zbyt małe (np. nie mniejsze niż 10 pikseli) ale jednocześnie muszą być możliwie
jednorodne spektralnie (homogeniczne pod względem wartości pikseli). W przypadku tych klas, gdzie
poszczególne obiekty znacząco różnią się wartościami odbicia (np. tereny uprawne) konieczne może
być wyznaczenie więcej niż jednego pola treningowego. Liczba wzorców nie może być mniejsza od
liczby klas pokrycia terenu.
12
Wyznaczanie pól treningowych
W programie ArcGIS pola treningowe można utworzyć na trzy sposoby, tworząc: plik ASCII, raster lub
wektor. W tym ćwiczeniu utworzymy pola treningowe tworząc nową warstwę wektorową. Zazwyczaj
źródłem informacji o polach treningowych są materiały inne niż zdjęcie satelitarne, które ma zostać
poddane klasyfikacji, my jednak w tym ćwiczeniu dla uproszczenia procedury wykorzystamy jako
podkład do digitalizacji zdjęcie satelitarne oraz dodatkowo dane CORINE Land Cover.
Przejdź do ArcCatalog i utwórz nowy plik shape: pola. Tworząc ten plik wybierz opcję ‘polygon’ oraz
importuj system odniesień przestrzennych z obrazu landrat.img: UTM 34N
Wróć do ArcMap i dodaj nowo utworzoną warstwę, do wektoryzacji będziesz używał narzędzi
dostępnych w pasku narzędzi Editor (więcej informacji o wektoryzacji znajdziesz w konspekcie: do
ćwiczenia 5).
Pola treningowe (ilość pól):
- obszary zabudowane (3),
- obszary rolnicze (3),
- lasy (3)
- wody (2)
(pamiętaj o dodaniu nowej kolumny: id_klas, w której wprowadzisz odpowiednie id w zależności od
klasy pokrycia terenu, której odpowiada dane pole, odpowiednio: obszary zabudowane (1), obszary
rolnicze (2), lasy (3), wody (4))
Sygnatury spektralne
Dla pól treningowych utwórz sygnatury wykorzystując narzędzie Create signatures (w
ArcToolsbox/Spatial analyst/Multivariate). W polu Input raster bands [kanały wejściowe rastra]
wprowadź obraz landsat, w polu Input raster or feature sample data [wejściowe rastrowe lub
wektorowe pola treningowe/próba] wektor z utworzonymi polami treningowymi: pola, w polu
sample field [pole pola treningowego/próby]
!
‘id_klas’, w polu Output signature file [wyjściowy plik
sygnatur – zapisywany w formacie *.gsg, który można odczytac np. w WordPad] wpisz nazwę
tworzonego pliku sygnatur; sygnatury, zaznacz opcję Compute covariance matrices [oblicz macierz
kowariancji
2
].
Po wygenerowaniu sygnatur warto sprawdzić jakie wartości jaskrawości są do nich przypisane oraz
jak te wartości jaskrawości dla sygnatur, dla poszczególnych klas pokrycia terenu rozkładają się
względem siebie. Można to zrobić na kilka sposobów
W celu porównania sygnatur otwórz, np. w programie Wordpad, plik sygnatury.gsg i odczytaj średnie
wartości jaskrawości dla każdego z pól treningowych w kanałach 2 oraz 3, zapisz je w arkuszu
programu Excel. Następnie, podobnie jak w pierwszej części ćwiczenia, utwórz wykres
przedstawiający zapisane wartości w obydwu analizowanych kanałach (przykład na rycinie poniżej).
2
więcej informacji o macierzy kowariancji znajdziesz m.in. na stronie:
http://pl.wikipedia.org/wiki/Macierz_kowariancji
13
Ocenę doboru sygnatur umożliwia również narzędzie Class Probability [Prawdopodobieństwo klasy]
(w ArcToolsbox/Spatial analyst/Multivariate). Pole Input raster bands [wejściowe kanały] = landsat,
Input signature file [wejściowy plik sygnatur] = sygnatury, Output multiband raster [wyjściowy
wielowarstwowy raster] = praw, pozostałe opcje pozostaw domyślne.
Dla każdej klasy tworzona jest, dzięki narzędziu Class Probability, warstwa przedstawiająca
prawdopodobieństwo (wartości od 0 do 100) przynależności każdego z pikseli do danej klasy. W ten
sposób można znaleźć obszary o niskim prawdopodobieństwie przynależności do danej klasy, które
mogą być np. błędnie sklasyfikowane lub są obszarami o mieszanym typie pokrycia terenu. Na bazie
map prawdopodobieństwa można np. stwierdzi czy potrzebne jest dodawanie nowych pól
treningowych lub korekta już istniejących pól.
Jeszcze jednym sposobem oceny sygnatur w ArcGIS jest dendrogram, który nie będzie omawiany na
tych zajęciach.
Algorytm klasyfikacyjny
Klasyfikacja w programie AcrGIS realizowana jest przez narzędzie, wykorzystujące algorytmy
największego podobieństwa/prawdopodobieństwa (maximum likelihood). Zaletą tej metody jest
uwzględnianie zmienności klas poprzez użycie macierzy kowariancji natomiast wadą uzależnienie od
rozkładu normalnego danych. Reguła ta oparta jest na wielkości prawdopodobieństwa, że piksel
należy do określonej klasy. Prawdopodobieństwo to wyliczane jest na bazie średnich wartości
sygnatur oraz macierzy kowariancji (czyli zmienności wartości w danej klasie). Metoda ta w ArcGIS
jest rozszerzona o regułę Bayesa (
http://pl.wikipedia.org/wiki/Twierdzenie_Bayesa
), która
14
wykorzystuje prawdopodobieństwa dla danych klas znane a priori. Piksel przypisywany jest do tej
klasy, gdzie wyliczone prawdopodobieństwo jest największe.
Przeprowadź klasyfikację wykorzystując narzędzie Maximum Likelihood Classification (w
ArcToolsbox/Spatial analyst/Multivariate). Pole Input raster bands [wejściowe kanały]= landsat,
Input signature file [wejściowy plik sygnatur] = sygnatury, Output classified raster [wyjściowy
sklasyfikowany raster] = klasyfikacja, pozostałe opcje pozostaw domyślne.
Po wyświetleniu mapy będącej wynikiem klasyfikacji zwróć uwagę, że widnieją na niej cztery klasy. Są
one oznaczone w legendzie nazwami poszczególnych sygnatur. Zastanów się, jakie błędy są widoczne
na mapie wynikowej? Porównaj wygenerowana mapę z mapą Corine. Sprawdź jak zostały
sklasyfikowane obszary, które miały niskie prawdopodobieństwo przynależności do poszczególnych
klas na warstwach wygenerowanych z wykorzystaniem narzędzia Class Probability.
Weryfikacja
Ostatnim etapem klasyfikacji jest jej weryfikacja. Klasyfikacja nigdy nie jest dokładna w stu
procentach. Istotne jest jednak, aby można było ocenić rzetelność klasyfikacji i oszacować, z jakim
błędem należy się liczyć. Ocena błędu (weryfikacja) klasyfikacji polega na porównaniu wyniku
klasyfikacji z danymi o pokryciu terenu zebranymi podczas niezależnych badań. Weryfikację
najczęściej przeprowadza się dla niewielkich powierzchni rozmieszczonych losowo na
klasyfikowanym obszarze. Dane o pokryciu terenu użyte do weryfikacji wyniku klasyfikacji powinny
być dokładniejsze od tych, których dostarcza zdjęcie satelitarne. Zwykle zbiera się je w czasie badań
terenowych, ale można w tym celu wykorzystać także zdjęcia lotnicze albo dokładne mapy.
Na zajęciach dokonasz tylko wstępnej wizualnej weryfikacji dla klasyfikacji w okolicach Kampusu UJ
wykorzystując jako dane referencyjne mapę topograficzną i ortofotomapę. Oceń na ile Twoja
klasyfikacja zgadza się z informacją o pokryciu terenu na danych referencyjnych. Porównaj również
otrzymany wynik z mapami prawdopodobieństwa – sprawdź jak zostały sklasyfikowane obszary o
najniższym prawdopodobieństwie.
Podsumowanie
Wykonując powyższe ćwiczenie poznałeś podstawowe właściwości zdjęć satelitarnych – cztery typy
rozdzielczość: przestrzenną, radiometryczna, spektralną i czasową – wizualną interpretację zdjęć
satelitarnych, podstawowe możliwości analizy zmian w czasie oraz klasyfikacji zdjęcia satelitarnego.
Literatura
Adamczyk J., Będkowski K., 2007, Metody cyfrowe w teledetekcji, Wydawnictwo SGGW
Jezioro P., Kozak J., 2004, Wprowadzenie do Systemów Informacji Geograficznej (Ćwiczenia), rozdz. 8.
Klasyfikacja zdjęcia satelitarnego, IGiGP, s. 57-66