setwd("J:\Henitz2006\Escritório\Contacto\Região Central\Shape") getwd() require(shapefiles) args(read.shapefile) RegiaoCentral <- read.shapefile("Pesquisa_Fin5_point") Border <- read.shapefile("Border_region") RegiaoCentral Border args(RegiaoCentral) require(geoR) require(MASS) regCenDF <- data.frame(RegiaoCentral$shp$shp[, 2:3],RegiaoCentral$dbf$dbf[5]) Bord <- data.frame(Border$shp$shp[[1]]$points[,1:2]) regCenDF RegCen <- as.geodata(regCenDF,coords.col=1:2,data.col=3) plot(RegCen) plot (RegCen,bord = Bord) help(likfit) ##Elias Krainski v <- variog(RegCen, max.dist=1200) plot(v) mod <- likfit(RegCen, ini.cov.pars=c(3E+05, 500)) lines(mod, col=2) apply(RegCen$coo, 2, range) loc <- expand.grid(seq(331520, 333160, leng=50), seq(7395512, 7396890, leng=50)) krig.kc <- krige.control(obj.mod=mod) Bord krig <- krige.conv(RegCen, locat=loc, border=Bord,krige=krig.kc) image(krig) ### salvando shape com resultados ### pontos da malha (loc) que pertencem o poligono (Border$...) require(splancs) id <- which(inout(loc, Bord)) ### salvando require(maptools) write.pointShape(as.matrix(loc[id, ]), data.frame(VUNIT=krig$pred), "Krig7") lido <- read.shape("Krig7") plot(lido)