R from SpatialPointsDataFrame to SpatialLines

Question

I have a SpatialPointsDataFrame loaded with

pst<-readOGR("/data_spatial/coast/","points_coast") 

And I would like to get SpatialLines in the output, I found somthing

 coord<-as.data.frame(coordinates(pst)) Slo1<-Line(coord) Sli1<-Lines(list(Slo1),ID="coastLine") coastline <- SpatialLines(list(Sli1)) class(coastline) 

it works, but when I try the plot (coastline), I have a line that should not be ...

Can anybody help me? The form file is here !

Answer 1

I looked at the shapefile. There is an id column, but if you are drawing data, it seems that the identifier is not ordered from north to south or something like that. Extra rows are created because the point order is not perfect, connecting the points that are next to each other in the table, but far from each other in terms of space. You can try to figure out the correct ordering of the data by calculating the distances between the points and then arranging them at a distance.

A workaround is to remove those lines whose length exceeds a certain distance, for example. 500 m. First find out where the distance between successive coordinates is greater than this distance: breaks . Then grab a subset of the coordinates between the two breaks and finally create Lines for that subset. As a result, you get a coastline consisting of several segments ( breaks-1 ) and without error.

 # read data library(rgdal) pst<-readOGR("/data_spatial/coast/","points_coast") coord<-as.data.frame(coordinates(pst)) colnames(coord) <- c('X','Y') # determine distance between consective coordinates linelength = LineLength(as.matrix(coord),sum=F) # 'id' of long lines, plus first and last item of dataset breaks = c(1,which(linelength>500),nrow(coord)) # check position of breaks breaks = c(1,which(linelength>500),nrow(coord)) # plot extent of coords and check breaks plot(coord,type='n') points(coord[breaks,], pch=16,cex=1) # create vector to be filled with lines of each subset ll <- vector("list", length(breaks)-1) for (i in 1: (length(breaks)-1)){ subcoord = coord[(breaks[i]+1):(breaks[i+1]),] # check if subset contains more than 2 coordinates if (nrow(subcoord) >= 2){ Slo1<-Line(subcoord) Sli1<-Lines(list(Slo1),ID=paste0('section',i)) ll[[i]] = Sli1 } } # remove any invalid lines nulls = which(unlist(lapply(ll,is.null))) ll = ll[-nulls] lin = SpatialLines(ll) # add result to plot lines(lin,col=2) # write shapefile df = data.frame(row.names=names(lin),id=1:length(names(lin))) lin2 = SpatialLinesDataFrame(sl=lin, data=df) proj4string(lin2) <- proj4string(pst) writeOGR(obj=lin2, layer='coastline', dsn='/data_spatial/coast', driver='ESRI Shapefile')