0

I have a data frame which contains point daily precipitation for 4 station for 2 years. I want to interpolate to 50m resoulution and write them in to 2 raster images. I used following code to achieve this...

library(ggplot2)
library(gstat)
library(maptools)
library(raster)
library(rgdal)
xcord<-c(100,200,300,400)
ycord<-c(100,200,300,400)
value1<-c(1,2,3,1)
value2<-c(2,5,7,3)
datas<-data.frame(xcord,ycord,value1,value2)
coordinates(datas) = ~xcord+ycord
mesh <- expand.grid(x=seq(0,500,50),y=seq(0,500,50))
coordinates(mesh)=~x+y
gridded(mesh) <- TRUE
oneidw = idw(value1~1,datas,mesh)
spplot(oneidw["var1.pred"], main = " inverse distance weighted    interpolations")

It worked. but i want to apply a loop to do it for another variable value2 (and so on...) without doing it manually. and i used this

    sym<-paste("value", 1:2,sep="")
    variable=as.vector(print(sym,quote=FALSE))
    for (i in 3:ncol(datas)){
     one<-idw((print(variable[i],quote=FALSE))~1,datas,mesh)
    }

but i got error too many spatial dimensions........ can anybody help me with this....

Frank
  • 66,179
  • 8
  • 96
  • 180
user3978632
  • 283
  • 4
  • 17
  • 1
    Can you make your example reproducible? Just construct a small example (with the `lapply` call) that demonstrates the error you're getting. – Roman Luštrik Mar 25 '15 at 11:23
  • I tried as you said with `list.idw <- lapply(1:2, function(n) idw(print(variable[n],quote=FALSE)~1, datas,mesh))` but still it didn't work and displays too many spatial dimesions......Please help me with this – user3978632 Apr 06 '15 at 03:34

1 Answers1

0

I'm not too familiar with spplot, but this worked for me using ggplot. 

library(ggplot2)
library(gstat)
library(sp)
library(maptools)
library(maps)
library(dplyr)
library(rgdal)

xcord<-c(100,200,300,400)
ycord<-c(100,200,300,400)
value1<-c(1,2,3,1)
value2<-c(2,5,7,3)
datas<-data.frame(xcord,ycord,value1,value2)
new_datas <- select(datas, xcord, ycord)
parse_by <- colnames(datas)[3:4] #change according to designated value columns 

for ( i in parse_by ) { 
  variable <- datas[i]
  new_datas2 <- cbind(new_datas, variable) #combine single variable col w/ coordinates
  colnames(new_datas2)[3] = "variable" #rename so that you can call to in idw formula
  coordinates(new_datas2) = ~xcord+ycord
  mesh <- expand.grid(x=seq(0,500,50),y=seq(0,500,50))
  coordinates(mesh)=~x+y
  gridded(mesh) = TRUE
  plot(mesh) #plot background so ggplot can use later
  points(new_datas2) #points for ggplot to use later
  one<-idw(formula = variable~1, locations = new_datas2, newdata = mesh) #idw formula
  one.output <- as.data.frame(one)
  names(one.output)[1:3] <- c("xcord", "ycord", "var1.pred") #rename for ggplot geom_tile
  ggplot() + geom_tile(data = one.output, alpha = 1, aes(x = xcord, y = ycord, fill = var1.pred)) + 
    geom_point(data = new_datas, aes(x = xcord, y = ycord)) + 
    labs(fill = "inverse distance weighted interpolations") 
  ggsave(paste(i,".png",sep="")) #save as .png according to variable name
}
Echelle Burns
  • 1