[R-sig-Geo] Spatial Interpolation of Regularly Gridded Data

[Miha Staut]

Dear Greg,

A while back I did a similar task with the function interp.surface() from the fileds package. It briefly does bilinear interpolation from a regular grid to arbitrary point. It is computationally quite inexpensive since the method of interpolation is quite simple. 

Hope I understood you correctly. 

All best
Miha


--- On Tue, 7/4/09, Greg King <gregspam at kooji.com> wrote:

> From: Greg King <gregspam at kooji.com>
> Subject: Re: [R-sig-Geo] Spatial Interpolation of Regularly Gridded Data
> To: Roger.Bivand at nhh.no
> Cc: r-sig-geo at stat.math.ethz.ch
> Date: Tuesday, 7 April, 2009, 10:32 PM
> Thanks for your reply Roger.  I
> little while back I got kriging working in
> the gstat package with my data, but found that for my
> purposes it was too
> computationally expensive (read slow).  Therefore I am
> looking for something
> a little quicker.  I now appreciate I cannot use
> interp due to it treating
> coordinates as planar.  I found a good tutorial for
> kriging in gstat at
> http://casoilresource.lawr.ucdavis.edu/drupal/node/442. 
> >From a quick play
> around I can't immediately see how functions like Tps in
> the field package
> meet my needs.  Could anyone provide directions to
> primers on the web?.
> Googling for R is quite tricky (even with rseek).
> 
> Essentially I want a fast simple method for interpolating a
> single point
> value from its surounding points.  It shouldn't be too
> difficult but I must
> admit I am struggling a little!
> 
> Greg
> 
> On Tue, Apr 7, 2009 at 8:26 AM, Roger Bivand <Roger.Bivand at nhh.no>
> wrote:
> 
> > 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.Bivand at nhh.no
> >
> >
> 
>     [[alternative HTML version deleted]]
> 
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>