Applications for CryoSat2 satellite magnetic data in studies of Earth's core field variations

We use 20 years of magnetic field measurements from the Oersted, CHAMP and Swarm satellite missions, supplemented by calibrated platform magnetometer data from the CryoSat2 satellite, to study time variations of the Earth's core field at satellite altitude and at the core-mantle boundary (CMB). From the satellite data we derive composite time series of the core field secular variation (SV) with 4month cadence, at 300 globally distributed Geomagnetic Virtual Observatories (GVO). GVO radial SV series display regional fluctuations with 5-10 years duration and amplitudes reaching 20 nT/yr, most notably at low latitudes over Indonesia (2014), over South America and the South Atlantic (2007, 2011 and 2014), and over the central Pacific (2017). Applying the Subtractive Optimally Localized Averages (SOLA) method, we map the SV at the CMB as a collection of locally averaged SV estimates. We demonstrate that using 2-year windows of CryoSat2 data, it is possible to reliably estimate the SV and its time derivative, the secular acceleration (SA), at the CMB, with a spatial resolution, corresponding to spherical harmonic degree 10. Along the CMB geographic equator, we find strong SA features under Indonesia from 2011-2014, under central America from 2015 to 2019, and sequences of SA with alternating sign under the Atlantic during 2004-2019. We find that data from CryoSat2 make a valuable contribution to the emerging picture of sub-decadal core field variations. Using 1 year windows of data from the Swarm satellites, it is possible to study SA changes at low latitudes on timescales down to 1 year, with spatial resolution corresponding to spherical harmonic degree 10. We find strong positive and negative SA features appearing side-by-side in the Pacific in 2017, and thereafter drift westward.


Introduction
The main part of the Earth's magnetic field is generated by motions in the electrical conducting liquid outer core, in a process known as the geodynamo. This magnetic field, termed the core field, exhibits both spatial and temporal changes over a broad range of scales. Magnetic measurements from satellites have increased the recovery of small-scale features of this field and revealed rapid changes in its temporal be-haviour (e.g. Alken et al., 2020a;Baerenzung et al., 2020;Finlay et al., 2020;Ropp et al., 2020;Sabaka et al., 2020). From ground and space magnetic observations, variations in the first and second time derivatives of the field, termed the secular variation (SV) and acceleration (SA), respectively, may now be resolved down to periods of about 1 to 2 years (Lesur et al., 2010(Lesur et al., , 2017Ropp et al., 2020).
Studies have revealed oscillating SA pulse-like field features at the core-mantle boundary (CMB) focused in the region around the geographical equator (Alken et al., 2020a;Chulliat and Maus, 2014;Chulliat et al., , 2015Sabaka et al., 2018). The interpretation and geophysical mechanisms responsible for driving such distinctive behaviour in the SA signal remains under debate (e.g. Aubert and Finlay, 2019;Buffett and Matsui , 2019;Gerick et al., 2020;Gillet, 2019), as is the connection to abrupt changes in the SV observed at ground observatories (Mandea et al., 2010). The secular acceleration must be characterized with care, paying attention to those spatial and temporal scales that are well resolved, as its observed spatial spectra at the core surface is blue, showing increasing power with spherical harmonic degree, and its observed temporal spectra seems to be rather flat, meaning that there could be important unresolved fast variations (Bouligand et al., 2016;Christensen et al., 2012;Gillet, 2019;Lesur et al., 2017). In this respect, assessing the limitations of the information that may be obtained from the measurements by analysing their resolving power becomes crucial when seeking to investigate the SA signal at small length scales and short time scales.
Magnetic field measurements from low Earth orbiting (LEO) satellites provide global observations of the field, which have proved important for mapping the spatial structures of the core field signal (Olsen and Stolle, 2012). Since the launch of the Danish Ørsted (1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014) satellite, the German CHAMP (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010) satellite and the European Swarm (2013-) satellite trio, satellites have provided high quality magnetic field measurements, which have enable investigations into the spatio-temporal variations of field. Unfortunately, the CHAMP mission ended in September 2010 and since reliable vector measurements from the Ørsted satellite extend only up to 2006, there is a gap from 2010 and 2014 in the satellite magnetic records (Finlay et al., 2016). However, other satellite missions, not dedicated to measuring the magnetic field, may offer a possibility to fill in this gap adding information about the field. In particular, the CryoSat-2 (2010-) mission, intended for measuring polar ice thickness, carries three platform magnetometers for navigational purposes. Calibrated CryoSat-2 measurements, where vector fluxgate magnetometer readings have been transformed into reliable magnetic field vector outputs, from August 2010 to December 2018 are now available , such that there are in total 20 years of continuous satellite measurements. With the data from CryoSat-2 now available, it is important to assess the quality and limitations of such calibrated platform magnetometer measurements and to test what contribution they can make to studying core field variations.
The standard approach for using satellite magnetometer data for core field studies is to construct spherical harmonic (SH) field models by least-squares inversion methods (see e.g. Langel , 1987). In such global models B-splines are often used to parameterize the model time-dependence. Typically this necessitates temporal regularization which modifies the time-dependence of the harmonics in a non-uniform manner.
An undesirable consequence is that for the higher harmonics, the first time derivative effectively becomes a time average over an increasingly long interval rather than an estimate of the instantaneous secular variation (Olsen et al., 2009). Cryosat-2 data has already been used in the construction of such time-dependent spherical harmonic field models by Alken et al. (2020a); Finlay et al. (2020); Kloss et al. (2021).
In an effort to access more detailed information on the spatial and temporal structure of the field, it is also of interest to look at alternative techniques for studying secular variation that can complement the traditional SH approach. In this paper we focus on two different local methods for studying core field variations as recorded in satellite measurements, with a particular focus on assessing the quality and resolving ability of CryoSat-2 magnetic data.
In a first assessment of the quality of the CryoSat-2 data and its ability to map the SV field at satellite altitude together with CHAMP and Swarm we use the Geomagnetic Virtual Observatory (GVO) technique developed by Mandea and Olsen (2006); Olsen and Mandea (2007). This technique involves computing time series of field estimates at specified target locations at satellite altitude, from satellite measurements taken nearby. We apply the processing algorithm recently developed to derive four-monthly Swarm GVO data series  and derive time series on a global grid of 300 GVOs. This network of GVOs allows field changes at satellite altitude to be investigated globally at fixed locations. The GVOs provide a useful compression of satellite magnetic measurements and are a convenient dataset for workers wishing to use constraints for studies of core dynamics. They provide an alternative to fitting core dynamics models to conventional spherical harmonic field models. They involve a series of independent local constraints rather than global constraints, and have well understood error covariances that can be assigned by methods similar to those used with ground observatories. GVOs have already been used by a number of groups for studies of core dynamics (e.g. Barrois et al., 2018;Domingos et al., 2019;Whaler and Beggan, 2015).
In the second part of this study we directly map the radial field SV at the coremantle boundary, using the technique of Subtractive Optimally Localized Averages (SOLA) that was adapted to geomagnetism by Hammer and Finlay (2019). The SOLA technique can be used to compute estimates of the radial field SV directly at the CMB, based on local spatial averages described by averaging kernels, of the SV field centred on target locations of interest, and time-averaged over chosen time windows. By collecting many individual SOLA estimates on a grid at the CMB, the SV field can be mapped on regional or global scales. An important foundation of the SOLA technique is that for noise-free data and a linear forward problem, any linear combination of the data provides a specific average of the true model. With noisy data, a variance is ascribed to this spatial average value such that a trade-off between resolution and variance arises (Oldenburg, 1984;Parker , 1994). The SOLA technique readily provides information on the resolution offered by the magnetic field observations, in the form of averaging kernels, and estimates of the variance of the locally averaged field. We compare SOLA-based maps of SV and SA estimates derived from CryoSat-2 and Swarm data, in order to asses the possibilities of the CryoSat-2 data for mapping these fields at the CMB. Demonstrating the usefulness of the CryoSat-2 data, we then take advantage of this data to map the time evolution of the SA field along the geographic equator at the CMB from 2001 to 2019.
We wish to emphasize that both the GVO and the SOLA methods can result in patterns of secular variation different from those seen in the CHAOS field model, despite the fact they use similar data selection schemes and the same magnetospheric field model. In the GVO and SOLA methods data close to a location of interest are used to effectively determine localized field or SV estimates. In contrast, estimation of the coefficients of truncated spherical harmonic expansions, in models such as CHAOS, is an inherently global approach that aims to find the best possible global model. In the GVO method only data from within a 700km radius of a target location is used to determine the local potential. In the SOLA method measurements far from the site of interest have essentially no influence on the estimated CMB SV because the data kernels, based on Green's functions for Laplace's equation under Neumann boundary conditions, have decreasing sensitivity far away from the target location.
Moreover, the CHAOS model involves temporal regularization whereby one minimizes global norms based on time derivatives of the CMB radial field, integrated over the entire timespan of the model. In the GVO method there is no temporal smoothing beyond the choice of 4 month data windows and use of annual differences to produce SV series. In the SOLA method we use one or two year windows to estimate the SV and then annual differences to estimate SA. The SOLA method involves a trade-off parameter specifying the balance between the spatial resolution and the variance of each local estimate, each estimate being time-averaged over a one or two year time window; there is no global regularization over longer time spans. Both methods therefore constitute a localized compression of information contained within the satellite data on the potential field near the location of interest within the specified time window. The GVOs and SOLA can thus give a different picture of the core field evolution compared to the spherical-harmonic based CHAOS model with its global temporal regularization, especially regarding rapid changes at short wavelengths scales when the data quality is high (see Section 4.2).
Section 2 describes the satellite measurements used in this study, including the CryoSat-2 platform magnetometer data, and how these have been selected and processed. Section 3 presents the GVO technique and results concerning global GVO time series. Section 4 presents the SOLA technique and results of applying this to estimate SV and SA at the CMB. Conclusions and perspectives are given in Section 5.

Data Selection
We use satellite vector magnetic field measurements from the Ørsted satellite between July 2000 and December 2005, from CHAMP taking L3 magnetic data between July 2000 and September 2010, and from the Swarm trio taking Level 1b MAG-L data, version 0505/0506, between January 2014 and April 2020. Most importantly for this study we make use of platform magnetometer data from the CryoSat-2 mission, taking calibrated vector measurements with a sampling rate of 4s from the FGM1 magnetometer dataset, version 3, between August 2010 to December 2018. This dataset has had extensive corrections applied for disturbances fields, and was calibrated using a reference field model -for full details see Olsen et al. (2020).
From the Ørsted, CHAMP and Swarm measurements we produced two data sets: dataset #1 used in the GVO application taking a 15 s subsampling of the vector field measurements from Ørsted, CHAMP and Swarm, while taking every 4th measurement from the CryoSat-2 dataset (i.e. a 16 s sub-sampling); dataset #2 used in the SOLA application taking a 5 s subsampling of the vector field measurements from Ørsted, CHAMP and Swarm, while taking every element of the CryoSat-2 dataset with its 4 s sampling rate. Field measurements having gross data outliers for which the vector field components deviated more than 500 nT from the CHAOS-7.2 internal field model predictions  were rejected. For both data sets we apply a dark geomagnetically quiet time selection criteria scheme. See table 1 below for full details of the selection requirements for the datasets used in the GVO and SOLA applications. In both cases we required the sun to be at least 10 • below the horizon, adding restrictions on the geomagnetic activity index K p and the change in the ring current index (see Olsen et al., 2014), as well as constraints on the merging electric field at magnetopause, and on the magnitudes of the B Y and B Z components of the interplanetary magnetic field (e.g. Finlay et al., 2020;Ritter et al., 2004). We used minute values of the IMF components and solar wind speed from the OMNI database, http://omniweb.gsfc.nasa.gov, computing two hourly means prior to the time of the considered datum (Finlay et al., 2016).
Since our focus is the core field we have applied corrections to the data for the lithospheric and external fields. For the lithospheric field model, we used the LCS-1 model (Olsen et al., 2017); the precise choice of lithospheric field is not crucial for studies of the SV and SA. For the solar-quiet ionospheric field and associated induced fields we used the CIY4 model (Sabaka et al., 2018). For the magnetospheric field and related induced fields we used the CHAOS-7 model . These models were chosen as they are well established and compatible with the data selection criteria described above.
< 500 nT from CHAOS < 500 nT from CHAOS LCS-1 Lithospheric field for n ∈ [14, 120] subtracted CHAOS-7.2 Magnetospheric (plus induced) fields subtracted CIY4 Ionospheric (plus induced) fields subtracted Table 1: Selection criteria, and model corrections applied for the GVO and SOLA datasets used in this study.
As noted in Table 1, for the GVO application we use the sums and differences of the magnetic field measurements. It has been shown by Olsen et al. (2015) and Sabaka et al. (2018), that taking differences of the satellite measurements alongtrack and east-west (between Swarm satellites A and C) helps the recovery of the small scale core field, as this reduces the impact of correlated errors caused by unmodelled large-scale external fields. Here we follow such an approach by taking differences of the measurements, but we also include along-track and east-west sums of the measurements in order to ensure sufficient constraint on the larger wavelengths of the field (Hammer , 2018;Sabaka et al., 2013). We denote the magnetic vector measurements by B k (r), where k is any of the three given vector component of the field, ∆d k and Σd k denote measurement differences and sums of this particular component, respectively. Here the along-track (AT) and East-West (EW) data differences are denoted by ∆d k = (∆d AT k , ∆d EW k ), and the data sums by Σd k = (Σd AT k , Σd EW k ). The along-track data differences are calculated using the 15 s differences ∆d AT k = [B k (r, t)−B k (r+δr, t+15s)]. The along-track sums were calculated as Σd AT k = [B k (r, t) + B k (r + δr, t + 15s)]/2. For Swarm, East-West differences were calculated as ∆d EW ] having an East-West orbit separation between the Swarm Alpha (SWA) and Charlie (SWC) satellites of ≈ 1.4 • corresponding to 155 km at the equator (Olsen et al., 2015). The East-West sums were calculated as Σd EW k = [B SWA k (r 1 , t 1 ) + B SWC k (r 2 , t 2 )]/2. For a particular orbit of Swarm Alpha the corresponding Swarm Charlie data were selected to be those closest in colatitude with the condition that |∆t| = |t 1 − t 2 | < 50s (Olsen et al., 2015).

Four monthly Core Field GVOs and Secular Variation Estimates
In the first application, we compute Geomagnetic Virtual Observatory time series derived from dataset #1. Of particular interest is the quality of the GVO series obtained from CryoSat-2 data compared with similar series obtained from CHAMP and Swarm data. The time series consist of estimates of the geocentric spherical polar vector components of the magnetic field at specified target points, referred to as GVOs (Mandea and Olsen, 2006;Olsen and Mandea, 2007). Here we use the same algorithm described in detail by Hammer et al. (2020) (see also http://www.spacecenter.dk/files/magnetic-models/ GVO/GVO_Product_Algorithm.pdf) to produce the Swarm GVO product, and derive global grids of 300 uniformly distributed GVO time series each having 4 month cadence. The GVOs are located in an approximately equal area grid based on the sphere computed using the algorithm of Leopardi (2006).
For each GVO in the grid we take data from within a cylinder of horizontal radius r cyl = 700 km. The GVO's have the spherical polar coordinates r GV O = (r, θ, φ), and are placed at fixed altitudes r = r a + h GV O where h GV O is the height above the Earth's mean spherical radius r a = 6371.2 km. For the CHAMP, CryoSat-2 and Swarm missions the GVO altitudes were chosen as h GV O = 370 km, h GV O = 727 km and h GV O = 490 km, respectively, such that the GVO's are located at approximately the mean orbital altitude for each mission during the time interval considered.
The input measurements of dataset #1 are provided in an Earth-Centred-Earth-Fixed (ECEF) coordinate frame by the spherical polar components B obs = (B r , B θ , B φ ). From the extracted vector field measurements surrounding each GVO, within a radius of 700 km and within a 4 month time window, a residual magnetic field, δB, is first computed by subtracting off estimates of the main field and noncore fields where the field estimates removed are: a) B M F the internal field for SH degrees n ∈ [1, 13] as given by IGRF-13 (Alken et al., 2020b), b) B lith the static internal field for SH degrees n ∈ [14, 185] as given by the LCS-1 model, (Olsen et al., 2017), c) B mag the magnetospheric and associated induced field as given by the CHAOS-7.2, model, , d) B iono the ionospheric and associated induced field as given by the CIY4 model, (Sabaka et al., 2018). Note that estimates of the main field from IGRF-13 (with linear time dependence over 5 year intervals) are subtracted here. At a later stage, main field estimates, again from IGRF-13, are added back for the specified GVO times and positions. Removal of a main field at this stage allows for a more effective pre-whitening of the data, such that Huber weights, used in the robust GVO estimation scheme, can be well determined. We emphasize that this approach still allows us to capture departures from the subtracted main field when required by the satellite data.
Of the remaining residual field, we are interested in the core field part of that signal. Although we have removed estimates of the external fields and their associated Earth-induced counterpart, contributions from non-core sources remain in the residual field. A particular concern is contamination from fields caused by rapidly varying ionospheric currents, for example polar electrojet currents at high latitudes in the E-layer, and possibly also signatures of F-layer currents at mid and low latitudes. Further work is needed on these aspects. We attempt to mitigate leakage of field-aligned currents by removing toroidal field estimates (Sabaka et al., 2010) obtained by performing an epoch-by-epoch spherical harmonic analysis performed on the global network of GVOs , see below for further details.
Next, the residual magnetic field data and their positions are transformed from the spherical system to a right-handed local topocentric Cartesian system (x, y, z) with origin at the GVO target location. At the GVO location and only at this location, x points towards geographic south, y points towards east and z points radially upwards . Assuming that the satellite measurements are made in a source free region, the residual magnetic field, δB, is a Laplacian potential field. In the local Cartesian coordinate system the magnetic scalar potential, V , can be expanded as a sum of polynomials of the form C abc x a y b z c (Backus et al., 1996). Here we use this expansion out to cubic terms V (x, y, z) = C 100 x + C 010 y + C 001 z + C 200 x 2 + C 020 y 2 + C 002 z 2 + C 110 xy + C 101 xz + C 011 yz + C 300 x 3 + C 030 y 3 + C 003 z 3 + C 210 x 2 y + C 201 x 2 z + C 120 y 2 x + C 021 y 2 z + C 102 z 2 x + C 012 z 2 y The forward problem linking the vector of GVO model coefficients, m = [C 100 , C 010 , ..., C 111 ] T , with the data vector d vec containing the residual field components, δB, can be written where G vec is a design matrix derived from the spatial derivatives of eq. (2). As noted above, instead of using residual vector field components to compute the potential, we use sums and differences of the residual field vector components such that the data vector is d = {∆d vec x , ∆d vec y , ∆d vec z , Σd vec x , Σd vec y , Σd vec z }, where ∆ and Σ denotes the differences and sums of the computed residual field as described in Section 2. The relevant design matrix is then constructed as To determine the GVO model coefficients, we use a robust iteratively reweighted least-squares inversion scheme, based on a diagonal weight matrix consisting of robust (Huber) weights for each entry in the data vector (e.g. Constable, 1988). We also include a down-weighting factor of 1/2 for the Swarm satellites Alpha and Charlie accounting for the fact that these two satellites fly side-by-side and therefore do not provide completely independent measurements. A minimum number of 30 data points were required for the computing the inversion. Using the resulting coefficients of the potential for a given GVO target location and time, derived from the associated sums and differences satellite data, a prediction for the mean residual field at the GVO target point and epoch can be computed as δB GV O (x, y, z) = −∇V (0, 0, 0) = −(C 100 , C 010 , C 001 ).
Moving back to the vector components in spherical polar coordinates, We then add back the IGRF-13 main field predictions for the given target point and epoch, The above procedure is then repeated for each epoch and each component to obtain time series of GVO estimates of the vector magnetic field at the GVO target locations.
The GVO method assumes that the residual field in eq.(1), is a potential field. However, because the satellite measurements are made in the ionospheric F-region, in-situ currents can cause non-potential fields to leak into the estimated potential (Olsen and Mandea, 2007). Therefore, in a final post-processing step, we carry out an epoch-by epoch spherical harmonic analysis of our 4-monthly GVOs, estimating external and toroidal field contributions to SH degree 13 or to degree 6 if fewer than 300 GVOs are available. For epochs having an insufficient number of GVOs available to ensure a stable solution, the external and toroidal coefficients were computed by a linear interpolation between nearby epochs. The external and toroidal field estimates, reaching a level of ±15nT at high latitudes, are then removed to obtain an estimated Core Field GVO time series . Secular variation at a given GVO location is computed using annual differences between values at time t + 6 months and at time t − 6 months. We have chosen to take annual differences, as this helps to avoid annual non-core signals that may persists in the GVO series despite our best efforts in reducing such contamination. Figure 1 presents the number of 4-monthly Core Field GVO estimates during the past 20 years. The maximum possible number of GVOs per epoch is 300. A strong dip in the number of GVOs is seen during 2002-2004 due to increased solar activity that meant there were fewer data meeting our selection criteria. As noted above, if fewer than 30 measurements are available within a GVO target cylinder during a given 4-month window, we were unable to reliably determine GVO estimates. The remaining epochs from 2004-2020 are well covered with only few epochs having less than 250 GVO's available. We find that using CryoSat-2 data, we are able to provide between 200 and 300 GVO estimates at all times during the gap between the end of the CHAMP mission and the start of the Swarm mission.  Table 2 presents the root-mean-square (rms) and mean of the residuals between the contributing satellite data (sums and differences) and GVO model predictions, summing over all GVO's for a given vector component and a given region (polar or nonpolar). Here we defined polar to be polweard of ±54 • geographic latitude. The polar rms values for both sums and differences are higher than the non-polar, and the CHAMP values are slightly higher than the Swarm values. The CryoSat-2 values are seen to be higher for all components but not unreasonable, given they are derived from platform magnetometer data. The non-polar rms values for all components are below 2 nT during both CHAMP and Swarm times. The CryoSat-2 GVO's rms values are as expected larger, and especially for the along-track differences, indicating that along-track correlated noise is less dominant, due to the presence of other noise sources in platform magnetometer data.   Table 2: GVO model rms misfit statistics between GVO estimates and the contributing data for the global grid of 300 GVO's during CHAMP, CryoSat-2 and Swarm. Here and ∆ represent data sums and data differences, respectively, split into the North-South (NS) and East-West (EW) components.

GVO Results: Global time series of Secular Variation from 2002 to 2020
The CryoSat-2 magnetometer data have been cleaned from known platform signals and calibrated as described by Olsen et al. (2020). This calibration relies on computing residuals with respect to a reference magnetic field which was taken from the CHAOS-6-x9 field model (Finlay et al., 2016). Here, we carried out an experiment to verify that the GVO secular variation signals obtain from CryoSat-2 data are independent of the main field model used for data calibration. To do this we computed CryoSat-2 GVO estimates using two different datasets, the first being the official dataset calibrated using the CHAOS-6-x9 field model and the second a test version of the CryoSat-2 data calibrated instead using IGRF-13 (Alken et al., 2020b). Having estimated global grids of GVO series for each dataset, to each series we fit cubic smoothing splines, with a knot spacing at every 4 month, and with the smoothing parameter determined using a GCV (generalized cross-validation) approach (Green and Silverman, 1993). Figure 2 presents SV series for the three field components at three example GVOs, with colatitude/longitude (36 • , 77 • ) (top), (72 • , 110 • ) (Center) and (84 • , −137 • ) (bottom). The SV time series derived from the CHAOS-6x9 calibrated dataset are shown with red dots (spline fit in the red line) and similar series from the IGRF-13 calibrated dataset are shown with blue dots (spline fit in the blue line). SV model predictions from the CHAOS-6-x9 model up to SH degree 16 (in green), and IGRF-13 (in black), are also shown for reference. There are clearly differences in the GVO SV estimates derived using the CHAOS-6x9 and IGRF-13 calibrated data. Since the CryoSat-2 calibration relies on computing residuals with respect to a reference field model (here either CHAOS or IGRF), we do expect such differences, especially since the IGRF model assumes a crude piecewise constant SV field. Interestingly, the IGRF calibrated SV are an rms of 1.2nT/yr for the horizontal components, and 3.9nT/yr for the radial component, closer to CHAOS-6x9 than IGRF-13, indicating that the Cryosat-2 data does possess an SV signal regardless of the model used to calibrate them.
We find that the CHAOS-6x9 and IGRF-13 calibrated CryoSat-2 GVO SV series show similar sub-decadal changes, neither of which match those seen in CHAOS-6x9. The similarity of the IGRF and CHAOS calibrated SV series give us confidence that the sub-decadal secular variation seen in CryoSat-2 data is independent of the field model used in its calibration. In addition, we note a strong acceleration in the radial SV component from 2014 to 2018 at the GVO (top plot) located in the Northwest Siberia, a change which is only partially observed in the IGRF-13 model.
In Table 3 we report the mean rms differences between the GVO SV estimates and GCV spline fits, separated by component (B r , B θ , B φ and also the intensity F ) and by polar and non-polar regions. These numbers give an indication of the scatter in the SV datasets and allow the quality of the GVO SV series obtained from CHAMP, Swarm and CryoSat-2 (both CHAOS-6x9 and IGRF-13 calibrated versions) to be compared. Similar results are seen for the CHAOS-6-x9 and IGRF-13 model calibrated GVO SV series, with rms differences between GCV fits and intensity SV data, taking the over all GVOs, of 3.5 nT/yr and 3.4 nT/yr, respectively. As expected, similar numbers for the CHAMP and especially the Swarm derived GVO's, are much smaller being 1.8 nT/yr and 0.9 nT/yr, respectively. This indicates the lower scatter in the CHAMP and Swarm GVO series.  Table 3: Mean of the rms differences (in nT/yr) between GVO SV series and GCV cubic spline fits. Results are shown for both the CHAOS-6-x9 and IGRF-13 calibrated CryoSat-2 datasets, and for GVO SV series derived using Swarm and CHAMP data.  Figure 2: Time series of annual differences of the spherical field components of CryoSat-2 GVO time series derived using CHAOS-6-x9 (red) and IGRF-13 (blue) calibration along with GCV spline fits at the three example grid locations : 26 (top), 100 (center) and 140 (bottom). Also shown are IGRF-13 (black) and CHAOS6-x9 predictions (green). differences of less than 1.0nT/yr in all three components of the SV series derived from the CHAOS-6-x9 and IGRF-13 calibrations. With the availability of CryoSat-2 magnetic data, it is now possible to use GVOs to study secular variation globally over the past twenty years. In order to illustrate the quality of information they can provide, in Figure 4 we present comparisons of GVO SV series, from CHAMP, CryoSat-2 and Swarm, with Revised Monthly Means (rmm) (Olsen et al., 2014) from high quality ground observatory SV time series from Kourou in South America (top plots), from Novosibirsk observatory in Siberia (middle plots) and from Honolulu ground observatory in the central Pacific (bottom plots). Each plot shows the spherical polar components of the annual differences of revised monthly means (black dots) computed from the ground observatory data (Olsen et al., 2014), and GVO time series derived from CHAMP (purple dots), CHAOS-6x9 calibrated CryoSat-2 (blue dots) and Swarm (red dots) data relocated to ground level. Here the GVO time series for each satellite mission have been mapped to a common altitude of 700km by subtracting the field difference as given by the CHAOS-7.2 model for SH degrees 1-20, between 700km altitude and GVO altitudes which are close to the mean orbital altitude for each mission. At all three locations the GVOs derived from CryoSat-2 data have more scatter -this is particularly noticeable in the θ-and φ-components for the examples shown, whereas the scatter in the radial component is closer to the level seen in the CHAMP and Swarm derived GVOs. Table 3 indicates that at non-polar latitudes the scatter in the r-and θ-components is generally similar, 3.4 nT/yr compared to 3.8 nT/yr. Notice that the scatter in the ground observatory rmm's is also enhanced in the horizontal components. The good agreement between the independently estimated CryoSat-2 and Swarm derived GVO's at overlapping epochs from 2014 to 2018, is particularly evident in the r-and φ-components. SV variations are coherent in both phase and amplitude between the CryoSat-2 GVO time series and the ground observatory records, thus confirming that the CryoSat-2 GVO's are able to track the same field changes as observed by ground observatory records on timescales of 1 year and longer.
Next, in Figure 5 we presents a global map of the SV of the radial field component over the past 18 years. Here, to ease visualization, GVO time series from CHAMP (covering 2002-2010), CHAOS-6x9 calibrated CryoSat-2 (covering 2010-2014) and Swarm (covering 2014-2020) have been mapped to a common altitude of 700 km, again using the CHAOS-7.2 field model, and combined into one composite time series. This allows for the investigation of global patterns of sub-decadal SV. We find that regions at low latitudes display strong sub-decadal variations, however, not simultaneously at all longitudes. For instance, we observe a change of slope in the radial SV field occurring over the south Atlantic region around 2007 and again in 2014, over Indonesia around 2014, and in the Pacific region centred in 2017. Some of these variations are characterised by distinct "Λ" and "V "-shaped behaviour occurring over time spans of 5-10 years and locally confined to specific regions, but otherwise reminiscence of often discussed geomagnetic jerks (Mandea et al., 2010). The availability of CryoSat-2 magnetic field data clearly plays a key role in permitting a continuous coverage from satellite-based time series without a gap between 2010 and 2013, thus allowing the study of global patterns in the time-varying core field secular variation over the past twenty years.  (Olsen et al., 2014), and 4-monthly Core Field GVOs derived from CHAMP (purple dots), CryoSat-2 (blue dots) and Swarm (red dots) data. Note the GVO estimates have been mapped from satellite altitude to ground in order to aid the comparison.

SOLA Method for local estimation of CMB radial field SV
We now move on to investigate the behaviour of the core field not at satellite altitude but down at the core-mantle boundary (CMB), on the edge of the region where it originates. To do this we use the Green's functions of the Neumann boundary value problem that links the magnetic field at satellite altitude to the radial field at the CMB. Following Hammer and Finlay (2019), a localized estimate, B r , of the radial magnetic field at a target location and time, (r 0 , t 0 ), at the CMB can be computed as a localized spatial average around the target location time-averaged over a specified interval. Because the CMB radial magnetic field is linearly related to the spherical polar components of the vector field at satellite altitude, we can write B r as a weighted linear combination of the satellite magnetic measurements with weights q n Gilbert, 1968, 1970;Hammer and Finlay, 2019) where d n are satellite magnetic measurements (n = 1, ..., N ) within a specified time window and q n are weighting coefficients to be determined. Here the data d n , at positions r n and times t n , are taken from dataset #2 and for simplicity we consider using only observations of the radial component of the field. Corrections for the lithospheric field for SH degrees n ∈ [14, 185] as given by the LCS-1 model (Olsen et al., 2017), for the magnetospheric and associated induced fields as given by the CHAOS-7.2 model , and for the ionospheric and associated induced fields as given by the CIY4 model, (Sabaka et al., 2018) are removed from the observations in a pre-processing step. The radial magnetic field measurements, d n (r n , t n ), are then related to the radial magnetic field B r (r , t n ), integrated over the CMB, by (Gubbins and Roberts, 1983) d n (r n , t n ) = S G r (r n , r )B r (r , t n )dS (6) where the surface element is dS = sinθ dθ dφ . The data kernel G r (r n , r ) is the radial derivative with respect to r, of the Green's functions for the exterior Neumann boundary value problem (e.g. Barton, 1989;Gubbins and Roberts, 1983) where h n = r /r n and r is the CMB radius, f n = R n /r n where R n = r 2 n + r 2 − 2r n r ζ n and ζ n = cos γ n = cos θ n cos θ +sin θ n sin θ cos(φ n − φ ), where γ n is the angular distance between a measurement at position (θ n , φ n ) and a position on the CMB (θ , φ ). The data kernel describes how a particular measurement samples the CMB radial field; radial magnetic measurements sample the CMB radial field most strongly directly below the measurement position. Regarding the time-dependence, we use a first order Taylor expansion around a reference time t 0 , such that The time difference, ∆t n = t n − t 0 , is computed with respect to the target time, t 0 . Inserting eq.(8) into eq.(5) we obtain where K(r 0 , t 0 , r ) andK(r 0 , t 0 , r ) are spatial averaging kernels for the CMB field and secular variation respectively, constructed from the weighting coefficients and the data kernels By varying the weight coefficients, q n , the shape of the averaging kernels change.
Notice that time differences ∆t n , between the measurement times and the target time, are effectively additional temporal weights applied to the kernel K in order to obtainK.
In order to obtain estimates of the secular variation of the radial field on the CMB, at the target location and time Ḃ r (r 0 , t 0 ), we minimize the following objective function (12) where λ is a trade-off parameter (units of [nT −1 ]), q is vector of the weighting coefficients, E is the data error covariance matrix which we define below andṪ is an SV target kernel that we choose to be a Fisher distribution on the sphere (Fisher , 1953) where γ 0 is the angular distance on the CMB between the target position (θ 0 , φ 0 ) and another position (θ , φ ). On the basis of tests carried out by Hammer and Finlay (2019), we set κ = 600 corresponding to a target kernel width of 15 • ; this is narrower than can be achieved forK with the available data, but it avoids excessive ringing associated with taking a Dirac delta function as the target kernel. When computing SOLA estimates for a given time window, we select a subset (n = 1, ..., N ) of the measurements. Using this data subset the data error covariance matrix E is defined as follows. Using all available measurements (m = 1, ..., M ), for each satellite mission within 2 degree bins of quasi-dipole (QD) latitude (Richmond , 1995), we first derived robust data error variances as a function of QD latitude where m are residuals with respect to predictions of the CHAOS-7.2 internal field model for SH degrees n ∈ [1, 13], µ are robust mean residuals within the considered bin and w m are Huber weights (e.g., Constable, 1988) for the data within each bin. Here we use QD coordinates, as this is appropriate for characterizing processes related to unmodelled ionospheric currents which we consider to be a likely source of contamination, especially at high latitudes. Figure 6 presents the resulting QDlatitude-dependent error estimates σ(θ QD ) for the radial field component used in this study, comparing the values for the Ørsted, CHAMP, CryoSat-2 and Swarm datasets. When computing SOLA estimates for a specified time window of e.g. 2yrs, we select a data subset of dataset #2. Using this data subset of N measurements, a data error covariance matrix E computed. Diagonal elements of this data error covariance matrix E are finally defined as E n,n = σ 2 (θ QD )/w n where w n are robust (Huber) weights determined a-priori for each datum (n = 1, ..., N ), based on their residual to CHAOS-7.2, in order to account for the expected long-tailed error distribution. Off-diagonal elements of E are set to zero.
In addition to minimizing eq.(12) we simultaneously impose the following constraint where the first term is in practise very small when estimating the SV, since it is minimized in the objective function. This constraint ensures that a valid averaging kernel is obtained.
Discretization of integrals over the CMB was carried out using Lebedev quadrature (Lebedev and Laikov , 1999) and the system of equations was solved for the coefficients, q n , using a Lagrange multiplier method, see Hammer and Finlay (2019) for further details.
Once SOLA estimates of the CMB radial field SV at the chosen target location and epoch are obtained, by minimizing eq. (12) subject eq. (16), we are able to easily appraise them based on (i) their averaging kernel width, which we define as the angular distance between the points at which the averaging kernel first reaches zero amplitude moving away from its maximum value, and (ii) the variance of the SOLA estimate which we computed aŝ By changing the parameter λ of eq.(12), a range of solutions can be computed which describes a trade-off between having an averaging kernel width as small as possible and the variance of the estimate being as small as possible (Parker , 1977). Below we discuss the effect of changing λ on our results. and Swarm data We begin by first comparing SOLA estimates for the CMB radial field SV obtained using separate data subsets from the Swarm and CryoSat-2 missions respectively. Firstly, Swarm and CryoSat-2 datasets are extracted from the main dataset #2 described in Section 2, and each covering the same two year time window from 2015.0 to 2017.0 Next, in order to obtain data subsets with suitable spatial and temporal coverage, we considered bins surrounding each point in an approximately equal-distance grid at satellite altitude of ≈ 2.5 • spacing, based on the partitioning algorithm of Leopardi (2006), and randomly sampled one datapoint from each bin, resetting the bins every two months. Data subsets spanning the full two year window from 2015.0 to 2017.0 were produced by accumulating these two monthly globally distributed subsets. The resulting Swarm and Cryosat-2 data subsets spanning 2015.0 to 2017.0 consisted of 62469 and 54685 radial field observations respectively.
In Figure 7 we compare maps collecting SOLA CMB radial field SV estimates centered on epoch 2016.0, derived using the CryoSat-2 and Swarm data subsets spanning 2015.0 to 2017.0. To ensure that we obtained SOLA estimates of comparable resolution, we first computed SOLA SV estimates using the Swarm data subset and taking λ = 3 × 10 −3 nT −1 . This resulted in well-behaved averaging kernels with widths ≈ 42 • . Next, we used these averaging kernels as the target kernels in order to derive similar estimates using the Swarm and CryoSat-2 data subsets, thus effectively seeking Swarm and CryoSat-2 SV estimates with the same spatial resolution as the target estimates. The maps in the top row of Figure 7 show the resulting global collections of SOLA SV estimates, obtained on a 1 • grid of target locations at the CMB, based on the CryoSat-2 (left plot) and S warm (right plot) data subsets, respectively.  and averaging kernel widths (bottom plots). Estimates are derived using CryoSat-2 (left plots) and Swarm (right plots) data, respectively.
Maps collecting the formal standard deviations for each SOLA estimate (derived using eq.17) and their averaging kernel widths are presented in middle and bottom rows of Figure 7. The standard deviations of the Swarm-based SOLA SV estimates are fairly homogeneous with values of 0.3−0.4µT/yr; those for CryoSat-2 SOLA are somewhat larger, being in the range ≈ 0.6 − 1.8µT/yr. We note that the CryoSat-2 error estimates are slightly lower at higher latitudes. The same is true for the Swarm map, but here the variations in errors estimates are in that case much smaller. Kernel widths in both cases are also fairly homogeneous except at auroral latitudes where distinct behaviour of the kernels are found related to the data error estimates having increased amplitude, as seen in Figure 6.
As seen from the kernel widths in the bottom row of Figure 7, very similar resolution has been obtained (by construction) in the Swarm and Cryosat-2 SV maps. The same field features are clearly identified in both maps. For instance, we notice high latitude SV patches in the northern hemisphere which have been associated with a high latitude jet of core flow (Livermore et al., 2017), and there is increased amplitude of SV over the hemisphere centered on the Atlantic in comparison with the Pacific hemisphere. This first test gives us confidence that the CryoSat-2 measurements can be used to reliably map SV features also at the CMB on a timescale of 2 years, and with a resolution down to 42 • degrees.
Next, we go further and investigate the second time derivative or secular acceleration (SA) of the radial field at the CMB, which is of great interest for investigating the dynamics of the core (e.g. Chi-Durán et al., 2020;Finlay et al., 2016). Again we first compare maps based on CryoSat-2 and Swarm data. To obtain SA estimates we initially use the accumulated change between SV estimates two years apart. In particular, in Figure 8 we show the SA in 2017 based on the difference between SOLA SV estimates in 2016.0 and 2018.0. Note, that this is not an instantaneous secular acceleration but a centered difference in SV estimates 2 years apart, each based on 2 years of data. To study the SA, we computed SOLA SV estimates from Swarm data taking λ = 1 × 10 −2 nT −1 , and used the resulting associated averaging kernels as the target kernels for the SOLA estimates from both Swarm and CryoSat-2 data. Figure 8 presents global grids of SOLA CMB radial field SA estimates centered on 2017.0, again with a 1 • spacing, derived from CryoSat-2 (left plot) and Swarm (right plot) data subsets. Here the map is centered on the Pacific region where there has been interesting SA activity during the past 6 years (Finlay et al., 2016). As for the SV maps, we find error estimates, computed assuming the contributing SV estimates have independent errors, are fairly homogeneous, with values ranging between 0.11− 0.27µT/yr 2 for the estimates derived using CryoSat-2 data and 0.06 − 0.08µT/yr 2 for the estimates derived using Swarm data. In both cases kernel widths are close to ≈ 42 • , except in the auroral region. Comparing the CyroSat-2 and Swarm-based SA maps, similar features can be observed. In particular, this is the case for the features seen under Asia and Indonesia. A distinctive feature reproduced in both maps is the sequence of intense patches of SA at low latitudes in a localized region below central America, having amplitudes of approximately 1.9 ± 0.3µT/yr 2 and 1.7 ± 0.1µT/yr 2 for the CryoSat-2 and Swarm maps, respectively. The location and amplitude of these features are similar in the two maps, confirming that CryoSat-2 data can be used to track such SA structures. Because the SOLA estimates are local averages, distant high latitude measurements, where ionospheric electrical currents may be prominent even during dark quiet times, will have little influence on such SV and SA estimates at low latitudes. Finally, we note a strong patch under the Bering Sea of amplitude approximately 1.4 ± 0.3µT/yr 2 and 1.0 ± 0.1µT/yr 2 for the CryoSat-2 and Swarm maps, respectively.  These initial investigations of the CMB radial field SV and SA using the SOLA technique indicate, as seen earlier in the GVO time series, that low latitude regions experience significant sub-decadal core field variations. We are therefore motivated to study the field time-dependence in the equatorial region in more detail. We do this in Figure 9 by constructing time-longitude (TL) plots of SOLA CMB radial field SA estimates along the geographic equator from 2002 to 2019, centered on the Pacifc. We again compute our SA estimates based on differences of SOLA SV estimates two years apart, each derived from 2 year data windows, and sliding the windows in 2 month steps. Here we use radial field data from the Ørsted, CHAMP, CryoSat-2 and Swarm satellites.
The top left plot in Figure 9 presents estimates based on averaging kernels obtained from Swarm data using λ = 3 × 10 −2 nT −1 . These associated error estimates between 0.02 − 0.08µT/yr 2 and kernel widths ≈ 50 • . Their resolution is lower that shown in Figure 8 and corresponds to approximately SH degree 8. For comparison the SA predicted by the CHAOS-7.2 model for SH degrees n ∈ [1, 8] are shown on the top right plot. Note here that the CHAOS-7.2 model makes use of uncalibrated vector magnetic data from CryoSat-2 between 2010 and 2014, and co-estimates magnetometer calibration parameters . The SOLA and CHAOS TL plots in the top row of Figure 9 show largely the same SA features, illustrating the convergence of the two techniques at long wavelengths of the SA and when constructing SOLA SA estimates from SV differences between consecutive 2 year time windows. For instance, the evolution of the SA features observed in Figure 8 under the central Americas, can be identified ranging from longitudes 240 • to 320 • centered on 2017. In addition, we find strong SA patches in the CryoSat-2 data around 2013 at longitudes 70 • to 160 • and 280 • to 320 • . Notice, that there seems to be a sign changing sequence occurring at longitudes 240 • to 320 • going from 2005 to 2019.
Next, we increase the spatial resolution by instead deriving SV estimates using λ = 1 × 10 −2 nT −1 , which leads to slightly larger error estimates in the range 0.1−0.5µT/yr 2 and kernel widths ≈ 42 • , i.e. similar to Figure 8. This is shown in the bottom left plot while the CHAOS-7.2 model predictions for SH degrees n ∈ [1, 10], matching approximately the kernel width, are shown on the bottom right plot. Although the SOLA TL-plot looks somewhat noisier than the CHAOS plot, similar coherent evolving structures having higher amplitudes can clearly be identified. The noiser appearance in the interval 2010-2014 likely indicates the limitations of the CryoSat-2 data, but they clearly provide useful information during this period.  Left column: SOLA radial field SA estimates derived from differences between consecutive SOLA SV estimates based on 2yr data windows moving in 2 month steps and using λ = 3 × 10 −2 nT −1 (top left), and λ = 1 × 10 −2 nT −1 (bottom left). Right column: SA predictions from the CHAOS-7.2 model for SH degrees up to 8 (top right) and 10 (bottom right); note that temporal regularization causes significant time-averaging of the high degree SA in the CHAOS model.
With data from the Swarm mission, it is possible to go further and also increase the temporal resolution of the SA by taking 1 year differences of SOLA SV estimates derived from 1 year data windows sliding in 1 month steps. The result of applying this procedure to obtain SA estimates on the geographic equator between 2015.0 and 2019.5 is shown in the left plot of Figure 10. These SOLA estimates have associated errors of 0.3 − 0.6µT/yr 2 and kernel widths ≈ 42 • . The right plot shows similar CHAOS-7.2 model predictions for SH degrees n ∈ [1, 10]. Both TL-plots shows the similar large scale features, for instance, the features under central America from longitudes 240 • to 320 • , which are elongated compared with Figure 9 due to the change in scale of the y-axis (time). However the SOLA results show significantly more time dependence, revealing features that were smoothed out by the temporal regularization of CHAOS-7.2. Changes of sign in the SA within about 1 year can be observed.  Particularly interesting is the appearance in the Pacific region around 2017, at longitudes 150 • to 220 • , of side-by-side positive and negative intense SA features, that have subsequently drifted westwards. This SA change coincides with the peak in the radial SV field observed in the Pacific region during Swarm time seen in the GVO map Figure 5. We note the presence of features in Figure 10 that appear to drift rapidly both eastwards and westwards, for example from 160 • East in 2015 to 220 • East in 2017. Such rapidly drifting behavior of low latitude SA patches is difficult to explain in terms of simple core flow advection processes. They may instead be a signature of wave propagation close to the core surface. A range of possible candidates for fast waves in the core have recently been described, some requiring only a strong magnetic field and rotation (Aubert and Finlay, 2019;Gerick et al., 2020) while others rely on the presence of a possible stratified layer at the top of the core (Buffett and Matsui , 2019). Though tempting, it may be dangerous to interpret such features that are at the limit of the present spatial resolution and temporal resolution (Gillet, 2019). It will be important to assess whether such features remain coherent in the future, as the resolution of the SA increases.

Conclusions
In this article we have studied global patterns and sub-decadal changes in geomagnetic secular variation during the past 20 years. We have shown there is now continuous coverage of magnetic field measurements from low Earth orbiting satellite missions during this period, provided one takes advantage of calibrated platform magnetometer data from the CryoSat-2 satellite. Using vector magnetic field measurements from Ørsted, CHAMP, Cryosat-2 and Swarm we have constructed Geomagnetic Virtual Observatory (GVO) time series that track sub-decadal changes in the core field secular variation at satellite altitude, and we have used the Subtractive Optimally Localized Averages (SOLA) technique to study the secular variation of the radial field and its time changes down at the core-mantle boundary. These are both local methods whereby field measurements in the vicinity of the location of interest are combined so as to estimate the field at that point; measurements far away from the site of interest, have little or no influence on the field estimates.
Using the GVO method, we derived composite time series of geomagnetic secular variation, spanning nearly 20 years, on a global grid of 300 GVOs. GVO time series derived from CHAOS-6-x9 and IGRF-13 calibrated CryoSat-2 data, show similar sub-decadal SV features, and comparable level of scatter. We found a scatter level for radial field SV of 3.5 nT/yr for the CryoSat-2 GVOs compared with 2.5 nT/yr for CHAMP and 1.0 nT/yr for Swarm. Comparing GVOs with overlapping epochs from 2014 to 2018, derived from CHAOS-6x9 calibrated CryoSat-2 data and Swarm data, we find similar sub-decadal SV changes, thus confirming the possibility of using CryoSat-2 for core field studies. In our 20 yr long composite record, we observe fluctuations in the radial SV field of up to 20 nT/yr in amplitude occurring at low latitudes over time periods of 5-10 years. For instance, we see a rapid change of slope over Indonesia around 2014, over the South America and South Atlantic region around 2007, 2011 and 2014, and in the central Pacific around 2017. Some of these events have previously been discussed in ground observatory records (Brown et al., 2013;Torta et al., 2015). They have the distinct "Λ" and "V "-shapes that are often associated with geomagnetic jerks, but as indicated in earlier investigations (Chulliat and Maus, 2014;Olsen andMandea, 2007, 2008) we find these events are localized rather than global features.
Using the SOLA technique, we mapped the radial field SV directly at the coremantle boundary using satellite data, computing spatially localized averages of the SV, time-averaged over specified windows. Taking differences between consecutive SV estimates we obtained estimates for the secular acceleration (SA) at the CMB. Using only CryoSat-2 measurements we are able to successfully map the SA at the CMB, down to spatial averaging widths of ≈ 42 • , corresponding approximately to SH degree 10 or length scales of 2500 km at the CMB. Comparing SV and also SA field maps at the CMB, derived from CryoSat-2 and Swarm measurements, the same features can be identified having similar amplitudes and latitude/longitude extend.
In time-longitude plots of SOLA-based radial field SA estimates along the geographic equator at the CMB we find strong SA features, with amplitudes of ±2.5µT/yr 2 , under Indonesia at longitudes 70 • to 160 • from 2011 to 2014 during CryoSat-2 time, under central America at longitudes 240 • to 320 • appearing from 2015 to 2019, and in sequences with alternating SA signs under South America and the South Atlantic region at longitudes 240 • to 360 • during 2004-2019. The imaged SA features around 2013 at longitudes 70 • to 160 • and 280 • 360 • demonstrates the usefulness of the CryoSat-2 measurements. Our results lend support to a sign changing sequence of SA, for observable length scales down to 2500 km in this study, at longitudes 240 • to 360 • from 2005 to 2019, that has been noticed in previous studies (Alken et al., 2020a;Chulliat et al., 2015). We have shown it is possible to increase the temporal resolution of SA estimates during Swarm era, compared to that seen for example in the CHAOS-7 model, by computing SA estimates from the differences of consecutive SOLA SV estimates derived using 1-yearly time windows. We find the similar coherent structures as seen in TL plots constructed using 2-yearly time window, but also see changes of sign in the SA within 1 year. In the central Pacific region at longitudes 150 • to 220 • we find strong positive and negative SA features appearing side-by-side in late 2017 drifting westwards until 2020. The results presented in Figure 9 and 10, demonstrate that estimates of core field SA different from those found in the CHAOS model can be obtained despite similar data selection and external field modelling schemes. This is due to the important role played by the model parameterization and regularization in the SA recovered in the CHAOS model, especially for the small length-scales and fast time changes which are at the limit of what can be reliably resolved from the data.
The rapid fluctuations of the core magnetic field described in this study are likely caused by time variations in the motions of the liquid metal outer core. In particular, changes in secular acceleration patterns at low latitudes provide constraints on the equatorial dynamics of the outer core (Aubert and Finlay, 2019;Kloss and Finlay, 2019). This is a topic of active research with various phenomenon recently proposed including equatorially-trapped MAC waves in a stratified layer close to the core surface (Buffett and Matsui , 2019;Chi-Durán et al., 2020;Gerick et al., 2020) or the equatorial focusing of hydrodynamic waves driven by turbulent convection deep within the core (Aubert and Finlay, 2019).
It is undeniable that much core dynamics occurs on time-scales either much longer, or much shorter, than can be can be resolved using the current available satellite and ground observations (Gillet, 2019). However, with 20 years of low-Earth orbit satellite measurements of the vector magnetic field, and with tools similar to those presented here, we are now able to probe and characterize core field changes with increasing spatial and temporal resolution. We have shown here that platform magnetometer data can help in this activity, provided they are appropriately calibrated.