[R-sig-Geo] Problem overlying Kriging grid with measurement points
Mauricio Zambrano
hzambran.newsgroups at gmail.com
Wed May 28 16:21:18 CEST 2008
Dear List,
I carried out a IDW interpolation and I got a plot of the
interpolations values over the grid that I defined.
However, for a better interpretation of the results, I would like to
plot, over the Kriging grid, the points that were used as data for the
interpolation, which are stored in a shapefile.
For doing this, I tried with:
spplot(pp.idw["var1.pred"], main="Daily Mean PP, [mm]", sp.layout =
list("SpatialPointsDataFrame", meteo_pt[,"PP_DAILY_MEAN_MM"] ) )
but I only get the plot of the interpolated value over the grid, and
the plot shows the following error message:
"error using packet 1"
"Missing 'data' argument, no omission value"
Could you give some hint about what I'm doing wrong ?.
Additional information is:
class(pp.idw)
[1] "SpatialPixelsDataFrame"
attr(,"package")
[1] "sp"
class(meteo_pt)
[1] "SpatialPointsDataFrame"
attr(,"package")
[1] "sp"
And complete tried procedure:
p4s <- CRS("+proj=utm +zone=30 +ellps=intl +units=m +no_defs")
meteo <- read.csv2("Meteorological_by_basin.csv")
# Selecting only those stations that have a NON "NA" value on the
field "PP_DAILY_MEAN_MM":
meteo_nna <- meteo[!is.na(meteo[,"PP_DAILY_MEAN_MM"]),]
# Settting the COORDINATES
coordinates(meteo_nna) <- ~EAST_ED50 + NORTH_ED50
# Projecting the coordinates of the meteorological stations into the
right system
proj4string(meteo_nna) = p4s
#proj4string(meteo_nna) <- CRS("+init=epsg:28992") ; # another way to
set up the coordinates
# Reading the boundary of the catchment
library(maptools)
catchment <- readShapePoly("only_catchment.shp", proj4string=p4s)
# Defining a sampling GRID of 1000mx1000m with regular spacing
# For avoiding a random grid (sampled randomly from the first cell),
and getting a fixed,
# and reproducable grid, it is neccessary to add the argument "offset
= c(0.5, 0.5)"
catchment.grid <- spsample(catchment, type="regular", cellsize=1000,
offset = c(0.5, 0.5))
# Makin possible that the grid can be used in the interpolations:
gridded(catchment.grid) <- TRUE
# Plotting the grid
spplot(catchment, sp.layout = list("sp.points", catchment.grid))
library(gstat)
# Interpolating with the INVERSE DISTANCE WEIGHTED
pp.idw <- krige(PP_DAILY_MEAN_MM~1, meteo_nna, catchment.grid)
# Plotting the interpolated values
spplot(pp.idw["var1.pred"], main="Daily Mean PP, [mm]")
# Reading the Meteorologicla Stations whitin the cachment
meteo_pt <- readShapePoints("meteo_by_basin.shp", proj4string=p4s)
#Attempting to get the overlay:
spplot(pp.idw["var1.pred"], main="Daily Mean PP, [mm]", sp.layout =
list("SpatialPointsDataFrame", meteo_pt[,"PP_DAILY_MEAN_MM"] ) )
Thanks in advance,
Mauricio
--
Linux user #454569 -- Ubuntu user #17469
More information about the R-sig-Geo
mailing list