[STATSGRASS] problems with kriging latlong areas

Carlos "Guâno" Grohmann carlos.grohmann at gmail.com
Sat Jan 6 12:03:14 EST 2007


First of all, Happy 2007 for everyone.

Well, gstat 0.9-35 came out, and I went to try it, since I want to do
some kriging interpolation on some latlong areas. I did exactly what I
was doing before with utm areas. The variogram modeling went well
(er.. not so well as I expected, see below for more on this), but when
it came to the kriging step, I get an erros msg about dynamic
memory... the area is quite small (about 15x15km, 167x167 rows/cols),
and I have 1.5Gb of ram, so maybe we have a bug here?

this is the hole thing:

> srtm<-readGDAL("/home/guano/artigos/DEM_krigagem/Los_Angeles/srtm_LA_90.tif");
/home/guano/artigos/DEM_krigagem/Los_Angeles/srtm_LA_90.tif has GDAL
driver GTiff
and has 167 rows and 167 columns
Closing GDAL dataset handle 0x191a9d0...  destroyed ... done.
>
>  summary(srtm)
Object of class SpatialGridDataFrame
Coordinates:
        min        max
x -118.77208 -118.63292
y   34.83958   34.97875
Is projected: FALSE
proj4string : [+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs]
Number of points: 2
Grid attributes:
 cellcentre.offset     cellsize cells.dim
x         -118.7717 0.0008333333       167
y           34.8400 0.0008333333       167
Data attributes:
    band1
 Min.   : 622
 1st Qu.:1135
 Median :1299
 Mean   :1312
 3rd Qu.:1478
 Max.   :1993
>
>
> offset<- c(srtm at grid@cellcentre.offset[1]+ (srtm at grid@cellsize[1]/3), srtm at grid@cellcentre.offset[2]+ (srtm at grid@cellsize[2]/3));
> size<- c((srtm at grid@cellsize[1]/3), (srtm at grid@cellsize[2]/3));
> dim<- c((srtm at grid@cells.dim[1]*3), (srtm at grid@cells.dim[2]*3));
> grd <- GridTopology(offset,size,dim);
> mask_SG <- SpatialGridDataFrame(grd, data=list(k=rep(1, dim[1]*dim[2])), proj4string=CRS(srtm at proj4string@projargs));
Closing GDAL dataset handle (nil)... done.
>
> coords<-coordinates(srtm);
> jcoords <- cbind(jitter(coords[,1]), jitter(coords[,2]));
> band1<-as.data.frame(srtm$band1);
> srtm2<-SpatialPointsDataFrame(jcoords,band1,proj4string=CRS(srtm at proj4string@projargs));
>
> variog1 <-variogram(srtm$band1~coords[,1]+coords[,2], loc=srtm, srtm);
>
> plot(variog1);
>
> plot(variog1, xlim=c(0,1), ylim=c(0,10000));
>
> vrg.eye<-(vgm(psill=9500, model="Gau", range=0.55, nugget=10));
>
> plot(variog1, model=vrg.eye, xlim=c(0,1), ylim=c(0,10000));
>
> OK_pred_n8 <- krige(srtm$band1 ~ 1, loc=srtm2, newdata=mask_SG, model=vrg.eye, nmax=8);
[using ordinary kriging]
Error in predict.gstat(g, newdata = newdata, block = block, nsim = nsim,  :
       gstat: out of dynamic memory
>


That's for kriging. about the variogram, if I just call variogram()
without specifing a width value (as in above), I get a normal
variogram. but if I call variogram() with a width=0.0004 (half cell,
as I do with utm regions), I get a cloud of points. Is this normal? I
don't think so. (see image in the link)

http://www.igc.usp.br/pessoais/guano/temp/variogram_cloud.jpg



Hope to hear from you guys

Carlos

-- 
+-----------------------------------------------------------+
              Carlos Henrique Grohmann - Guano
  Geologist M.Sc  - Doctorate Student at IGc-USP - Brazil
Linux User #89721  - carlos dot grohmann at gmail dot com
+-----------------------------------------------------------+
_________________
"Good morning, doctors. I have taken the liberty of removing Windows
95 from my hard drive."
--The winning entry in a "What were HAL's first words" contest judged
by 2001: A SPACE ODYSSEY creator Arthur C. Clarke




More information about the grass-stats mailing list