A better strategy for interpolating gravity and magnetic data

We present a new strategy for gravity and magnetic data interpolation and processing. Our method is based on the equivalent layer technique (EQL) and produces more accurate interpolations when compared with similar EQL methods. It also reduces the computation time and memory requirements, both of which have been severe limiting factors.

The equivalent layer technique (also known as equivalent source, radial basis functions, or Green’s functions interpolation) is used to predict the value of gravity and magnetic fields (or transformations thereof) at any point based on the data gathered on some observation points. It consists in estimating a source distribution that produces the same field as the one measured and using this estimate to predict new values. It generally outperforms other general-purpose 2D interpolators, like the minimum curvature or bi-harmonic splines, because it takes into account the height of measurements and the fact that these fields are harmonic functions. Nevertheless, defining a layout for the source distribution used by the EQL is not trivial and plays an important role in the quality of the predictions.

The most widely used source distributions are: (a) a regular grid of point sources and (b) one point source beneath each observation point. We propose a new source distribution: (c) divide the area into blocks, calculate the average location of observation points inside each block, and place one point source beneath each average location. This produces a smaller number of point sources in comparison with the other source distributions, effectively reducing the computational load. Traditionally, the source points are located: (i) all at the same depth or (ii) each source point at a constant relative depth beneath its corresponding observation point. Besides these two, we also considered (iii) a variable relative depth for each source point proportional to the median distance to its nearest neighbours. The combination of source distributions and depth configurations leads to seven different source layouts (the regular grid is only compatible with the constant depth configuration).

We have scored the performance of each configuration by interpolating synthetic ground and airborne gravity data, and comparing the interpolation against the true values of the model. The block-averaged source layout (c) with variable relative depth (iii) produces more accurate interpolation results (R² of 0.97 versus R² of 0.63 for the traditional grid layout) in less time than the alternatives (from 2 to 10 times faster on our test cases). These results are consistent between ground and airborne survey layouts. Our conclusions can be extrapolated to other applications of equivalent layers, such as upward continuation, reduction-to-the-pole, and derivative calculation. What is more, we expect that these optimizations can benefit similar spatial prediction problems beyond gravity and magnetic data.

The source code developed for this study is based on the EQL implementation available in Harmonica (fatiando.org/harmonica), an open-source Python library for modelling and processing gravity and magnetic data.