Skip to content

Create SpatialLinesDataFrame

6 messages · Michael Sumner, T.D. Rudolph, Roger Bivand

#
Using apply makes it difficult to set an appropriate ID for each Line,
so I've done it in a loop:

ranlines <- list()

for (irow in 1:nrow(temp)) {
	ranlines[[irow]] <- Lines(Line(rbind(as.numeric(temp[irow, c("x",
"y")]), as.numeric(temp[irow, c("X2", "Y2")]))),
		ID = as.character(irow))
}

sldf <- SpatialLinesDataFrame(SpatialLines(ranlines), temp, match.ID = FALSE)


I've also used named arguments as your indexes don't match the table
as you presented it, I presume that start is "x/y" and the end is
"X2/Y2"?

Then you can use this to write a shapefile to the current directory:

library(rgdal)
writeOGR(sldf, ".", "ranlines", "ESRI Shapefile")

I agree there needs to be more examples, the fact of multi-part lines,
polygons and points can be hard to grasp and the need for them makes
the structure of sp objects difficult to understand. It's not hard
once you do though. I posted this as I thought it might be of use:

https://stat.ethz.ch/pipermail/r-sig-geo/2010-February/007619.html

It's not a great example, but it shows how I explore these objects to
get an understanding. In case it's of use I also put this together
this morning, to work with the results of package "alphahull":

This example creates a (possibly multi-branched) line object and
writes to shapefile.

#R

library(alphahull)

library(rgdal)

# Uniform sample of size n=300 in the annulus B(c,0.5)\B(c,0.25),

# with c=(0.5,0.5).

n <- 300

theta<-runif(n,0,2*pi)

r<-sqrt(runif(n,0.25^2,0.5^2))

x<-cbind(0.5+r*cos(theta),0.5+r*sin(theta))



x <- do.call('cbind', locator())

# Value of alpha

alpha <- 0.1

# alpha-shape

x.shape <- ashape(x, alpha = alpha)

## construct a GIS-line object

l <- list()

for (i in 1:nrow(x.shape$edges)) {

 l[[i]] <-  Line(rbind(x.shape$edges[i, 3:4], x.shape$edges[i, 5:6]))

}

## one line object, one ID value

l <- list(Lines(l, as.character("1")))

## equivalent to a drawing with a single line object:

sldf <- SpatialLinesDataFrame(SpatialLines(l), data.frame(name =
"ashape"), match.ID = FALSE)



## write to desired OGR driver (see ogrDrivers() for available ones)

## this creates lines.shp/shx/dbf in current directory

writeOGR(sldf, ".", "lines", "ESRI Shapefile")


Cheers, Mike.





On Tue, Feb 16, 2010 at 11:42 AM, Tyler Dean Rudolph
<tylerdeanrudolph at gmail.com> wrote:
#
On Mon, 15 Feb 2010, Tyler Dean Rudolph wrote:

            
I have discussed with Achim Zeilis, who faces the same problem, the 
ethical appropriateness of mentioning our book on the Spatial Task View, 
and we've found until now that we shouldn't mention our own books. If 
opinion is that we should, and because they are in a number of university 
libraries, I can do that, and encourage Achim to do the same for 
Econometrics.

If you visit the book website (http://www.asdar-book.org), you will find 
links to code showing how to do this (chapter 2). In addition, under 
"Additional materials", you find links to course materials, especially to 
a course given at Imperial College in 2007, in which a good deal of these 
questions are covered (see "Representing Spatial Data").

Because the classes were designed as interface containers, most users do 
not need to know how to construct them from scratch, typically reading in 
data or coercing from other classes of objects. While the help pages could 
certainly be improved, we do provide a vignette that ships with the 
package, and where section 6 covers Line objects etc. See:

http://cran.r-project.org/web/packages/sp/vignettes/sp.pdf

or:

vignette("sp", package="sp")

in R.

I'm sure that Mike's answer resolves your immediate problem - this is more 
of an answer to your general question about where to look for information 
about sp classes and how to make and use them.

Hope this helps,

Roger

  
    
#
You need to put "flop" in a list of its own, as SpatialLines is a
collection class that may contain more than one "Lines" (think of
multi-branched line sequences that may not be connected, but exist as
a single entity/feature - a "Line" is a matrix of a coordinates that
are connected in sequence, a "Lines" is a collection of those, and a
SpatialLines is a collection of Lines - each Lines will have a row in
the attribute table). So SpatialLines(list(flop)) should work.

See in my code how I start with a generic empty list, then populate it
with many Lines:

ranlines <- list()

for (irow in 1:nrow(temp)) {
       ranlines[[irow]] <- Lines(Line
I use the rgdal package as it was my first introduction to read/write
for GIS data in R, and it is very general. It can be hard to install
as there are 3rd-party dependencies even for the base package (very
easy on Windows with a binary/zip package), though when you compile
from source you have greater control and the potental for adding new
drivers (formats) over and above the base ones. There may be important
differences between the way shapefiles are written but it's personal
preference and ease of use/installation as far as I can see.

cheers, Mike.

On Wed, Feb 17, 2010 at 10:28 AM, Tyler Dean Rudolph
<tylerdeanrudolph at gmail.com> wrote:
#
On Wed, 17 Feb 2010, Michael Sumner wrote:

            
I think that there is a misplacement too:

flop <- Lines(list(Line(rbind(as.numeric(temp[1,7:8]),
   as.numeric(temp[1,9:10])))), ID=row.names(temp))

looks better - Lines() takes a list of Line objects as its first argument.

Roger