Python c.d.n…

Witajcie!

Zapewne już przerobiliście kurs proponowany we wcześniejszym poście Pythona czas zacząć. Nie????? W sumie nie dziwię się po pracy są ciekawsze rzeczy do robienia niż siedzenie przed komputerem.

Jeśli nie to może być Wam ciężko ale zrobimy małą powtórkę w trakcie i damy radę 🙂

W tym poście wykonamy CLIP w ArcGIS przy pomocy Pythona, ale na początek trochę teorii i wprowadzenia.

Każdy choć trochę pracujący w programach GIS’owych używał narzędzi typu wytnij, przetnij, buforuj, połącz itp. Są to podstawowe narzędzia dostępne m.in. w ArcMap z menu podręcznego o nazwie „Geoprocessing” zlokalizowanego u góry okna programu.

Narzędzia te umożliwiają operację na plikach shp lub klasach obiektów znajdujących się w geobazie. Klasa obiektu to pewnego rodzaju kontenery przechowujące informacje. Tymi informacjami mogą być obiekty geometryczne (np:punkty, linie, poligony) oraz powiązane z nimi atrybuty.

Wyobraźmy sobie taką sytuację – mamy 3  pliki shp przedstawiające (nazwy plików oczywiście piszemy bez polskich znaków):

• gminy (zawiera tylko 63 gminy z obszary całego kraju),
• obiekty fizjograficzne (cała Polska)
• miejscowości (cała Polska)

W plikach tych zgromadzony jest wyżej podany zakres danych dla obszaru całej Polski.

Dostajemy zadanie aby dla wybranych 63 gmin przypisać miejscowości i obiekty fizjograficzne, w taki sposób aby w pliku miejscowości dla każdej pojedynczej miejscowości znalazła się nazwa gminy oraz dla każdego obiektu fizjograficznego znalazła się nazwa gminy.

Reasumując mamy taką sytuację w naszych danych źródłowych:

 

gminy.shp OBJECTID|Shape|jpt_nazwa|

obiekty_fizjograficzne.shp FID|Shape|nazwa_glowna|rodzaj_obi

miejscowosci.shp FID|Shape|naz_glowna|rodzaj_obi

Celowo ćwiczenie (Clipowanie) wykonujemy dla 2 plików dla lepszego zrozumienia, wyobrazić sobie możemy jednak sytuację, że takich plików mamy np:100 i co wtedy,  będziemy „wyklikiwać ” kolejno kroki w ArcMap. No jasne, że nie z pomocą rusza nam Python a ściślej mówiąc ArcPy (pakiet bibliotek w ArcGIS).

KROK 1.

ZADANIE: Wyciąć poligonem gminy.shp pliki obiekty_fizjograficzne.shp oraz miejscowosci.shp.

Jakie to proste Menu głowne -> Geoprocessing -> Clip i gotowe po co się męczyć. No i pewnie racja, nie ma co się porywać z motyką na słońce w przypadku wykonania klipa 2 razy. Można od biedy skorzystać z batch processingu. I tutaj odsyłam do postu „Batch processing w ArcGIS”.

Ale my zrobimy tą operację w Pythonie tak na wszelki wypadek jeśli kiedyś przyjdzie nam „zklipować” 100 albo więcej plików 🙂

a) Uruchamiamy PyScriptera.

I piszemy…

Na początku musimy zaimportować pakiet ArcPy oraz klasę odpowiadającą za zmienne środowiskowe.

A co robi ten ArcPy, a co to jest ta zmienna środowiskowa? Bardzo słuszne pytania. Pakiet ArcPy to podstawowy pakiet zawierający bardzo dużą liczbę funkcji pomocnych przy przetwarzaniu danych. Zmienne środowiskowe określają nam m.in. lokalizację danych, z których będziemy korzystać.

import arcpy
from arcpy import env

Tutaj określamy dokładną lokalizację przestrzeni roboczej, pamiętaj, że nie wystarczy skopiować linku do danych z exploratora, python czyta „/” a nie „\”, dlatego zamień backslasha na slasha.

# Zdefiniowanie przestrzeni roboczej
env.workspace = „E:/geoprofesja/Python/dane”

Dobrze to teraz musimy jeszcze zdefiniować dane wejściowe tzw. input. Nie musimy podawać pełnej ścieżki, ponieważ wcześniej określiliśmy przestrzeń roboczą dlatego interpreter będzie wiedział gdzie ma szukać nazw miejscowosci.shp oraz obiekty_fizjograficzne.shp

# Zdefiniowanie zmiannych
#Dane wejsciowe
dane_wejsciowe1 = „miejscowosci.shp”
dane_wejsciowe2 = „obiekty_fizjograficzne.shp”

Dobrze to piszemy dalej, teraz pora na zdefiniowanie warstwy, która posłuży do wycięcia. Zauważ, że interpreter poradzi sobie z tym, że nie podasz mu rozszerzenia *.shp. Znajdzie sobie warstwę o takiej nazwie, ponieważ korzysta ze zdefiniowanej wcześniej przestrzeni roboczej.

# Warstwa ograniczająca
warstwa_wycinajaca = „gminy”

Określamy teraz warstwy wynikowe.

#Dane wynikowe
wynik1 = „miejscowosciCLIP.shp”
wynik2 = „obiekty_fizjograficzneCLIP.shp”

I już ostatni krok teraz korzystamy z funkcji Clip_analysis. Zwróć uwagę, że po otwarciu nawiasu program podpowiada nam składnię, i tak najpierw piszemy INPUT poźniej CLIP_FEATURE, OUTPUT, TOLERANCE). Ostatnia opcja TOLERANCE umożliwia nam określenie zakresu tolerancji, ale domyślnie piszemy NONE nie określając tego parametru

# Wykonanie wycięcia
arcpy.Clip_analysis(dane_wejsciowe1, warstwa_wycinajaca, wynik1, None)
arcpy.Clip_analysis(dane_wejsciowe2, warstwa_wycinajaca, wynik2, None)

Zapisz skrypt i wczytaj go do ArcMap klikając w oknie Python prawym przyciskiem myszy i wybierając opcję „LOAD”

Po zaimportowaniu skryptu kliknij ENTER uruchamiając skrypt. Jeśli nie zrobiłeś błędów skrypt uruchomi się i wykona 2 clipy.

Ale zaraz, ktoś może powiedzieć, to jak będę mieć 100 klipów do wykonania to będę musiał pisać 100 razy arcpy.CLIP.. itd. to też się zejdzie dużo czasu. I masz rację, skoro tak uważasz to znaczy, że APETYT ROŚNIE, ROŚNIE, ROŚNIE…

W kolejnym poście zoptymalizujemy kod.