0

I'm having a problem exporting a density im from spatstat to a file format readable by ArcGIS. Here's my code

library(raster)
library(spatstat)
library(maptools)
library(sp)
# make a spatstat ppp with California boundary as window
ca <- readShapePoly("ca.shp")
o3 <- readShapePoints("o3.shp")
o3 <- as(o3, "SpatialPoints")
o3p <- as.ppp(o3)
o3p$window <- as.owin(ca)
# calculate density
d.o3p <- density.ppp(o3p)

which all works fine. But when I attempt to export to an ascii raster file

writeRaster(raster(d.o3p), filename="grid.asc", format="ascii", NAflag=-9999)

I get this problem

Error in .startAsciiWriting(x, filename, ...) : 
x has unequal horizontal and vertical resolutions. Such data cannot be  
stored in arc-ascii format

The data I'm using are for the state of California, and so the aspect ratio is not 1. So... how do I make the density im have square pixels?

1 Answers1

0

Do you have to use the ascii format for ArcGIS? It may be difficult to obtain perfectly square pixels with an arbitrary window. The pixel resolution is controlled by either the argument eps or the argument dimyx (not both), which are sent to as.mask and you can read more in that help file. Basically eps is the pixelsize, so the easiest is to set that to a given value. However, if the window size doesn't fit exactly with the window in one direction the pixels might be slightly non-square. E.g. this works perfectly:

X <- runifpoint(100, win = owin(c(0,9),c(0,10)))
d <- density.ppp(X, eps = 1)

But this gives non-square pixels:

d2 <- density.ppp(X, eps = 2)
Ege Rubak
  • 4,347
  • 1
  • 10
  • 18
  • Would upvote if I could (not enough reputation...) because this is helpful. Have also discovered since I posted that spatstat.options(npixel=c(nx,ny)) allows specification of a different image size to the default, which could enable writing a function to control this. –  Apr 03 '15 at 23:13