On Mon, 6 Apr 2009, Greg King wrote:
Hi,
I'm fairly new to R, but finding the experience a good one. However I am a
little overwhelmed by the number of packages and not sure I am necessarily
using the most appropriate ones.
Have you read the Spatial Task View on your nearest CRAN? There you should
find some information to help. Firstly, akima treats all coordinates as
planar, so don't use interp(). Both fields and gstat can interpolate using
geographical coordinates, gstat with IDW and kriging, fields with thin
plate splines (I believe using rdist.earth, but have not checked) and
kriging. You are trying to interpolate, so why not use a geostatistical
method? You could see whether the new automap package (running gstat
internally) can handle geographical coordinates - check by comparing the
output with runs with the proj4string set and unset.
Hope this helps,
Roger
Here is the background to what I am trying to achieve: I have a CSV file
which contains weather forecast data for latitude and longitude points (see
attached out.csv, the data I understand is on WGS84
http://www.ready.noaa.gov/faq/geodatums.html). The sample points are at
half degree intervals. My objective is to work out what the forecast data
is at any specific given latitude/longitude by interpolating data from the
0.5x0.5 degree grid. I am doing this for a number of different time points
using the following functions:
library(akima)
library(rgdal)
gks <- function(inLat,inLon,dframe,variab) {
wind.grid <- expand.grid(lat=inLat, lon=inLon)
coordinates(wind.grid) <- ~ lat+lon
proj4string(wind.grid) <- CRS("+init=epsg:4326")
pnt<-interpp(coordinates(dframe)[,1], coordinates(dframe)[,2],
z=as.data.frame(dframe)[,1],
coordinates(wind.grid)[,1],coordinates(wind.grid)[,2],linear=TRUE)
return(pnt$z)
}
interp_gk <- function(lat, lon) {
wind<-read.csv(file="/Users/greg/Out.csv",head=TRUE,sep=",",
colClasses=c("POSIXct","numeric","numeric","numeric","numeric"))
coordinates(wind)=~lat+lon
proj4string(wind) <- CRS("+init=epsg:4326")
times<-unique(wind$vt)
columns<-names(wind)[2:length(names(wind))]
dOut<-data.frame(dateTime=times[1])
for (i in 1:length(times)) {
dOut[i,"dateTime"]<-times[i]
for (j in 1:length(columns)) {
sset<-subset(wind, wind$vt==times[i], select=c(columns[j]))
dOut[i,columns[j]]<-gks(lat,lon,sset,columns[j])
}
}
dOut<-cbind(dOut, mag=sqrt(dOut$ugrd_0^2+dOut$vgrd_0^2))
return(dOut)
}
However, I have the following concerns:
- Should I really be using akima? The documentation states it is not for
use on regularly spaced grids - what are my alternatives?
- The interp funcrtion will not work for cubic spline "linear=FALSE"
interpolation (is it because my data is regularly gridded?). How can I
achieve cubic spline interpolation?
- Is my function really using the Cordinate reference system specified?
When I comment out the CRS lines, my functions return the same values?
Lots of questions I appreciate, but I am curious! It seems R can achieve
what I am trying to do... but I may just be missing some vital bits of
information.
Thanks,
Greg
--
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.biv...@nhh.no
_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
