Recovering external contribution from the monthly mean series of a given geomagnetic observatory ANNALS OF GEOPHYSICS, 59, 3, 2016, G0321; doi:10.4401/ag-6971 G0321 Recovering external contribution from the monthly mean series of a given geomagnetic observatory Bejo Duka1,*, Eni Duka2, Klaudio Peqini1 1 University of Tirana, Faculty of Natural Sciences, Department of Physics, Tirana, Albania 2 University of Zurich, Department of Informatics, Zurich, Switzerland ABSTRACT The differences between monthly mean values of the observed geomag- netic field and monthly values predicted by different models of the in- ternal geomagnetic field (named “model biases”) for the time period 2000-2015 at several geomagnetic observatories are analyzed. We no- tice that increasing the maximum degree of the model is not always fol- lowed by the decrease of such “model bias”. The time series of these “model biases” reduced by their average resulted to be approximately the same for all models and should represent the external (non-modeled) contribution to the observed geomagnetic field. These time series for dif- ferent observatories (close or away to each other) are compared and their power spectra are analyzed. Such spectra have common features like the annual and semi-annual variation with some possible sporadic cases of seasonal variation. 1. Introduction The magnetic field measured at any geomagnetic observatory is the result of several magnetic contribu- tions generated by various sources. More than 90% of the measured field is of internal origin and is mostly generated in the Earth’s outer fluid core. Also, of in- ternal origin is the magnetic field generated by the re- manent magnetization of the crust and by the induced magnetization of the crust produced from the core and external fields. The known contributions of external origin are magnetospheric and ionospheric fields, fields of ring currents and aligned field currents. The global models of the geomagnetic field and its secular variation (SV), like IGRF, WMM, EMM, POMME, CM, GUFM, MF, etc., are based on spherical harmonic analysis (SHA). The magnetic field potential (W) is expressed as the sum of spherical harmonic func- tions. The internal sources contribution to the radial variation of these functions is quite different from the external sources contribution: respectively (a/r)l+1 and (r/a)l, where a is the reference radius of the Earth, and l the spherical harmonic degree. The internal source magnetic field is presented by Br, Bi, B{ - orthogonal components in the geocentric reference frame or X (north), Y (east), Z (downwards) as: Using different data sets of ground or/and satel- lite measurements in disposal, different models define the Gauss coefficients (gml ; h m l ) of spherical harmonic expansion up to a maximum number of degrees L and their time derivative up to a maximum number of de- grees L’ < L. Generally, the denser and numerous the data used by the model, the greater is the maximum degree of harmonics (L) and the shorter is the minimum wave- length of harmonics i.e. more detailed is the model de- scription of the geomagnetic field, according to the formula [Backus et al. 1996]: (2) The separation of internal contributions into the core and the crustal (actually lithospheric) origin is . cos sin cos cos sin cos sin cos sin cos B l r a g m h m P B r a gQ m h m P B r a mgmg m mh m P 1 r m l l L l W l m l m l m l W m l l L l m l m l m l W m l l L l m l m l m 01 2 2 01 2 01 2 2 { { i { { i i { { i i = + + = + = - + { i == + + == + == Q R R R gQ Q Q Q Q Q l W l W V l W V V V V V V|| || || a L2.m r +Q / .1 2+ V Article history Received January 22, 2016; accepted June 24, 2016. Subject classification: Geomagnetic field model, Geomagnetic observatory, Internal geomagnetic field, Power spectrum, Spherical harmonics. (1) based on the spherical harmonics spatial power spec- trum [Lowes 1974, Langel and Estes 1982]. It is known that the spherical harmonics up to maximum degree L = 13 represent the core field, while the spherical har- monics of degree greater than 13 represent the Litho- sphere field and there is an overlap of contributions in between 12 and 14 [Cain et al. 1984]. In fact, lithosphere sources have some small contributions to the spherical harmonics of degrees < 13 and core sources have also some small contribution to the spherical harmonics of degrees > 13. Some models aim to describe the Earth’s magnetic field with high spatial and temporal resolution, for ex- ample CHAOS-5 model robustly determined spherical harmonics up to degree L = 90 for the lithosphere field, up to degree L = 20 for the core field and up to degree L’ = 16 for the time-varying core field [Olsen et al. 2014]. Temporal changes of the field of internal origin up to degree 16 are usually attributed to changes in the core field itself and temporal changes of the field of in- ternal origin produced by the induced part of the lith- osphere magnetization could dominate the core field signal beyond degree 22-24 [Hulot et al. 2009]. The time-varying lithosphere field dominates and conceals the time-varying core field beyond that critical degree (22-24), just as the permanent lithosphere field domi- nates and conceals the permanent of the core field be- yond degree 14 [Hulot et al. 2009]. It is more difficult to separate the crustal contri- bution into parts coming from remanent magnetiza- tion and from induced magnetization. [Lesur and Gubbins 2000] developed a new method for such sepa- ration, but its application for several observatories was unsuccessful. The so named “crustal bias” [Mandea and Langlais 2002] in an observatory is estimated by comparing the magnetic components measured in the observatory with those predicted by a geomagnetic model trun- cated to its nuclear part (i.e., up to degree and order of 13). These models are based only on the satellite data sets [Bloxham and Gubbins 1986], where the crustal sources contributions are negligible, therefore the dif- ferences between observed values and model values are considered as the signature of the crustal field of in- duced or remanent magnetization origin. Such “crustal biases” of annual means over 42 years period are well investigated [Verbanac et al. 2007a] for 46 European geomagnetic observatories by using models based on Magsat, Ørsted, CHAMP and SAC-C satellites. Also, the temporal evolution of the observatory monthly means (over 9 years period) “crustal biases” are ana- lyzed by Verbanac et al. [2015] using the G-field model [Lesur et al. 2015]. Even in this analysis, “crustal biases” are considered the differences of the magnetic compo- nents measured at the observatory with the values pre- dicted by a model obtained from satellite data only. Verbanac et al. [2007b] attempted to separate, in- terpret and explain the external field signals obtained by subtracting from the annual means the core field predicted by models truncated at the maximum degree 14. From the obtained external signal they removed magnetospheric contributions of the external field pre- dicted by the POMME-2.5 model and a Sq variation pre- dicted by the CM4 model. Averaging these residuals over representative European observatories, they found interesting results regarding the correlations of the residuals with external sources such as solar cycles. Here we will study a kind of “bias”, which is esti- mated as the difference of the monthly mean values of the geomagnetic field components registered at a given observatory with the respective values predicted by the geomagnetic field models based on ground and satel- lite data. The models are truncated at different maxi- mum degrees providing different “biases” that have mostly signatures of non-modeled local sources of the crust and of external sources. Different models provide different “biases” (named “model bias”) at the same ob- servatory. We will compare the results of three different models for different observatories that are geographically close to each other (WIK: latitude 48.265°N, longitude 16.318°E, altitude 0.4 km; NGK: latitude 52.072°N, lon- gitude 12.675°E, altitude 0.078 km; FUR: latitude 48.17°N, longitude 11.28°E, altitude 0.572 km) and two other observatories located far away from them (HER: latitude 34.43°S, longitude 19.23°E, altitude 0.026 m) (KAK: 36.23°N, longitude 140.18°E, altitude 0.036 km). We expect that the difference (“model bias”) be- tween the observed value and predicted value given by a model should decrease as the maximum degree of spherical harmonics of the model increases. But this does not happen always for all components at every place on Earth, as it will be shown in Section 2. The global geomagnetic models supply Gauss co- efficients and their secular variation for every 0.5 year, 1 year, 2.5 years, etc. Then the geomagnetic field at a given time is calculated from such models by interpo- lating the Gauss coefficient values in between time intervals, i.e. each model has limitations on time dis- tinction. When the monthly values of the internal ge- omagnetic field are predicted by a model, the time distinction of one month is considered. On the other side, the monthly values series of observatories are cal- culated by averaging the observatory minute registra- tions for all days of a month. Therefore the geomagnetic field variations with period less than a month (pulsa- tions, quite daily variations Sq, substorms/storms effects) DUKA ET AL. 2 3 are principally removed. By removing the averages of the differences be- tween the predicted monthly series and the respective monthly mean of observed value series, we can recover the un-modeled external source contributions to the ge- omagnetic field at the given observatory (see Section 3). In Section 4, we will compare the time series of un-modeled external contributions for different obser- vatories (NGK, FUR, WIK, KAK, HER) that are of the same length (2000-2015) and will try to find some com- mon features from their spectral analysis (Section 5). 2. Differences of different degree models We considered three different models that have a large maximum degree: CM4 (L = 65)/CM5 (L = 85), POMME-9 (L = 133), EMM2015 (L = 720). The comprehensive model, CM4 [Sabaka et al. 2004], covers the whole time interval from 1960 to 2002. It has been derived from quiet-time POGO, MAGSAT, Ørsted and CHAMP satellite data in combination with obser- vatory hourly means. The internal field is expanded to different number of maximum degree of spherical har- monics up to 65, resolving magnetic anomalies down to 611 km. The model also calculates primary and induced magnetospheric field, primary and induced ionospheric field. It offers calculation of the geomagnetic field by two steps. In each step, the user inputs the first and second maximum number of degrees. Using codes from (http:// core2.gsfc.nasa.gov/CM/codes.html), the monthly val- ues series of core field (L = 14) and lithospheric field (with different values of maximum degree, from L =15 to L = 65) at different observatory locations are generated. The new generation of CM models, CM5 model [Sabaka et al. 2015] used an extended data set and defined the Gauss coefficients of the core field up to maximum de- gree of 20 for the period 2000-2015 and the Gauss co- efficients of the crust field up to maximum degree of 100. Using cubic B-spline, the Gauss coefficients of SV with 6 month knot spacing are saved. Using Matlab code from (http://www.spacecenter.dk/files/magnetic-models/ CHAOS-5), the monthly series (2000-2015) of core and lithospheric field are generated. POMME-9 (http://geomag.org/models/POMME- 9.html) is an internal field model representing the geo- magnetic field in the region from the Earth’s surface to an altitude of a couple of thousand kilometers [Maus et al. 2004, 2006, 2010]. The time variations of the inter- nal field are given by a piece-wise linear representation of the spherical harmonic (Gauss) coefficients of the magnetic potential. POMME-9 was produced from CHAMP satellite vector magnetic measurements from July 2000 up to September 2010, Ørsted satellite total field measure- ments from January 2010 to June 2014 and Swarm satel- lite vector magnetic measurements from December 2013 to January 2015. It has the same parameterization of the magnetospheric field as POMME-6, POMME-7 and POMME-8. POMME-9 extends the maximum de- gree of spherical harmonics to 133, resolving magnetic anomalies down to 300 km. The Gauss coefficients of every year (from 2000 to 2015) are embedded into the header (.h) files for the core and crust field. The pro- gram calculates the core geomagnetic field every day from the starting year 2000, while the crust field coeffi- cients are static from degree 16 to 133. Enhanced magnetic model (EMM2015) based on the ellipsoidal harmonic representation of Earth’s lithos- pheric magnetic field [Maus 2010] represents the main field with maximum degree L =15 and the lithospheric field up to L = 720. It was compiled using a vast amount of data from satellite, marine, aeromagnetic and ground magnetic surveys. It also includes data from the Euro- pean Space Agency’s Swarm satellite mission. This model (https://www.ngdc.noaa.gov/geomag/EMM/EMM SurveySPH.shtml) covers the 2010-2020 period resolv- ing magnetic anomalies down to 56 km. NON-MODELED EXTERNAL CONTRIBUTION Figure 1. Comparison of short monthly series of geomagnetic field observed at NGK observatory and predicted by different models (CM4, CM5, POMME-9, EMM2015). We tried to find a common time interval for all considered models. All considered models have a short common time interval 2000-2002. We firstly compare the monthly mean series of the NGK geomagnetic ob- servatory with respective monthly values series pre- dicted by these models for that short time period (Figure 1). One can see that series predicted by CM4, CM5, EMM2015 represent a smooth variation, while the series predicted by POMME-9 model have some ir- regular variations much like the series of observations. Even the monthly values series received by averaging the daily values series generated by POMME-9 show the same behavior. It seems that this model fits better to the observed monthly values series than other mod- els. The CM5 model represents a greater baseline jumps from the observed monthly values series. In Figure 2, the “model biases” are presented, i.e. the differences of monthly means of X, Y, Z compo- nents observed at Niemegk observatory (NGK) during the common time interval 2000-2015 with respective values predicted by: a) EMM2015 model of maximum degree of 720 (dashed line), b) POMME-9 model of maximum degree of 133 (solid line) and c) CM5 model of maximum degree of 85 (dotted line). Comparing the results of different models, one can see that for the X component, the POMME-9 model fits the observations better (see also Table 1, where the val- ues of averaged “bias” and their errors for different mod- els and different maximum degree are presented). While for the Y and Z components the EMM2015 model (re- spective averaged “bias”: −1.9 nT and −7.2 nT) fits bet- ter with the observations than the POMME-9 model (respective averaged “bias”: 4.9 nT and −56.1 nT). The CM5 model fits the observations worse than the two other models for all three components. It seems that CM5 model gives worse results than CM4 model. See for example Figure 4, where the results of CM5 at WIK observatory are presented, and Table 2, where the av- eraged “biases” of CM4 model at NGK for a long pe- riod (1960-2002) are presented. The same comparison for two observatories, Furstenfeldbruck (FUR) and Wien Kobenzl (WIK) ob- servatories are presented respectively in Figure 3 and Fig- ure 4. As it is seen in these figures, the POMME-9 and EMM2015 models fit the observed data better than CM5 model. Also all models provide the same time variations. In order to study the impact of the maximum de- gree number of the model on the fitting quality between DUKA ET AL. 4 Figure 2. Differences between monthly mean values observed at Niemegk observatory and respective values predicted by EMM2015 model (dashed line), POMME-9 model (solid line) and CM5 model (dotted line). NGK 2000-2015 EMM2015 model averaged “bias” (in nT) POMME-9 model averaged “bias” (in nT) CM5 model averaged “bias” (in nT) L = 720 L = 520 L = 320 L = 120 L = 133 L = 95 L = 65 L = 45 L = 100 L = 85 L = 60 L = 40 X -24.5± 0.6 -21.7±0.6 -20.5±0.6 -3.6±0.6 1.5±0.7 -20.2±0.7 -43.5±0.7 -43.8±0.7 -24.5±0.6 -21.7±0.6 -20.5±0.6 -3.6±0.6 Y -1.9±0.2 12.1±0.2 5.6±0.2 28.1±0.2 4.9±0.3 25.9±0.3 4.0±0.3 -3.5±0.3 -1.9±0.2 12.1±0.2 5.6±0.2 28.1±0.2 Z -7.2±0.5 17.2±0.5 7.1±0.5 -32.9±0.5 -56.1±0.6 -34.9±0.6 1.4±0.6 -33.8±0.6 -7.2±0.5 17.3±0.5 7.1±0.5 -32.9±0.5 Table 1. Averaged “bias” of different models at NGK observatory (2000-2015 period). 5 the observed and predicted series, the monthly values series are generated by EMM2015, POMME-9, CM4 and CM5 models with different numbers of maximum de- grees at different observatories. Some results are pre- sented in the plots of Figures 5, 6, 7, and in Table 1. In Figure 5, one can see that decreasing the maxi- mum number of degrees from L = 133 to L = 95, or L = 65 for the POMME-9 model worsened the fitting qual- ity for X component, but not for Z component. In Figure 6 (EMM2015 model), one can see that changing the maximum number of degrees from 720 to 520 gives better fitting for the X component (see also Table 1), while for Y and Z components the fitting qual- ity is poorer. When the maximum number of degrees is decreased to 120 we see the fitting is worsening for all components. In Figure 7 (CM5 model), one can see that fitting behaves differently for different components when the maximum number of degrees (L) is changed. The dif- ference series of Y component is always decreased when L is changed from 100 to 40, while difference series of NON-MODELED EXTERNAL CONTRIBUTION NGK 1960-2002 CM4 model averaged “bias” (in nT) L = 65 L = 45 L = 35 L = 25 X -62.25 ± 0.43 -58.83 ± 0.43 -64.1 ± 0.43 -50.6 ± 0.43 Y 17.1 ± 0.15 4.84 ± 0.15 -8.75 ± 0.15 1.53 ± 0.15 Z 8.57 ± 0.25 -16.7 ± 0.25 -36.55 ± 0.25 -54.9 ± 0.25 Figure 3. Differences between monthly mean values observed at FUR observatory and respective values predicted by EMM2015 model (dashed line), POMME-9 model (solid line) and CM5 model (dotted line). Table 2. Averaged “bias” of CM4 model at NGK observatory (1960- 2002 period). Figure 4. Differences between monthly mean values observed at WIK observatory and respective values predicted by EMM2015 model (dashed line), POMME-9 model (solid line) and CM5 model (dotted line). DUKA ET AL. 6 Figure 5. X, Y, Z - component differences between monthly mean values at NGK observatory and respective values predicted by POMME-9 model with maximum degrees L: a) 133 (dash-dotted line), b) 95 (solid line), c) 65 (dotted line), d) 45 (dashed line). Figure 6. X, Y, Z - component differences between monthly mean values at NGK observatory and respective values predicted by EMM2015 model with maximum degrees L: a) 720 (dash-dotted line), b) 520 (solid line), c) 320 (dotted line), d) 120 (dashed line). Figure 7. X, Y, Z - component differences between monthly mean values at NGK observatory and respective values predicted by CM5 model with maximum degrees L: a) 100 (dash-dotted line), b) 85 (solid line), c) 60 (dotted line), d) 40 (dashed line). 7 X and Z components are considerable increased when L changes from 100 to 85 then shift a little when L changes from 85 to 60 and from 60 to 40. Regarding the CM4 model, as its validity period is 1960-2002 we have applied this model for the period 1987-2002 to WIK observatory and for the period 1960- 2002 to NGK observatory. As can be seen in Figure 8 and in Table 2, there is a tendency of decreasing of the average “bias” by increasing of the maximum number of degrees from 25 to 65 for the Z component at NGK observatory, that is not so clear for X and Y components Almost the same happened for WIK observatory (the results are not shown here). 3. Recovering external contribution We have considered the averaged “biases” as con- tributions of the crustal sources: remanent magnetiza- tion and induced magnetization from the internal sources that are not modeled by the given model. Therefore, the residuals of subtracting such average from the model “bias” should represent the contribution of non mod- eled external sources at a given place (observatory). The time series of differences between the geo- magnetic field values observed at an observatory (that are the superposition of magnetic fields from internal and external sources) and respective values predicted by a model of the internal geomagnetic field, are: where i = 1, 2 …N indexes the series of months, k in- dexes the different models. The averages of these series are: where the first term is the averaged of the internal part not modeled by the k-th model and the second one is the average of external contribution to the observations. Subtracting these averages from the series of dif- ferences Dxi, we have: Considering the un-modeled part of internal con- tribution is almost unchanging during the considered time interval (16 years), then: and These series of residuals do not depend on the k- index, i.e. the series are principally the same for all models and should represent the non-constant external contri- bution to the given observatory. These series are calcu- lated as residuals of average subtraction from differences between monthly mean values and respective values pre- dicted by different models. These series for three different models EMM2015 (L = 720), POMME-9 (L = 133) and CM5 (L = 85), at NGK observatory are plotted in Figure 9. The graphs are almost the same for all models and they coincide for each component X, Y, Z. That means such series have the time variation contributions of the exter- nal sources at the NGK observatory for the period 2000- . int int int x xD averagage k x x exext x x x exext ( ) ( ) ( ) ( ) mod mod mod mod i k i k i obs i obs i k un k obs Dl xD = - + -= + + - Dl r rQ Q Q Q Q QV V V V V V int int intx x x( ) ( )mod modi obs i k un k, + rQ Q QV V V .x x exext x exext( )modi k i obs obsD = -l rQ QV V ,V int int int x N x exext N x x exext averaragage k N x 1 1 1 ( ) ( )mod mod i obs i obs i i k i obs i un k = + - - = +r r Q Q Q Q Q Q V V V V V ,V | | | ,Vint int intxD x x x x xexext ( ) ( ) ( ) mod mod mod i k i obs i k i obs i obs i k xD = - = = + -Q Q Q QV V V ,V NON-MODELED EXTERNAL CONTRIBUTION Figure 8. X, Y, Z - component differences between monthly mean values at NGK observatory and respective values predicted by CM4 model with maximum degrees L: a) 65 (solid line), b) 45 (dashed line), c) 35 (dash-dotted line), d) 25 (dotted line). (2) (3) 2015. While the averages x– obs(ext) provide the averaged “bias” at NGK observatory that is not modeled by the given model. Their values are presented in Table 1. The same calculations are carried out for different observatories. When plotted for the same observatory, the series of residuals calculated by different models should collapse together. Here we present the results for FUR nearby NGK observatory (see Figure 10) and two other observatories far away from NGK: HER (see Figure 11) and KAK (see Figure 12). In all figures one can see that after removing the averages, the “model bias” is the same for all models. Again, we noticed that averaged “bias” (non modeled static contribution) for CM5 model is greater than for other models (see Table 3, where the averaged “bias” for different models and different components are presented). 4. Comparison of external contributions from different geomagnetic observatories Plotting the external contributions, calculated as described in Section 3, at three nearby observatories NGK, FUR, and WIK, one can see (Figure 13) that the external contributions in these observatories are almost the same. This can be seen in all components. The ex- ternal contribution to the eastward component (Y), which is supposed to be the least affected by the exter- nal fields [Mandea et al. 2010], looks like noise. In Figure 14, the external contributions at three observatories (NGK, KAK and HER), located far away each other, are plotted. It can be seen that the differ- ences between such contributions, at the three obser- vatories that are more than 100° in latitude or longitude apart, are small ones, especially for the X component. DUKA ET AL. 8 Figure 9. The external contribution to the NGK observatory, calculated as shown in the text for three models: POMME-9 (solid line), EMM2015 (dashed line) and CM5 (dotted line). Figure 10. The external contribution to the FUR observatory, calculated as shown in the text for three models: POMME-9 (solid line), EMM2015 (dashed line) and CM5 (dotted line). 9 While the time series for the Z components are similar, but the series of north observatories (NGK, KAK) are inverted in phase from the series of the south observa- tory (HER). This reflects the sign change of the Z com- ponent from North to South hemisphere. Regarding the Y component, one can notice the expressed noise-like character of these series and the phase inversion can be noticed more between eastern observatories (KAK) and less eastern observatories (NGK and HER). 5. Spectral analyses of the external contribution The external contributions to the monthly mean series of the geomagnetic field observed at different observatories should contain signatures from the same external sources. In order to investigate any common contributors in the series plotted in Figures 9 to 14, we have applied the same spectral analysis method (FFT) on all of these series. The frequency contents (until the Nyquist frequency) of power spectrum of NON-MODELED EXTERNAL CONTRIBUTION Observatory 2000-2015 POMME-9 model averaged “bias” (in nT) EMM2015 model averaged “bias” (in nT) CM5 model averaged “bias” (in nT) X Y Z X Y Z X Y Z HER -15.67±0.76 -8.58±0.35 -34.68±0.37 -19.51±0.77 -18.23±0.33 -24.97±0.35 83.75±0.80 31.78±0.22 -42.29±0.33 KAK 16.9±0.85 9.81±0.27 -79.01±0.31 -1.81±0.87 -1.24±0.14 -44.52±0.27 300.67±0.86 37.68±0.18 -503.78±0.33 Table 3. Averaged “bias” at HER and KAK observatories. Figure 11. The external contribution to the HER observatory, calculated as shown in the text for three models: POMME-9 (solid line), EMM2015 (dashed line) and CM5 (dotted line). Figure 12. The external contribution to the KAK observatory, calculated as shown in the text for three models: POMME-9 (solid line), EMM2015 (dashed line) and CM5 (dotted line). these time series (16 years long) with 1 month time sampling are shown in Figures 15, 16 and 17, respec- tively: Figure 15 for POMME-9 model, Figure 16 for EMM2015 model and Figure 17 for CM5 model. The resolution of these power spectra, i.e. the ability to discriminate the smaller spectral frequency feature of the signal, is about 0.0052 month-1. In these plots, one can see several peaks that are common for all series, lightly shifted for different components, different observatories and different models. The first one being in the frequency interval 0.0039 ÷ 0.0055 month-1 (period interval 21.37 ÷ 15.15 years), should belong to the time series length (16 years). The second peak is in frequency interval 0.02 ÷ 0.025 month-1 (period interval 4.16 ÷ 3.2 years) and the third one is in the frequency interval 0.049 ÷ 0.051 month-1 (period interval 1.7 ÷ 1.6 years). These two last peaks are not so clear. However there are two distinct and clear peaks in the intervals: 0.082 ÷ 0.086 month-1 (period in- terval 1.01 ÷ 0.97 years) and 0.165 ÷ 0.168 month-1 (pe- riod interval 6.06 ÷ 5.96 months). The amplitudes of these peaks are different for dif- ferent components and different observatories. In case of POMME-9 model, the greatest peak is reached at KAK observatory for X and Y components and at NGK observatory for the Z component at the first frequency and the feeble maximum is reached at HER observatory for the Z component. The reason is that the absolute values of X and Y components at KAK observatory are greater than those at other observatories. It can be noticed that the frequency content of the Y component for the POMME-9 model (only for this model) represents a series of small peaks that correspond to the frequencies: 0.3555 month-1 (period 2.82 months); DUKA ET AL. 10 Figure 13. The external contribution (by POMME-9 model) to: a) NGK observatory (solid line), b) FUR observatory (dashed line), c) WIK observatory (dotted line). Figure 14. The external contribution (by POMME-9 model) to: a) NGK observatory (solid line), b) HER observatory (dashed line), c) KAK observatory (dotted line). 11 NON-MODELED EXTERNAL CONTRIBUTION Figure 15. Power spectrum of X. Y, Z - components of external contribution (POMME-9 model) at different observatories: NGK (dotted line), FUR (dashed line), KAK (dash-dotted line) and HER (solid line). Figure 16. Power spectrum of X, Y, Z - components of external contribution (EMM model) at different observatories: NGK (dotted line), FUR (dashed line) KAK (dash-dotted line) and HER (solid line). Figure 17. Power spectrum of X, Y, Z - components of external contribution (CM5 model) at different observatories: NGK (dotted line), FUR (dashed line) KAK (dash-dotted line) and HER (solid line). 0.4375 month-1 (period 2.3 months); 0.4766 month-1 (period 2 months). Maybe these frequencies can indi- cate some external contributions with seasonal origin. We think that such distinction of the Y component of POMME-9 model can be explained by the way the time series of this model is generated. The model is more accurate if the present state of the magnetos- phere is provided as an input in the form of magnetic indices [Lühr and Maus 2010]. Unfortunately, these in- dices are presently not available in real-time, so during real-time evaluation of the model we have to set these indices to zero. Then the time series (especially of Y component) of this model should bear some additional magnetosphere contribution. Similar results are found in the cases of the EMM2015 and CM5 models. The most visible peaks are located in the following frequency intervals: 0.0039 ÷ 0.0059 month-1 (period 21.37 ÷ 14.12 years), 0.0215 ÷ 0.0234 month-1 (period 3.88 ÷ 3.56 years), 0.0332 ÷ 0.0352 month-1 (2.51 ÷ 2.37 years), 0.0488 ÷ 0.0508 month-1 (1.71 ÷ 1.64 years), 0.0820 ÷ 0.0840 month-1 (1.02 ÷ 0.99 years) and 0.1602 ÷ 0.1621 month-1 (6.24 ÷ 6.17 months) for the X component; 0.0039 ÷ 0.0059 month-1 (21.37 ÷ 14.12 years), 0.0332 ÷ 0.0352 month-1 (2.51 ÷ 2.37 years), 0.0488÷ 0.0508 month-1 (1.71 ÷ 1.64 years), 0.0781 ÷ 0.0801 month-1 (1.07 ÷ 1.04 years) and 0.1602 ÷ 0.1621 month-1 (6.24 ÷ 6.17 months) for the Y component; 0.0039 ÷ 0.0059 month-1 (period 21.37 ÷ 14.12 years), 0.0215 ÷ 0.0234 month-1 (period 3.88 ÷ 3.56 years) and 0.0820 ÷ 0.0840 month-1 (1.02 ÷ 0.99 years) for the Z component (Figures 16, 17). Interestingly none of the contributions with sea- sonal origin are visible in the power spectrum of the Y component in both EMM2015 and CM5 models. How- ever the Y component in the case of the CM5 model has a spiky power spectrum, especially in the high fre- quency part, although there are not any distinguishable peaks. This fact indicates the contributions from exter- nal sources of different time scales. A major problem in the spectra of finite length time series is that they suffer from “spectral leakage”. While the infinite-length signal has its power concen- trated exactly at the discrete frequencies, the truncated signal has a continuum of power “leaked” around these discrete frequencies. This “spectral leakage” is more ev- ident when data series are short as our monthly time series are (192 month long). Thompson’s multitaper method (MTM) provides an improved PSD estimate [Thomson 1982] by using a bank of optimal bandpass filters. These optimal filters, i.e. with minimal leakage, are derived from a set of sequences known as discrete prolate spheroidal sequences (DPSSs, also known as Slepian sequences) [Percival et al. 1993]. Another prob- lem when working with truncated time series is vari- ance. It is shown that while reducing spectral leakage one can increase variance and viceversa. Thus an important issue is the appropriate choice of parameters so that a trade-off between spectral leakage and variance is achieved [Percival et al. 1993]. Wardinski and Mandea [2006] used this method in order to investigate the geographic dependence of the annual and semi-annual variations. They use the nota- tion of previous workers, defining K, the number of Slepian tapers, and p which is basically the number of points in the bandwidth W of the time series in the fre- quency domain. The relation between the two param- eters is K = 2p − 1. It is proven that only 2p − 1 tapers are good enough for reducing spectral leakage, while av- eraging their results [Percival et al. 1993] also helps re- ducing variance. For further details we refer the reader to formulas (1)-(3) of Wardinski and Mandea [2006]. The MTM method with parameters described above is provided by “pmtm” Matlab function. For large values of time series length N, one can construct the sequence of Slepian tapers for all the values of p usually used (1, 2, 2.5, 3, 3.5, 4). In our case (N = 192), we have used p = 3.5, that means K = 6. So we use 6 Slepian ta- pers which is considered accurate enough for short time series like what we use in this study [Park et al. 1987]. Apart from the power spectrum density function of frequency: PSD ( f ), the “pmtm” Matlab function calculates the confidence interval. According to our cal- culations, for a 95% confidence level, the true spectral density function lies within the interval from 0.52 to 2.3 of calculated PSD ( f ) function. Decreasing the confi- dence level one can reduce this interval, but the form of PSD ( f ) function does not change. We present here the results for the time series of three different observatories (POMME-9 model; see Figure 18). Similar plots are obtained for EMM and CM5 models that are not shown here. Beside the greatest peak corresponding to the series length (192 months), there are two clear and wide peaks centered in the frequen- cies: 0.08 and 0.16 month-1 corresponding to annual (12.5 months) and semiannual (6.25 months) variations. In case of the Y component PSD (only for POMME-9 model) there appear to be three peaks centered in the frequen- cies: 0.36, 0. 44 and 0.48 month-1, corresponding respec- tively to the periods: 2.78, 2.27 and 2.08 months. 6. Conclusions Different models of the geomagnetic field of in- ternal origin use different maximum degrees of the spherical harmonic expansion of the geomagnetic field components. Comparing the differences between monthly mean values of the geomagnetic field at a given place on DUKA ET AL. 12 13 the Earth and the respective values predicted by different models, we noticed that these differences have almost the same time dependency for the period considered (2000- 2015). These time dependencies are shifted from each other and this shift can be removed by subtracting from these differences their averages. The resulting residuals representing the external contributions for the different models (EMM2015, POMME-9, CM5) at different obser- vatories (WIK, NGK, KAK, HER) show some common features. They are almost in the same amplitude at all observatories (order of 20 nT for the X and Z compo- nents, order of 5 nT for the Y component) are almost the same at the nearly observatories (NGK, FUR, WIK). In particular the X component varies almost the same in all observatories. The time series of the Z component at ob- servatories located in different hemispheres show inver- sion of the phase differences, whilst the time series of the Y component show much more noisy-like behavior. The spectral analyses of these time series for dif- ferent models, different observatories and different com- ponents evidenced some common frequency contents, with maximum of power spectrum for periods of 16 years, 1 year, and 0.5 year. These latest variations show the presence of an annual non-ionospheric and a sea- sonal modulation of Sq [Wardinski and Mandea 2006]. By varying the maximum degree for a given model of the internal field, the changes on the predictions at given place are seen. We noticed that increasing the maximum degree of the model does not go along with the best fitting of predicted values to the observed val- ues of any geomagnetic field component at any place, i.e. the un-modeled internal or external contributions are not always decreased by increasing the maximum degree. However, we noticed that the shape of the plot- ted time series is not affected by the change of maxi- mum degree as it is always shifted a little bit up or down according to the considered component, model and maximum degree. By removing this shift, we can easily detect the external contribution in a given place that is not modeled by the models of internal field with dif- ferent maximum degree. It is easy to notice that the time dependence is practically the same for all series. We have not studied the correlation of these un- modeled external signals with DST-index or for other indices, such as F10.7, suggested by Wardinski and Holme [2011], that would help in specifying their ex- ternal sources. Removal of these signals enhances the resolution of fine-scale detail in secular variation; this is useful in considering the phenomenology of geomag- netic jerks [Wardinski and Holme 2011]. Acknowledgements. The results presented in this paper rely on the data collected at magnetic observatories: Koblenz (WIK) - Austria, Niemegk (NGK) and Furstenfeldbruck - Germany, Kakioka (KAK) - Japan and Hermanus (HER) - South Africa. We thank: Zen- tralanstalt für Meteorologie und Geodynamik (Vienna) for supply- ing data to us, the national institutes that supported NGK, KAK and HER observatories and INTERMAGNET for promoting high stan- dards of magnetic observatory practice (www.intermagnet.org). References Backus, G., R.L. Parker and C. Constable (1996). Foun- dations of Geomagnetism, Cambridge University Press, Cambridge, UK. Bloxham, J., and D. Gubbins (1986). Geomagnetic field analysis-IV. Testing the frozen-flux hypothesis, Geo- phys. J. R. Astr. Soc., 84, 139-152; doi:10.1111/j.1365- 246X.1986.tb04349.x. Cain, J., D.R. Schmitz and L. Muth (1984). Small-Scale Features in the Earth’s Magnetic Field Observed by Magsat, J. Geophys. Res., 89 (b2), 1070-1076. Hulot, G., N. Olsen, E. Thebault and K. Hemant (2009). NON-MODELED EXTERNAL CONTRIBUTION Figure 18. Power spectrum density (PSD in nT2/frequency units) calculated by PMTM method for external contributions series of POMME- 9 model at different observatories: NGK (solid line), FUR (dashed line), HER (dotted line), KAK (dash-dotted line). Crustal concealing of small-scale core-field secular variation, Geophys. J. Int., 177, 361-366; doi:10.1111/ j.1365-246X.2009.04119.x. Langel, R.A., and R.H. Estes (1982). A geomagnetic field spectrum, Geophys. Res. Lett., 9 (4), 250-253; doi:10.1029/GL009i004p00250. Lesur, V., and D. Gubbins (2000). Using geomagnetic secular variation to separate remanent and induced sources of the crustal magnetic field, Geophys. J. Int., 142, 889-897. Lesur, V., K. Whaler and I. Wardinski (2015). Are geo- magnetic data consistent with stably stratified flow at the core-mantle boundary?, Geophys. J. Int., 201 (2), 929-946; http://doi.org/10.1093/gji/ggv031. Lowes, F.J. (1974). Spatial power spectrum of the main geomagnetic field, and extrapolation to the core, Geophys. J. R. Astr. Soc., 36, 717-730. Lühr, H., and S. Maus (2010). Solar cycle dependence of quiet-time magnetospheric currents and a model of their near-Earth magnetic fields, Earth Planets Space, 62 (10), 843-848. Mandea, M., and B. Langlais (2002). Observatory crustal magnetic biases during MAGSAT and Ørsted satellite missions, Geophys. Res. Lett., 29 (15), ORS 4-1 - ORS 4-4; doi:10.1029/2001GL013693. Mandea, M., R. Holme, A. Pais, K. Pinheiro, A. Jackson and G. Verbanac (2010). Geomagnetic Jerks: Rapid Core Field Variations and Core Dynamics, Space Sci. Rev., 155, 147-175. Maus, S., H. Lühr, G. Balasis, M. Rother and M. Man- dea (2004). Introducing POMME, Potsdam Mag- netic Model of the Earth; http://www.geomag.us/ info/Smaus/Doc/CHAMP-2-Pomme.pdf. Maus, S., M. Rother, C. Stolle, W. Mai, S.-C. Choi, H. Lühr, D. Cooke and C. Roth (2006). Third genera- tion of the Potsdam Magnetic Model of the Earth (POMME), Geochem. Geophys. Geosyst., 7, Q07008; doi:10.1029GC001269. Maus, S. (2010). An ellipsoidal harmonic representation of Earth’s lithospheric magnetic field to degree and order 720, Geochem. Geophys. Geosyst., 11, Q06015; doi:10.1029/2010GC003026. Maus, S., C. Manoj, J. Rauberg, I. Michaelis and H. Lühr (2010). NOAA/NGDC candidate models for the 11th generation International Geomagnetic Reference Field and the concurrent release of the 6th genera- tion POMME magnetic model, Earth Planets Space, 62 (10), 729-735; doi:10.5047/eps.2010.07.006. Olsen, N., L. Hermann, C.C. Finlay, T.J. Sabaka, I. Michaelis, J. Rauberg and L. Tøffner-Clausen (2014). The CHAOS-4 geomagnetic field model, Geophys. J. Int., 197, 815-827. Park, J., C.R. Lindberg and F.L. Vernon (1987). Multi- taper spectral analysis of high-frequency seismo- grams, J. Geophys. Res., 92, 12675-12684. Percival, D.B., and A.T. Walden (1993). Spectral Analy- sis for Physical Applications. Multi-Taper and Con- ventional Univariate Techniques, Cambridge University Press, Cambridge, UK. Sabaka, T.J., N. Olsen and M. Purucker (2004). Extend- ing comprehensive models of the Earth’s magnetic field with Ørsted and CHAMP data, Geophys. J. Int., 159 (2), 521-547. Sabaka, T.J., N. Olsen, H.T. Robert and A. Kuvshinov (2015). CM5, a pre-Swarm comprehensive geomag- netic field model derived from over 12 yr of CHAMP, Ørsted, SAC-C and observatory data, Geophys. J. Int., 200 (3), 1596-1626. Thomson, J.D. (1982). Spectrum estimation and har- monic analysis, IEEE Proc., 70 (9), 1055-1096. Verbanac, G., M. Korte and M. Mandea (2007a). On long-term trends in European geomagnetic obser- vatory biases, Earth Planets Space, 59, 685-695. Verbanac, G., H. Lühr, M. Rother, M. Korte and M. Mandea (2007b). Contributions of the external field to the observatory annual means and a proposal for their corrections, Earth Planets Space, 59, 251-257. Verbanac, G., M. Mandea, M. Bandi and S. Subaši (2015). Magnetic observatories: biases over CHAMP satellite mission, Solid Earth, 6, 775-781; doi:10.519 4/se-6-775-2015. Wardinski, I., and M. Mandea (2006). Annual and semi- annual variations of the geomagnetic field compo- nents analysed by the multi-taper method, Earth Planets Space, 58, 785-791. Wardinski, I., and R. Holme (2011). Signal from noise in geomagnetic field modelling: denoising data for sec- ular variation studies, Geophys. J. Int., 185, 653-662; doi:10.1111/j.1365-246X.2011.04988.x. * Corresponding author: Bejo Duka, University of Tirana, Faculty of Natural Sciences, Department of Physics, Tirana, Albania; email: bejo.duka@fshn.edu.al. © 2016 by the Istituto Nazionale di Geofisica e Vulcanologia. All rights reserved. DUKA ET AL. 14 << /ASCII85EncodePages false /AllowTransparency false /AutoPositionEPSFiles false /AutoRotatePages /None /Binding /Left /CalGrayProfile (Dot Gain 20%) /CalRGBProfile (sRGB IEC61966-2.1) /CalCMYKProfile (U.S. Web Coated \050SWOP\051 v2) /sRGBProfile (sRGB IEC61966-2.1) /CannotEmbedFontPolicy /Warning /CompatibilityLevel 1.3 /CompressObjects /Tags /CompressPages true /ConvertImagesToIndexed true /PassThroughJPEGImages true /CreateJobTicket false /DefaultRenderingIntent /Default /DetectBlends true /DetectCurves 0.1000 /ColorConversionStrategy /LeaveColorUnchanged /DoThumbnails false /EmbedAllFonts true /EmbedOpenType false /ParseICCProfilesInComments true /EmbedJobOptions true /DSCReportingLevel 0 /EmitDSCWarnings false /EndPage -1 /ImageMemory 1048576 /LockDistillerParams true /MaxSubsetPct 100 /Optimize false /OPM 1 /ParseDSCComments true /ParseDSCCommentsForDocInfo true /PreserveCopyPage true /PreserveDICMYKValues true /PreserveEPSInfo true /PreserveFlatness true /PreserveHalftoneInfo false /PreserveOPIComments false /PreserveOverprintSettings true /StartPage 1 /SubsetFonts true /TransferFunctionInfo /Apply /UCRandBGInfo /Preserve /UsePrologue false /ColorSettingsFile (None) /AlwaysEmbed [ true /AndaleMono /Apple-Chancery /Arial-Black /Arial-BoldItalicMT /Arial-BoldMT /Arial-ItalicMT /ArialMT /CapitalsRegular /Charcoal /Chicago /ComicSansMS /ComicSansMS-Bold /Courier /Courier-Bold /CourierNewPS-BoldItalicMT /CourierNewPS-BoldMT /CourierNewPS-ItalicMT /CourierNewPSMT /GadgetRegular /Geneva /Georgia /Georgia-Bold /Georgia-BoldItalic /Georgia-Italic /Helvetica /Helvetica-Bold /HelveticaInserat-Roman /HoeflerText-Black /HoeflerText-BlackItalic /HoeflerText-Italic /HoeflerText-Ornaments /HoeflerText-Regular /Impact /Monaco /NewYork /Palatino-Bold /Palatino-BoldItalic /Palatino-Italic /Palatino-Roman /SandRegular /Skia-Regular /Symbol /TechnoRegular /TextileRegular /Times-Bold /Times-BoldItalic /Times-Italic /Times-Roman /TimesNewRomanPS-BoldItalicMT /TimesNewRomanPS-BoldMT /TimesNewRomanPS-ItalicMT /TimesNewRomanPSMT /Trebuchet-BoldItalic /TrebuchetMS /TrebuchetMS-Bold /TrebuchetMS-Italic /Verdana /Verdana-Bold /Verdana-BoldItalic /Verdana-Italic /Webdings ] /NeverEmbed [ true ] /AntiAliasColorImages false /CropColorImages true /ColorImageMinResolution 150 /ColorImageMinResolutionPolicy /OK /DownsampleColorImages true /ColorImageDownsampleType /Bicubic /ColorImageResolution 300 /ColorImageDepth -1 /ColorImageMinDownsampleDepth 1 /ColorImageDownsampleThreshold 1.10000 /EncodeColorImages true /ColorImageFilter /DCTEncode /AutoFilterColorImages true /ColorImageAutoFilterStrategy /JPEG /ColorACSImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /ColorImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /JPEG2000ColorACSImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /JPEG2000ColorImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /AntiAliasGrayImages false /CropGrayImages true /GrayImageMinResolution 150 /GrayImageMinResolutionPolicy /OK /DownsampleGrayImages true /GrayImageDownsampleType /Bicubic /GrayImageResolution 300 /GrayImageDepth -1 /GrayImageMinDownsampleDepth 2 /GrayImageDownsampleThreshold 1.10000 /EncodeGrayImages true /GrayImageFilter /DCTEncode /AutoFilterGrayImages true /GrayImageAutoFilterStrategy /JPEG /GrayACSImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /GrayImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /JPEG2000GrayACSImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /JPEG2000GrayImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /AntiAliasMonoImages false /CropMonoImages true /MonoImageMinResolution 1200 /MonoImageMinResolutionPolicy /OK /DownsampleMonoImages true /MonoImageDownsampleType /Bicubic /MonoImageResolution 1200 /MonoImageDepth -1 /MonoImageDownsampleThreshold 1.08250 /EncodeMonoImages true /MonoImageFilter /CCITTFaxEncode /MonoImageDict << /K -1 >> /AllowPSXObjects false /CheckCompliance [ /None ] /PDFX1aCheck false /PDFX3Check false /PDFXCompliantPDFOnly false /PDFXNoTrimBoxError true /PDFXTrimBoxToMediaBoxOffset [ 0.00000 0.00000 0.00000 0.00000 ] /PDFXSetBleedBoxToMediaBox true /PDFXBleedBoxToTrimBoxOffset [ 0.00000 0.00000 0.00000 0.00000 ] /PDFXOutputIntentProfile (None) /PDFXOutputConditionIdentifier () /PDFXOutputCondition () /PDFXRegistryName (http://www.color.org) /PDFXTrapped /Unknown /CreateJDFFile false /SyntheticBoldness 1.000000 /Description << /ENU (Use these settings to create PDF documents with higher image resolution for high quality pre-press printing. The PDF documents can be opened with Acrobat and Reader 5.0 and later. These settings require font embedding.) /JPN /FRA /DEU /PTB /DAN /NLD /ESP /SUO /NOR /SVE /KOR /CHS /CHT /ITA >> >> setdistillerparams << /HWResolution [2400 2400] /PageSize [595.000 842.000] >> setpagedevice