diff --git a/notebooks/04_ortofotomapa.Rmd b/notebooks/04_ortofotomapa.Rmd index c286dcb..f946455 100644 --- a/notebooks/04_ortofotomapa.Rmd +++ b/notebooks/04_ortofotomapa.Rmd @@ -54,44 +54,41 @@ bliskiej podczerwieni (NIR) ## Wyszukiwanie -Centroid (środek geometryczny) Lublina znajduje się na 22,56° długości -geograficznej (X) i 51,22° szerokości geograficznej (Y). W celu wygenerowania tego -punktu, możemy uprzednio stworzyć macierz, gdzie w wierszach znajdą się punkty -(w naszym przypadku tylko jeden), a w kolumnach kolejno współrzędne X i Y. Następnie -należy dokonać konwersji do obiektu wektorowego używając funkcji `vect()` -i definiując układ współrzędnych. +W katalogu `dane` znajduje się warstwa poligonowa `Lublin.gpkg`, która reprezentuje +granicę Lublina w układzie `EPSG:2180`. ```{r} -wspolrzedne = matrix(c(22.56, 51.22), ncol = 2) -centroid = vect(wspolrzedne, type = "points", crs = "EPSG:4326") -centroid +lublin = vect("../dane/Lublin.gpkg") +lublin ``` -Dokonajmy również konwersji ze współrzędnych geograficznych na układ metryczny -`EPSG:2180`. +Funkcja `perim()` pozwala obliczyć długość granicy, natomiast funkcja `expanse()` +wskazuje powierzchnię obszaru. ```{r} -centroid = project(centroid, "EPSG:2180") +perim(lublin) / 1000 # wynik w km +expanse(lublin, unit = "km") ``` -W kolejnym kroku stwórzmy bufor, który częściowo pokryje obszar miasta. +Kolejny etap dotyczy wyszukania dostępnych danych dla naszego obszaru zainteresowania. +Dla uproszczenie wyznaczmy centroid. Jego współrzędne można odczytać używając +funkcji `crds()`. ```{r} -bufor = buffer(centroid, width = 1000) +punkt = centroids(lublin) +crds(punkt) ``` -Stworzyliśmy bufor o szerokości 1 km, a teraz przygotujmy prostą wizualizację. - ```{r} -plot(bufor, main = "Bufor Lublin") -plot(centroid, add = TRUE) + +plot(lublin, main = "Lublin") +plot(punkt, col = "blue", add = TRUE) ``` -Kolejny etap, po określeniu obszaru zainteresowania, związany jest z wyszukaniem -dostępnych danych. W tym celu wykorzystamy funkcję `ortho_request()`. +Teraz wyszukajmy dostępne dane za pomocą funkcji `ortho_request()`. ```{r} -dane = ortho_request(bufor) +dane = ortho_request(punkt) ``` Możemy wyświetlić część otrzymanej ramki danych lub alternatywnie przeglądać @@ -120,10 +117,10 @@ dane[order(-dane$year), ] Jako przykład pobierzmy dwie kompozycje tego samego obszaru wykonane w naturalnych barwach i z kanałem bliskiej podczerwieni z 2022 r. -(ID: `76746_1143331_M-34-34-A-c-2-3` i `76745_1143332_M-34-34-A-c-2-3`). +(ID: `76746_1143325_M-34-34-A-c-1-4` i `76745_1143326_M-34-34-A-c-1-4`). ```{r} -id = c("76746_1143331_M-34-34-A-c-2-3", "76745_1143332_M-34-34-A-c-2-3") +id = c("76746_1143325_M-34-34-A-c-1-4", "76745_1143326_M-34-34-A-c-1-4") dane_sel = dane[dane$filename %in% id, ] ``` @@ -174,4 +171,5 @@ b) do jednego wirtualnego pliku *.vrt* używając funkcji `vrt()` Sprawdź zajmowaną ilość miejsca przez te pliki na dysku wykorzystując funkcję `file.size()` (wynik zwracany jest w bajtach). Sprawdź również zawartość pliku *.vrt* (czym on jest w rzeczywistości?). Następnie, zmniejsz rozdzielczość -mozaiki do 10 m i zapisz wynik. +mozaiki do 10 m i zapisz wynik. Jak zmieniła się jakość w stosunku do +oryginalnego zdjęcia? diff --git a/notebooks/04_ortofotomapa.html b/notebooks/04_ortofotomapa.html index 667d336..5bbdfcd 100644 --- a/notebooks/04_ortofotomapa.html +++ b/notebooks/04_ortofotomapa.html @@ -1531,51 +1531,52 @@

Ortofotomapa

Wyszukiwanie

-

Centroid (środek geometryczny) Lublina znajduje się na 22,56° -długości geograficznej (X) i 51,22° szerokości geograficznej (Y). W celu -wygenerowania tego punktu, możemy uprzednio stworzyć macierz, gdzie w -wierszach znajdą się punkty (w naszym przypadku tylko jeden), a w -kolumnach kolejno współrzędne X i Y. Następnie należy dokonać konwersji -do obiektu wektorowego używając funkcji vect() i definiując -układ współrzędnych.

-
wspolrzedne = matrix(c(22.56, 51.22), ncol = 2)
-centroid = vect(wspolrzedne, type = "points", crs = "EPSG:4326")
-centroid
+

W katalogu dane znajduje się warstwa poligonowa +Lublin.gpkg, która reprezentuje granicę Lublina w układzie +EPSG:2180.

+
lublin = vect("../dane/Lublin.gpkg")
+lublin
##  class       : SpatVector 
-##  geometry    : points 
+##  geometry    : polygons 
 ##  dimensions  : 1, 0  (geometries, attributes)
-##  extent      : 22.56, 22.56, 51.22, 51.22  (xmin, xmax, ymin, ymax)
-##  coord. ref. : lon/lat WGS 84 (EPSG:4326)
-

Dokonajmy również konwersji ze współrzędnych geograficznych na układ -metryczny EPSG:2180.

-
centroid = project(centroid, "EPSG:2180")
-

W kolejnym kroku stwórzmy bufor, który częściowo pokryje obszar -miasta.

-
bufor = buffer(centroid, width = 1000)
-

Stworzyliśmy bufor o szerokości 1 km, a teraz przygotujmy prostą -wizualizację.

-
plot(bufor, main = "Bufor Lublin")
-plot(centroid, add = TRUE)
-

-

Kolejny etap, po określeniu obszaru zainteresowania, związany jest z -wyszukaniem dostępnych danych. W tym celu wykorzystamy funkcję +## extent : 741124, 756250.7, 369483.5, 387356.7 (xmin, xmax, ymin, ymax) +## source : Lublin.gpkg +## coord. ref. : ETRF2000-PL / CS92 (EPSG:2180) +

Funkcja perim() pozwala obliczyć długość granicy, +natomiast funkcja expanse() wskazuje powierzchnię +obszaru.

+
perim(lublin) / 1000 # wynik w km
+
## [1] 81.4903
+
expanse(lublin, unit = "km")
+
## [1] 147.4631
+

Kolejny etap dotyczy wyszukania dostępnych danych dla naszego obszaru +zainteresowania. Dla uproszczenie wyznaczmy centroid. Jego współrzędne +można odczytać używając funkcji crds().

+
punkt = centroids(lublin)
+crds(punkt)
+
##             x      y
+## [1,] 748243.3 379235
+
plot(lublin, main = "Lublin")
+plot(punkt, col = "blue", add = TRUE)
+

+

Teraz wyszukajmy dostępne dane za pomocą funkcji ortho_request().

-
dane = ortho_request(bufor)
+
dane = ortho_request(punkt)

Możemy wyświetlić część otrzymanej ramki danych lub alternatywnie przeglądać całość używając funkcji View().

# wyświetl 10 pierwszych wierszy i 6 pierwszych kolumn
 dane[1:10, 1:6]
##            sheetID year resolution composition  sensor        CRS
-## 1  M-34-34-A-c-2-3 2022       0.25         CIR Digital    PL-1992
-## 2  M-34-34-A-c-2-3 2020       0.25         RGB Digital    PL-1992
-## 3  M-34-34-A-c-2-3 2019       0.10         RGB Digital    PL-1992
-## 4  M-34-34-A-c-2-3 2018       0.10         RGB Digital    PL-1992
-## 5      8.151.08.08 2015       0.25         RGB Digital PL-2000:S8
-## 6      8.151.08.13 2015       0.25         RGB Digital PL-2000:S8
-## 7      8.151.08.18 2015       0.25         RGB Digital PL-2000:S8
-## 8  M-34-34-A-c-2-3 2015       0.25         CIR Digital    PL-1992
-## 9  M-34-34-A-c-2-3 2015       0.25         RGB Digital    PL-1992
-## 10 M-34-34-A-c-2-3 2013       0.10         RGB Digital    PL-1992
+## 1 M-34-34-A-c-1 2003 0.25 B/W Analog PL-1992 +## 2 M-34-34-A-c-1-4 2022 0.25 CIR Digital PL-1992 +## 3 M-34-34-A-c-1-4 2020 0.25 RGB Digital PL-1992 +## 4 M-34-34-A-c-1-4 2019 0.10 RGB Digital PL-1992 +## 5 M-34-34-A-c-1-4 2018 0.10 RGB Digital PL-1992 +## 6 8.151.08.12 2015 0.25 RGB Digital PL-2000:S8 +## 7 M-34-34-A-c-1-4 2015 0.25 RGB Digital PL-1992 +## 8 M-34-34-A-c-1-4 2015 0.25 CIR Digital PL-1992 +## 9 M-34-34-A-c-1-4 2013 0.10 RGB Digital PL-1992 +## 10 8.151.08.12 2013 0.10 RGB Digital PL-2000:S8

Standardowo dane możemy filtrować z uwzględnieniem zadanych parametrów.

dane[dane$year > 2016, ]
@@ -1588,9 +1589,9 @@ 
Wyszukiwanie

 
Pobieranie

 
Jako przykład pobierzmy dwie kompozycje tego samego obszaru wykonane
 w naturalnych barwach i z kanałem bliskiej podczerwieni z 2022 r. (ID:
-76746_1143331_M-34-34-A-c-2-3 i
-76745_1143332_M-34-34-A-c-2-3).

-
id = c("76746_1143331_M-34-34-A-c-2-3", "76745_1143332_M-34-34-A-c-2-3")
+76746_1143325_M-34-34-A-c-1-4 i
+76745_1143326_M-34-34-A-c-1-4).

+
id = c("76746_1143325_M-34-34-A-c-1-4", "76745_1143326_M-34-34-A-c-1-4")
 dane_sel = dane[dane$filename %in% id, ]

 
Po selekcji potrzebnych danych, można je pobrać wykorzystując funkcję
 tile_download(). Możliwe jest również wskazanie katalogu,
@@ -1607,17 +1608,17 @@ 
Pobieranie

 do plików (full.names = TRUE).

 
pliki = list.files("dane", pattern = "\\.tif$", full.names = TRUE)
 pliki

-
## [1] "dane/76745_1143332_M-34-34-A-c-2-3.tif"
-## [2] "dane/76746_1143331_M-34-34-A-c-2-3.tif"

+
## [1] "dane/76745_1143326_M-34-34-A-c-1-4.tif"
+## [2] "dane/76746_1143325_M-34-34-A-c-1-4.tif"

 
W ostatnim kroku możemy kolejno wczytać rastry i je wyświetlić.

 
# kompozycja z bliską podczerwienią
 r1 = rast(pliki[1])
 plot(r1)

-

+

 
# kompozycja w naturalnych barwach
 r2 = rast(pliki[2])
 plot(r2)

-

+
@@ -1634,7 +1635,8 @@

Zadanie

funkcję file.size() (wynik zwracany jest w bajtach). Sprawdź również zawartość pliku .vrt (czym on jest w rzeczywistości?). Następnie, zmniejsz rozdzielczość mozaiki do 10 m i -zapisz wynik.

+zapisz wynik. Jak zmieniła się jakość w stosunku do oryginalnego +zdjęcia?