13

I was trying to make a raster image from an irregularly spaced point database. The data looks like-

> head(s100_ras)
         x       y         z
1 267573.9 2633781 213.29545
2 262224.4 2633781  69.78261
3 263742.7 2633781  51.21951
4 259328.4 2633781 301.98413
5 264109.8 2633781 141.72414
6 255094.8 2633781  88.90244

I want these 'z' values within a mesh which I created by

# Create a fine mesh grid
my_mesh=expand.grid(seq(min(s100_ras$Y),max(s100_ras$Y),l=100),
                    seq(min(s100_ras$X),max(s100_ras$X),l=100))

I also want the z-values to be assigned as 'NA' for those mesh points that are outside the data points. The points over the mesh looks like this: https://drive.google.com/file/d/0B6GUNg-8d30vYzlwTkhvaHBFTnc/edit?usp=sharing when I plot

plot(my_mesh)
points(s100_ras$Y, s100_ras$X, pch="*", col='blue')

The problem is that I'm not sure how to build on this, the following steps don't work because my mesh grid and data points are not of the same scale!!

library(rgdal)
library(raster)
xyz<-cbind(my_mesh, s100_ras)
r <- rasterFromXYZ(xyz)
image(r)

If I try to make a raster just by using the data points (without any mesh), R throws an error since my data is irregularly spaced!

library(sp)
s100_ras <- data.frame(expand.grid(x = s100_ras$Y, y = s100_ras$X), 
                       z = as.vector(s100_ras$mean))
coordinates(s100_ras) <- ~x+y
proj4string(s100_ras) <- CRS("+proj=utm +zone=46 +datum=WGS84")
gridded(s100_ras) = TRUE

suggested tolerance minimum: 0.916421 
Error in points2grid(points, tolerance, round) : 
  dimension 1 : coordinate intervals are not constant

Moreover, I was trying to play with 'rasterize' function (for irregular grids) of 'raster package', but couldn't get a way with it :(. I know how to interpolate and make a regular grid, but for the sake of originality, I want to AVOID interpolation. Is it possible to make a raster of irregularly spaced data points without idw or kriging methods?

  • The problem with irregular spaced grids is that the algorithm fails if points lie too close/far together. A (non-optimal) workaround: Why not take the minimal distance between cells and create a rectangular vector grid. Then join the average values of the points to that grid and rasterize it. – Curlew Nov 30 '13 at 23:20
  • I had the same problem - the solution was to use SpatialPixelsDataFrame with the suggested tolerance argument (0.916421 in your case). – Tomas Nov 13 '14 at 10:48
18

I assume you want your irregular point data on a regular raster. In that case, rasterize should work, and the examples in ?rasterize show how. Here is something based on your data

s100 <- matrix(c(267573.9, 2633781, 213.29545, 262224.4, 2633781, 69.78261, 263742.7, 2633781, 51.21951, 259328.4, 2633781, 301.98413, 264109.8, 2633781, 141.72414, 255094.8, 2633781, 88.90244),  ncol=3,  byrow=TRUE)
colnames(s100) <- c('X', 'Y', 'Z')

library(raster)
# set up an 'empty' raster, here via an extent object derived from your data
e <- extent(s100[,1:2])
e <- e + 1000 # add this as all y's are the same

r <- raster(e, ncol=10, nrow=2)
# or r <- raster(xmn=, xmx=,  ...

# you need to provide a function 'fun' for when there are multiple points per cell
x <- rasterize(s100[, 1:2], r, s100[,3], fun=mean)
plot(x)
  • excellent solution to the problem@Robert – ToNoY Dec 1 '13 at 16:16
  • Sorry can you clarify why you add e <- e + 1000 ? – mmann1123 Sep 10 '15 at 3:00
  • @mman1123 That is just to make things work with these odd example data. All y coordinates are the same and hence the extent is zero in the y direction, so I add 1000 --- entirely arbitrary --- to be able to make a raster from the extent. – Robert Hijmans Sep 10 '15 at 3:12
  • Is there a pythonic solution for this? – raaj 2 days ago
5

This worked for me - the solution was to use SpatialPixelsDataFrame with the suggested tolerance argument (0.916421 in your case):

points <- SpatialPoints(s100_ras[,c('x','y')], s100_ras[,c('z')])
pixels <- SpatialPixelsDataFrame(points, tolerance = 0.916421, points@data)
raster <- raster(pixels[,'z'])

though, due to the high tolerance value, the raster doesn't fit very well to the original points. Could be fit much better.

  • 1
    The first line of code has a closing bracket missing? – Antti Aug 21 at 8:57
  • @Antti thanks, corrected! – Tomas Aug 21 at 14:16

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.