[GRASSLIST:10001] Re: R, kriging and grass6

[Kirk R. Wythers]

Thanks Roger. That was a huge help. I appreciate your taking the time.

Kirk
On Jan 25, 2006, at 4:43 AM, Roger Bivand wrote:

> On Mon, 23 Jan 2006, Kirk R. Wythers wrote:
>
>> Roger,
>>
>> Thanks for putting this together. I do have a very basic question
>> about getting the vector points from a GRASS vector to R. Can you
>> roll this description back a couple of steps and start with a GRASS
>> vector? Is this example based on the maas data?
>
> Yes, this uses the data in the sp package for R:
>
> library(sp)
> library(maptools)
> data(meuse)
> coordinates(meuse) = c("x", "y")
> writePointsShape(meuse, "meuse")
> data(meuse.grid)
> coordinates(meuse.grid) = c("x", "y")
> meuse.grid_sp <- SpatialPixelsDataFrame(meuse.grid, data=as 
> (meuse.grid,
>   "data.frame"))
> writeAsciiGrid(meuse.grid_sp, "soil", fn="soil.aai")
> meuse.grid_sp$ffreq1 <- as.integer(meuse.grid_sp$ffreq)
> writeAsciiGrid(meuse.grid_sp, "ffreq1", fn="ffreq1.aai")
>
> then in GRASS in an arbitrary location:
>
> r.in.gdal input=soil.aai output=soil location=meuse
>
> and moving to location meuse (more or less):
>
> g.mapset PERMANENT
> g.proj -c proj4='+proj=stere +lat_0=52.15616055555555
>   +lon_0=5.38763888888889 +k=0.999908 +x_0=155000 +y_0=463000
>   +ellps=bessel +units=m'
> # watch line breaks
> g.mapset rsb
> # some fiddling to get the projection data right
> # then some more to get ewres and nsres to 40
> r.in.gdal -o input=ffreq1.aai output=ffreq1
> v.in.ogr -o dsn=. layer=meuse output=meuse
>
> then R and the kriging code.
>
> Hope this helps,
>
> Roger