2. Einführung in Leaflet
Vorwort zu working directory, Daten und Packages
Die Daten für dieses Tutorial werden über Github zur Verfügung gestellt. Wenn Sie den Code mit eigenen Daten reproduzieren wollten, ersetzen Sie die URL mit ihrem lokalen Dateipfad.
Definieren Sie in diesem Fall Ihren “working directory”, ihr Arbeitsverzeichnis.
setwd("YOUR/FILEPATH")
# this will not work. if you want to use it, define a path on your computer
Dabei handelt es sich um den Dateispeicherort, wo Ihre Daten und Skripte gespeichert sind. Weitere Infos finden Sie hier.
Der Code für das Speichern von Daten und Ergebnissen ist in diesen Tutorials deaktiviert. Ersetzen Sie zum Speichern “YOUR/FILEPATH” mit ihrem gewünschten Dateipfad.
Wir haben den Code zum Laden von packages so geschrieben, dass diese installiert werden, sofern diees noch nicht auf Ihrem PC abgespeichert sind. Anschließend werden diese direkt geladen. Sofern nicht anderweitig spezifiziert werden die packages an einer standardmäßigen Stelle auf Ihrem Computer abgespeichert.
Was ist Leaflet
Leaflet ist eine beliebte JavaScript Bibliothek zur Erstellung von dynamischen, interaktiven webbasierten Karten (s. https://leafletjs.com/).
Typisches Ablaufschema
- Erstelle eine sog. “map widget” mit der Funktion
leaflet(). - Füge dem “map widget” neue Ebenenen aus Geodaten hinzu indem die Funktionen
addTiles,addMarkers,addPolygonsverwendet werden. - Schritt 2 kann beliebig häufig wiederholt werden.
- Führe das fertige “map widget” aus um die erstellte Karte zu visualisieren.
Voraussetzungen
Für die nachfolgen Schritte müssen zunächst folgende R-Pakete installiert werden:
- leaflet - als führende Open-Source-JavaScript-Bibliothek für mobile, interaktive Karten.
- rgdal - für den Import von Shapefiles nach R als Spatial data frame.
- spdplyr - zur Manipulation der Attributdaten innerhalb des Spatial data frame.
- sf - stellt Methoden und Klassen für die Bearbeitung von Geodaten (Vektordaten) zur Verfügung.
- sp - Standardpaket für Vektordaten bevor das sf Paket existierte.
- geojsonio - für die Konvertierung des Spatial data frame in das GeoJSON-Format und die Speicherung im Dateisystem.
- rmapshaper - zur Manipulation von Geometrieinformationen (Polygon, Linie, Marker) einer GeoJSON-Datei.
- geojson - Allgemeine nützliche Funktionen für GeoJSON-Dateien.
- raster - stellt Methoden und Klassen für die Bearbeitung von Rasterdaten zur Verfügung.
if (!require(leaflet)){install.packages("leaflet"); library(leaflet)}
if (!require(rgdal)){install.packages("rgdal"); library(rgdal)}
if (!require(spdplyr)) {install.packages("spdplyr"); library(spdplyr)}
if (!require(sp)){install.packages("sp"); library(sp)}
if (!require(sf)) {install.packages("sf"); library(sf)}
if (!require(geojsonio)) {install.packages("geojsonio"); library(geojsonio)}
if (!require(rmapshaper)) {install.packages("rmapshaper"); library(rmapshaper)}
if (!require(geojson)) {install.packages("geojson"); library(geojson)}
if (!require(raster)) {install.packages("raster"); library(raster)}
Beispiel
m <- leaflet() %>% ## Add empty leaflet frame (map widget)
addTiles() %>% # Add default OpenStreetMap map tiles
addMarkers(lng=7.84505099, lat=47.9934184, popup="Werthmannstraße 4, CIP-Pool")
m # Print the map
Map widget anpassen (Map methods)
Für das map widget können spzifische Parameter vergeben werden, z.B. kann die Zoomstufe angegeben werden.
m <- leaflet(options = leafletOptions(minZoom = 0, maxZoom = 18)) %>%
addTiles() %>% # Add default OpenStreetMap map tiles
addMarkers(lng=7.84505099, lat=47.9934184, popup="Werthmannstraße 4, CIP-Pool")
m # Print the map
Auch die initiale Position und Zoomstufe der Hintergrundkarte kann verändert werden. Gewünschte Koordinaten können unter https://www.mapcoordinates.net/de gefunden werden.
m <- leaflet(options = leafletOptions(minZoom = 0, maxZoom = 18)) %>%
setView(lng=8.9, lat=49.0, zoom = 12) %>%
addTiles() %>% # Add default OpenStreetMap map tiles
addMarkers(lng=7.84505099, lat=47.9934184, popup="Werthmannstraße 4, CIP-Pool")
m # Print the map
Auch die Bounding Box, d.h. der Kartenausschnitt kann bestimmt werden:
m <- leaflet(options = leafletOptions(minZoom = 0, maxZoom = 18)) %>%
fitBounds(7.65,47.90,7.95,48.08) %>%
addTiles() %>% # Add default OpenStreetMap map tiles
addMarkers(lng=7.84505099, lat=47.9934184, popup="Werthmannstraße 4, CIP-Pool")
m # Print the map
Mit dem “Data Object” und Vektordaten arbeiten
Eine einfache Matrix mit lng/lat Werten laden
Es kann ein gewöhnliches Datenfenster (xls, csv, txt, etc.) mit lng/lat Spalten (Variablen) als Data Object genutzt werden. Merke: die Projektion sollte bekannt sein und in EPSG: 43261 vorliegen.
CSV-Daten können als Data Object beispielsweise folgendermaßen in R eingeladen werden:
# lade die Datei aus Github
urlfile <- "https://raw.githubusercontent.com/RafHo/teaching/master/geodaten_kommunizieren/datasource/lng_lat1.txt"
df <- read.csv(urlfile, sep = ",", header = TRUE, encoding = "UTF-8")
Nun können die Werte aus den lng/lat Spalten mittels Leaflet visualisiert werden; ggf. müssen sie manuell definiert werden:
#lng/lat Daten einladen
m <- leaflet(df) %>%
addTiles() %>% # Add default OpenStreetMap map tiles
addMarkers()
m # Print the map
Eine Shapefile laden (Point)
Um Shapefiles in R zu importieren, können Sie die Funktion st_read() aus dem Paket’sf’ verwenden.
#load data from Github
download.file("https://raw.githubusercontent.com/RafHo/teaching/master/geodaten_kommunizieren/datasource/lng_lat.zip", destfile = "lng_lat.zip")
# unzip data
unzip("lng_lat.zip")
#read data in with the sf package
shp_point <- st_read("lng_lat.shp")
## Reading layer `lng_lat' from data source
## `V:\github\nils_website\blogdown\content\courses\geoinformationen_kommunizieren\lng_lat.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 466 features and 5 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: 7.791445 ymin: 47.976 xmax: 7.862785 ymax: 48.02543
## Geodetic CRS: WGS 84
#shp_point <- readOGR("lng_lat.shp", verbose = FALSE)
#plot(shp_point)
Nach dem Import sehen Sie, dass ein Spatial Data Frame erstellt wurde. Wenn Ihre Shapefiles Polygone (Flächen) enthalten, ist es nun ein sog. SpatialPolygonsDataFrame erhalten, wenn es Linien waren, ist es nun ein sog. SpatialLinesDataFrame. In unserem Fall ist es ein sog. SpatialPointDataFrame
Der Attributteil der ersten 10 Datenreihen kann eingesehen werden, indem Sie folgendes ausführen:
# Show the first 10 rows
head(shp_point, 10)
## Simple feature collection with 10 features and 5 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: 7.818589 ymin: 47.976 xmax: 7.858387 ymax: 48.02543
## Geodetic CRS: WGS 84
## position latitude longitude timestamp
## 1 Freiburg, Tullastrasse 48.02543 7.847819 2017-01-13 08:00:00
## 2 Freiburg, Vauban 47.97601 7.824701 2017-01-13 08:00:00
## 3 Freiburg, Sundgauallee 30 48.00735 7.818589 2017-01-13 08:00:00
## 4 Freiburg, Vauban 47.97601 7.824701 2017-01-13 08:00:00
## 5 Freiburg, Barbarastraße 8 48.00373 7.839775 2017-01-13 08:00:00
## 6 Freiburg, Okenstraße 48.01042 7.858387 2017-01-13 08:00:00
## 7 Freiburg, Engelbergerstrasse 41 47.99721 7.833161 2017-01-13 08:00:00
## 8 Freiburg, Okenstraße 48.01042 7.858387 2017-01-15 08:00:00
## 9 Freiburg, Barbarastraße 8 48.00373 7.839775 2017-01-13 08:00:00
## 10 Freiburg, Barbarastraße 8 48.00373 7.839775 2017-01-13 08:00:00
## node_label geometry
## 1 shortterm precipitation POINT (7.847819 48.02543)
## 2 cloud cover POINT (7.824701 47.97601)
## 3 temperature level POINT (7.818589 48.00735)
## 4 shortterm precipitation POINT (7.824701 47.97601)
## 5 shortterm precipitation POINT (7.839775 48.00373)
## 6 cloud cover POINT (7.858387 48.01042)
## 7 cloud cover POINT (7.833161 47.99721)
## 8 cloud cover POINT (7.858387 48.01042)
## 9 temperature level POINT (7.839775 48.00373)
## 10 shortterm precipitation POINT (7.839775 48.00373)
Nun können die Werte aus den lng/lat Spalten mittels Leaflet visualisiert werden. Standardmäßig sucht die Funktion leaflet() nach korrekt beschrifteten lng/lat Spalten; ggf. müssen sie manuell definiert werden:
#lng/lat Daten einladen
m <- leaflet(shp_point) %>%
addTiles() %>% # Add default OpenStreetMap map tiles
addMarkers()
m # Print the map
Eine Shapefile laden (Polygon) und das CRS transformieren
Shapefiles können als Data Object beispielsweise folgendermaßen in R eingeladen und gefiltert werden:
# load from Github
download.file("https://raw.githubusercontent.com/RafHo/teaching/master/geodaten_kommunizieren/datasource/Administrative_Boundaries_2016.zip", destfile = "Administrative_Boundaries_2016.zip")
# unzip shapefile
unzip("Administrative_Boundaries_2016.zip")
shp_poly <- st_read("Administrative_Boundaries_2016.shp")
## Reading layer `Administrative_Boundaries_2016' from data source
## `V:\github\nils_website\blogdown\content\courses\geoinformationen_kommunizieren\Administrative_Boundaries_2016.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 1857 features and 20 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: 336371.5 ymin: 5214597 xmax: 492299.8 ymax: 5466220
## Projected CRS: ETRS89 / UTM zone 32N
#shp_poly <- readOGR("Administrative_Boundaries_2016.shp", encoding = "UTF-8", use_iconv = FALSE)
# vorsicht, Encoding
shp_poly <- subset(shp_poly, shp_poly$name %in% c(
"Freiburg im Breisgau","Emmendingen", "Vrstetten", "Denzlingen", "Reute",
"Gundelfingen", "Glottertal", "Heuweiler"))
Leaflet kann nur mit nicht projezierten latitude/longitude Koordinaten umgehen. Entpsrechend ist das CRS ggf. anzupassen. Das Ziel CRS sollte immmer EPSG:4326 sein. Der Proj4-String hierfür lautet: “+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs”.
Die Transformation kann entweder mit einer GIS-Software (QGIS, ArcGIS, …) durchgeführt werden oder direkt mit R und der Funktion spTransform sp aus dem sp package.
shp_poly <- st_transform(shp_poly, CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))
# or with `spTransform` from the sp package
#shp_poly <- spTransform(shp_poly, CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"))
Mit der korrrekten Projektion können die Polygone nun visualisiert werden:
leaflet(shp_poly) %>%
addTiles() %>%
addPolygons()
Die Visualisierung kann um weitere Optionen ergänzt werden. Dies ist Gegenstand der kommenden Sitzungen.
leaflet(shp_poly) %>%
addTiles() %>%
addPolygons( color = "#444444", weight = 1, smoothFactor = 0.5,
opacity = 1.0, fillOpacity = 0.5,
fillColor = ~colorQuantile("YlOrRd", pop13)(pop13),
highlightOptions = highlightOptions(color = "white", weight = 2,
bringToFront = TRUE))
Die Attribute eines Shapefile verändern
Um zu zeigen, wie einfach die Attributstruktur einer Shapefile verändert werden kann, beziehen wir einen neuen Beispieldatensatz. Über die Website https://gadm.org kann auf eine offene Geodatenbank zugregriffen werden, die Verwaltungsgrenzen der Welt (z.B. Länder, Provinzen, etc.) beinhalten. In verschiedenen Formaten können dei Daten für die nicht-kommerzielle Nutzung verwendet werden. Im folgenden Beispiel werden die Verwaltungsgrenzen von Deutschland verwendet: https://biogeo.ucdavis.edu/data/gadm3.6/shp/gadm36_DEU_shp.zip
# load from Github
download.file("https://raw.githubusercontent.com/RafHo/teaching/master/geodaten_kommunizieren/datasource/gadm_36_BW_2.zip", destfile = "gadm_36_BW_2.zip")
# unzip shapefile
unzip("gadm_36_BW_2.zip")
bw <- st_read("gadm_36_BW_2.shp")
## Reading layer `gadm_36_BW_2' from data source
## `V:\github\nils_website\blogdown\content\courses\geoinformationen_kommunizieren\gadm_36_BW_2.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 45 features and 13 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 7.51159 ymin: 47.5024 xmax: 10.49575 ymax: 49.79148
## Geodetic CRS: WGS 84
#bw <- readOGR("gadm_36_BW_2.shp", verbose=FALSE)
Anhand des R-Paket spdplyr kann ein spatial data frame in das einfache Tabellenformat tibble überführt werden, sodass alle Funktionen des R-Pakets dplyr darauf angwendet werden können (s. https://dplyr.tidyverse.org/. Der Vorteil dabei ist, dass weiterhin alle dazugehörigen Geometrien im Datenfenster erthalten bleiben.
Es kann sinnvoll sein die zugörigen Datensätze entsprechend einer bestimmten Raumeinheit zu filtern, sodass z.B. nur die administrativen Raumeinheiten im neuen Datensatz enthalten sind:
bw_fr <- filter(bw, NAME_2 == "Freiburg im Breisgau")
bw_fr
## Simple feature collection with 1 feature and 13 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 7.662455 ymin: 47.90339 xmax: 7.930314 ymax: 48.07111
## Geodetic CRS: WGS 84
## GID_0 NAME_0 GID_1 NAME_1 NL_NAME_1 GID_2
## 1 DEU Germany DEU.1_1 Baden-Württemberg <NA> DEU.1.12_1
## NAME_2 VARNAME_2 NL_NAME_2 TYPE_2 ENGTYPE_2 CC_2 HASC_2
## 1 Freiburg im Breisgau <NA> <NA> Stadtkreis District 08311 DE.BW.FB
## geometry
## 1 MULTIPOLYGON (((7.677702 47...
Wenn nur Attributnamen in den Daten umbenannt werden sollen, kann der Befehl rename wie unten beschrieben verwendet werden:
bw_bl <- rename(bw, Bundesland = NAME_1)
bw_bl
## Simple feature collection with 45 features and 13 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 7.51159 ymin: 47.5024 xmax: 10.49575 ymax: 49.79148
## Geodetic CRS: WGS 84
## First 10 features:
## GID_0 NAME_0 GID_1 Bundesland NL_NAME_1 GID_2
## 1 DEU Germany DEU.1_1 Baden-Württemberg <NA> DEU.1.1_1
## 2 DEU Germany DEU.1_1 Baden-Württemberg <NA> DEU.1.2_1
## 3 DEU Germany DEU.1_1 Baden-Württemberg <NA> DEU.1.3_1
## 4 DEU Germany DEU.1_1 Baden-Württemberg <NA> DEU.1.4_1
## 5 DEU Germany DEU.1_1 Baden-Württemberg <NA> DEU.1.5_1
## 6 DEU Germany DEU.1_1 Baden-Württemberg <NA> DEU.1.6_1
## 7 DEU Germany DEU.1_1 Baden-Württemberg <NA> DEU.1.7_1
## 8 DEU Germany DEU.1_1 Baden-Württemberg <NA> DEU.1.8_1
## 9 DEU Germany DEU.1_1 Baden-Württemberg <NA> DEU.1.9_1
## 10 DEU Germany DEU.1_1 Baden-Württemberg <NA> DEU.1.10_1
## NAME_2 VARNAME_2 NL_NAME_2 TYPE_2 ENGTYPE_2 CC_2
## 1 Alb-Donau-Kreis <NA> <NA> Landkreis District 08425
## 2 Baden-Baden <NA> <NA> Stadtkreis District 08211
## 3 Biberach <NA> <NA> Landkreis District 08426
## 4 Böblingen <NA> <NA> Landkreis District 08115
## 5 Bodensee <NA> <NA> Water body Water body <NA>
## 6 Bodenseekreis <NA> <NA> Landkreis District 08435
## 7 Breisgau-Hochschwarzwald <NA> <NA> Landkreis District 08315
## 8 Calw <NA> <NA> Landkreis District 08235
## 9 Emmendingen <NA> <NA> Landkreis District 08316
## 10 Enzkreis <NA> <NA> Landkreis District 08236
## HASC_2 geometry
## 1 DE.BW.AD MULTIPOLYGON (((10.00971 48...
## 2 DE.BW.BB MULTIPOLYGON (((8.127194 48...
## 3 DE.BW.BR MULTIPOLYGON (((9.692888 47...
## 4 DE.BW.BL MULTIPOLYGON (((8.846692 48...
## 5 <NA> MULTIPOLYGON (((9.555235 47...
## 6 DE.BW.BD MULTIPOLYGON (((9.6394 47.6...
## 7 DE.BW.BH MULTIPOLYGON (((8.191485 47...
## 8 DE.BW.CL MULTIPOLYGON (((8.600422 48...
## 9 DE.BW.EM MULTIPOLYGON (((8.125374 48...
## 10 DE.BW.EZ MULTIPOLYGON (((8.550905 48...
Wenn nur bestimmte Spalten aus dem importierten spatial data frame von Interesse sind, können sie mit dem Befehl select ausgewählte werden:
bw_select <- dplyr::select(bw, NAME_0, NAME_1, NAME_2, TYPE_2, ENGTYPE_2)
bw_select
## Simple feature collection with 45 features and 5 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 7.51159 ymin: 47.5024 xmax: 10.49575 ymax: 49.79148
## Geodetic CRS: WGS 84
## First 10 features:
## NAME_0 NAME_1 NAME_2 TYPE_2 ENGTYPE_2
## 1 Germany Baden-Württemberg Alb-Donau-Kreis Landkreis District
## 2 Germany Baden-Württemberg Baden-Baden Stadtkreis District
## 3 Germany Baden-Württemberg Biberach Landkreis District
## 4 Germany Baden-Württemberg Böblingen Landkreis District
## 5 Germany Baden-Württemberg Bodensee Water body Water body
## 6 Germany Baden-Württemberg Bodenseekreis Landkreis District
## 7 Germany Baden-Württemberg Breisgau-Hochschwarzwald Landkreis District
## 8 Germany Baden-Württemberg Calw Landkreis District
## 9 Germany Baden-Württemberg Emmendingen Landkreis District
## 10 Germany Baden-Württemberg Enzkreis Landkreis District
## geometry
## 1 MULTIPOLYGON (((10.00971 48...
## 2 MULTIPOLYGON (((8.127194 48...
## 3 MULTIPOLYGON (((9.692888 47...
## 4 MULTIPOLYGON (((8.846692 48...
## 5 MULTIPOLYGON (((9.555235 47...
## 6 MULTIPOLYGON (((9.6394 47.6...
## 7 MULTIPOLYGON (((8.191485 47...
## 8 MULTIPOLYGON (((8.600422 48...
## 9 MULTIPOLYGON (((8.125374 48...
## 10 MULTIPOLYGON (((8.550905 48...
Eine Shapefile in das GeoJSON-Format überführen
GeoJSON hat sich zunehmend zu einem Standardformat zur Speicherung von Geoinformationen entwickelt, insbesondere für Webanwendungen. Anwendungen wie Leaflet vereinfachen es, Geodaten im GeoJSON-Format mit Java Script zu visualisieren.
Ausführliche Informationen über das GeoJSON-Format werden in einem Artikel von Tom MacWright in leicht verständlicher Weise erklärt:
Die meisten heute öffentlich zugänglichen Geodaten liegen dennoch meist im sog. Shapefile-Format vor. Es definiert ebenso einen Standard für die Speicherung der Geometrien Polygone, Linien und Punkte sowie deren Attributinformationen.
Eine Konvertierung von Shapefiles in das GeoJSON-Format ist in R dank des Pakets rgdal leicht zu bewerkstelligen. Im folgenden wird gezeigt, wie GeoJSON-Dateien aus Shapefiles erzeugt und die Daten angepasst werden können, sodass sie der Visualisierung mittels Leaflet dienen können.
Es kann Anwendungsfälle geben bei denen es sinnvoll ist Geometrien für Visualisierungen zu vereinfachen, um die Datengröße zu reduzieren und gleichzeitig die Integrität der Geometrie zu erhalten.
Bevor die Geometriedaten bearbeitet werden, soll der spatial data frame zunächst in das GeoJSON-Format konvertiert werden. Es ist selbstverständlich auch möglich die Änderung von Geometrien im spatial data frame vorzunehmen, aber dank der R-Pakete geojsonio und rmapshaper lassen sie sich sehr simple im GeoJSON-Format vereinfachen.
Mittels der Funktion geojson_json kann ein spatial data frame in das GeoJSON-Format konvertiert werden:
bw_json <- geojsonio::geojson_json(bw_select)
GeoJSON-Dateien können oftmals größer sein, da sie eine Menge Informationen enthalten, um die Polygone oder Linien zu zeichnen. Deshalb ist es ratsam sie zum Laden, Verarbeiten und Zeichnen auf einer Karte zu vereinfachen. Mittels der Funktion ms_simpify aus dem Paket rmapshaper können Detailinformationen von Geometrien für Visualisierungszwecke entfernt werden.
if (!require(rmapshaper)){install.packages("rmapshaper"); library(rmapshaper)}
bw_json_simplified <- ms_simplify(bw_json)
Sie können explizit festlegen, wie viel Sie vereinfachen möchten, indem Sie das Argument keep verwenden.
bw_json_simplified <- ms_simplify(bw_json, keep = 0.1)
Geometrien löschen und ausschneiden
Mit den Funktionen ms_clip und ms_erase aus dem Paket rmapshaper kann ein Ausschnit defniert werden der erhalten bleiben soll bzw. es kann ein Bereich defniert werden der gelöscht werden soll. Hier ist ein Beispiel mit der Funktion ms_clip:
#Calculate a bounding box
bbox_fr <- c(7.65, 47.90, 7.95, 48.08)
bw_json_clipped <- ms_clip(bw_json_simplified, bbox = bbox_fr)
Wenn Koordinatenwerte unbekannt sind, kann ihre bounding Box mit der Funktion lawn_bbox_polygon aus dem Paket lawn wie unten beschrieben unmittelbar in R-Studio visualisiert werden:
#load package "lawn"
#if (!require(lawn)) {install.packages("lawn"); library(lawn)}
#Visualize a bounding box
#view(lawn_bbox_polygon(bbox_fr))
Das Paket lawn bietet eine R-Schnittstelle für’turf.js’, eine Javascrip-Bibliothek, die die Analyse und Transformation von GeoJSON-Daten erleichtert. Es hält viele nützliche Analysewerkzeuge bereit.
GeoJSON exportieren
Fertig berarbeitet GeoJSON-Daten können mittels der Funktion geojson_write vom Paket geojsonio gespeichert bzw. exportiert werden:
# speichert im working directory sofern nicht näher spezifiziert
geojson_write(bw_json_clipped, file = "bw_json_clipped.geojson")
Mit dem “RasterLayerObject” und Rasterdaten arbeiten
Zweidimensionale RasterLayerObjects, die mittels des Pakets raster eingeladen wurden, können mit der Funktion addRasterImage in Bilder umgewandelt, die einer Leaflet-Karte bzw. dem sog. “map widget hinzugefügt werden können.
Am Beispiel eines Jahresrasters der Anzahl der Heißen Tage für Deutschland vom DWD wird gezeigt, wie ein RasterLayerObject eingeladen wird
urlfile2 <- "https://raw.githubusercontent.com/RafHo/teaching/master/geodaten_kommunizieren/datasource/TADXCDGE30_17_2018_01.asc" # get url from github repositury
# lade das raster
dwd <- raster(urlfile2)
CRS defnieren und
Es ist wichtig, dass das RasterLayerObject mit einem geeigneten Koordinatenbezugssystem versehen ist. Viele Rasterdateien enthalten diese Informationen, aber einige nicht. D.h. das RasterLayerObject muss ggf. mit der entsprechenden Information (d.h. dem korrekten Proj4-String) versehen werden, z.B.:
crs(dwd) <- sp::CRS("+proj=tmerc +lat_0=0 +lon_0=9 +k=1 +x_0=3500000 +y_0=0 +ellps=bessel +datum=potsdam +units=m +no_defs")
Bildgröße reduzieren
Die Funktion addRasterImage arbeitet so, dass sie jede Zelle in eine RGB-Farbe kodiert, um ein PNG-Bild zu erzeugen. Dieses Bild wird dann in das Map-Widget eingebettet.
Sobal ein RasterLayerObject in das Map-Widget einbettet wird, vergrößert sich die Größe der erzeugten HTML-Datei proportional. Um unzulässige Downloadzeiten und Speicherverbrauch zu vermeiden, soll das erzeugte PNG die nicht größer als 4MB sein. Größere Rasterebenen können mit den Funktionen raster::resample oder raster::aggregate verringert werden (s. hierzu die Funktionshilfen).
Beispiel: Ein RasterLayerObject in Leaflet visualisieren
Um die Rasterebene als Bild darzustellen, muss jeder Zellwert in eine RGB-Farbe umgewandelt werden. Sie können die Farbskala mit dem Argument colors angeben, das eine Vielzahl von Farbspezifikationen akzeptiert:
- Der Name einer Farbpalette von Color Brewer 2. Wenn kein Farbargument angegeben wird, ist “Spectral” der Standard.
- Ein Vektor, der die geordnete Liste der Farben darstellt, die den Daten zugeordnet werden sollen. Es kann jede Farbspezifikation verwendet werden, die von
grDevices::col2rgbakzeptiert wird. - Eine Farbskalierungsfunktion, z.B.
colors <- colorBin("Greens", Domain = NULL, bins = 5, na.color = "transparent").
pal <- colorNumeric(c("#fcf0dc", "#cc9500", "#940a0a"), values(dwd),
na.color = "transparent")
Die Geovisualisierung des RasterLayerObjects kann mit Leaflet umgesetzt werden:
leaflet() %>% addTiles() %>%
addRasterImage(dwd, colors = pal, opacity = 0.75) %>%
addLegend(pal = pal, values = values(dwd),
title = "Heiße Tage")
Die Verwendung andere Projektionen ist derzeit nicht empfohlen. Hierin liegt eine häufige Ursache für Fehlermeldungen↩︎