[R-sig-Geo] reboxing with package=GSIF

[Tom Roche]

summary: how to use GSIF::make.3Dgrid for general-purpose "reboxing," or
3D interpolation? Particularly, can one use varying vertical bounds
rather than the same "standard depths" for both input and output?

details:

Apologies if the following questions reveal a profound ignorance of
spatiality or geostatistics (which I'm learning slowly and "on-the-fly")
and if I'm missing something simple regarding GSIF; unfortunately the
tutorial

http://gsif.r-forge.r-project.org/tutorial_eberg.php

seems to require more geostatistical background than I currently have.

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/016578.HTML
> I've produced some 3D interpolation examples in the experimental
> package GSIF (under construction). You can actually create some kind
> of voxels via the SpatialPixelsDataFrame class

> http://gsif.r-forge.r-project.org/make.3Dgrid-method.html

> and then predict at those locations.

From

http://gsif.r-forge.r-project.org/make.3Dgrid-method.html

it appears one can provide as input

* an unprojected SpatialPixelsDataFrame or RasterBrick

* a PROJ.4 string defining the output CRS

* a "list of standard depths"

and get an output list of SpatialPixelsDataFrame's (one per standard
depth). From the description, it seems like GSIF::make.3Dgrid is doing a
2D regridding at each standard depth, such that the heights of the input
box/voxels must match that of the output box/voxels. Is that correct?

If not, how can one rebox from an input with one set of layer heights
(or vertical bounds) to an output with another set of layer heights?
Can `gdalwarp` do that, or is it (like, IIUC, raster::projectRaster)
restricted to 2D?

Your assistance is appreciated! Tom Roche <Tom_Roche at pobox.com>