插值到栅格

Question

我有一个包含 LAT、LON 和温度数据的三列数据集 我想生成一个光栅图像，根据 24 个数据记录器数据点预测景观的温度。数据集可以在这里访问：DATA

这是我到目前为止所尝试的：

    #Lets try to interpolate the data onto a raster 
library (raster)
library (gstat)
library (sp)
#Temp and XY data
temp<-read.csv ('test_temp.csv')
#create a blank raster to the extent of the system
r<- raster (nrows=300, ncols=100, xmn=-84.95, xmx=-84.936, ymn=45.7, ymx=45.74)

#build a prediction model temp ~ LAT*LON
loc<-temp [,c(3,5,6)]
loc<-na.omit (loc)
#TPS model 
tps<-Tps(loc, loc$TemperatureC)
#gstat model
mod1<-gstat (data=temp, formula=TemperatureC ~ 1, locations=loc )
summary (mod1)

r2<-interpolate(r, model=tps)
    Error in scale.default(x, xc, xs) : 
  length of 'center' must equal the number of columns of 'x'

r2<-interpolate(r, model=mod1)
Error in bbox(dataLst[[1]]$data) : object not a >= 2-column array

Ultimatley 我想创建一系列插值数据的栅格，以显示一天中不同时间的温度变化。关于如何做到这一点的任何想法？

Answer 1

结合使用代码和示例，我想出了一些对我来说看起来正确的东西：

library(ggplot2) # start needed libraries
library(gstat)
library(sp)
library(maptools)

my_location<-c(lon= -84.94, lat = 45.72)
basemap <- get_map(location = my_location, zoom = 14, maptype = "hybrid")
plot (basemap)
#load data 
field_data<-read.csv ("test_temp.csv")
test_temp <- field_data # duplicate air temp. data file
test_temp$x <- test_temp$LON # define x & y as longitude and latitude
test_temp$y <- test_temp$LAT
test_temp<-na.omit(test_temp)
coordinates(test_temp) = ~x + y #set spatial coordinates to create a Spatial object

plot(test_temp)
x.range <- as.numeric(c(-84.94, -84.935)) # map extent
y.range <- as.numeric(c(45.705, 45.735))
grd <- expand.grid(x = seq(from = x.range[1], to = x.range[2], by = 0.001), 
                   y = seq(from = y.range[1], to = y.range[2], by = 0.001)) # expand points to grid
coordinates(grd) <- ~x + y 
gridded(grd) <- TRUE
plot(grd, cex = 1.5)
points(test_temp, pch = 1, col = "red", cex = 1)
idw <- idw(formula = TemperatureC ~ 1, locations = test_temp, newdata = grd) # apply idw model for the data
idw.output = as.data.frame(idw)
names(idw.output)[1:3] <- c("long", "lat", "var1.pred") # give names to the modelled variables
# plot results
ggplot() +
  geom_tile(data = idw.output, alpha=0.5, aes(x = long, y = lat, fill=var1.pred, 0)) +
  geom_point(data=field_data, aes(x=LON, y=LAT), shape=21, colour="red") + 
  scale_fill_gradient(low="cyan", high="orange") + theme_bw() 

ggmap(basemap, extent = "device") + 
  geom_tile(data = idw.output, alpha=0.5, aes(x = long, y = lat, fill=var1.pred, 0)) +
  geom_point(data=field_data, aes(x=LON, y=LAT), shape=16, size=.5, colour="red") + 
  scale_fill_gradient(low="cyan", high="red") + theme_bw()