[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