What you describe doing here is substantially different from what you
previously described. Also what you describe doing must be substantially
different, yet again, from what you actually did, since errors would have
been thrown by what you describe doing.
Finally, what you describe doing makes little if any sense.
See comments in line below.
On 19/01/12 13:38, Jose Funes wrote:
> Thanks Rolf for your reply.
>
> I am attaching some data. I am a bit confused and I dont know were i
> am making the mistakes. I followed p.50 from "Analysing spatial point
> patterns in R" workshop notes. I get the window from the code below:
>
> It creates the window and the matrix with no error message
> image to matrix
> >sudan1<-readGDAL("sudanbndy")
> >im.sudan<- as.im <http://as.im/>(sudan1)
> >Vsudan<- im.sudan$v
> >Wsudan<- levelset(Vsudan,21,"==")
This cannot have been what you did, since would have thrown an
error. The function levelset() takes an object of class "im" as its
(first) argument --- ***NOT*** a matrix.
You may have done
Wsudan <- levelset(im.sudan,21,"==")
> >matrix.sudan<- as.matrix(Wsudan)
This makes no sense at all; Wsudan is an object of class "owin" and
there is no ".owin" method for as.matrix(). Consequently
as.matrix.default()
gets called here, and yields a mess, resulting in the error you got.
(It is probably true that as.matrix.default() should have better
error traps,
but that is not really relevant here.)
>
> Revising the outputs:
> The binary image mask seems fine
> >Wsudan
> window: binary image mask
> 290 x 251 pixel array (ny, nx)
> enclosing rectangle: [0.5, 251.5] x [0.5, 290.5] units
>
> But when I do check the output for the matrix, it does seems to be a
> problem
> >matrix.sudan
> window: Error in switch(x$type, rectangle = { : EXPR must be a length
> 1 vector
As explained above, matrix.sudan is not a matrix; it's a mess.
You *may* have wanted to do
matrix.sudan <- Wsudan$m
which extracts the matrix of TRUE/FALSE pixel values from Wsudan.
>
> Then I capture xrange and yrange as follow,
> >xrange<-Wsudan$xrange
> >yrange<-Wsudan$yrange
>
> ###########Text file#############
> >sudan98<- read.csv("sudan98.txt")
Sigh! There is no such file. The file you sent me is called:
sudan1998.txt
So I did:
sudan98<- read.csv("sudan1998.txt")
> >sudan98sa<- ppp(sudan98$Longitude, sudan98$Latitude,
> owin(mask=matrix.sudan,xrange,yrange))
Well this obviously won't work for *your* "matrix.sudan" since it
isn't a matrix at all. Moreover, even with a sensible value for
"matrix.sudan"
it still won't work because the third argument to ppp() is "..."
--- which consists
of "arguments passed to owin to create the window, if window is
missing".
(RTFM.)
So you need to do:
> sudan98sa<- ppp(sudan98$Longitude, sudan98$Latitude,
window=owin(mask=matrix.sudan,xrange,yrange))
With this amendment and with my alternative proposal for the value of
"matrix.sudan", i.e Wsudan$m, it does "work".
But what on earth is the point of playing ring-around-the-rosie with
matrix.sudan anyway? The object Wsudan is already a window (of type
"mask") and appears to be the window you want. Why are you diddling
about and re-constructing it? Note:
U <- owin(mask=matrix.sudan,xrange,yrange)) # With *my* version
of matrix.sudan.
all.equal(Wsudan,U)
TRUE
So you are just reinventing a rather boring wheel.
Moreover, all the fooling around you are doing with "levelset()" is
a complete
waste of time since "im.sudan" has only *one* value, namely 21, to
start with.
So you are not imposing any restriction by extracting the level 21
level set.
You could simply have done:
V <- as.owin(im.sudan)
Note:
all.equal(V,Wsudan)
TRUE
It seems to me that all you want is to obtain a window from
"sudanbndy".
In which case all that you need to do is:
sudan1 <- readGDAL("sudanbndy")
W <- as.owin(sudan1)
sudan98 <- read.csv("sudan1998.txt")
sudan98sa <- ppp(sudan98$Longitude,sudan98$Latitude,window=W)
Finally, note that there is a great deal of duplication in the
points of the pattern
you are playing with --- there are only 59 unique points (and two
of these points
fall outside of the window obtained from sudan1, presumably due to
an overly
coarse pixellation of the map of Sudan)
You can get rid of the duplication by using unique() (which has a
method for
class "ppp"). You could try replacing the points which fall
outside of the
Sudan window by the centres of the pixels inside the window which are
nearest to these two points.
But mostly you need to understand what you are doing rather than
blindly
typing R commands and hoping. Also when sending inquiries to this
list,
show the commands that you actually used; don't make them up.
cheers,
Rolf Turner
[[alternative HTML version deleted]]
_______________________________________________
R-sig-Geo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
