A harmonic spline magnetic main field model for Southern Africa combining ground and satellite data to describe the evolution of the South Atlantic Anomaly in this region between 2005 and 2010

CHAMP satellite and relatively densely spaced ground-based data measured over southern Africa between 2005 and 2010 are combined in a new regional, harmonic spline based, core field model. This new SACFM-2 model is compared to the regional SARM model, which is based only on satellite data, and the global CHAOS-6 model. The results agree well in the vertical (Z) component, with somewhat larger differences in the horizontal components. The Z component and total intensity F are used to investigate the evolution of the South Atlantic Anomaly in this region. The computed maps of main field of the Z component and total intensity F show a steady decrease in the field over the years during the study period, indicating the evolution of the South Atlantic Anomaly over southern Africa and suggesting an increase in the area of this feature.


Introduction
The largest part of the magnetic field measured at the Earth's surface can be ascribed to a self-exciting geodynamo process that takes place inside the outer core. This outer core is fluid and consists mostly of iron and nickel that are good electrical conductors. The main field is generated by electrical currents flowing in the outer core (e.g. Hollerbach 1996). It undergoes longterm changes with decrease and increase in the field at irregular and unpredictable rates. This phenomenon of a long-term time variation of geomagnetic field is also known as secular variation, while abrupt changes in the secular variation trend are referred to as geomagnetic jerks (Courtillot et al. 1978). There are additional sources of the Earth's magnetic field. The lithospheric field results from the rocks magnetisation in the crust and upper mantle. Electrical currents in the ionosphere and magnetosphere, resulting from interactions between solar wind and the near-Earth space environment, add a significant contribution of external sources.
Monitoring the field evolution by direct observations has revealed that the magnetic dipole of the Earth oriented along the Earth's rotational axis has decreased in the last 175 years and it has become now 9% weaker than it was in 1840 (e.g. . The decrease in the field is attributed to the changes taking place in the molten iron outer core. In a recent study,  found that the majority of the dipole decay is caused by meridional flux advection since 1840 with magnetic diffusion making a small contribution. The chosen area of study, southern Africa, lies in a region where magnetic field intensity has decreased rapidly and this region covers the area from southern Africa over the South Atlantic Ocean to South America. (e.g. Badhwar 1997). Southern Africa borders a region known as the South Atlantic Anomaly (SAA) where the field has already weakened by approximately 30% compared to other regions located at similar latitudes. The study of the long-term variation of the geomagnetic field over southern Africa can contribute to monitoring the evolution of the South Atlantic Anomaly (SAA) (Fig. 1). The SAA is essential in explaining how the dipole decay rate may be controlled by a planetary-scale gyre in the molten outer core . Some geophysicists see this anomaly as an indicator of forthcoming geomagnetic transition, such as a dipole excursion or reversal (Pavon-Carrasco and DeSantis 2015). The weakening in field intensity can have a negative effect on the magnetic field shielding of the Earth against high energetic particles and cosmic rays that can penetrate into the magnetosphere and atmosphere, creating space weather hazards for technology such as, e.g. communication satellites (e.g. Heirtzler et al. 2002).
An investigation into a combination of both satellite and ground-based data to develop a Southern Africa Core Field Magnetic Model (SACFM-2) is conducted using a harmonic spline technique based on data recorded by the CHAMP satellite between 2005 and 2010, ground-based data collected at 38 geomagnetic repeat stations in southern Africa and observatory data recorded at four INTER-MAGNET observatories (Hermanus, Hartebeesthoek, Tsumeb and Keetmanshoop) during the same time period. Regional field studies of the long-term time variation of the geomagnetic field in southern Africa using the harmonic spline technique were conducted before using only either ground-based data ) or satellite data (Nahayo et al. 2015). The method is summarised in section "Methods", and the utilised data are described in section "Data selection and processing". The new model and its fit to the different data types is presented and compared to previous models in section "Results" with the help of the selected points over the area of study in Table 1, and further discussed in section "Discussion". Finally, the evolution of the SAA in southern Africa between 2005 and 2010 is discussed based on the new model, before concluding with an outlook to future monitoring of the magnetic field and its evolution in the region (Section "Conclusions").

Methods
The aim of this paper is to use a regional main field model, derived from CHAMP satellite and ground-based magnetic data, to investigate the evolution of the SAA over southern Africa. This regional model is derived using harmonic splines, following the approach recently applied by Geese et al. ( , 2011 to southern African ground data. The harmonic splines were first introduced by Shure et al. (1982) for global magnetic field modelling using Magsat data, but they are also appropriate for regional modelling. The harmonic spline functions meet Fig. 1 The long-term decrease in total field intensity over southern Africa. Evolution of total field intensity from 1940 to 2010 at Hermanus (HER), Hartebeesthook (HBK) and Tsumeb (TSU) magnetic observatories in southern Africa. There is a decrease of 20% of total field intensity at Hermanus between 1940 and 2010 the Laplace's equations requirements; they allow a potential field approach combining the individual magnetic field components in a physical significant way (Geese et al. 2011). The technique allows the successful integration of both satellite and ground-based observations at different altitudes above the Earth's surface. We summarise the method starting with the spatial aspect and introducing at the end a time dependency in the model. According to Lesur (2006) and , the Earth's magnetic field can be written as: where the functions F Lr j (ϑ, ϕ, r), F Lϑ j (ϑ, ϕ, r) and F Lϕ j (ϑ, ϕ, r) are defined as follows: where L l,m denotes the double sum L l=1 l m=−l , ϑ, ϕ and r are the geocentric colatitude, longitude and radius, a is the Earth's reference radius (6371.2 km), Y m l (ϑ, ϕ) are the usual Schmidt semi-normalised spherical harmonic functions of degree l and order m, L is the maximum degree of expansion of internal sources and f l = 2l+1 4�(l+1) 4 l 2 . The parameter f l defines the shape of the localised functions.
Considering the above description of the magnetic field representation, and provided that the functions are defined at the observation locations, the fitting of threecomponent data as a potential field is possible, so long as the data contain spherical harmonic degrees less than or equal to L. Equation (1) is not time dependent and it determines a linear system of equations where α r j , α ϑ j and α ϕ j are unknowns, the data are magnetic field measurements and the matrix is expressed by combinations of the functions F Lr j , F Lϑ j and F Lϕ j . The solution of this linear system can be found directly. For time representation, each of the coefficients α r j , α ϑ j and α ϕ j is expanded on a basis of B-splines, i.e piecewise polynomials between spline knots: and α ϑ j and α ϕ j are expanded in a similar way. No regularisation was applied, i.e. the spatial smoothness of the model is determined by degree L and the temporal variability by the number of splines. To suppress the contributions of the lithospheric field to the derived core field model, the maximum degree of expansion of internal sources L was set to 12. The degree 12, with the minimum wavelength of 4000 km, is small enough to retain mostly the field contributions from the core. For the time representation, we use order six B-splines with 17 spline knots spaced at 1 year intervals between 2005.0 and 2011.0. A linear system is built from Eqs. (1) and (5) and using a least squares method, we solve for the unknown model parameters β r k j , β ϑ k j and β ϕ k j as defined in Eq. (5).

Vector satellite data
Data selection for magnetically quiet times are commonly done in order to eliminate external field influences as best possible from core field models when using satellite data (e.g. Lesur et al. 2011). We performed data selection on CHAMP vector magnetic data of 1 Hz sampling interval recorded between approximately 300 and 380 km of altitude. We selected night data recorded between 2005 and 2010 during geomagnetic quiet time over the southern African region covering the area between 16 • S and 36 • S in latitude and 10 • E and 34 • E in longitude. We considered only quiet time data recorded when a Dst index was between − 20 and 0 nT during local night times between 20:00 and 06:00. The time interval was chosen to eliminate the influence of the solar quiet day variation. The external field source dominating at mid-latitudes is the ring current. The proxy for the strength of the ring current, Dst index (e.g. Fares Saba et al. 1997), was therefore exclusively used in data selection. Before setting the data selection criterion of Dst values to be between − 20 and 0 nT (Nahayo et al. 2015) a trade-off was made between getting quiet time data and a good data coverage. The model input data were obtained by creating 525 data bins of 0.2 • × 0.2 • on a grid of 1 • of latitude and 1 • of longitude. Additionally 120 data bins on a grid of 2 • of latitude and 2 • of longitude and offset by 0.5 • were used to increase the number of data points. The middle of each bin is named a data centre in the following. Six hundred and forty five (645) satellite data centres were created (Fig. 2). We averaged data recorded at the same epoch (within 1 day) in each bin, which provided values at the same geographic location at different epochs and different altitudes. A data centre thus can be considered similar to a geomagnetic repeat station, with the differences that the data might be measured at different altitudes between 300 and 380 km and that repeat intervals are not annual. The altitudes of the values in the same bin recorded at the same epoch were found to be very close (with an average scatter value of less than 10 m) so that no altitude correction was necessary before averaging. The altitude of the averaged value was taken as the average of altitudes of individual values in the bin recorded at the same epoch. Considering all used bins, the average scatter values between data values in the same bin recorded at the same epoch are 14.5, 5.9 and 3.8 nT for X, Y and Z components, respectively.

Ground data
The hourly mean values derived from 1-min data from the four continuous recording magnetic observatories Hermanus (HER), Tsumeb (TSU), Hartebeesthoek (HBK) and Keetmanshop (KMH) (also available at http://www. intermagnet.org) were selected using the same criteria as for satellite data selection. Monthly mean values were computed as averages of the obtained hourly values to get the input data for the model.
The field survey data collected at 38 geomagnetic repeat stations from South Africa, Namibia and Botswana were also used as input data (also available at http://www.geomag.bgs.ac.uk/data_service/data/surveydata.shtml). These data were collected by running a LEMI-008 fluxgate vector magnetometer over night and taking absolute observations in the evening and morning using the standard absolute equipment of a DIfluxgate magnetometer (http://www.bartington.com/ mag-01h-declinometer-inclinometer-system.html) and an Overhauser magnetometer (http://www.gemsys.ca/ technology/ground/). The absolute observations provide baseline values for the variometer values, which are averaged over the night hours between 20:00 and 06:00 LT to obtain the main field values at this geomagnetic repeat station for this particular night when the survey is conducted (Korte et al. 2007). The repeat stations thus provide approximately annually spaced values for quiet night times. The location of all ground stations is indicated in Fig. 2.

Model input data
The data used to produce the model fall in three categories, i.e. the monthly mean values from observatory data, the spot values from geomagnetic repeat stations and the averaged values from satellite data centres ( Table 2). The selection of weights for the different data categories (see Table 2) was chosen assuming that the external field (noise) has cancelled out better in the ground data, in particular the time-averaged observatory data, than the satellite data and taking into account the different numbers of data with strong dominance of the satellite data. Remaining external field contributions to data were further removed using the comprehensive model phase 4 (CM4) (Sabaka et al. 2004;Onovughe 2016). The computed CM4 external field values, primary and induced magnetospheric and ionospheric contributions, were removed from all data before the averaging process. The crustal field was removed from both satellite and ground-based data before the averaging process using the Potsdam Magnetic Model of the Earth (POMME-9) (Maus et al. 2006). This model includes in its representation the field contribution from crustal magnetisation up to spherical harmonic degree and order 133 (300 km resolution) and the maximum expansion degree for the core field, that was not removed, was set to 15.

Results
The SACFM-2 model is based on harmonic spline fitting of vector satellite and ground-based data recorded between January 2005 and December 2010. The rootmean-square (rms) values of the difference between input data and model values are 10.7, 6.3 and 4.7 nT for the X, Y and Z components, respectively. Spline-based models often suffer from end effects and we therefore consider SACFM-2 to be reliable from 2005.5 to 2010.5.
A comparative evaluation of SACFM-2 and CHAOS-6 ) at 0.8 km altitude (average altitude for ground-based data) and for all three components X, Y and Z over the whole area of study is presented in Tables 3 and 4 for main field and secular variation, respectively. CHAOS-6 is a global spherical harmonic model spanning the time interval 1999.0-2016.5. For the time of interest here it is based on monthly means from 160 ground observatories and CHAMP and Oersted satellite quiet time vector and scalar data. Differences in secular variation between the two models are generally small on average except towards the end of the time interval in X and Z. There is a systematic, temporally nearly constant difference between the two models of the same order of magnitude in all three components. Figure 3 shows comparative maps of all three field components from the two models for epoch 2008.0. Figure 4 presents the time evolution of the main field of X, Y and Z components for both models SACFM-2 and CHAOS-6 at four magnetic observatories in the region (HER, HBK, TSU and KMH).
Three groups of input data, the satellite data, the observatory data and the repeat station data, and one group of independent data of three reference ground points were used to evaluate the SACFM-2 model in comparison with the regional model SARM (Nahayo et al. 2015) and global model CHAOS-6. The independent data at satellite level were also used to compare the performance of both SACFM-2 and CHAOS-6 models (Tables 5 and 6). SARM is built by the same modelling technique but based on satellite data only, and the spatial data coverage is less dense (data centres are on a grid of 2 • × 2 • of latitude and longitude) with the maximum degree of expansion of internal sources set to 12. The same maximum degree 12 was used in getting synthesised data from CHAOS-6. Quiet time external field correction and lithospheric field elimination had not been performed on the SARM data set. The results of rms values for this comparison are given in Tables 7 and 8. The rms fit to the input data of SACFM-2 is very close, in general below 12 nT, in all components and for all three data types. For the satellite data, the misfit of SARM is two to three times larger, and the misfit of CHAOS-6 is clearly much larger, in particular for X and Z. For the observatory data the fit of SACFM-2 is much better than that of SARM or CHAOS-6 in all Table 2 The SACFM-2 model input data in three catego-

SAT OBS GRS
Number of data centres 645 4 35 Number of data points 6150 290 210 Weight of data 25 100 50 components (less than 7 nT rms misfit in SACFM-2 compared to 43 to 168 nT in the other two models). CHAOS-6 fits this data type clearly closer than SARM in X and Z, with some apparent time dependence, and the fits to Y are of the same order for both these models. The misfit to the repeat station data are of similar magnitude for SARM and CHAOS-6, slightly lower for X in the latter, and notably lower in all components for SACFM-2. In Figs. 5, 6, 7 the secular variation from data, SACFM-2 and CHAOS-6 is shown for a couple of widely distributed locations. The secular variation values for each component were calculated from monthly observatory and model values as follows: where the unit of t is months. SACFM-2 shows notably higher temporal variability than CHAOS-6, most pronounced in the X component. Figures 8 and 9 show the evolution of the Z component and total intensity F over the time interval 2005.5 to 2010.5 as predicted by the SACFM-2 model, indicating a weakening of the field strength over most of the studied area during that time.

Discussion
The comparative evaluation of SACFM-2 using regional model SARM based only on satellite data and global model CHAOS-6 firstly shows a strong dependence of the different models on their input data (Tables 7 and 8). This is not surprising and has to be kept in mind when using and interpreting models. The fact that the satellite data with external and lithospheric field elimination used here are fit less well by SARM than SACFM-2, which used satellite data without these corrections, indicates that SARM contains more external field leakage and / or lithospheric field influence and clearly is superseded by the new SACFM-2. External field signals in general tend to be stronger than lithospheric field signatures at satellite altitude, so the influence of the latter might be negligible, though. This is supported by the observation that SARM shows the least fit to all ground data among the three models.
Interpreting the differences between the regional SACFM-2 and the global CHAOS-6 model (Tables 6, 7, 8 and Fig. 3) is less straightforward. Except for basis functions, details of quiet time data selection and data processing two further differences between these two models are that CHAOS-6 uses less ground information as it does not include repeat station data, and that CHAOS-6 should better minimise large-scale magnetospheric field leakage by taking this contribution into account separately from the core field in its parameterisation. The main field description of SACFM-2 might be influenced by local lithospheric field residuals at the ground stations that have not fully been eliminated in the data processing (Fig. 4). This would not influence the secular variation of SACFM-2.
The observed model differences are mostly largescale in all components (Fig. 3), yet do not show the geometry of ring current dominated magnetospheric

Table 3 RMS values of the differences between SACFM-2 and CHAOS-6 main field models between 2005 and 2010 at 0.8 km altitude
These differences were computed using 684 data points over the modelled area  field, which could be assumed to leak more strongly into SACFM-2 than CHAOS-6. On the other hand, the scale of the differences argues against local lithospheric field influence in SACFM-2, although it seems surprising that spherical harmonic degree up to 133 as used for the lithospheric field elimination from the data should be enough to remove local crustal fields to residuals of less than 10 nT in ground data (Tables 7  and 8). Future expansions of the modelling method to co-estimate crustal biases for the ground stations (e.g. Geese et al. 2011;Talarn et al. 2017) might shed light on the influence of the lithospheric field on the regional model. Trends of increasing residuals to observatory data in X and Z from 2005 to 2008 in CHAOS-6 not seen in SACFM-2 might indicate a more accurate secular variation description by the regional model (Table 7). Figures 5 and 7 suggest that this might result from a slight lack of temporal variability in CHAOS-6 in the southern African region, where secular variation and its gradients are strong in global comparison. However, a similar trend is only seen in the X component for the repeat station data (Table 8).
The comparative evaluation of the SACFM-2 using independent data at four groups of reference points at satellite altitude suggest that the regional model SACFM-2 performs better in all three components compared to Table 5 The presentation of satellite data and computed data from the SACFM-2 and CHAOS-6 main field models for the Z component at four reference group points (G1, G2, G3 and G4) at selected epochs between 2005 and 2010  Table 6 The RMS values of the differences between the SACFM-2 and CHAOS-6 main field models and measured satellite data at four reference group points (G1, G2, G3 and G4) for selected epochs in Table 5 Component SACFM-2-satellite data CHAOS-6-satellite data CHAOS-6 ( Table 6). The noise in X and Y might be the source of the deterioration of the SACFM-2 model at ground level in these components. The noise signal also amplifies with downward continuation. The inclusion of more ground data seems to have improved the model at the Earth's surface in many regions compared to downward predictions from a model from satellite data only, and helps to track the rapid secular variation observed in the southern African region. The vertical field and total intensity maps in Figs. 8 and 9, respectively, show a continuing decrease in field strength over southern Africa from 2005.5 to 2010.5, contributing to a deepening of the SAA. The decrease in the field in magnitude ranges between 0 to 200 nT and 0 to 250 nT in the Z component and total intensity F, respectively, with faster decrease in the southwest corner of the studied area. However, there is a very small area in the south-east corner that shows the increase in the field strength. The secular variation curves in Fig. 7 confirm that the decrease ranges between 0 and 55 nT/year, with faster changes in the south (HER, SWP) than the north (TSU, MAU). Strong secular variation gradients manifest themselves in the observed patterns. In the northern part of the region, the rate of change in the Z component decreases from east to west (see MAU and TSU), while it decreases Table 7 The RMS values of the differences between satellite and observatory input data and the three models, SARM,

SACFM-2 and CHAOS-6 for the period 2005-2010
In 2005 and 2010, only the input satellite or observatory data recorded after 2005.5 and before 2010.5 were used, respectively from west to east in the southern part (see KMH and HBK or SWP and HER). The rate of change from west to east increases in the X component (Fig. 5) and does not show a clear trend in the Y component (Fig. 6). The magnitude of the Y component is decreasing in the northern part and south western part of the region (see TSU, MAU and SWP) and increasing in the south eastern part (see HER and HBK). Both the global and the regional model (CHAOS-6 and SACFM-2) identify the 2007 geomagnetic jerk in Y and Z components (Figs. 6 and 7). There is also an indication of another geomagnetic jerk during 2009 of opposite sign to the 2007 jerk that is indicated by some of the SACFM-2 Z component curves (SWP, KMH, TSU and MAU). This confirms results by Chulliat et al. (2015) and Kotze and Korte (2016), who used only ground-based field survey and observatory data.

Conclusions
We have presented a new regional magnetic field model for the southern African region, named SACFM-2, which for the first time incorporates both ground and satellite data. The ground data coverage is improved compared to global models by using data from 38 repeat stations distributed through South Africa, Namibia and Botswana. The use of both CHAMP satellite and ground-based data has shown the possibility to develop accurate regional models that can be used to characterise the variation of the geomagnetic field over a small area. A comparison of SACFM-2 to the global CHAOS-6 model showed similar results within the area of observation, revealing highly nonlinear time-varying characteristics of the geomagnetic field in southern Africa between 2005 and 2010. SACFM-2 supersedes the previous SARM regional model, which was based only on satellite data without elimination of external and lithospheric field influences during data processing. From comparisons, we conclude that SARM suffered more strongly than SACFM-2 from external field leakage, while the global CHAOS-6 model might slightly underestimate secular variation in this area.
Like other models, SACFM-2 shows a steady decrease in the Earth's magnetic field over southern Africa, thus confirming the deepening of the SAA in this region.
The comparative evaluation of these two models with the help of measured satellite and ground-based data shows that the regional model SACFM-2 correlates better with observations than CHAOS-6 for the Z main field applications. The SACFM-2 fits all component data better at satellite level (Table 8), but noise apparently gets amplified by the downward continuation in X and Y in regions not well-constrained by ground data.
Given the high magnetic field horizontal gradient in Southern Africa, an increase in the number of data centres for better coverage might improve the model further. Moreover, an improved approach for data selection and noise removal is required to improve the horizontal Table 8 The RMS values of the differences between geomagnetic repeat stations input data (GRS) and the three models, SARM, SACFM-2 and CHAOS-6 for the period

2005-2009
The year 2010 is not included because the ground vector survey was carried out in the second half of the year, outside the interval where SACFM-2 is considered to be reliable   components X and Y for downward continuation, and an extension of the modelling method to co-estimate crustal biases of the ground data might avoid lithospheric field leakage into the model. A careful monitoring of the SAA evolution is envisaged using high-quality geomagnetic field observations by European Space Agency Swarm satellite constellation (Olsen et al. 2014), providing denser data coverage than a single satellite mission and an additional third observation level. Future models will expand our current knowledge about the magnetic anomaly region, its geographic expansion towards inside the African continent and how fast the field intensity decreases.