2015-03-25 87 views
1

我正在嘗試在R中使用Gstat軟件包進行通用剋剋裏格處理。我有一個腳本,我被幫助,但現在我卡住了,不能要求原始來源的援助。 問題是我無法更改cokriged數據的輸出分辨率。我想將插入的地圖導入ArcMap,點對柵格以非常低的分辨率離開我。R:Gstat通用剋剋裏格分辨率

我的腳本如下:

library(raster) 
library(gstat) 
library(sp) 
library(rgdal) 
library(FitAR) 

加載我的數據集,即containes座標和採樣值:

kova<-read.table("katvus_point_modif3.txt",sep=" ",header=T) 
coordinates(kova)=~POINT_X+POINT_Y 

在同一座標先前加載的深度值,這是我的協:

Sygavus<-read.table("sygavus_point_cokrig.txt",sep=" ",header=T) 
coordinates(Sygavus)=~POINT_X+POINT_Y 

overlay <- over(kova,Sygavus) 
kova$Sygavus <- overlay$Sygavus 

這應該是爲我的插值設置邊界,該文件是從ArcMap中導出的Shape文件:

border <- shapefile("area_2014.shp") 
projection(kova)=projection(border) 

據說這是爲了創造協同克里金法網格和RES =應該讓我指定我所要的輸出是,但不管我用什麼號碼什麼的分辨率輸出不會改變

grid <- spsample(border,type="regular",res=25) 

我刪除overlaping點:

zero <- zerodist(kova) 
kova <- kova[-zero[,2],]    

我在深度協光柵文件加載。這是從ArcMap導出爲ascii格式的深度柵格:

depth <- raster("htp_depth_covar.asc") 
projection(depth)=projection(border) 

overlay <- extract(depth,kova) 
kova$depth <- overlay 

我刪除了na!從疊加式的深度值的值(這些值應該是相同的,在各自座標上以前加載的深度協表,但是如果我離開的那部分,腳本停止運行

kova <- kova[!is.na(kova$depth),] 


kova.gstat <- gstat(id="Kova",formula=kova~depth,data=kova) 
kova.gstat <- gstat(kova.gstat,id="Sygavus",formula=Sygavus~depth,data=kova) 

var.kova <- variogram(kova.gstat) 
plot(var.kova) 

kova.gstat <- gstat(kova.gstat,id=c("Kova","Sygavus"),model=vgm(psill=cov(kova$kova,kova$Sygavus),model="Mat",range=12000,nugget=0)) 
kova.gstat <- fit.lmc(var.kova,kova.gstat,model=vgm(psill=cov(kova$kova,kova$Sygavus),model="Mat",range=12000,nugget=0)) 

plot(var.kova,kova.gstat$model) 

overlay <- extract(depth,grid) 
grid <- as.data.frame(grid) 
grid$depth <- overlay 
coordinates(grid)=~x1+x2 
projection(grid)=projection(border) 

krige <- predict.gstat(kova.gstat,grid) 

spplot(krige,c("Kova.pred")) 

write.table(krige, "kova.raster1.ck.csv", sep=";", dec=",", row.names=F) 

在任何幫助瞭解gstat cokriging和腳本整體將不勝感激!

回答

0

因爲你不提供一個可重複的例子,我只能猜測,但我認爲spsample忽略res=25參數。改爲嘗試n=1000,然後增加該值以獲得更高的分辨率。