S U M M A R YWe show the numerical applicability of a multiresolution method based on harmonic splines on the 3-D ball which allows the regularized recovery of the harmonic part of the Earth's mass density distribution out of different types of gravity data, for example, different radial derivatives of the potential, at various positions which need not be located on a common sphere. This approximated harmonic density can be combined with its orthogonal anharmonic complement, for example, determined out of the splitting function of free oscillations, to an approximation of the whole mass density function. The applicability of the presented tool is demonstrated by several test calculations based on simulated gravity values derived from EGM96. The method yields a multiresolution in the sense that the localization of the constructed spline basis functions can be increased which yields in combination with more data a higher resolution of the resulting spline. Moreover, we show that a locally improved data situation allows a highly resolved recovery in this particular area in combination with a coarse approximation elsewhere which is an essential advantage of this method, for example, compared to polynomial approximation. and the references therein). From the mathematical point of view, the inverse gravimetric problem is based on a Fredholm integral equation of the first kind involving Newton's law of gravitation. Thus, the relation between the gravitational potential V of the Earth and the density distribution ρ is given byIn reality, the gravitational potential V can be determined outside the Earth whereas the density ρ represents the unknown part. This means that we have to consider the inversion of the equation above. Problems of such kind are called inverse problems. By Hadamard's classification, that is the uniqueness, the existence and the stability of the solution, we divide inverse problems in ill-posed and well-posed problems. It is a well-known fact that the inverse gravimetric problem is ill-posed because each of these three criteria can be invalid.For finding an approximation for the solution of such an ill-posed inverse problem several different methods were developed. The classical approach, for example, is a method using a truncated singular value decomposition (see, e.g. Škorvanek 1981;Tscherning & Sünkel 1981;Tscherning & Strykowski 1988). This method has several well-known disadvantages like the non-localizing character of the used spherical C 2008 The Authors