[STATSGRASS] problem with kriging in gstat( cooordinates) (long)

Jarek Jasiewicz jarekj at amu.edu.pl
Thu Aug 31 15:44:12 EDT 2006


I try to krige some data, I do following grass wiki tutorial (that with 
gaiciture)

1) I read data points as spatial grid data frame and g region
2) I make SpatialGridDataFrame according instruciton I I recived spatial 
data frame
3) I compared the created data frame which I read from grass (grid) soon 
after importing points the geographic coordinates seems look same

4) I made variogram and fit the model - all OK, model is fitted well

5) finally I try kriging:

krig=krige(lenght~1, locations=kr, newdata=sg, model=vgf)

where: lenght~1 the modeled data kr-spatialpointdataframe, 
sg-spatialgriddataframe generated in R, vgf - the fitted variogram

I I recive the error:

Error in predict.gstat(g, newdata = newdata, block = block, nsim = nsim,  :
    data items in gstat object and newdata have different coordinate 
reference systems

In fact, ONLY bounding box of grid and points deffer a little - is a 
litlle smaller (see at the end of the message) because points not cover 
all region

So, I try to change region definition to be exaclty the same with point 
layer and point bounding box read to R

I have recived the same error

I thing something is missing, but I completly don't know what!

I thing also (but only thing) If my problems are not connected with 
region problems mentioned (also by me) in grass 6.1

regards
Jarek Jasiewicz

PS1. I have no problem with tutorial data (meuse) according man pages, 
but I can't repead this on my data with spatial reference

PS2. I tried to krige the same data with sgeostat but I recive error 
that matrix (534x667) is too big ....hmmm... it looks like that sgeostat 
is for tutorial data only

PS3. here are strucutures of data frames used by me

kr - SpatialPointsDataFrame Readed from Grass (6.1)

Formal class 'SpatialPointsDataFrame' [package "sp"] with 5 slots
  ..@ data       :Formal class 'AttributeList' [package "sp"] with 1 slots
  .. .. ..@ att:List of 14
  .. .. .. ..$ cat   : num [1:492] 11 40 43 49 140 148 175 181 199 205 ...
// (JJ........ lot of data I cut it )
  .. .. .. ..$ lenght: num [1:492] 512 434 302 247 302 ...
  ..@ coords.nrs : num(0)
  ..@ coords     : num [1:492, 1:2] 3568825 3567206 3571253 3568405 
3571613 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : NULL
  .. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
  ..@ bbox       : num [1:2, 1:2] 3560820 5635012 3575900 5654902
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
  .. .. ..$ : chr [1:2] "min" "max"
  ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slots
  .. .. ..@ projargs: chr " +proj=stere +lat_0=51.67083333333333 
+lon_0=16.67222222222222 +k=0.999800 +x_0=3703000 +y_0=5627000 
+ellps=krass +units=m +no_"| __truncated__

sg - SpatialGridDataFrame created in R:

Formal class 'SpatialGridDataFrame' [package "sp"] with 6 slots
  ..@ data       :Formal class 'AttributeList' [package "sp"] with 1 slots
  .. .. ..@ att:List of 1
  .. .. .. ..$ k: num [1:356178] 2 2 2 2 2 2 2 2 2 2 ...
  ..@ grid       :Formal class 'GridTopology' [package "sp"] with 3 slots
  .. .. ..@ cellcentre.offset: num [1:2] 3560010 5635012
  .. .. ..@ cellsize         : num [1:2] 30 30
  .. .. ..@ cells.dim        : int [1:2] 534 667
  ..@ grid.index : int(0)
  ..@ coords     : num [1:2, 1:2] 3560010 3575990 5635012 5654992
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : NULL
  .. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
  ..@ bbox       : num [1:2, 1:2] 3559995 5634997 3576005 5655007
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
  .. .. ..$ : chr [1:2] "min" "max"
  ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slots
  .. .. ..@ projargs: chr " +proj=stere +lat_0=51.67083333333333 
+lon_0=16.67222222222222 +k=0.999800 +x_0=3703000 +y_0=5627000 +no_defs 
+a=6378245 +rf=29"| __truncated__

and for comparing: SpatialGridDataFrame readed from Grass

Formal class 'SpatialGridDataFrame' [package "sp"] with 6 slots
  ..@ data       :Formal class 'AttributeList' [package "sp"] with 1 slots
  .. .. ..@ att:List of 1
  .. .. .. ..$ dem.30: num [1:356178] 71.0 70.9 70.7 70.5 70.3 ...
  ..@ grid       :Formal class 'GridTopology' [package "sp"] with 3 slots
  .. .. ..@ cellcentre.offset: num [1:2] 3560010 5635012
  .. .. ..@ cellsize         : num [1:2] 30 30
  .. .. ..@ cells.dim        : int [1:2] 534 667
  ..@ grid.index : int(0)
  ..@ coords     : num [1:2, 1:2] 3560010 3575990 5635012 5654980
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : NULL
  .. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
  ..@ bbox       : num [1:2, 1:2] 3559995 5634997 3576005 5654995
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
  .. .. ..$ : chr [1:2] "min" "max"
  ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slots
  .. .. ..@ projargs: chr " +proj=stere +lat_0=51.67083333333333 
+lon_0=16.67222222222222 +k=0.999800 +x_0=3703000 +y_0=5627000 +no_defs 
+a=6378245 +rf=29"| __truncated__

So I thing spatial grid is genreated well

 




More information about the grass-stats mailing list