1

I am trying to convert a raster object to an .im object for use with a point process model in the spatstat package in R. I begin by creating the raster from a tiff file using the raster() package. No problems there. I then proceed by cropping the raster according to a given extent. Again, no problem there. I then specify a spatial window (owin) defined using the same extent. Still no problems. When I then proceed to the final step of converting the raster to the im object using as.im(), the function runs and the new im object is created, but it has somehow lost the pixel information that was contained in the original raster such that each pixel now has the same value in the im object.

The date file used is here: https://www.dropbox.com/s/n67djm3n0tfa6sx/MGVF_2001_30_arc_sec.tif?dl=0

I also have a pre-cropped file, which is much smaller in size, here: https://www.dropbox.com/s/co2rsb4vbgn5kog/mgvf.2001.africa.tif?dl=0

And the R code is as follows:

library(raster)
library(spatstat)

# First set the geographic extent we'll be using
e <- extent(-20, 60, -40, 35)

# Then read in the Maximum Green Vegetation Fraction tiff and crop it
mgvf <- raster("MGVF_2001_30_arc_sec.tif")
mgvf.2001.africa <- crop(mgvf, e)

# Now let's create a window for in spatstat
SP.win <- as(e, "SpatialPolygons")
W <- as(SP.win, "owin")

# Finally, we create the .im object
mgvf.img <- as.im(X = "mgvf.2001.africa", W = W)

# Notice, there are no errors thrown. However, compare the plots below and see the loss of information:
plot(mgvf.2001.africa)
plot(mgvf.img)

Incidentally, I have tried the above as shown as well as trying to replace the NAs in the raster prior to converting to im. The result is the same.

  • Example raster is too big, could you plz crop or resize it to have several MBs maximum? – SS_Rebelious Sep 25 '14 at 18:38
  • Thanks @SS_Rebelious, here is the link for the cropped tiff: dropbox.com/s/co2rsb4vbgn5kog/mgvf.2001.africa.tif?dl=0 – thegowda Sep 25 '14 at 20:00
  • It's because you as.im a character string. Converting to im is easy enough, but I don't want to write it just now. I can dig it up later. We should get the method added to raster. – mdsumner Sep 25 '14 at 20:43
  • 1
    Actually I thought it already was there: gis.stackexchange.com/a/73357/1204 need to check. – mdsumner Sep 25 '14 at 20:47
  • Thanks @mdsumner; that method appears to apply to a raster brick, however. Whereas I am having the problem with a single raster. – thegowda Sep 25 '14 at 21:00
1

If you would check help for as.im() you would notice that raster data type is not supported. We should thank @mdsumner for bringing up my old question (I do not remember the occasion that forced me to ask it - suppose it was just some testing and I wouldn't recall that question by myself). I used there a geostatsp package that provides needed functionality.

In your case it is just:

library(geostatsp)
mgvf.img <- as.im(mgvf.2001.africa)

Notice that if you don't want to use geostatsp (for no good reason), as.im() will work with matrices. So you will need to as.matrix(mgvf.2001.africa) and then rotate it (several times probably - I'm too sleepy to test it) with

rotate = function(mat) t(mat[nrow(mat):1,,drop=FALSE]) # mat - matrix

and then pass rotated matrix to as.im().

0

I just ran across this post while trying to convert a raster to an im object.

I'm sure the accepted post was correct at the time, but as of now geostatsp has deprecated as.im(). The function is now available in maptools as as.im.RasterLayer()

https://www.rdocumentation.org/packages/geostatsp/versions/1.5.4/topics/asImRaster

0

Something changed in the im class which depreciated a previous function of mine for coercing raster to im. Unfortunately, it is now a little bit more complicated than just rotating the matrix.

Here is a function will coerce a raster into an im class object.

raster.as.im <- function(im) {
  r <- raster::res(im)
  orig <- sp::bbox(im)[, 1] + 0.5 * r
  dm <- dim(im)[2:1]
    xx <- unname(orig[1] + cumsum(c(0, rep(r[1], dm[1] - 1))))
    yy <- unname(orig[2] + cumsum(c(0, rep(r[2], dm[2] - 1))))
  return(spatstat::im(matrix(raster::values(im), ncol = dm[1], 
         nrow = dm[2], byrow = TRUE)[dm[2]:1, ], 
         xcol = xx, yrow = yy))
}

Now, test the function

library(raster)
library(spatstat)
r <- raster(system.file("external/test.grd", package="raster"))
  plot(r)

r.im <- raster.as.im(r)
  class(r.im)
  plot(r.im)

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.