Hello experts,
Here's another novice to R. Hope you can help me out. Apologies in
advance if questions have been asked before (in which case, could you
tell me the thread).
My dataset is a table with daily Precipitation for the month January in
2000, measured at more than 300 weather stations. Each row is a weather
station and there are 31 columns (or should I say attribute), each day
it's own column.
I use a for-loop to go through the columns (the days) and each loop
performs an idw and autokrige interpolation, plots images and saves them
as jpg to hard disk.
Here's the script:
Pd200001 <- read.table(paste(wordir, "/Input/", sep=""), sep=("\t"),
header = TRUE)
coordinates(Pd200001) = ~RD_LOCATION_X + RD_LOCATION_Y
nl.grd <- read.asciigrid(paste(wordir, "/Input/nl_clip_1km.asc",
sep=("")), as.image=FALSE, plot.image=TRUE)
for (i in 1:31) {
kol <- paste(names(Pd200001[i+2]), "~1", sep="") # first two columns
are id and stationnr
# INVERSE DISTANCE
Pd.idw <- krige(as.formula(kol), Pd200001, nl.grd)
image(Pd.idw["var1.pred"], col = topo.colors(20))
title(paste("Pd_idw", 20000100+i, sep = ""))
dev.copy(jpeg, file = paste( wordir, "/img/Pd_idw_", 20000100+i, ".jpg",
sep = ""), height = 17, width = 17, units = "cm", quality = 100,
res=150, bg="white") ; dev.off()
rm(Pd.idw)
# ORDINARY KRIGING
Pd.kri <- autoKrige(as.formula(kol), Pd200001, nl.grd)
plot(Pd.kri)
dev.copy(jpeg, file = paste( wordir, "/img/Pd_kri_", 20000100+i, ".jpg",
sep = ""), height = 17, width = 17, units = "cm", quality = 100,
res=150, bg="white") ; dev.off()
rm(Pd.idw)
}
The loop works and I get images in files but I also get Warning
Messages, both after interpolation and after plot().
Pd.idw <- krige(as.formula(kol), Pd200001, nl.grd)
[inverse distance weighted interpolation]
Warning message:
In points2grid(points, tolerance, round, fuzz.tol) :
grid has empty column/rows in dimension 2
spplot(Pd.idw, "var1.pred")
Warning message:
In data.row.names(row.names, rowsi, i) :
some row.names duplicated: 2,3,4,5,6,7,8,9,10,
Question 1
How "bad" are the Warning Messages, how should I prevent them
Question 2
I use 2 different plotting techniques, which makes the images hard to
compare. At this stage I do not need real map data for further analysis,
just some jpegs to quickly compare the results. How can I plot the idw's
and de krige results the same way? They are of different class.
Question 3
automapPlot(plot_data, zcol, col.regions, ...)
The ... is intriguing: "Arguments that are passed on to spplot." This is
also the case for many R functions, which makes R very flexible, but for
a beginner this seems endlessly layered. For instance, I would like to
create my own colour ramp, in either image(), plot() or spplot(). But I
can't figure it out. Is it also possible to reverse the standard ramps
heat.colors, topo.colors etc?
Question 4
After the jpg is written to HD, I remove the objects Pd.idw and Pd.kri.
Need I do this, or is the object overwritten anyway, at the start of the
next loop? Also, is there another method to deal with the results, for
instance can I create an object Pd.idw.many which contains all of the 31
interpolated maps? This would be convenient for further R processing.
Many thanks in advance!
[[alternative HTML version deleted]]
_______________________________________________
R-sig-Geo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
