Skip to content

Regridding

miturbide edited this page Oct 6, 2016 · 4 revisions

Regridding/interpolation

Function interpGrid performs the interpolation of gridded datasets into a new user-defined grid using bilinear weights or nearest-neighbour methods. In the following, we illustrate the different options that are supported using gridded data from the NCEP/NCAR Reanalysis 1 dataset. Note the the same operation can be done on multimember grids. Optionally, parallelization is supported for bilinear interpolation.

In general, the preferred choice for interpolation is nearest-neighbours, being the default method. If a bilinear interpolation approach is chosen, there are two different algorithm implementations, that can be selected using the argument bilin.method:

  • The "akima" choice (the default when bilinear method is in use) uses akima::interp. Faster, but does not support missing data in the input grid.
  • The "fields" option uses the fields::interp.surface.grid algorithm implementation. Much slower than akima, but able to handle missing data.

Grid data:

data(NCEP_Iberia_tas)

plotClimatology(climatology(NCEP_Iberia_tas, list(FUN = mean, na.rm = T)), 
                backdrop.theme = "countries", 
                scales = list(draw = T),
                main = attr(NCEP_Iberia_tas$Variable, "longname"))

# Bilinear interpolation to a smaller domain centered in Spain using a 0.5 degree resolution
# in both X and Y axes
GtoG <- interpGrid(NCEP_Iberia_tas, 
                           new.coordinates = list(x = c(-10,5,.5), y = c(36,44,.5)),
                           method = "bilinear")

## [2016-02-18 10:47:31] Performing bilinear interpolation... may take a while
## [2016-02-18 10:47:31] Done

plotClimatology(climatology(GtoG, list(FUN = mean, na.rm = T)), 
                backdrop.theme = "countries", 
                scales = list(draw = T),
                main = attr(NCEP_Iberia_tas$Variable, "longname"))

Note that new attributes "interpolation", "resX" and "resY" indicate that the original data have been interpolated:

attributes(GtoG$xyCoords)

## $names
## [1] "x" "y"
## 
## $projection
## [1] "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0"
## 
## $resX
## [1] 0.5
## 
## $resY
## [1] 0.5
## 
## $interpolation
## [1] "bilinear"

We could use method = "nearest" to perform the nearest neighbour interpolation in all the cases of this worked example, for instance:

GtoG <- interpGrid(NCEP_Iberia_tas, 
                           new.coordinates = list(x = c(-10,5,.5), y = c(36,44,.5)),
                           method = "nearest")

## [2016-02-18 10:47:32] Performing nearest interpolation... may take a while
## [2016-02-18 10:47:44] Done

plotClimatology(climatology(GtoG, list(FUN = mean, na.rm = T)), 
                backdrop.theme = "countries", 
                scales = list(draw = T),
                main = attr(NCEP_Iberia_tas$Variable, "longname"))

The new.coordinatesargument can be filled with an existing grid of other loaded data. We can use function getGridto extract the coordinates:

data(tp_forecast)
par(mfrow = c(1,3))
plotMeanGrid(tp_forecast)
plotMeanGrid(NCEP_Iberia_tas)

# Bilinear interpolation to a smaller domain centered in Spain using a 0.5 degree resolution
# in both X and Y axes
new.coordinates <- getGrid(tp_forecast)
GtoG.new <- interpGrid(NCEP_Iberia_tas, new.coordinates = new.coordinates,
                            method = "bilinear")

## [2016-02-18 10:47:45] Performing bilinear interpolation... may take a while
## [2016-02-18 10:47:45] Done

The warning messages indicate that the bounding box of the new.coordinates definition is outside the boundaries of the input grid bounding box.

plotMeanGrid(GtoG.new)
par(mfrow = c(1,1))


<-- Home page of the Wiki