Billet

Cartographie avec R (suite)

Billet publié le 15/12/2009

Je cherche à donner à voir, par des points sur une carte, la localisation d’églises (ou de boulangeries, ou de sex-shops, ou de lobbyistes…) en région parisienne. Il est possible de créer un “mashup” avec google maps, ou une carte dans google earth, mais cela ne donne pas de jolis fichiers PDF utilisables dans une publication scientifique qui se respecte. Imaginons que je dispose des données “Longitude / Latitude” des églises.

Il me faut un fond de carte. On trouve une carte de la France (avec les frontières administratives) sur “cloudmade” : http://downloads.cloudmade.com/europe/france. Il faut télécharger le fichier : “france.shapefiles.zip”
On trouve aussi, ailleurs, une carte des principales rues, routes, autoroutes… d’Île de France : http://download.geofabrik.de/osm/europe/france/ : il faut télécharger la carte de l’Île de France : ile-de-france.shp.zip

Ces cartes “open source” proviennent du projet OpenStreetMap : il y a des erreurs, des morceaux non complets, des manques. Mais à notre échelle, cela suffira. Les fichiers téléchargés sont des “shapefiles”. Ils consistent en 4 fichiers différents : un fichier .prj qui contient des informations concernant la projection, puis trois autres fichiers contenant les données elles-mêmes (un fichier dbf, un fichier shp et un fichier shx).

Ouvrons maintenant R.

library(maptools) #charge le package "maptools"
france<-readShapeLines(
"Desktop/france/france_administrative.shp",
proj4string=CRS("+proj=longlat")
)

l’instruction précédente demande à R de charger les informations de la carte de France dans “france”.

summary(france) # donne la structure de "france" 

On constate que dans cette “Data frame” il est indiqué, par “ADMIN_LEVE” le type de frontière administrative: 8 pour les communes, 6 pour les départements.

routesidf<-readShapeLines(
"Desktop/ile-de-france/roads.shp",
proj4string=CRS("+proj=longlat")
)
summary(routesidf)

permet de constater que le type de route est indiqué par “type” : “primary”, “secondary”, “residential”…

les fichiers peuvent être longs à se charger : ce sont des objets très lourds et il serait préférable de demander à ne charger qu’une petite partie des fichiers (par exemple les routes principales et pas tous les chemins communaux). Mais je ne sais pas le faire… pas encore du moins.

plot(france,xlim=c(2.35,2.45),ylim=c(48.87,48.97),lty=3)

donne l’image suivante. Seul un regard averti y discernera le nord de Paris et une partie de la Seine-Saint-Denis :
Paris-Nord
Rendons cette carte un peu plus lisible :

plot(france[france$ADMIN_LEVE==6,],add=TRUE,lwd=2)
plot(routeidf[routeidf$type=="primary",],add=TRUE,lwd=2,col="lightgray")
plot(routeidf[routeidf$type=="secondary",],add=TRUE,lwd=2,col="lightgray")

Paris-Nord2
J’ai ajouté les routes principales (de type “primary” et “secondary”), j’ai indiqué certaines des frontières départementales par un trait noir. Je vais maintenant ajouter mes églises, qui sont dans l’objet “coordeglises” : X indiquant la longitude et Y la latitude. :

points(coordeglises$X,coordeglises$Y,pch=20,cex=2,col="red")

Paris-Nord3

Il me semble pouvoir remarquer que mes églises s’installent assez souvent à proximité de ces grandes routes, voire même à proximité du croisement de deux de ces grandes routes.

Note : Mis à part le bel iMac sur lequel j’ai réalisé ces cartes, tout le reste fut “gratuit”. Open Source ou non. Seashore, R, OpenStreetMaps… et l’indispensable géocodage offert par google….

7 commentaires

Un commentaire par Frédéric Dejean (23/12/2009 à 21:23)

et dire que j’ai placé mes point un par un en utilisant un fond de carte communal…. Et tu arrives à exporter tes cartes sous adobe pour les retravailler?

Un commentaire par Baptiste Coulmont (23/12/2009 à 22:03)

>Frederic : oui, R génère des PDF qu’Illustrator peut lire.

Un commentaire par Lamireau (08/01/2010 à 16:54)

Je réagis un peu tard, mais WOAW !.. Ce travail cartographique, présenté avec une telle clarté, donne envie de se plonger dans les méandres des SIG “libres”. C’est en tout cas l’une de mes résolutions pour 2010.

Un commentaire par Baptiste Coulmont (09/01/2010 à 10:27)

>Lamireau : “méandres”, le terme est bon… il y a un peu d’apprentissage à faire et des informations à chercher dans des forums…

Un commentaire par Ouro-Bang'na wahabou (10/02/2010 à 22:11)

Bonsoir!
J’essaie de procéder comme vous l’avez indiqué. Cependant cela ne marche pas. Je reçois un message d’erreur.

> togo<-readShapeLines("Desktop/togo/togo_administrative.shp", proj4string=CRS("+proj=longlat"))

Erreur dans getinfo.shape(filen) : Error opening SHP file

J'ai téléchargé dezipé et laisser sur le bureau mon dossier togo.chapefiles avant de passer au R.

Que faire?

Cordialement,
Wahabou

Un commentaire par Baptiste Coulmont (11/02/2010 à 11:11)

> Wahabou : A mon avis, votre problème vient simplement du fait que “Desktop” ne désigne pas votre bureau. Vous devez peut-être écrire quelque chose du type “C://WINDOWS/Bureau…”

Un commentaire par Regis (23/07/2012 à 16:45)

Bonjour,

En découvrant votre site j’ai eu l’envie de tester votre tutorial sur un shapefile de la belgique. En précisant correctement le chemin d’accès au.shp ,j’ai le même message d’erreur que Wahabou.
Que faire?
Cordialement