2014-07-09 1 views
4

J'essaie d'obtenir des quartiers de Varsovie et les dessiner sur google map. En utilisant ce code, où 2536107 est le code de la relation pour OpenStreetMap seul quartier de Varsovie, me donne presque ce que je veux, mais avec quelques bugs. Il y a un contour général mais aussi des lignes entre les points qui ne devraient pas être connectés. Qu'est-ce que je fais mal?Tracer OpenStreetMap avec ggmap

map <- get_googlemap('warsaw', zoom =10) 
warszawa <- get_osm(relation(2536107), full = T) 
warszawa.sp <- as_sp(warszawa, what='lines') 
warsawfort <- fortify(warszawa.sp) 

mapa_polski <- ggmap(map, extent='device', legend="bottomleft") 
warsawfort2 <- geom_polygon(aes(x = long, y = lat), 
       data = warsawfort, fill="blue", colour="black", 
       alpha=0.0, size = 0.3) 

base <- mapa_polski + warsawfort2 
base 

Edit: je me suis dit, il doit être en quelque sorte lié à l'ordre de tracer tous les points/ligne, mais je ne sais pas comment résoudre ce problème.

Répondre

9

Il y a une façon de générer votre carte sans utiliser les paquets externes: ne pas utiliser osmar ...

This link, à l'excellent site Mapzen, fournit un ensemble de shapefile des zones administratives en Pologne. Si vous téléchargez et décompressez-le, vous verrez un ensemble shapfile appelé warsaw.osm-admin.*. Il s'agit d'un fichier de formes polygonal de tous les districts de Pologne, indexé convenablement par osm_id (!!). Le code ci-dessous suppose que vous avez téléchargé le fichier et décompressé dans le "répertoire avec vos fichiers de formes".

library(ggmap) 
library(ggplot2) 
library(rgdal) 

setwd(" <directory with your shapefiles> ") 
pol <- readOGR(dsn=".",layer="warsaw.osm-admin") 
spp <- pol[pol$osm_id==-2536107,] 
wgs.84 <- "+proj=longlat +datum=WGS84" 
spp <- spTransform(spp,CRS(wgs.84)) 

map <- get_googlemap('warsaw', zoom =10) 
spp.df <- fortify(spp) 

ggmap(map, extent='device', legend="bottomleft") + 
    geom_polygon(data = spp.df, aes(x = long, y=lat, group=group), 
       fill="blue", alpha=0.2) + 
    geom_path(data=spp.df, aes(x=long, y=lat, group=group), 
      color="gray50", size=0.3) 

Deux nuances: (1) Les ID osm sont stockés sous forme de nombres négatifs, de sorte que vous devez utiliser, par exemple,

spp <- pol[pol$osm_id==-2536107,] 

pour extraire le district concerné, et (2) le shapefile n'est pas projeté en WGS84 (long/lat). Nous devons donc reprojeter en utilisant:

spp <- spTransform(spp,CRS(wgs.84)) 

La raison osmar ne fonctionne pas est que les chemins sont dans le mauvais ordre. Votre warszawa.sp est un SpatialLinesDataframe, composé d'un ensemble de chemins (12 dans votre cas), constitués chacun d'un ensemble de segments de ligne. Lorsque vous utilisez fortify(...) à ce sujet, ggplot essaie de les combiner en une seule séquence de points. Mais puisque les chemins ne sont pas dans un ordre convexe, ggplot essaie, par exemple, de relier un chemin qui se termine au nord-est, à un chemin qui commence au sud-ouest. C'est pourquoi vous obtenez toutes les lignes supplémentaires. Vous pouvez le voir en coloriant les segments:

xx=coordinates(warszawa.sp) 
colors=rainbow(11) 
plot(t(bbox(warszawa.sp))) 
lapply(1:11,function(i)lines(xx[[i]][[1]],col=colors[i],lwd=2)) 

Les couleurs sont dans l'ordre "rainbow" (rouge, orange, jaune, vert, etc.). Clairement, les lignes ne sont pas dans cet ordre.

EDIT Réponse au commentaire de @ ako.

Il existe un moyen de "réparer" l'objet SpatialLines, mais ce n'est pas trivial. La fonction gPolygonize(...) du package rgeos prend une liste de et la convertit en objet SpatialPolygons, qui peut être utilisé dans ggplot avec fortify(...). Un énorme problème (que je ne comprends pas, franchement), c'est que l'objet warszaw.sp de l'OP a 12 lignes, dont deux semblent être des doublons - cela fait gPolygonize(...) échouer. Donc, si vous créez une liste SpatialLines avec seulement les 11 premiers chemins, vous pouvez convertir warszawa.sp en un polygone. Ce n'est pas général cependant, car je ne peux pas prédire comment ou si cela fonctionnerait avec d'autres objets convertis à partir de osm. Voici le code, qui mène à la même carte que ci-dessus.

library(rgeos) 
coords <- coordinates(warszawa.sp) 
sll <- lapply(coords[1:11],function(x) SpatialLines(list(Lines(list(Line(x[[1]])),ID=1)))) 
spp <- gPolygonize(sll) 
spp.df <- fortify(spp) 
ggmap(map, extent='device', legend="bottomleft") + 
    geom_polygon(data = spp.df, aes(x = long, y=lat, group=group), 
       fill="blue", alpha=0.2) + 
    geom_path(data=spp.df, aes(x=long, y=lat, group=group), 
      color="gray50", size=0.3) 
+0

+1 pour expliquer l'erreur si clairement. Si le fichier alternatif n'était pas disponible, existe-t-il un moyen trivial d'utiliser la topologie pour réordonner par adjacence, c'est-à-dire; réparer plutôt que de remplacer? – ako

+0

L'homme, votre or! Je me suis rendu compte qu'il s'agissait de lignes tracées dans le mauvais ordre mais je ne me suis jamais souvenu que quelque part il y avait une solution aussi simple que ce fichier de Varsovie (je pensais que ce serait trop beau pour un ensemble de données particulier) - Maintenant tout ce que je dois faire est d'importer warsaw.osm-admin.shp et tout est pris en charge! Un grand merci pour vous homme! –

+1

@ako Oui, il y en a mais ce n'est pas trivial. Voir mes modifications à la fin. – jlhoward

2

Je ne suis pas sûr que ce soit un blocage général - je peux reproduire votre exemple et voir le problème. Ma première pensée était que vous ne fournissiez pas group = id, qui est généralement utilisé pour les polygones avec beaucoup de lignes, mais vous avez des lignes, donc cela ne devrait pas être nécessaire.

La seule façon de l'afficher correctement était de changer vos lignes en script polygone off. Qgis' line to polygon n'a pas obtenu ce "droit", obtenir un grand trou de beignet, alors j'ai utilisé ArcMap, qui a produit un polygone complet. Si c'est un hors fonction qui peut fonctionner pour votre flux de travail. Les chances sont ce n'est pas. Dans ce cas, RGDAL peut peut-être transformer des lignes en polygones, en supposant qu'il s'agisse bien d'un problème général. En lisant le fichier de formes polygone et en le fortifiant, votre code s'est exécuté sans problème. enter image description here