[R-sig-Geo] extract point values from contour lines
Barry Rowlingson
b.rowlingson at lancaster.ac.uk
Wed Sep 4 13:06:39 CEST 2013
On Wed, Sep 4, 2013 at 11:12 AM, Barry Rowlingson
<b.rowlingson at lancaster.ac.uk> wrote:
> One approach might be to consider the points of the contour line as
> samples from the underlying data and do a smoothing such as kriging or
> inverse distance weighting, but I suspect that having a lot of samples
> with the same height at a range of distances is going to mess up the
> variogram somewhat...
x <- 10*1:nrow(volcano)
y <- 10*1:ncol(volcano)
image(x,y,volcano)
lines = contourLines(x, y, volcano)
Following on from that, the Kriging approach doesn't look so bad:
require(plyr)
# convert output of contourLines to a spatial points data frame with
height in 'h':
xyhd=ldply(l,function(e){data.frame(x=e$x,y=e$y,h=e$level)})
coordinates(xyhd)=~x+y
require(automap)
a = autoKrige(h~1,xyhd)
plot(a)
that shows you how the points of the contours have been interpolated
over the space. If we want to just get predictions at locations we
construct a new points data set and feed that in:
pnew = data.frame(x=runif(100,100,800),y=runif(100,50,500))
coordinates(pnew)=~x+y
anew = autoKrige(h~1,xyhd,pnew)
plot(anew)
and you get point estimates and kriging standard errors - these are in
anew$krige_output as a spatial points data frame.
Now I suppose I should compare those estimates with direct samples
from the original volcano DEM to see how well we've done.... Left as
an exercise...
There seems to be a bit of literature on interpolating from contour
lines via various other methods too.
Barry
More information about the R-sig-Geo
mailing list