[R-sig-Geo] How should i change the shaplefile of a polygon to be used as a window in ppp function of spatstat?--still wrong?

[Roger Bivand]

On Fri, 1 Jun 2007, zhijie zhang wrote:

> Dear Prof .Roger,
>   There is still a little problem with your method, the following are the
> programs:

Thanks for attaching a small data file, makes it easier to give advice.

> 
> a<-readShapePoly("d:/deleting/a/guichi2.shp")
> W <-as(as(a, "SpatialPolygons"), "owin")  # polygons window
> plot(W)
> 
> cases<-coordinates(readShapePoints('d:/deleting/a/cases.shp'))  #points
> cases<-data.frame(cases)
> names(cases)<-c("x","y")
> cases[1:2,]
> 
> plot(W)
> points(cases$x,cases$y)  #no errors,83 points were inside the window
> 
> W<-as.owin(W)
> points<-ppp(cases$x,cases$y,win=W)  *#????????????
> **#error*:83 points were rejected as lying outside the specified window in:
> ppp(cases$x, cases$y, win = W)

In this case, the argument in the call to ppp() must be given by its full 
name:

library(maptools)
library(spstat)
library(spatstat)
a <- readShapePoly("guichi2.shp")
W <- as(as(a, "SpatialPolygons"), "owin")
cases <- coordinates(readShapePoints("cases.shp"))
cases
colnames(cases)<-c("x","y")
p2 <- ppp(cases[,1], cases[,2], window=W)
plot(p2)
points(cases, pch=3)

works, as does:

p1 <- as.ppp(cases, W=W)
plot(p1)
points(cases, pch=3)

Naming the argument "win" after the "..." was not being seen as a partial 
match of window, but being passed out in the "..." list.

Roger



> *In fact,* 83 points were inside the specified window,what's wrong?
> (This packages:
> spatstat     mgcv maptools  foreign       sp
> "1.11-6" "1.3-24" "0.6-12" "0.8-20" "0.9-13")
>  BTW, i put the data in the attachment so that you can check it.
> Thanks.
> 
> On 6/1/07, Roger Bivand <Roger.Bivand at nhh.no> wrote:
> >
> > On Fri, 1 Jun 2007, zhijie zhang wrote:
> >
> > > Dear friends,
> > >   I have a polygon, which is the shapefile format. I want to use it as
> > the
> > > window of points in ppp(x,y,window) of spatstat package, how should i
> > do?
> > > The following programs can't work,:
> > > a<-readShapePoly("d:/deleting/a/a.shp")  # read the boundary data
> > > a2<-Polygon(coordinates(a), hole=as.logical(NA))
> > > ppp(x,y,a2)
> >
> > Using the shapefile shipped with maptools:
> >
> > library(maptools)
> > library(spatstat)
> > xx <- readShapePoly(system.file("shapes/sids.shp", package="maptools")[1])
> > xx1 <- xx[1,] # and taking just the first county
> > plot(xx1)
> > W <- as(as(xx1, "SpatialPolygons"), "owin")
> > plot(W)
> > points(rpoispp(1000, win=W))
> >
> > (This for:
> >
> > other attached packages:
> > spatstat     mgcv maptools  foreign       sp
> > "1.11-5" "1.3-24" "0.6-12" "0.8-20" "0.9-13")
> >
> > Hope this helps,
> >
> > Roger
> >
> >
> > >  THANKS.
> > >
> >
> > --
> > Roger Bivand
> > Economic Geography Section, Department of Economics, Norwegian School of
> > Economics and Business Administration, Helleveien 30, N-5045 Bergen,
> > Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
> > e-mail: Roger.Bivand at nhh.no
> >
> >
> 
> 
> 

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no