Carlos,
Thank you! That helps me as well. I got to the point following Roger's example using my own
dataset, but:
OK_pred <- krige(IMMERSIONE~ 1,locations=giaciture,newdata=mask_SG,model=efitted)
returns a DataFrame, not a SpatialGridDataFrame! So, writeRast6sp() failed for me. I could
use R's image()
command to see the kriged surface in R, but had not seen how to get it into GRASS.
Interestingly, my equivalent
of the line:
cvgm <- variogram(IMMERSIONE~1,locations=giaciture,width=400,cutoff=4000)
would not work. I had to use:
cvgm <- variogram(giaciture$IMMERSIONE~1,locations=giaciture,width=400,cutoff=4000)
Regards,
Tom
----- Original Message -----
From: "Carlos \"Guâno\" Grohmann" <carlos.grohmann@gmail.com>
Date: Wednesday, February 22, 2006 5:00 pm
Subject: [GRASSLIST:10527] Re: help with R and GRASS6
Well, I managed to get the residuals of the trend-surface analysis
in a
Spatial data Frame.Here's what I did:
library(spgrass6) ; library(spatial);
srtm <- getSites6sp("toposrtm");
coords<-coordinates(srtm);
class(srtm);
G <- gmeta6();
grd <- GridTopology(cellcentre.offset=c(G$west+(G$ewres/2),
G$south+(G$nsres/2)), cellsize=c(G$ewres, G$nsres),
cells.dim=c(G$cols,G$rows));
mask_SG <- SpatialGridDataFrame(grd, data=list(k=rep(1,
G$cols*G$rows)),proj4string=CRS(G$proj4));
class(mask_SG);trend <- surf.ls(1, x=coords[,1], y=coords[,2], z=srtm$cat);
resid1<-residuals(trend);
resid2<-as.data.frame(resid1);
resid.trend<-SpatialPointsDataFrame(coords,resid2);
resid.img <- SpatialGridDataFrame(grd, resid.trend);Now I'll start playing with variograms and kriggind.
Thanks all!
Carlos
On 2/21/06, Roger Bivand <Roger.Bivand@nhh.no> wrote:
>
> Hi Carlos,
>
> On Tue, 21 Feb 2006, Carlos "Guâno" Grohmann wrote:
>
> > Hello list
> >
> > I've away from R a while, and now I'm struggling with it...
> > I read Bivand's artile in GRASS-News3, but I dindn't understant
exactly> how
> > the data is stored in R.
> >
>
> Please see a note in R News in November 2005 by Edzer Pebesma and
me.>
> > If I want to make a trend-surface analysis with surf.ls(), I
need the
> x,y
> > and z coordinates of each point (or pixel), but I don't know
how to get
> > these values
>
> The object returned by getSites6sp() is a SpatialPointsDataFrame.
You can
> get the coordinates by saying:
>
> pointlayer <- getSites6sp("mypoints")
> coords <- coordinates(pointlayer)
>
> and assuming z is called z:
>
> myTSA <- surf.ls(2, x=coords[,1], y=coords[,2], z=pointlayer$z)
>
> but you can also do TSA straight away in gstat:
>
> myTSA_with_gstat <- gstat(id="tr2", formula = z ~ 1, data =
pointlayer,> degree=2)
>
> (the most recent gstat just takes the location coordinates from the
> SpatialPointsDataFrame object.)
>
> There was a thread about this recently, I think copied to:
>
>
> http://grass.gdf-
hannover.de/twiki/pub/GRASS/WebHome/krigingusingRGRASS6.txt>
> where a grid is created from the GRASS region settings to predict
to.>
> Hope this helps,
>
> Roger
>
> >
> > the whole idea is to make TSA, then variogram fitting in the
residuals,> then
> > krigging (maybe with gstat or geoR, still haven't decided)
> >
> > thanks
> >
> > --
> > +-----------------------------------------------------------+
> > 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
> >
>
> --
> Roger Bivand
> Economic Geography Section, Department of Economics, Norwegian
School of
> Economics and Business Administration, Helleveien 30, N-5045 Bergen,
> Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
> e-mail: Roger.Bivand@nhh.no
>
>--
+-----------------------------------------------------------+
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