[R-sig-Geo] Reply to: Spatial interpolation using thin plate spline [fields] through time

Robert Way rway024 at uottawa.ca
Thu Jul 9 00:36:48 CEST 2015


This is a reply to a prior post I had written here:
https://stat.ethz.ch/pipermail/r-sig-geo/2014-December/022082.html

The solution to the problem mentioned is to create a data frame which has
the to-be-interpolated measurements organized in columns per time step
corresponding to the independent x-variables. This code will allow an
individual to reuse the same independent variables (x position, y position,
elevation) with changing measurements at each time step in a loop. The loop
uses an input DEM to interpolate and output a .geotiff at each time step.

Below the rough code required for this type of thin plate spline
implementation is described. A cleaned up and fully commented version will
be made available as a supplemental to the following publication when it is
accepted:

Way, R.G., Lewkowicz, A.G. and Bonnaventure, P.P. Creating high-resolution
gridded air temperature and degree-day maps for the Labrador-Ungava region
of northern Canada. To be submitted to International Journal of Climatology.


### Load packages (not all required) ###
library(fields)
library(raster)
library(spatial.tools)
library(gdalUtils)
library(rgdal)
library(gstat)

### Use all computer resources ###
sfQuickInit()

### Load Independent X-Data ###
xdata <- file.choose()
x <- read.csv(xdata, TRUE, ",")

### Load Dependent Y-Data organized by time step ###
ydata <- file.choose()
y <- read.csv(ydata, TRUE, ",")

### Load elevation data ###
elev <- raster("elev.tif")

### Loop ###
### Make sure to set up name and appropriate file path ###
### This takes a while - Status Bar may be a good add-on ###

for (i in 1:(ncol(y))){

fit <- Tps(x, y[,i])
result <- interpolate(elev,fit, xyOnly=FALSE)
myfile <- file.path("...",
paste0("Name", "_", i, ".tif"))
writeRaster(result, myfile, overwrite=TRUE)
}

sfQuickStop()

-- 



*Robert WayB.A., M.Sc., PhD StudentUniversity of Ottawa*

	[[alternative HTML version deleted]]



More information about the R-sig-Geo mailing list