Having done "regridding" (aka 2D interpolation), I'd like to know how to do "reboxing" (aka 3D interpolation) in R. Can anyone point me to docs explaining how to do this, or example code that does this? What I mean, why I ask: I'm a student attempting to incorporate N2O emissions inventories (EIs), initial conditions (ICs), and boundary conditions (BCs) from relatively coarse-scaled global inventories (with horizontal grids on the order of degrees lon-lat) into a relatively fine-scaled regional atmospheric model (12-km LCC horizontal). The emissions are all from surface, hence effectively 2D, and I now know how to "regrid," or interpolate from one 2D grid to another. However the model and the N2O IC/BCs are 3D: the cells are boxes, not just grids. And of course the model (on the one hand) and the N2O IC/BCs have different vertical layering (i.e., numbers of layers, and top height, or maximum vertical extent), as well as different horizontal grids. I'd like to know: is there R code available that will "rebox," i.e., interpolate vertically as well as horizontally? thanks in advance, and feel free to forward, Tom Roche <Tom_Roche at pobox.com>
reboxing, or 3D interpolation?
3 messages · Tom Roche, Edzer Pebesma
Package gstat provides 3D interpolation: library(gstat) demo(gstat3D) in this case, if you'd leave out the variogram model, inverse distance is used for interpolation. An typical issue is, independent whether the third dimension represents z (vertical space) or time, that variability, or spatial correlation, in x/y differs from that in z for the same distance unit (where "same" is only meaningful when z is vertical space). You can then rescale z, or use anisotropic variogram models.
On 11/01/2012 07:44 AM, Tom Roche wrote:
Having done "regridding" (aka 2D interpolation), I'd like to know how to do "reboxing" (aka 3D interpolation) in R. Can anyone point me to docs explaining how to do this, or example code that does this? What I mean, why I ask: I'm a student attempting to incorporate N2O emissions inventories (EIs), initial conditions (ICs), and boundary conditions (BCs) from relatively coarse-scaled global inventories (with horizontal grids on the order of degrees lon-lat) into a relatively fine-scaled regional atmospheric model (12-km LCC horizontal). The emissions are all from surface, hence effectively 2D, and I now know how to "regrid," or interpolate from one 2D grid to another. However the model and the N2O IC/BCs are 3D: the cells are boxes, not just grids. And of course the model (on the one hand) and the N2O IC/BCs have different vertical layering (i.e., numbers of layers, and top height, or maximum vertical extent), as well as different horizontal grids. I'd like to know: is there R code available that will "rebox," i.e., interpolate vertically as well as horizontally? thanks in advance, and feel free to forward, Tom Roche <Tom_Roche at pobox.com>
_______________________________________________ R-sig-Geo mailing list R-sig-Geo at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Edzer Pebesma Institute for Geoinformatics (ifgi), University of M?nster Weseler Stra?e 253, 48151 M?nster, Germany. Phone: +49 251 8333081, Fax: +49 251 8339763 http://ifgi.uni-muenster.de http://www.52north.org/geostatistics e.pebesma at wwu.de
3 days later
Apologies if the following reveals a profound lack of understanding of spatiality or geostatistics (esp kriging), which I'm learning slowly and "on-the-fly": summary: for "reboxing" from a regular global unprojected grid to a regular local projected grid (both 3D), should one use * gstat::krige? If so, what model? * gstat::idw? If so, what formula? or something completely different? details: https://stat.ethz.ch/pipermail/r-sig-geo/2012-November/016573.html
[I seek to] incorporate N2O emissions inventories (EIs), initial conditions (ICs), and boundary conditions (BCs) from relatively coarse-scaled global inventories (with horizontal grids on the order of degrees lon-lat) into a [finer-scaled] regional atmospheric model (12-km LCC horizontal). The [EIs] are [surface/2D] and I now know how to "regrid," or interpolate from one 2D grid to another.
However the model and [its] IC/BCs are 3D[, and] have [both] different vertical layering (i.e., [number and height of layers]) [and] different horizontal grids. I'd like to know: is there R code available that will "rebox," i.e., interpolate vertically as well as horizontally?
To further specify: my global IC/BC estimates mass concentration for each box/voxel in a regular, unprojected grid 1.875? x 2.5? x 56 levels. I want to "rebox" from global IC/BCs to model-ready IC/BCs, where the latter must estimate mass concentration over a regular grid that horizontally projects to LCC over North America at 12-km resolution, and which has 34 vertical levels. (Both the vertical and horizontal extents of the model-ready IC/BC are subsets of the global extents.) I can compute the sizes of the boxes, so I can compute the mass in each input box, and recompute the output concentrations, presuming I can learn how to rebox the masses. https://stat.ethz.ch/pipermail/r-sig-geo/2012-November/016574.html (rearranged)
[A] typical issue is [whether] variability, or spatial correlation, in x/y differs from that in z for the same distance unit[.]
I believe the vertical and horizontal correlations of N2O concentrations differ strongly: * vertical distribution of N2O is governed by atmospheric processes (advection and diffusion) very different from the mostly biological pedospheric processes that govern the horizontal distribution of its emissions (though the two sets of processes are coupled via precipitation and temperature) * N2O in the troposphere (which unfortunately is all my group can now model) has a long life (scale ~= century) relative to the temporal scale of the processes that produce/emit it (which vary diurnally and seasonally, with soil moisture and temperature) That being said, N2O datasets are quite appallingly sparse, given its current significance as a depleter of stratospheric O3 and as a GHG.
Package gstat provides 3D interpolation:
library(gstat) demo(gstat3D)
in this case, if you'd leave out the variogram model, inverse distance is used for interpolation.
I note R session
demo(gstat3D)
...
# $Id: gstat3D.R,v 1.5 2007-02-23 13:34:07 edzer Exp $ # simple demo of 3D interpolation of 50 points with random normal # values, randomly located in the unit cube n <- 50 data3D <- data.frame(x=runif(n), y=runif(n), z=runif(n), v=rnorm(n)) coordinates(data3D)=~x+y+z range1D <- seq(from=0, to=1, length=20) grid3D <- expand.grid(x=range1D, y=range1D, z=range1D) gridded(grid3D)=~x+y+z res3D <- krige(formula=v ~ 1, data3D, grid3D, model=vgm(1, "Exp", .2))
So for my usecase, 1 'coordinates(data3D)' would be the coordinates of the centers of the boxes defined by the global/unprojected grid 2 'data3D' would be the pairs of coordinates and mass concentrations for the global/unprojected grid 3 'grid3D' would be the coordinates of the centers of the local/projected gridboxes 4 'res3D' would be the pairs of coordinates and mass concentrations for the global/unprojected grid correct? (Though I note that, in the demo
names(res3D at data)
[1] "var1.pred" "var1.var" and I don't know why there are 2 data columns rather than 1.) If so, I can create data3D and get its coordinates from my input grid, and obtain grid3D by projecting coordinates(data3D) using R package=M3. But to get res3D from grid3D, for my usecase, should I use * gstat::krige? If so, what model? (or no model?) * gstat::idw? If so, what formula? * or something else? Your assistance is appreciated! Tom Roche <Tom_Roche at pobox.com>