PBSmapping to SpatialPolygonsDF

Hi,

I am also interested in this convertion function. In fact, yesterday
night I wrote some code to convert from PolySet to SP, but Roger's code
is much better than mine... 

If nobody is already doing it, can I take Roger's code, add a few checks
and make a function to be included in, say, maptools? 

Best regards,

Virgilio

El lun, 26-03-2007 a las 16:21 -0400, Andrew Niccolai escribi?:

Roger!!! Thank you very much.  I needed to adjust your code only to account
for the fact that my tesselated object did not have a SID column and also I
needed to add a step to "repeat" the first set of XY cords for each polygon
to close the loop.  Other than that, it worked like a charm out of the gate.

Thank you.

Here are my changes to Roger's code in case others have data sets w/o closed
loops.

##--- convert from Polyset in PBSmapping to SPDF
## example(PolySet2SpatialPolygons ---written by R.Bivand)
res0 <- split(cb.tess, cb.tess$PID)
for(i in 1:length(res0)){
res0[[i]] <- rbind(res0[[i]],res0[[i]][1,])
}
outPolygons <- vector(mode="list", length=length(res0)) 
for (i in seq(along=outPolygons)) {
outPolygons[[i]] <- Polygons(lapply(res0[[i]], function(x)
Polygon(cbind(res0[[i]]$X, res0[[i]]$Y))), ID=as.character(i))
} 
outSP <- SpatialPolygons(outPolygons)
plot(outSP)

Andrew Niccolai
Doctoral Candidate
Yale School of Forestry
203-432-5144
 

-----Original Message-----
From: Roger Bivand [mailto:Roger.Bivand at nhh.no] 
Sent: Monday, March 26, 2007 3:34 PM
To: Andrew Niccolai
Cc: r-sig-geo at stat.math.ethz.ch
Subject: Re: PBSmapping to SpatialPolygonsDF

On Mon, 26 Mar 2007, Andrew Niccolai wrote:

Is there a pre-existing function or a clever way to pull a PBSmapping

object

of class "Polyset", "data.frame" into a SpatialPolygonsDataFrame?

I found the SpatialPolygons2PolySet function and was hoping that there is

a

way to go the other direction.

I have created a voronoi tessellation for a set of 627 points and then
clipped that by the convex hull of those points.  The object looks like
this:

class(cb.tess)

[1] "PolySet"    "data.frame"

names(cb.tess)

[1] "PID" "POS" "X"   "Y"  

head(cb.tess)

  PID POS        X        Y
1   1   1 945062.2 860290.8
2   1   2 945073.5 860273.0
3   1   3 945062.7 860272.6
4   1   4 945062.5 860272.7
5   1   5 945055.9 860277.4
6   1   6 945057.8 860291.4

tail(cb.tess)

  PID POS        X        Y
1 626   9 944972.3 859813.1
2 627   1 945009.7 859794.3
3 627   2 944993.8 859769.0
4 627   3 944985.3 859763.3
5 627   4 944990.9 859804.2
6 627   5 945008.2 859804.6

I would like to slip back into the world of SPDF if possible.

Thanks in advance for help, suggestions or code.

Not yet a function, but if you could try this and polish it a bit, I think 
it's workable:

example(SpatialPolygons2PolySet)
# gives the input PolySet nor_coast_poly_PS
res0 <- split(nor_coast_poly_PS, nor_coast_poly_PS$PID)
res1 <- lapply(res0, function(x) split(x, x$SID))
outPolygons <- vector(mode="list", length=length(res0))
for (i in seq(along=outPolygons)) {
   outPolygons[[i]] <- Polygons(lapply(res1[[i]], function(x) 
   Polygon(cbind(x$X, x$Y))), ID=as.character(i))
}
outSP <- SpatialPolygons(outPolygons)
plot(outSP)

# get the projection/zone from attr(nor_coast_poly_PS, "projection") and
# attr(nor_coast_poly_PS, "zone") if UTM, use in constructing
# the SpatialPolygons object. Now the ID= in Polygons() is just the loop 
# index - should be the PID values if they are the key value for adding
# in a data frame.

Hope this helps,

Roger

cheers

Andrew Niccolai
Doctoral Candidate
Yale School of Forestry
203-432-5144