Change of support: details

Figure A.1: Block discretization: locations $s_i$ (+) and weights $w_i$ using (a) 4 $\times$ 4 Gaussian quadrature, (b) 4 $\times$ 4 regular

(a)

(b)

The average of a function $f(\cdot)$ over a block (or line or volume) $B$,

$$f(B)=\vert B\vert^{-1} \int_B f(s)ds$$

with $\vert B\vert$ the block area (or lenght or volume) is approximated by:

$$f(B) \approx \sum_{i=1}^{N} w_i f(s_i)$$

with $\sum_{i=1}^{N} w_i = 1$, and $s_i$ the points that discretise the block $B$ and where $w_i$ are the weights for each point $s_i$.

Rectangular blocks

For rectangular blocks, gstat calculates block averages (for semivariances, (generalised) covariances, coordinate polynomials or inverse distance weighted interpolations) by using Gauss quadrature with 4 points in each block dimension (Fig. A.1, [1]). Regular discretization ($w_i = N^{-1}$) is obtained by setting nblockdiscr to the number of points in each direction (section 4.4). Blocks are always centred at prediction locations.

Block-to-block covariances $C(B_1 , B_2 )= \vert B_1\vert^{-1} \vert B_2\vert^{-1} \int_{B_1} \int_{B_2} C(x) dudv$ are calculated as

$$C(B_1 , B_2 ) \approx \sum_{i=1}^{N_1} \sum_{j=1}^{N_2} w_i w_j C(s_i , s_j )$$

with $s_i$ discretizing $B_1$ and $s_j$ discretizing $B_2$. For block kriging, block-to-block (generalised) covariances $C ( B_0 , B_0 )$ are calculated only once per variogram model. For Gaussian simulation of block averages though, this double sum is recalculated for each pair of (simulated) block averages in a kriging neighbourhood. This takes a while.

Non-rectangular blocks

Mean values for arbitrary shaped, e.g. non-rectangular `blocks' are obtained when the area:'file', ... ; command is set (using the data syntax) to specify the points that describe the shape. This area should, depending on the dimensions, be centred around the location (0), (0,0), or (0,0,0) in order to obtain predictions centred around the prediction locations. Weights $w_i$ of points discretizing the area are set to $N^{-1}$, or, if the V field is defined for the area, they are set to the values of that field.

If area is specified but no prediction locations are specified, then the area average of the (points discretizing) area itself will be calculated, and written to output. In this case, the area should not (necessarily) be centred around zero, but should be the actual area for which area-average predictions are required.

Base functions block averages

When base functions are used for the trend, gstat assumes that the user-defined base function values at the prediction location are block average values: gstat cannot average user-defined base functions over the prediction block (it does so for coordinate polynomial base-functions).