n,z,x,y
上のプロットデータをRで取り込んで内挿格子データを出力する「press.R」
library(gstat)
d <- read.csv("press.csv")
dd <- d[c(3,4,2)]
coordinates(dd) <- ~ x+y
bubble(dd, zcol='z', fill=FALSE, do.sqrt=FALSE, maxsize=2)
x.range <- as.integer(range(dd@coords[,1]))
y.range <- as.integer(range(dd@coords[,2]))
grd <- expand.grid(x=seq(from=0, to=600, by=10), y=seq(from=0, to=580, by=10) )
coordinates(grd) <- ~ x+y
gridded(grd) <- TRUE
## plot(grd, cex=0.5)
## points(dd, pch=1, col='red', cex=0.7)
g <- gstat(id="press", formula=z ~ 1, data=dd)
plot(variogram(g, map=TRUE, cutoff=400, width=20), threshold=10)
v <- variogram(g, alpha=c(0,45,90,135))
v.fit <- fit.variogram(v, model=vgm(model='Lin' , anis=c(0, 0.5)))
plot(v, model=v.fit, as.table=TRUE)
g <- gstat(g, id="press", model=v.fit )
p <- predict(g, model=v.fit, newdata=grd)
pts <- list("sp.points", dd, pch = 4, col = "black", cex=0.5)
spplot(p, zcol="press.pred", col.regions=terrain.colors(20), cuts=19, sp.layout=list(pts), contour=TRUE, labels=FALSE, pretty=TRUE, col='brown', main='OK Prediction')
n <- 1
for(i in 1:59){
cat("[")
for(j in 1:61){
cat(p@data$press.pred[n],",",sep="")
n <- n + 1
}
cat("],")
}
