Basis functions for the consistent and accurate representation of surface mass loading

Inversion of geodetic site displacement data to infer surface mass loads has previously been demonstrated using a spherical harmonic representation of the load. This method suffers from the continent-rich, ocean-poor distribution of geodetic data, coupled with the predominance of the continental load (water storage and atmospheric pressure) compared with the ocean bottom pressure (including the inverse barometer response). Finer-scale inversion becomes unstable due to the rapidly increasing number of parameters which are poorly constrained by the data geometry. Several approaches have previously been tried to mitigate this, including the adoption of constraints over the oceanic domain derived from ocean circulation models, the use of smoothness constraints for the oceanic load, and the incorporation of GRACE gravity field data. However, these methods do not provide appropriate treatment of mass conservation and of the ocean's equilibrium-tide response to the total gravitational field. Instead, we propose a modified set of basis functions as an alternative to standard spherical harmonics. Our basis functions allow variability of the load over continental regions, but impose global mass conservation and equilibrium tidal behaviour of the oceans. We test our basis functions first for the efficiency of fitting to realistic modelled surface loads, and then for accuracy of the estimates of the inferred load compared with the known model load, using synthetic geodetic displacements with real GPS network geometry. Compared to standard spherical harmonics, our basis functions yield a better fit to the model loads over the period 1997-2005, for an equivalent number of parameters, and provide a more accurate and stable fit using the synthetic geodetic displacements. In particular, recovery of the low-degree coefficients is greatly improved. Using a nine-parameter fit we are able to model 58 per cent of the variance in the synthetic degree-1 zonal coefficient time-series, 38-41 per cent of the degree-1 non-zonal coefficients, and 80 per cent of the degree-2 zonal coefficient. An equivalent spherical harmonic estimate truncated at degree 2 is able to model the degree-1 zonal coefficient similarly (56 per cent of variance), but only models 59 per cent of the degree-2 zonal coefficient variance and is unable to model the degree-1 non-zonal coefficients.