Czytasz plik KML do R?

42

Pracuję z dużymi plikami .kml (do 10 Gb) i potrzebuję wydajnego sposobu na ich odczytanie do R. Do tej pory przekształcałem je w pliki shapefile za pomocą QGIS, a następnie z powrotem w R z readShapePoly i readOGR (ten ostatni , nawiasem mówiąc, jest ~ 1000 szybszy niż poprzedni). Idealnie chciałbym wyciąć etap pośredni QGIS, ponieważ jest on uciążliwy i powolny.

Jak odczytywać bezpośrednio pliki .kml?

I zobaczyć można również dokonać z readOGR . Niestety nie widzę sposobu wdrożenia działającego przykładu (po długim przygotowaniu pliku .kml:) xx <- readOGR(paste(td, "cities.kml", sep="/"), "cities"). Wydaje się, że „miasta” to tutaj nazwy obiektów przestrzennych.

Roger Bivand przyznaje, że „odkrycie tej nazwy nie jest oczywiste, ponieważ sterownik KML w OGR potrzebuje go do uzyskania dostępu do pliku. Jedną z możliwości jest:

system(paste("ogrinfo", paste(td, "cities.kml", sep="/")), intern=TRUE)

Ale to też nie działa dla mnie. Oto testowy plik .kml do wypróbowania. Z tym w moim katalogu roboczym readOGR("x.kml", "id")generuje ten komunikat o błędzie:

Error in ogrInfo(dsn = dsn, layer = layer, encoding = encoding, use_iconv = use_iconv) : 
  Cannot open layer . 

I system(paste("ogrinfo", "x.kml"), intern=TRUE)generuje:

[1] "Had to open data source read-only."   "INFO: Open of `x.kml'"               
[3] "      using driver `KML' successful." "1: x (3D Polygon)"  

, którego po prostu nie rozumiem.

Czy getKMLcoordinates{maptools} byłoby prawidłową alternatywą?

Próbowałem także:

tkml <- getKMLcoordinates(kmlfile="x.kml", ignoreAltitude=T)
head(tkml[[1]])
tkml <- SpatialPolygons(tkml, 
                        proj4string=CRS("+init=epsg:3857"))

Współrzędne są generowane poprawnie, ale moja próba przekształcenia ich z powrotem w obiekt wieloboku nie powiodła się z następującym komunikatem:

Error in SpatialPolygons(tkml, proj4string = CRS("+init=epsg:3857")) : 
  cannot get a slot ("area") from an object of type "double"
RobinLovelace
źródło
1
Możesz uzyskać warstwy w kml za pomocą funkcji rGdal ogrListLayers.
Mario Becerra

Odpowiedzi:

37

Aby odczytać plik KML za pomocą sterownika OGR, podaj mu nazwę pliku i nazwę warstwy.

Komentarz Rogera jest taki, że nazwa warstwy jest ukryta w pliku KML, i jeśli nie wiesz, jak KML został utworzony, nie możesz wywnioskować nazwy warstwy na podstawie nazwy pliku KML.

Patrząc na twój przykładowy KML, widzę:

<?xml version="1.0" encoding="utf-8" ?>
<kml xmlns="http://www.opengis.net/kml/2.2">
<Document><Folder><name>x</name>
<Schema name="x" id="x">

Co mówi mi, że nazwa warstwy to x: nie id, a więc:

> foo = readOGR("/tmp/x.kml", "x")
OGR data source with driver: KML 
Source: "/tmp/x.kml", layer: "x"
with 1 features and 2 fields
Feature type: wkbPolygon with 2 dimensions

działa ładnie.

Teraz można spróbować uzyskać nazwę przez parsowania KML jako XML przy użyciu parsera XML R, można też może spróbować go czytać w R jako plik tekstowy, aż znajdziesz plakietkę.

Inne podejście polega na uruchomieniu programu ogrinfo z linii poleceń, który wyrzuca nazwy warstw pliku KML:

$ ogrinfo /tmp/x.kml 
Had to open data source read-only.
INFO: Open of `/tmp/x.kml'
      using driver `KML' successful.
1: x (Polygon)

tutaj pokazano, że istnieje warstwa wielokąta o nazwie x.

Spacedman
źródło
Dzięki za odpowiedź Spaced - od razu rozwiązałem problem. Takie jasne wyjaśnienie sprawia, że ​​uwielbiam wymianę stosów! Jedno pytanie „punkt bonusowy”: czy mogę użyć tego samego polecenia do odczytania podzbioru danych (np. Pierwszego miliona wielokątów)? W przeciwnym razie będzie szukać podziału dużych kmls za pomocą zewnętrznego programu.
RobinLovelace 16.04.13
2
KML będący XML nie jest tak naprawdę przeznaczony do losowego dostępu. Prawdziwym rozwiązaniem jest umieszczenie danych przestrzennych w przestrzennej bazie danych i uzyskanie pewnych wskaźników przestrzennych dotyczących prędkości. Sprawdź PostGIS.
Spacedman
OK, dobry plan - powiedziałem klientowi, że PostGIS jest sposobem na osiągnięcie tak dużych zbiorów danych i jestem przekonany, że jest to właściwa opcja dla rzeczy, które chce robić. Dobra wymówka, abym się tego nauczył!
RobinLovelace 16.04.13
Istnieje również przestrzenne rozszerzenie sqlite , baz danych plikowych, które nie wymagałyby instalowania usługi i wymagałyby mniejszej konfiguracji niż PostGIS.
Frank
dziwnie systemw badania potrzebne path.expandna ~na ogrinfodo pracy, mimo że pracował dobrze na nierozszerzonym ścieżki w linii poleceń (MacOS; Sys.which('ogrinfo')i which ogrinfowrócił te same ścieżki)
MichaelChirico
5

Jeśli chcesz zrobić alternatywny sposób za pomocą maptool, powinno to działać:

tkml <- getKMLcoordinates(kmlfile="yourkml.kml", ignoreAltitude=T)
#make polygon
p1 = Polygon(tkml)
#make Polygon class
p2 = Polygons(list(p1), ID = "drivetime")
#make spatial polygons class
p3= SpatialPolygons(list(p2),proj4string=CRS("+init=epsg:4326"))

Kluczem tutaj jest to, że musisz przejść kilka kroków, aby stworzyć klasę wielokąta przestrzennego.

Widziany
źródło
cześć @Seen, próbowałem twojego podejścia, ale wydaje się, że nie działa? Mam błąd: błąd w wielokącie (tkml): współrzędne muszą być macierzą dwukolumnową> head (tkml) [[1]] [1] -87.88141 30.49800 i mam go jako listę .. czy uważasz, że jego konwersja jest dobra lista współrzędnych do matrycy? tahnks!
maycca
1

Nie wiem, czy to nadal stanowi problem dla kogokolwiek innego, ale przez jakiś czas biegałem w kółko. Co w końcu dla mnie zadziałało, poniżej. Używa XMLpakietu, aby dostać się xmlValuedo właściwego węzła. Musiałem ustawić layerparametr readOGRna nazwę jednego z folderów w pliku kml. Gdy layerustawię parametr na plik kml, otrzymam ten sam błąd, co RobinLovelace opisany powyżej.

Poniżej pokazano wiele linii kodu, które pokazują tylko, jak zobaczyć różne poziomy węzłów dokumentu kml. Myślę, że będzie to nieco inne w zależności od źródła kml. Ale powinieneś być w stanie użyć tej samej logiki, aby określić poprawną wartość parametru.

Ponadto, stworzyłem listę plików KML więc może być łatwo przekształcony funkcji, które można umieścić w lapply- do.callpary. Może to następnie pobrać dane z długiej listy plików kml. Lub wiele podfolderów w jednym pliku kml, jak się wydaje, readOGRnie radzi sobie z wieloma podfolderami w pliku kml.

library(rgdal); library(XML)

# SET WORKING DIRECTORY FIRST!!
dir <- getwd()

kmlfilelist <- list.files(dir, pattern =".kml$", full.names=TRUE, recursive=FALSE)

doc0 <- xmlTreeParse(kmlfilelist[2], useInternal = TRUE)
rootNode0 <- xmlRoot(doc0)
rootName0 <- xmlName(rootNode0)
element1Name0 <- names(rootNode0)

nodeNames <- names(rootNode0[1][[1]])

# entire rootNode - kml Document level
rootNode0[[1]]

# 1st element of rootNode - kml file name
rootNode0[[1]][[1]] 

# 2nd element of rootNode - kml Style Map 
rootNode0[[1]][[2]] 

# 3rd element of rootNode - Style
rootNode0[[1]][[3]]

# 4th element of rootNode - Style
rootNode0[[1]][[4]] 

# 5th element of rootNode - kml Folder with data in it.
rootNode0[[1]][[5]] 

# 5th element 1st subelement of rootNode - kml Folder name with data in it. 
#  What to set readOGR() layer parameter to.
rootNode0[[1]][[5]][[1]] 

kmlfoldername <- xmlValue(rootNode0[[1]][[5]][[1]]) # Folder name to set = layer.

readOGR(dsn=kmlfilelist[2], layer =  kmlfoldername)
Glina
źródło
0

Nie wiem, czy powinienem zmodyfikować poprzednią odpowiedź. Być może, ale dotyczy to niektórych kwestii, których nie ma w tej odpowiedzi, więc postanowiłem to zostawić.

W każdym razie poniższy kod działa dla mnie dobrze. Wyszukuje wszystkie xmlNodes w pliku kml o nazwie „Folder”, a następnie ustawia layerparametr readOGRna to xmlValue. Testowany na katalogu roboczym z około 6 oddzielnymi plikami kml. Dane wyjściowe to lista importowanych obiektów SpatialDataFrames. Każdy SpatialDataFrame może być łatwo podzbiór z listy.

Nadal nie adresuje plików kml z wieloma węzłami folderów. Ale tę funkcję można łatwo dodać za pomocą innej applyfunkcji zagnieżdżonej .

library(rgdal); library(XML)

# SET WORKING DIRECTORY FIRST!!
dir <- getwd()

kmlfilelist <- list.files(dir, pattern =".kml$", full.names=TRUE, recursive=FALSE)

ImportKml <- function (kmlfile) {
  doc0 <- xmlTreeParse(kmlfile, useInternal = TRUE)
  rootNode0 <- xmlRoot(doc0)
  rootName0 <- xmlName(rootNode0)
  element1Name0 <- names(rootNode0)

  kmlNodeNames <- unname(names(rootNode0[1][[1]]))
  kmlFolderNodeNum <- which(kmlNodeNames == "Folder")
  kmlFolderNodeName <- xmlValue(rootNode0[[1]][[kmlFolderNodeNum]][[1]])

  kmlIn <- readOGR(dsn=kmlfile, layer = kmlFolderNodeName)
}
ImportedKmls <- lapply(kmlfilelist, ImportKml)
Glina
źródło