Créer un polygone de coque convexe à partir de points et enregistrer sous shapefile

besoin d'aide concernant un problème de conversion en R.

j'ai calculé le coque convexe d'un nuage de points. Je voudrais, à partir des points formant la coque convexe, construire un objet polygone et enregistrer en tant que shapefile qui peut être lu par un logiciel SIG (ArcMap ou similaire).

Mon code ressemble à ceci:

gps <- read.csv(f)  ##reads the lat-long coordinates  file 
x <- gps$LONGITUDE  ##tells R which columns is which
y <- gps$LATITUDE
z<-chull(x,y)       ##calculates the convex hull --this is just a list of x-y points, N vertex 
dfHull <-cbind(x[z],y[z])  ##the convex hull expressed as a list of selected x-y points
plot(dfHull)     ##this plots the vertex of the polygon, just a check
lines(dfhull)    ##plots the polygon in screen

##generate polygon shapefile, from dfHull, and save it externally as a shapefile ???

le fichier source ne contient que des coordonnées lat-long, E. g:

52.73336     N  0.365974
52.7332  N  0.366051
52.73289     N  0.36636
52.73297     N  0.366258
52.73298     N  0.366243
52.733   N  0.366112
52.73308     N  0.365942
52.73317     N  0.365881
52.73321     N  0.36593
52.73328     N  0.365942
52.73352     N  0.36579
52.73362     N  0.365678
52.73391     N  0.365536
52.7373  N  0.36543
52.73289     N  0.36728

je sachez qu'il y a des paquets (rgdal,maptools,..) pour aider avec ceux-ci, mais je suis très peu familier avec les choses spatiales. Tout ce dont j'ai besoin, c'est de générer l'objet polygon et de le sauvegarder sous forme de shapefile.

Toute aide appréciée. Merci d'avance, dev.

15
demandé sur Lennert 2014-09-01 17:03:18

1 réponses

Voici un exemple simple pour créer un SpatialPolygonsDataFrame, qui peut être sauvegardé comme un shapefile avec rgdal::writeOGR():

set.seed(1)
dat <- matrix(stats::rnorm(2000), ncol = 2)
ch <- chull(dat)
coords <- dat[c(ch, ch[1]), ]  # closed polygon

plot(dat, pch=19)
lines(coords, col="red")

convex hull

library("sp")
library("rgdal")

sp_poly <- SpatialPolygons(list(Polygons(list(Polygon(coords)), ID=1)))
# set coordinate reference system with SpatialPolygons(..., proj4string=CRS(...))
# e.g. CRS("+proj=longlat +datum=WGS84")
sp_poly_df <- SpatialPolygonsDataFrame(sp_poly, data=data.frame(ID=1))
writeOGR(sp_poly_df, "chull", layer="chull", driver="ESRI Shapefile")
22
répondu rcs 2014-09-01 14:37:04