FULL PAPER

Open Access

# Global maps of the magnetic thickness and magnetization of the Earth's lithosphere

<figure><img src="image_1.png" /></figure>

CrossMark

Foteini Vervelidou<fnref n="1" />* and Erwan Thébault²

**Abstract**

We have constructed global maps of the large-scale magnetic thickness and magnetization of Earth's lithosphere. Deriving such large-scale maps based on lithospheric magnetic field measurements faces the challenge of the masking effect of the core field. In this study, the maps were obtained through analyses in the spectral domain by means of a new regional spatial power spectrum based on the Revised Spherical Cap Harmonic Analysis (R-SCHA) formalism. A series of regional spectral analyses were conducted covering the entire Earth. The R-SCHA surface power spectrum for each region was estimated using the NGDC-720 spherical harmonic (SH) model of the lithospheric magnetic field, which is based on satellite, aeromagnetic, and marine measurements. These observational regional spectra were fitted to a recently proposed statistical expression of the power spectrum of Earth's lithospheric magnetic field, whose free parameters include the thickness and magnetization of the magnetic sources. The resulting global magnetic thickness map is compared to other crustal and magnetic thickness maps based upon different geophysical data. We conclude that the large-scale magnetic thickness of the lithosphere is on average confined to a layer that does not exceed the Moho.

**Keywords:** Lithospheric magnetic field; Curie isotherm; Moho discontinuity; Magnetization; Spectral analyses; R-SCHA; WDMAM; CHAMP satellite

**Background**

The magnetic field on the Earth's surface results from the superposition of various sources, both of internal and external to the Earth origin. The surface field is dominated by the field generated within Earth's outer core (see, e.g., Hulot et al. 2015); however, there are important contributions from the magnetized rocks of Earth's lithosphere (for reviews on the lithospheric magnetic field see, e.g., Purucker and Whaler 2015; Thébault et al. 2010). The core and the lithospheric magnetic field overlap in the spectral domain. The core field dominates the power spectrum over spherical harmonic (SH) degrees 1 to 15 (Langel and Estes (1982)) approximately, and it is only at larger degrees that the lithospheric field becomes the primary contributor.

As a consequence of the dominant role of the core field, the large-scale lithospheric magnetic field corresponding

to SH degrees 1 to 15 cannot be recovered through magnetic field measurements alone. To overcome this constraint, studies inferring the thickness and magnetization of the magnetic sources of the lithospheric field have relied on a priori information (e.g., Purucker et al. 2002; Purucker and Whaler 2004). Specifically, these studies assume that over SH degrees 1 to 15 the Moho, a seismic discontinuity that separates the crust from the mantle, is also a magnetic boundary. This approximation is supported by studies that have shown the mantle to be mainly depleted in magnetic minerals (e.g., Wasilewski and Mayhew 1982). However, recent studies, e.g., Ferré et al. 2013; 2014, Friedman et al. 2014; Martin-Hernandez et al. 2014, propose that there is evidence in some areas for a magnetized upper mantle. This suggests the lower boundary of the magnetic thickness does not necessarily coincide with the Moho. In this study, the term magnetic thickness refers to the part of the lithosphere that contributes to the observable magnetic signal, as opposed to the seismically defined crustal thickness.

<footnote>*Correspondence: foteini@gfz-potsdam.de</footnote>

<footnote n="1">¹Helmholtz Centre Potsdam-GFZ German Research Centre for Geosciences, Section 2.3 Earth's magnetic field, Telegrafenberg, 14473 Potsdam, Germany
Full list of author information is available at the end of the article</footnote>
New ways of recovering Earth's large-scale magnetic thickness and magnetization, that do not rely on the assumption that Moho is a magnetic boundary, would provide some independent assessment. This in turn could contribute to studies on the composition and the thermal state of Earth's lithosphere (e.g., Fox Maule et al. 2005), since the minerals are only magnetized down to the Curie isotherm.

We propose such a methodology in this paper. We choose to work in the spectral domain and rely upon statistical assumptions about the nature of lithospheric field sources, along the lines of studies seeking either to model the large scale lithospheric field or to infer information about magnetic thickness and magnetization. Such studies have been carried out either on a global scale by means of Spherical Harmonic Analysis (SHA) (e.g., Jackson 1990, 1994; Korte et al. 2002; O'Brien et al. 1999; Voorhies 1998; Voorhies et al. 2002; Voorhies 2008) or on local scales using planar geometry (e.g., Bouligand et al. 2009; Li et al. 2013; Maus et al. 1997). However, existing global SH statistical expressions do not represent accurately recent lithospheric field models, and planar geometry is not well suited for large-scale studies (e.g., Langel and Hinze 1998, chapter 5.3.2 or Vervelidou 2013, Figure 5.27b). Recently, Thébault and Vervelidou 2015 proposed a statistical SH expression for the power spectrum of Earth's lithospheric magnetic field that improves the fit to recent lithospheric field models while also being in agreement with global forward models of continental and oceanic magnetization.

In this study, we make use of this expression to construct global maps of Earth's large-scale magnetic thickness and magnetization. In section "Surface spherical cap harmonic power spectrum", we use the R-SCHA methodology (Thébault et al. 2004, 2006, Thébault 2008) and derive an expression for the surface spherical cap harmonic power spectrum, which we show can be directly related to any SH power spectrum. In section "A statistical expression for the Earth's lithospheric magnetic field", we present briefly the statistical expression of Thébault and Vervelidou 2015. In section "Synthetic analyses", we combine this expression with synthetic spherical cap spectra and demonstrate the ability to provide magnetic thickness and magnetization estimates through a misfit analysis. In section "Global models of the Earth's mean magnetic thickness and magnetization", we calculate the observational spherical cap spectra based on the satellite-, aeromagnetic- and marine-based NGDC-720 lithospheric field model (Maus 2010), and we present the magnetic thickness and magnetization maps derived for the Earth. We note that a similar approach, based on the statistical expression of Voorhies et al. 2002 and the spatiospectral localization techniques of Dahlen and Simons 2008; Wieczorek and Simons 2007,

was followed by Lewis and Simons 2012 for inferring the crustal magnetic thickness of Mars. Finally, we transform our magnetic thickness map into spherical harmonics and compare it to the SH power spectra of other crustal and magnetic thickness models.

## Methods

### Surface spherical cap harmonic power spectrum

#### Definition

The spherical cap harmonic analysis introduced by Haines 1985 was revised by Thébault et al. 2004, 2006 to account for the upward/downward continuation of the magnetic field with altitude. The magnetic potential described by the Revised-Spherical Cap Harmonic Analysis (R-SCHA) is the solution of a boundary value problem with continuity conditions on each boundary of the considered domain. The magnetic field under consideration is assumed to be measured in source-free regions so that the magnetic potential solves the Laplace's equation. The domain, $\Omega$, is a spherical cone resulting from the intersection of two spherical caps defined by two concentric spheres of radius $\alpha$ and $\beta$ and an infinite cone of half angle $\theta_0$, whose apex lies at the Earth's center. This general problem takes the formal expression

$$ \nabla^2 V(r, \theta, \phi) = 0, \text{ within the volume } \Omega, \qquad (1a) $$

$$ \delta_i V(r, \theta, \phi) + \gamma_i \frac{\partial V(r, \theta, \phi)}{\partial n_i} = G_i, \text{ on each surface } \partial \Omega_i, \qquad (1b) $$

where $n_i$ is the normal vector of the surface $\partial \Omega_i$ pointing outward; $\delta_i$ and $\gamma_i$ are two arbitrary constants; and $r, \theta, \phi$ are radius, colatitude and longitude, respectively, in the reference frame of the spherical cap. $\partial \Omega_{\theta_0}$ refers to the surface that forms the lateral boundary of the cap and $\partial \Omega_{\alpha}$ ($\partial \Omega_{\beta}$) to the surface that forms the lower (upper) boundary of the cap. The boundary value problem of Eq. 1 can be split into two sub-problems so as to reveal homogeneous boundary conditions on each of the three boundaries of the cone:

$$ \nabla^2 V_1 = 0 \text{ in } \Omega, \\ \delta_{\theta}^1 V_1|_{\theta_0} + \gamma_{\theta}^1 \frac{\partial V_1}{\partial \theta} \Big|_{\theta_0} = 0 \text{ on } \partial \Omega_{\theta_0}, \qquad (2) $$

$$ \delta_r^1 V_1|_{\alpha, \beta} + \gamma_r^1 \frac{\partial V_1}{\partial r} \Big|_{\alpha, \beta} = G_{\alpha, \beta} \text{ on } \partial \Omega_{\alpha} \text{ and } \partial \Omega_{\beta}, $$
and

$$\nabla^2 V_2 = 0 \text{ in } \Omega,$$

$$\delta_{\theta}^{2} V_{2}|_{\theta_{0}} + \gamma_{\theta}^{2} \left. \frac{\partial V_{2}}{\partial \theta} \right|_{\theta_{0}} = G_{\theta_{0}} \text{ on } \partial \Omega_{\theta_{0}}, \quad (3)$$

$$\delta_r^2 V_2|_{\alpha, \beta} + \gamma_r^2 \left. \frac{\partial V_2}{\partial r} \right|_{\alpha, \beta} = 0 \text{ on } \partial \Omega_{\alpha} \text{ and } \partial \Omega_{\beta}.$$

Now, the solution of the boundary value problem of Eq. 1 can be written as the superposition of potentials $V_1$ and $V_2$, such that

$$V(r, \theta, \phi) = V_1(r, \theta, \phi) + V_2(r, \theta, \phi), \quad (4)$$

with $V_1$ and $V_2$ being the solutions of the boundary values problems 2 and 3.

Although there is an infinite number of possible boundary conditions, here we examine two that give rise to solutions previously considered for geophysical applications (Thébault et al. 2004, 2006).

### Neumann-Neumann Boundary Value Problem

One problem involves Neumann boundary conditions on each of the boundaries of the conical domain. This Neumann-Neumann Boundary Value Problem (BVP) is defined by setting $\delta_{\theta}^{1} = \delta_{\theta}^{2} = \delta_{r}^{1} = \delta_{r}^{2} = 0$ and $\gamma_{\theta}^{1} = \gamma_{r}^{2} = \gamma_{\theta}^{2} = \gamma_{r}^{1} = 1$ in Eqs. 2 and 3. The formal solution involves the so-called Legendre solution

$$V_1(r, \theta, \phi) = \alpha \sum_{k \ge 1} \sum_{-k \le m \le k} \left[ \left(\frac{\alpha}{r}\right)^{n_k+1} G_{n_k}^{i,m} + \left(\frac{r}{\alpha}\right)^{n_k} G_{n_k}^{e,m} \right] \beta_{n_k}^m(\theta, \phi), \quad (5)$$

where

$$\beta_{n_k}^m(\theta, \phi) = \cos m\phi P_{n_k}^m(\cos\theta), \quad 0 \le m \le k \quad (6)$$

$\beta_{n_k}^m(\theta, \phi) = \sin|m|\phi P_{n_k}^{|m|}(\cos\theta), \quad -k \le m \le -1, \quad (7)$
with $P_{n_k}^{|m|}(\cos\theta)$ the generalized Legendre functions and $\{G_{n_k}^{i,m}, G_{n_k}^{e,m}\}$ the harmonic coefficients of order $m$ and real degree $n_k$, and the so-called Mehler solution

$$V_2(r, \theta, \phi) = \alpha \sum_{p \ge 0} \sum_{m \ge 0} G_p^m \gamma_p^m(r, \theta, \phi), \quad (8)$$

with $\gamma_p^m$ the Mehler basis functions and $G_p^m$ the harmonic coefficients. For the purpose of this work, the exact expression of the Mehler basis functions is not needed but we refer the reader to Thébault et al. 2004 for a detailed expression. These basis functions have a decreasing contribution as a function of colatitude $\theta$ from the lateral boundaries of the cap to its axis of revolution.

### Dirichlet-Neumann Boundary Value Problem

The second boundary values problem differs from the previous one as it considers Dirichlet instead of Neumann

boundary conditions on the lateral boundary of the cap, $\partial\Omega_{\theta_0}$. This problem is obtained by setting $\delta_{\theta}^{1} = \delta_{\theta}^{2} = \gamma_{r}^{1} = \gamma_{\theta}^{2} = 1$ and $\gamma_{\theta}^{1} = \gamma_{\theta}^{2} = \delta_{r}^{1} = \delta_{r}^{2} = 0$ in Eqs. 2 and 3. With these boundary conditions, the potential $V_1$ has formally the same form as the one given in Eq. 5 but is numerically different. The degrees $n_k$ are now the roots of $P_{n_k}^m(\cos\theta_0) = 0$ rather than that of $\left.\frac{dP_{n_k}^m}{d\theta}\right|_{\theta_0} = 0$, as is the case for the Neumann boundary condition. The solution $V_2$ is identical to Eq. 8.

From the above description of the magnetic potential, we seek the definition of a spherical cap power spectrum that could be readily compared with a SH power spectrum. In spherical harmonics, the spatial power spectrum is defined as the average square magnetic field contribution per SH degree $n$ over a spherical surface (Lowes 1974). Therefore, we restrict ourselves to a spherical cap surface. Over a spherical cap surface, $\vartheta\Omega_{\rho}$ at radius $\alpha \le \rho \le \beta$, the square magnetic field $E_{\vartheta\Omega_{\rho}}$ is given by

$$E_{\vartheta\Omega_{\rho}} = \|\nabla V_1 + \nabla V_2\|_{\vartheta\Omega_{\rho}}^2 = \|\nabla V_1\|_{\vartheta\Omega_{\rho}}^2 + \|\nabla V_2\|_{\vartheta\Omega_{\rho}}^2 + 2\langle\nabla V_1, \nabla V_2\rangle_{\vartheta\Omega_{\rho}}, \quad (9)$$

where $\langle \cdot, \cdot \rangle$ stands for the inner product and $\|\cdot\|$ for the norm, both defined on the spherical cap's surface. In this expression, the gradients of the potentials $V_1$ and $V_2$ are not orthogonal on a spherical cap surface when they are solutions of the Neumann-Neumann or the Dirichlet-Neumann boundary conditions (see, e.g., Thébault et al. 2004 or Thébault et al. 2006). As a result, the contributions of the two basis functions cannot be separated from each other. However, if we assume that the contribution of the Mehler basis functions is zero, we can benefit from the orthogonality properties of the Legendre basis functions over the surface $\vartheta\Omega_{\rho}$ of a spherical cap (see for instance Haines 1985). Based on this assumption, the $E_{\vartheta\Omega_{\rho}}$ can be approximated by

$$E_{\vartheta\Omega_{\rho}} = \|\nabla V_1 + \nabla V_2\|_{\vartheta\Omega_{\rho}}^2 \approx \|\nabla V_1\|_{\vartheta\Omega_{\rho}}^2 = \sum_{m=0}^{\infty} (1 + \delta_{m,0}) \pi \rho^2 \sum_{k=0}^{\infty} \left\{ (n_k + 1)(2n_k + 1) \left(\frac{a}{\rho}\right)^{2n_k+4} \left[ (G_{n_k}^{i,m})^2 + (H_{n_k}^{i,m})^2 \right] + n_k(2n_k + 1) \left(\frac{\rho}{a}\right)^{2n_k-2} \left[ (G_{n_k}^{e,m})^2 + (H_{n_k}^{e,m})^2 \right] \right\} \|P_{n_k}^m\|^2, \quad (10)$$

where, according to the common convention, coefficients $H_{n_k}^{i,m}$ ($H_{n_k}^{e,m}$) stand for the coefficients $G_{n_k}^{i,m}$ ($G_{n_k}^{e,m}$) of Eq. 5 with negative order $m$. By further normalizing the expression (10) by $S_{\vartheta\Omega_{\rho}} = 2\pi\rho^2(1 - \cos\theta_0)$, the area of the surface $\vartheta\Omega_{\rho}$, we obtain the estimation of the mean square magnetic field.
## Relationship to a spherical harmonic power spectrum

Contrary to the case of the Spherical Harmonic Analysis (SHA), where the whole sphere is under consideration, in the R-SCHA only part of the sphere is under consideration. This has direct consequences for the spectral resolution of the R-SCHA. It is known that the higher the concentration of a function in the space domain, the poorer its spectral resolution (e.g., Wieczorek and Simons 2005). In the case of SHA, the spectral resolution, i.e. the distance between two consecutive terms of its power spectrum, is equal to one ($n_{i+1} - n_i = 1$). In R-SCHA, the spectral resolution depends on the cap's half-angle, $\theta_0$, and is equal to $\pi/\theta_0$, where $\theta_0$ is given in radians. The larger the half-angle $\theta_0$, the smaller the distance between two consecutive terms of its power spectrum becomes and therefore the better the spectral resolution. For the limit case where $\theta_0 = \pi$, we recover the spectral resolution of the SHA.

The above considerations are taken into account for the numerical implementation of Eq. 10. Each term of the double sum of Eq. 10 corresponds to a different degree $n_k$, where $n_k$ is a real number. For this spectrum to be useful for geophysical interpretations, each of its terms should be related to a characteristic length scale. For this, it is convenient to group together terms that correspond to degrees that lie close to one another. The size of these bins is dictated by the available spectral resolution, therefore they are $\pi/\theta_0$ large. Each bin can be represented by the mean of its degrees, $\bar{n}$. We can therefore write $E_{\vartheta\Omega_\rho}$ of Eq. 10 as

$$E_{\vartheta\Omega_\rho} = \sum_{\bar{n}} E_{\bar{\vartheta}\Omega_\rho}^{\bar{n}}, \qquad (11)$$

where $E_{\bar{\vartheta}\Omega_\rho}^{\bar{n}}$ is the sum of the terms of Eq. 10 whose degrees $n_k$ belong to the same bin. We can then relate each term $E_{\bar{\vartheta}\Omega_\rho}^{\bar{n}}$ to a specific wavelength, $\bar{\lambda}$, in the same way as in SHA (see Backus et al. 1996, page 103):

$$\bar{\lambda} \approx \frac{2\pi\rho}{\bar{n} + \frac{1}{2}}, \qquad (12)$$

with $\rho$ the radius of the spherical surface where we estimate the wavelength.

According to Eq. 11, the terms of the R-SCHA power spectrum correspond to bins of harmonic degrees, rather than to a single harmonic degree. Consequently they carry a part of the total power that is $\pi/\theta_0$ times larger than that of the terms of a SHA power spectrum. For this reason, a direct comparison between SHA and R-SCHA power spectra is straightforward, provided that the R-SCHA power spectra are scaled by the factor of $\pi/\theta_0$. Defining $E^n$ the terms of a SH power spectrum that has the same total energy as a spherical cap power spectrum $E_{\vartheta\Omega_\rho}$, it follows that

$$E_{\bar{\vartheta}\Omega_\rho}^{\bar{n}} = \frac{\pi}{\theta_0} E^n, \qquad (13)$$

where $\theta_0$ is given in radians. This expression shows that any regional and global power spectra, experimental or theoretical, despite their differing spectral resolution, can be readily compared with each other.

## The implicit assumption

The assumption that the contribution of the Mehler basis functions to Eq. 9 is negligible, puts us in a specific framework. The Mehler basis functions are essential for ensuring the upward/downward continuation of the magnetic field (or conversely for processing data measured at various altitudes). Moreover, they are required when the magnetic field has large-scale features and has therefore non zero components on the lateral boundaries of the spherical cap (Thébault et al. 2004). Prior to interpreting Eqs. 10 or 11 as a complete power spectrum, and not as representing the energy of one particular set of basis functions, some investigations are required on the consequences of this assumption.

Considering the Dirichlet-Neumann BVP, Vervelidou 2013 shows that the Mehler basis functions carry a negligible part of the energy in the volume of the spherical cap when the wavelengths considered are smaller than the cap's aperture. This "rule of thumb" has long been known by the users of the spherical cap analysis initially defined by Haines 1985. Under this condition, we can safely neglect them also over a specific spherical cap's surface. The situation is even more conveniently illustrated for the Neumann-Neumann BVP. For this BVP, neglecting the Mehler functions over the surface is equivalent to assuming $\left. \frac{\partial V(r,\theta,\phi)}{\partial\theta} \right|_{\theta=\theta_0} = 0$. Applying the divergence free theorem in its integral form leads to the equality

$$\int_{\vartheta\Omega_\alpha} \frac{\partial V}{\partial r} ds = \int_{\vartheta\Omega_\beta} \frac{\partial V}{\partial r} ds. \qquad (14)$$

This equality requires that the radial component of the field at the lower and upper boundary of the spherical cap averages to the same value. Due to the decrease of an internal potential field with increasing altitude, this requirement can only hold true if the radial component of the field averages to zero over every spherical cap surface $\vartheta\Omega_\rho$, with $a \le \rho \le \beta$. This, again, is correct to a good approximation when the magnetic field under study has wavelengths smaller than the cap's aperture. Possible applications of the power spectrum therefore include studies of lithospheric magnetic fields, whose large scales are currently not considered in magnetic field modelling due to the masking effect of the core field. Such conditions are common in spectral analyses that seek to avoid aliasing and leakage, in particular when Fourier transform is applied to magnetic anomaly data, e.g., for upward/downward continuation purposes (see, e.g., Blakely 1996).
Consequently, Eq. 10 can be interpreted as a spatial power spectrum under the condition that the signal under study averages to zero within the spherical cap. The fulfilment of this condition should be investigated case-by-case. In this respect, we note that Hwang and Chen 1997 proposed a similar expression to Eq. 10 for working with sea-level data. However, they did not discuss the implicit assumptions behind neglecting the Mehler basis functions and therefore the applicability of this expression to their context. They applied the expression to a case where the signal does not average out within the cap, thus constructing aliased power spectra.

Once the conditions are fulfilled to neglect the contributions of the Mehler basis functions, it is wise to decide whether the spectral analysis should be performed with Legendre functions computed with the Neumann or with the Dirichlet boundary condition. Since the Neumann boundary condition has the advantage to converge faster (see, e.g., Thébault et al. 2004), it should be preferred.

However, we do not advise to perform regional modelling of magnetic field measurements available at different altitudes using only the Legendre basis functions (or to apply the SCHA as proposed originally by Haines 1985). Both the Legendre and the Mehler basis functions are needed for an accurate representation of the field in the spatial domain. Such a representation becomes possible by solving the Dirichlet-Neumann BVP (Thébault et al. 2006). This BVP is advantageous in terms of the inverse problem because the gradient of its basis functions are orthogonal in the volume of the spherical cap. In this case, the Legendre basis functions are computed with the Dirichlet boundary condition, and they can be used to calculate the regional power spectrum.

In this paper, we solve the Neumann boundary condition for the Legendre basis functions, making sure the signal under study is close to zero on average. We do not, however, interpret the model in the physical space. For this, the complete solution, based both on the Legendre and the Mehler basis functions would be necessary.

## A statistical expression for the Earth's lithospheric magnetic field

Several statistical expressions describing the power spectrum of the lithospheric magnetic field can be used for comparison to regional power spectra. Jackson 1990 for instance, proposes a form based on a random distribution of magnetization uncorrelated in space, whereas Jackson 1994 includes small-scale spatial correlation. Voorhies 1998; Voorhies et al. 2002; Voorhies 2008 explore a large variety of solutions by considering different combinations of uncorrelated and correlated sources. Recently, Thébault and Vervelidou 2015 proposed a statistical form

that depends on the mean magnetization and thickness of the magnetic lithosphere and a power law index related to the SH power spectrum of the magnetic susceptibility. This form does not take into account the effect of remanent magnetization to the power spectrum. This assumption is justified by the fact that the global magnetic field power spectra of induced and remanent magnetization are very similar (see for instance Voorhies 1998; Thébault and Vervelidou 2015). In addition, this power spectrum shows a good agreement with recent lithospheric field models over a wide range of SH degrees (from 16 to about 600) (see Vervelidou 2013, Figure 4.8 and Thébault and Vervelidou 2015, Figure 4). Moreover, its predictions for the mean magnetic thickness and magnetization compare well with those estimated by independent models (see Thébault and Vervelidou 2015, Fig. 4). This statistical power spectrum takes the form (Thébault and Vervelidou 2015, their Equation 26)

$$E^{SHA} = \sum_{n=1}^{\infty} E_{n}^{SHA} = \sum_{n=1}^{\infty} \frac{1}{2}(n+1)[\mu_{0}MF_{n}(\epsilon)]^{2}n^{-\gamma}C_{n} \quad (15)$$

with

$$F_{n}(\epsilon) = \frac{(1 - (1 - \epsilon/a)^{n-1})}{n-1} \quad (16)$$

and

$$C_n = \frac{n(20n^3 + 8n^2 - 13n + 3)}{3(2n + 3)(2n + 1)(2n - 1)} \quad (17)$$

and where a is the Earth's reference radius, equal to 6371.2 km, μ₀ is the vacuum permeability in H/m, M is the mean apparent induced magnetization in A/m, ε is the magnetic thickness and γ is the power law of the susceptibility power spectrum. As discussed in section “Relationship to a spherical harmonic power spectrum”, this form, when scaled by the factor π/θ₀, is directly comparable to any spherical cap power spectrum (see Eq. 13), provided the factor 1/2 in Eq. 15, which comes from averaging the magnetic field intensity of the axial dipole over the sphere, is replaced by the corresponding regional factor. This regional factor accounts for the increase of the Earth's main magnetic field with latitude and is given by (see for example Love and Constable 2003, their Equation 48):

$$F_{\text{latitude}} = 1 + 3 \sin(\text{lat})^2, \quad (18)$$

where lat stands here for the latitude of the center of the cap.

Therefore, the statistical expression to which the regional power spectrum of a spherical cap centred at lat
and having a half-angle of $\theta_0$ can be compared with, is written as

$$E^{R-SCHA} = \frac{\pi}{\theta_0} \frac{1}{F_{\text{latitude}}} (n+1) [\mu_0 MF_n(\epsilon)]^2 n^{-\gamma} C_n, \quad (19)$$

where $F_n(\epsilon)$, $C_n$ and $F_{\text{latitude}}$ are given by Eqs. 16, 17, and 18, respectively.

By a misfit analysis between the regional power spectrum, as described by Eqs. 10 and 11 and the statistical expression of Eq. 19 it is possible to carry out regional analyses and estimate for any spherical cap the mean magnetization and the mean magnetic thickness. In the following, all regional spectra, observational or statistical, have been downscaled by the factor $\frac{\pi}{\theta_0}$ (according to Eq. 13).

## Results

### Synthetic analyses

We assess our methodology with a set of synthetic numerical tests. We first test the compatibility between the global and the regional estimators. We draw a random set of the SH coefficients whose power spectrum coincides with the statistical expression of Eq. 15 for an arbitrary mean magnetic thickness equal to 40 km, a mean magnetization equal to 1 A/m, and a susceptibility power law equal to 1.36. We consider coefficients of SH degree 16 to 100. Using these coefficients, we create a global dataset of synthetic vector magnetic field values at Earth's reference radius. We then distribute 2000 spherical caps centred on the nodes of a global equal area grid. The half-angle of the spherical caps is set to $\theta_0 = 15$ so that the caps overlap. Using an aperture of 30 for each cap allows lithospheric field wavelengths as large as $\approx 3000$ km (SH degree 12) to be represented. Considering that lithospheric magnetic field data and models do not include larger wavelength

due to the masking by the core field, this cut-off value is sufficient to minimize the values of the mean magnetic field. Therefore, spectral leakage is avoided, and the contributions from the Mehler functions are kept negligible within each cap, as required by Eq. 14. We then estimate by least-squares the regional coefficients representing the synthetic magnetic field values over each spherical cap. Finally, we compute the regional power spectra with Eq. 11 up to an arbitrarily large degree, $\bar{n} = 480$. Figure 1 shows that despite an arbitrarily large expansion, the regional power spectra are all concentrated within the waveband covered by the synthetic data, with only limited spectral leakage from and to neighbouring harmonics. The observed fluctuations among the regional power spectra depicts the expected natural regional variability, for which the global power spectrum is only expected to represent its mean. Indeed, the mean of the 2000 regional spectra (shown in blue) is nearly identical to the global SHA power spectrum (shown in red) for degrees 16 to 100. This analysis shows that the spherical cap power spectrum can be meaningfully compared to a global power spectrum.

We carry out a second synthetic analysis to test the ability of the regional power spectra to provide regional estimates of the mean magnetic thickness and magnetization by misfit to the statistical global power spectrum. We now generate a synthetic global dataset for SH degrees 16 to 600. Again, the magnetic thickness is set to 40 km, the magnetization to 1 A/m and the susceptibility power law to 1.36. For each spherical cap, we verify that the magnetic field averages to zero, as described in section “The implicit assumption”. This is illustrated by Fig. 2, where the ratio between the mean and the rms value of the Z component is displayed for each of the 2000 caps. We see that for 98 % of the caps the ratio is less than $\pm 3 \%$, while the ratio does not exceed $\pm 6 \%$ for any cap. We then invert this data

<figure><img src="image_2.png" /><figcaption>**Fig. 1** The result of the first synthetic test. The 2000 regional spectra (grey curves), their average (blue curve) and the SHA power spectrum for a mean magnetic thickness of 40 km, a mean magnetization of 1 A/m and susceptibility power law of 1.36 (red curve) (see text for details)</figcaption></figure>
<figure><img src="image_3.png" /><figcaption>**Fig. 2** The normalised histogram of the mean value of component Z of the second synthetic test. It corresponds to the test case with 40 km magnetic thickness (see text for details). The average of the Z component is normalised in respect to its root mean square value</figcaption></figure>

set for its regional Legendre coefficients and calculate the corresponding regional power spectra. Finally, we conduct a misfit analysis between each regional spectrum and the statistical expression of Eq. 19. The free parameters of the statistical expression, $\epsilon, M$ and $\gamma$ were chosen to minimize the following misfit function (Voorhies et al. 2002, their Equation 14a)

$$ \chi^2 = \sum_n \left[ (\ln(R_n) - \ln(E\{R_n\}))^2 \right], \qquad (20) $$

where $R_n$ is the regional power spectrum of the synthetic data and $E\{R_n\}$ its expected value as given by Eq. 19.

The co-estimation of the three unknown parameters leads to some uncertainty as they are not easy to separate (see also Bouligand et al. 2009; Maus et al. 1997). Therefore, we set the susceptibility power law to 1.36 for the whole globe (in accordance to the preferred estimated mean global value of Thébault and Vervelidou 2015) and compute the value of $\chi^2$ for a variety of values for the parameters $\epsilon$ and $M$. The magnetic thickness $\epsilon$ is allowed to vary between 0 and 80 km and the magnetization $M$ between 0 A/m and 2 A/m. This investigation provides us with the 2D function $\chi^2(\epsilon, M)$, whose minimum corresponds to the maximum likelihood values for the mean magnetic thickness and the mean magnetization.

The mean value of the 2000 regional mean magnetic thickness estimations is 39.3 km, with a standard deviation due to regional variability of 4.3 km. The mean value for the mean magnetization is 0.94 A/m with a standard deviation of 0.2 A/m. The expected values lie within the uncertainty of the estimated values, despite a small bias in the estimated maximum likelihood.

We further test the accuracy of the algorithm for two other global synthetic data sets. For these two cases, the magnetic thickness are 10 and 100 km, and the other

parameters are left unchanged. As before, we verify that the Z component averages to zero. The mean estimated values are 10.8 ± 2.3 km and 0.91 ± 0.23 A/m, and 92.2 ± 9.9 km and 0.95 ± 0.2 A/m, for the 10 km and the 100 km test cases, respectively. Again, despite the small bias, we observe that the error bars span the expected values. Figure 3 shows the 2000 regional spectra, in grey, together with their mean (blue curve) and the realization of the statistical expression of Eq. 19 based on the estimated values (red curve), for all three test cases. We see that for all cases the blue and the red curves follow closely each other.

From this figure we also see that the maximum of the power spectrum moves to higher degrees as the magnetic thickness reduces (see also Thébault and Vervelidou 2015). Moreover, the maximum of the spectrum becomes smoother and spreads over a larger bandwidth. We estimate the 10 km magnetic thickness to be the lower limit of our resolution. Observational spectra of higher harmonic degree would be necessary to detect thinner magnetic layers as discussed by Thébault and Vervelidou 2015.

## Global models of the Earth's mean magnetic thickness and magnetization

### Data

We follow the same methodology as in section “Synthetic analyses” to infer the large-scale mean magnetic thickness and magnetization for Earth. To achieve this, we need global magnetic field measurements to a high spatial resolution. Many aeromagnetic, marine, and satellite magnetic measurements have been collected worldwide and merged together into the World Digital Magnetic Anomaly Map (WDMAM) (see Khorhonen et al. 2007 for its first version and Dyment et al. 2015a, 2015b for the recently published version 2.0, which can be downloaded
<figure><img src="image_4.png" /><figcaption>**Fig. 3** The spectra of the second synthetic test. From top to bottom: The 2000 regional spectra (grey curves), their average (blue curve) and the realization of the statistical expression based on the values estimated through the misfit analysis (red curve), for the 10, 40, and 100 km test case</figcaption></figure>
at www.wdmam.org). The WDMAM results from vast international efforts to represent globally lithospheric magnetic field structures at a wide range of spatial scales. Three different grids were built as candidates for the first edition of the WDMAM (Hemant et al. 2007; Hamoudi et al. 2007; Maus et al. 2007). Differences between these candidate grids, which rely on different processing and merging techniques, show that the WDMAM, although self-consistent, is imperfect. An attempt to improve, this grid was proposed by Maus et al. 2009 who constructed the EMAG2 grid using updated aeromagnetic and marine compilations and adding a priori information in the oceanic domain. We choose not to use the WDMAM or the EMAG2 scalar grid to compute the regional power spectra, as the computation of regional models from scalar data only proved to be unstable near the equatorial regions. This is a result of the Backus effect and of the fact that the radial component of Earth's main magnetic field cancels out at the magnetic equator. Therefore, we use the NGDC-720 model (Maus 2010), built upon the EMAG2 grid, for which this problem was addressed by a global regularization. The NGDC-720 model is a hybrid construction that uses the MF6 lithospheric magnetic field model (Maus et al. 2008) derived from CHAMP satellite (Reigber et al. 2002) measurements for SH degrees 16 to 130 and then the EMAG2 grid up to SH degree 720. This model therefore offers a maximum horizontal spatial resolution of about 50 km. We refer to Thébault and Vervelidou 2015 for a discussion on the limitations of the EMAG2 grid and the NGDC-720 model for carrying out spectral analyses.

Using the NGDC-720 Gauss coefficients, we calculate the vector magnetic field over an equal area grid. We verify, as before, that the average of the Z component over all spherical caps is close to zero. All the parameters involved

in the inversion of the NGDC-720 model and the subsequent misfit analysis are kept equal to those used for the synthetic tests. The histogram of the residuals between the NGDC-720 values and the R-SCHA estimations is shown in Fig. 4. The distribution of the residuals is centered on zero and the residual mean square for the 2000 caps lie between 4.9 and 51.3 nT. Figure 5 shows the 2000 regional spectra of the NGDC-720 model (grey lines), together with their average (red line) and the NGDC-720 SHA power spectrum (blue line). Despite the expected regional variability, we see that the average of the regional spectra follows closely the global NGDC-720 power spectrum, thus giving us confidence in the relevance of the computation.

### A quality index for the global maps

The spectral content of the NGDC model is characterised by a high regional variability due to gaps in the current aeromagnetic and marine coverage of the Earth but also due to differences in the acquisition and subsequent processing of the various implicated datasets. This has a direct impact on the reliability of the regional estimates, which results in a high regional variability. Therefore, we compute a quality index for our maps. This quality estimator, $QI$, has been chosen to be (after Voorhies et al. 2002, see their Eq. 14b)

$$QI [\%] = \frac{100}{e^{\frac{\chi}{\sqrt{N}}}}, \qquad (21)$$

where $\chi$ is given by Eq. 20 and $N$ are the number of the terms of the regional spectrum.

When $QI = 100\%$, the fit is perfect whereas $QI = 50\%$, for instance, indicates that the terms of the statistical spectrum are half or twice as large than the corresponding terms of the observational spectrum. Figure 6 shows the

<figure><img src="image_5.png" /><figcaption>Fig. 4 The histogram of the residuals of the inversion of the NGDC-720 model</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Fig. 5** Results of the spectral analysis of the NGDC-720 model. The 2000 regional spectra of the NGDC-720 model (grey curves), their average (red curve) and the NGDC-720 SHA power spectrum (blue curve)</figcaption></figure>

geographical distribution of *QI* estimated for each spherical cap. For most of the caps (80%) *QI* ≥ 67%, but there are also some caps (19%) for which 50% < *QI* < 67%. Lower *QI* characterises only 1% of the caps, and these correspond only to regions covered in Fig. 6 by white dots, which are regions over which EMAG2 is void of

near-surface data. Over these regions, the values of the NGDC-720 model are based on the satellite-based model MF6 (Maus et al. 2008), (see Maus 2010, paragraph 3.1). We note, however, that not all regions covered by data are based on near-surface measurements, as gridding and extrapolation techniques were used to produce the

<figure><img src="image_7.png" /><figcaption>**Fig. 6** The map of the quality index *QI*. Large values of *QI* indicate reliable magnetic thickness and magnetization estimates. The regional spectra of the four caps, whose centre is marked with a yellow star, are shown in Fig. 7. The regional spectra of the four caps, whose centre is marked with a pink star, are shown in Fig. 8. The white dots correspond to regions over which EMAG2 is void of near-surface data (see text for details)</figcaption></figure>
EMAG2 grid (Maus et al. 2009), especially over the oceans.

The largest *QI* values correlate well with regions where aeromagnetic coverage is dense and of high resolution, i.e., southern Africa, North America, and Australia. Figure 7 shows a spectrum over each of these regions, together with a spectrum of a cap over the Atlantic ocean, also characterised by a large *QI*. The exact location of the center of these caps is marked with a yellow star on Fig. 6. We see that these spectra agree well with the statistical expression. The spectra are presented in a descending order in respect to the estimated mean magnetic thickness. The spectrum of the cap over North America presents a pick at SH degree lower than 100, while the picks of the other spectra lie at larger SH degrees. The spectrum over southern Africa shows a maximum at a lower SH degree than the spectrum over Australia but is also less energetic. This leads to a larger mean magnetic thickness and a lower mean magnetization estimate for the cap over southern Africa than for the cap over Australia. We see here an example of how the magnetic thickness and the magnetization information can be disentangled.

The pink stars on Fig. 6 show the center of four caps for which $QI <= 60\%$. The corresponding spectra are shown in Fig. 8.

The regional spectrum over Antarctica, a region where EMAG2 has a sparse aeromagnetic coverage (see Fig. 7 of Maus et al. 2009), presents an abrupt drop of energy after degree 130. This pattern is similar to the cut-off illustrated in Fig. 1. We attribute the drop in power to the lack of near-surface data. The immediate consequence of this so-called spectral gap is that the large degrees are artificially low. This drop off causes an underestimation of the mean magnetization. At the same time, the comparatively large energy at low degrees leads to an overestimation of the magnetic thickness. Consequently, estimations from such regions should not be considered.

Another example of a poor fit is around the North Pole. Figure 8 shows the spectrum of a cap located over central Greenland. We see again a drop of energy for degrees larger than 130. However, Greenland has good aeromagnetic coverage (see again Fig. 7 of Maus et al. 2009), suggesting this feature could be genuine. This would indicate the statistical expression for the power spectrum is locally not in agreement with the statistical properties of the magnetic sources. Alternatively it could indicate a poor spectral content of the respective aeromagnetic grids. Interestingly, this spectrum is also characterised by a drop of energy between degrees 50 and 100. This, again, could be a genuine feature. But, it could also be related to how the satellite data were corrected for the external fields along the satellite orbits or how the underlying MF6 model was regularized (see, e.g., Lesur and Maus 2006; Maus et al. 2008; Thébault et al. 2012). In this respect,

we note that over the same cap satellite-based magnetic field models more recent than MF6 (Maus et al. 2008) e.g., MF7, GRIMM (Lesur et al. 2013), CM5 (Sabaka et al. 2015) are more energetic over the same bandwidth, although still less energetic than predicted by the statistical expression.

Another continental example of mismatch is found in Russia, over the broader area around the Kursk magnetic anomaly. For this major magnetic anomaly (see, e.g., Taylor et al. 2003), the underlying hypotheses of the statistical expression may not hold.

A final example of a spectrum with low *QI* is for a cap located over the Eastern Pacific Ocean. This spectrum is in clear disagreement with the form given by the statistical expression. Further caps, mainly located over oceanic regions, exhibit this behaviour. Testing different values for the power law $\gamma$ does not improve the fit. This mismatch could result from the extrapolation scheme used to fill in oceanic areas void of direct marine measurements (see Maus et al. 2009).

The above examples illustrate the spectral variability of the NGDC-720 model, which depends on a multitude of factors. Lack of aeromagnetic and marine data, incomplete data processing, and/or genuine geophysical features are all possible sources of incompatibility between the observational spectra and the statistical expression. Clearer information about the underlying datasets would help to discern the origin of the mismatch on a case-by-case basis. In any case, for the purposes of our study, we find that the *QI* index suitably accounts for various possible sources of incompatibility, and we use it as an indicator for the degree of reliability of the magnetic thickness and magnetization estimates.

## The maps

Figures 9 and 10 show maps of the estimated mean magnetic thickness and mean magnetization, respectively, where mean stands for the average value over the area of the spherical cap. These maps should be interpreted in tandem with Fig. 6. On the magnetization map, brown dots indicate the center of the caps for which a thickness of less than 10 km was estimated (see section "Synthetic analyses" regarding the depth resolution limit).

The point-to-point product of these two maps corresponds to the vertically integrated magnetization (VIM) and is given in Fig. 11a. VIM is a quantity calculated in studies based solely on satellite data, as the effect of magnetic thickness and magnetization on satellite measurements are not fully separable (see, e.g., section 6.4 in Thébault et al. 2010). While near-surface data gaps do not allow our magnetic thickness and magnetization estimates to be globally reliable, satellite data are globally available. We therefore expect to be able to estimate the VIM on a global scale based only on the satellite part
<figure><img src="image_8.png" /><figcaption>**Fig. 7** Examples of spectra for caps with a high *Ql*. The regional spectra (blue curves), together with their best fit realization of the statistical expression (red curves), of the four caps marked with a yellow star on Fig. 6. From the top to the bottom: The cap over North America (*Ql* = 73.3, *ϵ* = 53 km, *M* = 0.4 A/m), Southern Africa (*Ql* = 79, *ϵ* = 37 km, *M* = 0.3 A/m), Australia (*Ql* = 75.7, *ϵ* = 18 km, *M* = 0.7 A/m) and Atlantic ocean (*Ql* = 82.3, *ϵ* = 17 km, *M* = 0.4 A/m). See the text for a discussion</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>**Fig. 8** Examples of spectra for caps with a low *Ql*. The regional spectra (blue curves), together with their best fit realization of the statistical expression (red curves), of the four caps marked with a pink star on Fig. 6. From the top to the bottom: The cap over Antarctica (*Ql* = 32.9), Greenland (*Ql* = 60.6), Pacific ocean (*Ql* = 57.4) and Russia (*Ql* = 58.6). See the text for a discussion</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>**Fig. 9** The mean magnetic thickness map. The map should be read together with the quality index map given at Fig. 6</figcaption></figure>

of NGDC-720. We follow the same misfit procedure as before, but this time, we fit the statistical expression only to the degrees between 16 and 130 of the regional spectra. The point-to-point product of the resulting magnetic thickness and magnetization is shown in Fig. 11b. Interestingly, both maps of Fig. 11 are very similar, indicating that spectral analyses of satellite data can lead to VIM estimates.

We then transform the magnetic thickness map into spherical harmonics. Figure 12 shows the resulting SH power spectrum (blue line), together with an uncertainty estimate (light blue shaded curve). This uncertainty estimate results from the combination of two types of errors. The first type is the 95% confidence interval estimated numerically by the misfit analysis. The second type corresponds to the upper and lower boundary of the estimated magnetic thickness when $\gamma = 1.29$ and $\gamma = 1.5$ are considered instead of $\gamma = 1.36$. These two $\gamma$ values correspond to the 95% confidence interval found for $\gamma$ on a global scale by (Thébault and Vervelidou 2015).

In Fig. 12, we see that the spectrum carries energy from SH degree 0 to 36. The maximum SH degree that we can recover is dictated by the spatial distribution of the spherical caps over Earth's surface. The more caps, the denser the spatial sampling of the magnetic thickness and

therefore the higher the recovered maximum degree. For our distribution of 2000 caps, the maximum SH degree is equal to 36. Whether the thickness map has energy or not up to that degree depends on the cap's half-angle, since it is this parameter that defines the actual spectral content of our map. For the extreme case of $\theta_0 = \pi$, the magnetic thickness would be a constant for the whole Earth and so only SH degree 0 would carry energy. For $\theta_0 = 15$, which is our case, the thickness map carries energy up to degree 36. Even more interestingly, it also carries energy over SH degrees 1–16, which correspond to wavelengths that are absent from the known lithospheric magnetic field due to the masking by the core field.

## Discussion

The most prominent feature of the derived magnetic thickness map (Fig. 9) is the dichotomy between the oceans and the continents. The oceans are characterised by a smaller magnetic thickness than the continents. This is in agreement with our current understanding about crustal thickness (see, e.g., the CRUST1.0 model, Pasyanos et al. 2014). This dichotomy is visible in the map as its spectral content starts from very low SH degrees.

In Fig. 12, we show the power spectra of two other models: (1) the seismic crustal model CRUST1.0
<figure><img src="image_11.png" /><figcaption>**Fig. 10** The mean magnetization map. The brown dots correspond to the centres of the caps, for which a magnetic thickness of less than 10 km has been estimated (see text for details). The map should be read together with the quality index map given at Fig. 6</figcaption></figure>

(Pasyanos et al. 2014, http://igppweb.ucsd.edu/~gabi/crust1.html) (magenta line) and a recent version (pers. comm. M. Purucker) of the magnetic thickness map based on the methodology of Purucker et al. 2002, which was applied over Antarctica by FoxMaule et al. 2005 and on a global scale by FoxMaule et al. 2009 (cyan line). The second model is based on the 3SMAC seismic model (Nataf and Ricard 1996) for the SH degrees up to 15 and on the lithospheric magnetic field model MF7 (see Maus et al. 2008 for the latest published version of the MF series) for SH degrees larger than 15. We see that our model follows very closely the spectra of the other two models from SH degree 11 up to degree 18. It is, however, somewhat more energetic than the other models from SH degree 1 up to SH degree 10. This probably results from the overestimation of the magnetic thickness over large areas such as Antarctica, Africa, and the Arctic (see section “A quality index for the global maps” for an explanation). A greater data coverage will allow more realistic estimates for the magnetic thickness, and we believe it will lead to a decrease in the energy of our power spectrum over these wavelengths.

From SH degree 19 and onwards, the spectrum of our model follows closely the spectrum of Fox Maule et al. 2009, while both of them are a little more energetic than

the spectrum of CRUST1.0. This could indicate the depth of Moho is actually shallower than the depth of the Curie isotherm over this waveband, and that parts of the upper mantle are magnetized. However, both the spectrum of the hybrid model of Fox Maule et al. 2009 and of the crustal thickness model CRUST1.0 lie within the errors of the spectrum of our model, indicating that these two models are statistically equally energetic to ours.

Comparing our magnetic thickness and magnetization maps, we observe differences in their main features. There are three particularly notable examples. Australia's magnetic thickness pattern is almost reversed with respect to its magnetization pattern. The low magnetization feature over China and parts of Indonesia does not correspond to a particular feature on the magnetic thickness map. The same is true for the low magnetization feature that covers almost all of South America. Interestingly, the low magnetization pattern over South America has also been observed by Purucker and Whaler 2004 (see also Purucker and Whaler 2015). Due to their methodology (they transform the MF3 lithospheric field model (Maus et al. 2006) to a minimum magnetization map (see Whaler and Langel 1996), while considering the magnetic thickness to be constant), their magnetization map is expected to resemble more closely our VIM map rather than our
<figure><img src="image_12.png" /></figure>

magnetization map. Indeed, comparing their magnetiza-
tion map (see Fig. 10, Purucker and Whaler 2015) with
our VIM map (Fig. 11), we find common features that
include the signatures of the Bangui and Kursk magnetic
anomalies, but also, e.g., an anomaly over central Aus-
tralia and a dichotomy in the signals of North and South
America. Figure 11b also shows an apparent dichotomy
between West and East Antarctica. This feature could
have important implications for inferring the thermal
state of the crust in Antarctica and is in close agreement
with previous works (see, for instance Fig. 1 in Fox Maule
et al. 2005).

Conclusions

In this study, we present the surface spherical cap har-
monic power spectrum, based on the R-SCHA modelling
method (Thébault et al. 2006). It can be directly inter-
preted in terms of wavelengths, in the same way as the
SHA power spectrum. Such a regional spectrum can have
various applications that go beyond magnetic field studies,
as long as the signal under study averages to zero over
the spherical cap surface. It can also be applied on planets
other than Earth (see for example Lewis and Simons 2012
for an application of a different regional spectrum (Dahlen
and Simons 2008; Wieczorek and Simons 2007) for Mars).

<figure><img src="image_13.png" /><figcaption>**Fig. 12** The SHA power spectrum of our magnetic thickness model. The SHA power spectrum of our magnetic thickness model (**blue line**) together with its uncertainty estimation (**light blue shaded curve**), of the crustal layer of the model CRUST1.0 (Pasyanos et al. 2014) (**magenta line**) and of a recent version of the hybrid model of (Fox Maule et al. 2009) (**cyan line**)</figcaption></figure>
Here, we use this regional power spectrum to calcu-
late global maps of the large wavelengths of the mean
magnetic thickness and mean magnetization of Earth. We
do so by combining spherical cap observational spec-
tra with a recently proposed statistical expression for the
power spectrum of Earth's lithospheric magnetic field
(Thébault and Vervelidou 2015). This expression depends
on the mean magnetic thickness, the mean magnetiza-
tion and a power law of the susceptibility power spec-
trum. The regional observational spectra are obtained
by inversion of the NGDC-720 lithospheric field model
(Maus 2010) over a series of spherical caps, covering the
whole Earth.

The NGDC-720 model relies on information from satel-
lite, aeromagnetic, and marine magnetic measurements.
The resulting regional power spectra therefore cover a
broad range of wavelengths. The first outcome of our
study is that by exploiting simultaneously satellite and
near surface information, separate estimates of the mean
magnetization and mean magnetic thickness are provided.

However, our estimates are not equally reliable glob-
ally. For this reason, we accompany them with a quality
index, called *QI* (see Fig. 6). This variability stems mainly
from incomplete, and inhomogeneous in terms of quality
and resolution, near-surface data coverage. For regions,
however, with dense coverage of good quality near-surface
data, a low *QI* can act as an indicator that the regional
properties of the magnetic sources are different from
the assumptions underlying the statistical expression of
Eq. 15, e.g., the assumption of purely induced magneti-
zation. Consequently, a second outcome of our study, as
illustrated through Figs. 7 and 8, is that the R-SCHA,
and its power spectrum offer a way to study and evaluate
regionally the spectral content of data compilations and
models over a large bandwidth. This can be very useful
for regional magnetic studies. The WDMAM project (see
Dyment et al. 2015a, 2015b for its recently published ver-
sion 2.0, which can be downloaded at www.wdmam.org)
and regional efforts to enrich the available near-surface
data sets (e.g., Gaina et al. 2011; Golynsky et al. 2013),
provide currently a favourable framework for such studies.

Maps of magnetic thickness and magnetization can be
valuable in understanding further the temperature, com-
position, and structure of Earth's lithosphere. Efforts to
produce such maps are hindered because of the masking
effect of the core field on the large scales of the litho-
spheric magnetic field. A third outcome of our study, as
illustrated by Fig. 12, is the feasibility of estimating the
large scales of the magnetic thickness and magnetization
using magnetic data. Keeping in mind the uncertainties
involved, we conclude from Fig. 12 that the large scales
of the magnetic thickness are on average confined to a
layer that does not exceed the Moho. Updates to our
maps (Figs. 9 and 10), and consequently to the SH power

spectrum shown in Fig. 12, are especially welcome, as
the worldwide compilation of magnetic data gets richer,
both in the spatial and the spectral domain. Both the
SWARM satellite mission (see, e.g., Olsen et al. 2013) and
the continuous updating of the WDMAM move in this
direction.

Competing interests

The authors declare that they have no competing interests.

Authors' contributions

FV and ET designed the methodology and analysed the results. FV carried out the numerical analyses and drafted the manuscript. ET contributed to the drafting of the manuscript. Both authors read and approved the final manuscript.

Acknowledgements

The authors would like to thank three anonymous reviewers, the handling editor Prof. Toshi Yamazaki and the editor in chief, Prof. Yasuo Ogawa, for their constructive comments that helped ameliorate the manuscript. V. Lesur for useful discussions and M. Brown for an early review of the manuscript. FV was partly funded by Région Ile-de-France and partly by the DFG SPP1488 PlanetMag. ET was partly funded by the Centre National des Etudes Spatiales (CNES) within the context of the project of the "Travaux préparatoires et exploitation de la mission Swarm" and partly by INSU through the "Programme National de Planétologie". Figures 6, 9, 10 and 11 were drawn with the Generic Mapping Tools (Wessel and Smith 2001).

Author details

Received: 30 June 2015 Accepted: 14 September 2015

Published online: 26 October 2015

References

Backus G, Parker R, Constable CG (1996) Foundations of geomagnetism. Cambridge University Press, Cambridge

Bouligand C, Glen JMG, Blakely RJ (2009) Mapping Curie temperature depth in the western United States with a fractal model for crustal magnetization. J Geophys Res 114. (B11). doi: 10.1029/2009JB006494

Blakely RJ (1996) Potential theory in gravity and magnetic applications. Cambridge University Press, Cambridge

Dahlen FA, Simons FJ (2008) Spectral estimation on a sphere in geophysics and cosmology. Geophys J Int 174(3):774-807

Dyment J, Lesur V, Hamoudi M, Choi Y, Thébault E, Catalan M, The WDMAM Task Force, The WDMAM Evaluators, The WDMAM Data Providers: World Digital Magnetic Anomaly Map version 2.0 (2015a) The WDMAM Task Force: Dyment J (chair), Catalan M (co-chair), de Santis A, Hamoudi M, Ishihara T, Korhonen J, Lesur V, Litvinova T, Luis J, Meyer B, Milligan P, Nakanishi M, Okuma S, Pilkington M, Purucker M, Ravat D, Thébault E. (alphabetical order). The WDMAM Evaluators: Gaina C, Luis J, Maus S, Meyer B, Nakanishi M, Purucker M, Quesnel Y, Saltus R, Taylor P. (alphabetical order). Avalaible at http://www.wdmam.org

Dyment J, Choi Y, Hamoudi M, Lesur V, Thébault E (2015b) Global equivalent magnetization of the oceanic lithosphere. Earth Planet Sci Lett 430:54–65. doi:10.1016/j.epsl.2015.08.002

Ferré EC, Friedman SA, Martín-Hernández F, Feinberg JM, Conder JA, Ionov DA (2013) The magnetism of mantle xenoliths and potential implications for sub-moho magnetic sources. Geophys Res Lett 40(1):105–110

Ferré EC, Friedman SA, Martín-Hernández F, Feinberg JM, Till JL, Ionov DA, Conder JA (2014) Eight good reasons why the uppermost mantle could be magnetic. Tectonophysics 624:3–14

Friedman SA, Feinberg JM, Ferré EC, Demory F, Martín-Hernández F, Conder JA, Rochette P (2014) Craton vs. rift uppermost mantle contributions to magnetic anomalies in the United States interior. Tectonophysics 624:15–23

<footnote n="1">¹Helmholtz Centre Potsdam-GFZ German Research Centre for Geosciences, Section 2.3 Earth's magnetic field, Telegrafenberg, 14473 Potsdam, Germany.</footnote>

<footnote n="2">²UMR CNRS 6112, University of Nantes, Laboratoire de Planétologie et de Géodynamique, 2, rue de la Houssinière BP 92208, 44322 Nantes Cedex 3, France.</footnote>
Fox Maule C, Purucker ME, Olsen N, Mosegaard K (2005) Heat flux anomalies in Antarctica revealed by satellite magnetic data. Science 309(5733):464–467

Fox Maule C, Purucker ME, Olsen N (2009) Inferred magnetic crustal thickness and geothermal heat flux from crustal magnetic field models. Technical Report 9, Danish Climate Center

Gaina C, Werner SC, Saltus R, Maus S (2011) Circum-Arctic mapping project: new magnetic and gravity anomaly maps of the Arctic. Geol Soc Lond Mem 35(1):39–48

Golynsky A, Bell R, Blankenship D, Damaske D, Ferraccioli F, Finn C, Group AW (2013) Air and shipborne magnetic surveys of the Antarctic into the 21st century. Tectonophysics 585:3–12

Haines GV (1985) Spherical cap harmonic analysis. J Geophys Res 90(B3):2583–2591

Hamoudi M, Thébault E, Lesur V, Mandea M (2007) GeoForschungsZentrum Anomaly Magnetic Map (GAMMA): a candidate model for the World Digital Magnetic Anomaly Map. Geochem Geophys Geosyst 8(6). doi:10.1029/2007GC001638

Hemant K, Thébault E, Mandea M, Ravat D, Maus S (2007) Magnetic anomaly map of the world: merging satellite, airborne, marine and ground-based magnetic data sets. Earth Planet Sci Lett 260:56–71

Hulot G, Sabaka TJ, Olsen N, Fournier A (2015) 5.02 - the present and future geomagnetic field. In: Schubert G (ed). Treatise on Geophysics (Second Edition). Elsevier, Oxford. pp 33–78

Hwang C, Chen SK (1997) Fully normalized spherical cap harmonics: application to the analysis of sea-level data from TOPEX/POSEIDON and ERS-1. Geophys J Int 129(2):450–460

Jackson A (1990) Accounting for crustal magnetization in models of the core magnetic field. Geophys J Int 103(3):657–673

Jackson A (1994) Statistical treatment of crustal magnetization. Geophys J Int 119(3):991–998

Khorhonen JK, Fairhead JD, Hamoudi M, Hemant K, Lesur V, Mandea M, Maus S, Purucker ME, Ravat D, Sazonova T, Thébault E (2007) Magnetic Anomaly Map of the World-Carte des Anomalies Magnétiques du Monde. Scale: 1:50,000,000, 1st edition

Korte M, Constable CG, Parker RL (2002) Revised magnetic power spectrum of the oceanic crust. J Geophys Res: Solid Earth 107(B9). doi:10.1029/2001JB001389

Langel RA, Estes RH (1982) A geomagnetic field spectrum. Geophys Res Lett 9(4):250–253

Langel RA, Hinze WJ (1998) The Magnetic Fiel of the Earth's Lithosphere: the Satellite Perspective. Cambridge University Press, Cambridge

Lesur V, Maus S (2006) A global lithospheric magnetic field model with reduced noise level in the polar regions. Geophys Res Lett 33(13). doi:10.1029/2006GL025826

Lesur V, Rother M, Vervelidou F, Hamoudi M, Thébault E (2013) Post-processing scheme for modelling the lithospheric magnetic field. Solid Earth 4:105–118

Lewis KW, Simons FJ (2012) Local spectral variability and the origin of the martian crustal magnetic field. Geophys Res Lett 39(18)

Li C, Jian W, Jian L, Tingting W (2013) Thermal evolution of the North Atlantic lithosphere: new constraints from magnetic anomaly inversion with a fractal magnetization model. Geochem Geophys Geosyst 14(12):5078–5105

Lowes FJ (1974) Spatial power spectrum of the main geomagnetic field, and extrapolation to the core. Geophys J Int 36(3):717–730

Love JA, Constable CG (2003) Gaussian statistics for palaeomagnetic vectors. Geophys J Int 152(3):515–565

Martín-Hernández F, Ferré EC, Friedman SA (2014) Remanent magnetization in fresh xenoliths derived from combined demagnetization experiments: magnetic mineralogy, origin and implications for mantle sources of magnetic anomalies/rectangular harmonic analysis revisited. Tectonophysics 624:24–31

Maus S, Gordon D, Fairhead D (1997) Curie-temperature depth estimation using a self-similar magnetization model. Geophys J Int 129(1):163–168

Maus S, Rother M, Hemant K, Stolle C, Lühr H, Kuvshinov A, Olsen N (2006) Earth's lithospheric magnetic field determined to spherical harmonic degree 90 from CHAMP satellite measurements. Geophys J Int 164(2):319–330

Maus S, Sazonova T, Hemant K, Fairhead JD, Ravat D (2007) National Geophysical Data Center candidate for the World Digital Magnetic Anomaly Map. Geochem Geophys Geosyst 8(6). doi:10.1029/2007GC001643

Maus S, Yin F, Lühr H, Manoj C, Rother M, Rauberg J, Michaelis I, Stolle C, Müller RD (2008) Resolution of direction of oceanic magnetic lineations by the sixth-generation lithospheric magnetic field model from CHAMP satellite magnetic measurements. Geochem Geophys Geosyst 9(7). doi:10.1029/2008GC001949

Maus S, Barckhausen U, Berkenbosch H, Bournas N, Brozena J, Childers V, Dostaler F, Fairhead JD, Finn C, Von Frese RRB, Gaina C, Golynsky S, Kucks R, Lühr H, Milligan P, Mogren S, Müller RD, Olesen O, Pilkington M, Saltus R, Schreckenberger B, Thébault E, Caratori Tontini F (2009) A 2-arc min resolution Earth Magnetic Anomaly Grid compiled from satellite, airborne, and marine magnetic measurements. Geochem Geophys Geosyst 10(8). doi:10.1029/2009GC002471

Maus S (2010) An ellipsoidal harmonic representation of Earth's lithospheric magnetic field to degree and order 720. Geochem Geophys Geosyst 11(6). doi:10.1029/2010GC003026

Nataf HC, Ricard Y (1996) 3SMAC: an a priori tomographic model of the upper mantle based on geophysical modeling. Phys Earth Planet Inter 1:101–122

O'Brien MS, Parker RL, Constable CG (1999) Magnetic power spectrum of the ocean crust on large scales. J Geophys Res 104(B12):29189–29201

Olsen N, Friis-Christensen E, Floberghagen R, Alken P, Beggan CD, Chulliat A, Doornbos E, Teixeira da Encarnação J, Hamilton B, Hulot G, van den IJssel J, Kuvshinov A, Lesur V, Lühr H, Macmillan S, Maus S, Noja M, Olsen PEH, Park J, Plank G, Püthe C, Rauberg J, Ritter P, Rother M, Sabaka TJ, Schachtschneider R, Sirol O, Stolle C, Thébault E, Thomson AWP (2013) The Swarm Satellite Constellation Application and Research Facility (SCARF) and Swarm data products. Earth Planets Space 65(11):1189–1200

Purucker ME, Whaler KA (2015) 5.06 - crustal magnetism. In: Schubert G (ed). Treatise on Geophysics (Second Edition). Elsevier, Oxford. pp 185–218

Pasyanos E, Masters TG, Laske G, Ma Z (2014) LITHO1.0: An updated crust and lithospheric model of the Earth. J Geophys Res: Solid Earth 119(3):2153–2173

Purucker ME, Langlais B, Olsen N, Hulot G, Mandea M (2002) The southern edge of cratonic North America: Evidence from new satellite magnetometer observations. Geophys Res Lett 29(15). doi:10.1029/2001GL013645

Purucker ME, Whaler KA (2004) Recognizing and interpreting the longest wavelength lithospheric magnetic fields obscured by overlap with the core field. Eos Trans Am Geophys Union 85(47). Fall Meeting Supplement Abstract GP31A-0821

Reigber C, Lühr H, Schwintzer P (2002) CHAMP mission status. Adva Space Res 30(2):129–134

Sabaka TJ, Olsen N, Tyler RH, Kuvshinov A (2015) CM5, a pre-Swarm comprehensive geomagnetic field model derived from over 12 yr of CHAMP, Ørsted, SAC-C and observatory data. Geophys J Int 200(3):1596–1626

Taylor PT, Frawley JJ, Kim HR, von Frese R, Kim JW (2003) Comparing Magsat, Ørsted and CHAMP crustal magnetic anomaly data over the Kursk Magnetic Anomaly, Russia. First CHAMP Mission Results for Gravity Magnetic Atmospheric Studies. Springer, Berlin Heidelberg. pp 302–308

Thébault E, Schott JJ, Mandea M, Hoffbeck JP (2004) A new proposal for spherical cap harmonic modelling. Geophys J Int 159:83–103

Thébault E, Schott JJ, Mandea M (2006) Revised spherical cap harmonic analysis (R-SCHA): Validation and properties. J Geophys Res 111(B1). doi:10.1029/2005JB003836

Thébault E (2008) A proposal for regional modelling at the Earth's surface, R-SCHA2D. Geophys J Int 174(1):118–134. doi:10.1111/j.1365-246X.2008.03823.x

Thébault E, Purucker ME, Whaler KA, Langlais B, Sabaka TJ (2010) The magnetic field of the Earth's lithosphere. Space Sci Rev 155(1):95–127

Thébault E, Vervelidou F, Lesur V, Hamoudi M (2012) The along-track satellite analysis in planetary magnetism. Geophys J Int 188(3):891–907

Thébault E, Vervelidou F (2015) A statistical spatial power spectrum of the Earth's lithospheric magnetic field. Geophys J Int 201(2):605–620

Voorhies C (1998) Elementary theoretical forms for the spatial power spectrum of Earth's crustal magnetic field. Technical Report 208608, National Aeronautics and Space Administration, Goddard Space Flight Center

Voorhies CV, Sabaka TJ, Purucker ME (2002) On magnetic spectra of Earth and Mars. J Geophys Res: Planets 107(E6). doi:10.1029/2001JE001534

Voorhies CV (2008) Thickness of the magnetic crust of Mars. J Geophys Res: Planets 113(E4). doi:10.1029/2007JE002928

Vervelidou F (2013) Contribution à la modélisation et à l'interprétation multi-échelle du champ magnétique de la lithosphère terrestre. PhD thesis, Université Paris Diderot - Paris 7
Wasilewski P, Mayhew MA (1982) Crustal xenolith magnetic properties and long wavelength anomaly source requirements. Geophys Res Lett 9(4):329-332

Wessel P, Smith WH (2001) The Generic Mapping Tools. [2013-01-01]
http://gtm.soest.hawaii.edu

Whaler KA, Langel RA (1996) Minimal crustal magnetizations from satellite data. Phys Earth Planetary Interiors 98(3):303-319

Wieczorek MA, Simons FJ (2005) Localized spectral analysis on the sphere. Geophys J Int 162(3):655-675

Wieczorek MA, Simons FJ (2007) Minimum-variance spectral analysis on the sphere. J Fourier Anal Appl 13(6):665-692

**Submit your manuscript to a SpringerOpen® journal and benefit from:**

► Convenient online submission

► Rigorous peer review

► Immediate publication on acceptance

► Open access: articles freely available online

► High visibility within the field

► Retaining the copyright to your article

Submit your next manuscript at ► springeropen.com