1. Introduction

There is a clear need to measure and model the values of the main geomagnetic field and its secular variation at global and regional scales. In this paper we present two geomagnetic reference field models for the region comprising Albania, south-east Italy and surrounding sea areas from 1990 to 2010 with a prediction to 2012 using data from magnetic repeat station surveys and satellites. The repeat stations are the points of a magnetic network where the three-component magnetic field is periodically measured at intervals of some years (e.g. Newitt et al., 1996). The region under investigation is however not uniformly covered by ground magnetic measurements. The complete Italian territory on the other hand is comprehensively covered by magnetic repeat station measurements whereas the Albanian territory has always been poorly surveyed. Thanks to the short distance between these two countries, during the last 20 years there have been continuous collaboration between the “Istituto Nazionale di Geofisica e Vulcanologia” (INGV; “Istituto Nazionale di Geofisica”, ING till 2001), which is the Italian Institution responsible of the magnetic monitoring in Italy, and the Academy of Science of Albania, Centre of Geochemistry and Geophysics (CGG) and Physics Department of Tirana University (PDTU), i.e. the Albanian Institutions that usually perform magnetic measurements in Albania. This collaboration has leaded to coordinate specific campaigns of geomagnetic surveys over a designed repeat station network in Albania. As mentioned in a previous paper (Duka et al., 2004), INGV provided some personnel and instruments during the performance of measurement campaigns and afterwards with the competence in the modelling procedure.

In order to develop our regional reference field model, in this paper we make extensive use of the datasets of magnetic measurements in Albania from 1990 till present, together with recent results from the Italian repeat station network and from the ∅rsted and CHAMP satellite missions. First we compute a regional model by means of SCHA (Haines, 1985) with polynomial time dependency. Then, we present a new technique (RS-SCHA) which is a simplification of the former in the radial variation function. Both the techniques provide an optimal representation of the geomagnetic field over the area of investigation, which is an improvement on the International Geomagnetic Reference Field (IGRF-11). Since most of the recent analysed data were not used to form IGRF-11, this work represents also an independent test of the global model.

2. Data

The data used in this study have been collected during different Italian (vector and scalar measurements from 1990.5 to 2010.0) and Albanian (vector and scalar measurements from 1990.0 to 2009.7) magnetic repeat station surveys (Duka et al., 2004). In particular, 12 stations of the Albanian network have been taken into account in addition to 6 stations located in southern Italy. Because of the poor coverage of this area with geomagnetic data, especially at the borders of this cap, we calculated from IGRF-11 field model the corresponding synthetic Cartesian X, Y, Z components at sea level in 19 fixed points of this area, at years 1990, 1995, 2000 and 2005 (circles in Fig. 1(a)). Locations of all ground geomagnetic measurements and selected ∅rsted and CHAMP satellite data are shown in Fig. 1(a), 1(b) and 1(c), respectively.

Fig. 1.
figure 1

(a) Location of Italian and Albanian ground magnetic data with 18 vector and scalar magnetic field data points (stars) and 28 scalar magnetic field data points (squares) from 1990 to 2010; 19 synthetic vector data points from IGRF-11 (circles) at 5-year intervals from 1990 to 2005: (b) and (c) are the location of ∅rsted (circles) and CHAMP (circles) satellite data used for the SCHA, respectively, from 1999.4 to 2009.4.

The older dataset of magnetic measurements used here is the total intensity F measured by means of a proton magnetometer around 1990.0 by the former Geophysical Enterprise of Tirana and covers all Albanian territory with an array of 28 sites (Duka and Bushati, 1991). In September 1994 (epoch 1994.75), in the framework of a joint project between the CGG, PDTU and the former ING a new vector magnetic survey measuring F, inclination I and declination D (Chiappini et al., 1997, 1999) covered the Albanian territory. The measurements were taken by using a Geometrics proton magnetometer and a Bartington fluxgate theodolite, together with a gyro-theodolite for the absolute determination of the geographical azimuth. In August 2003 (epoch 2003.6; Duka et al., 2004) another scalar field campaign was conducted for the total intensity F only using an Overhauser effect magnetometer, repeating the scalar measurements at 10 of the 28 sites of the previous F survey. At each site data were recorded for 1–2 hours in order to obtain a satisfactory accumulation of total intensity data.

Another three component survey in the Albanian territory was made in September 2004 when a total of 12 (with the addition of Berat station) vector measurements were undertaken. In order to remove the effect of the diurnal variation, for the full duration of the campaign, a temporary station with a fluxgate variometer (30 sec. sampling) was installed not far from Tirana for the full duration of the campaign to reduce all values to the closest night time.

During September 2009 (epoch 2009.7) INGV carried out the last campaign of absolute measurements at all the 2004 points, measuring F, I and D magnetic elements. These data were reduced using data from L’Aquila Observatory because of some malfunctions in an installed temporary station in Tirana (Albania).

To derive a more representative and accurate model, we analyse together with the measurements in Albania, also those in Italian territory, in particular we considered 6 locations from the Italian magnetic repeat stations network placed in South-East Italy (Dominici et al., 2007). For a better temporal behaviour of the model and to improve the stability of the inversion, we have synthesized X, Y and Z components at the limits of the considered temporal interval, i.e. 1990.0 and 2010.0, for the 12 Albanian stations reducing the real vector measurements available for 1994.75 and 2009.7 to the closest extreme epoch. The corresponding temporal reduction was made applying the secular variation predicted by an updated version of the Italian geomagnetic reference field model (ITGRF; De Santis et al., 2003), a model that has been demonstrated to predict the temporal change of the magnetic field in this area better than global models, such as the IGRF-11 field model. With the aim of temporal stability, we added also L’Aquila observatory annual means from 1990.5 to 2010.0.

To overcome the non-uniform distribution of the data in this region (especially for the sea area), and since our models will take into account the proper altitude variation of the field, together with the ground data we also included a set of total intensity magnetic field measurements from ∅rsted and CHAMP satellite missions, selected between 1999.4 and 2009.4 during low external magnetic activity characterized by K p ≤ 3 and |Dst| ≤ 10 nT magnetic indices.

The different datasets were weighted according to the reciprocal of the variance of the total error associated with the measurements (Duka et al., 2006). Each variance was computed as the sum of the error σm intrinsically involved in the measurement process (instrumental errors, etc.) and that related to the lithospheric field, σ1. We assumed that the ground data have σ1 = 50 nT and σm = 10 nT, whereas no lithospheric contribution was assigned to satellite data, and their σm was set equal to that of ground data.

In the next section we will apply our regional techniques to the data set obtained from the original observations after removing the main field predicted by IGRF-10 (Macmillan and Maus, 2005) and IGRF-11 (http://www.ngdc.noaa.gov/IAGA/vmod/igrf.html) up to degree 13 at the central epoch (2000.0).

3. Spherical Cap Harmonic Analysis Model and a Radially Simplified Version

Modelling the geomagnetic field over the whole terrestrial sphere is usually approached using spherical harmonic analysis (SHA). When we are interested only in details of a certain area with typical lengths shorter than a given size (say, one thousand km, like the area under investigation), we must resort to a local analysis. The first regional model we propose is based on the Spherical Cap Harmonic Analysis (SCHA; Haines, 1985). The SCHA is a powerful analytical tool for modelling Laplacian potential and the corresponding field components over a spherical cap, overcoming most of the problems (e.g., non-orthogonality and huge number of model coefficients) that arise when global spherical harmonic models are applied to restricted areas. Introduced for the first time in 1985 by Haines, SCHA has evolved progressively during the recent years (e.g. De Santis, 1991, 1992) finally arriving to the Revised SCHA (R-SCHA; Thébault et al., 2004, 2006).

The solution of Laplace’s equation in spherical coordinates (r,θ,ϕ) for the magnetic potential due to internal sources (thus the subscript “int” for K) over a spherical cap can be written as an expansion of non integer-degree spherical harmonics:

((1))

where the polynomial time dependency is included: this allows to easily extrapolate the model forward in time (we propose here a prediction up to 2012); and are the spherical cap harmonic coefficients that determine the model; are the associated Legendre functions that satisfy appropriate boundary conditions (null potential or co-latitudinal derivative at the border of the cap) and have integer order m and generally non-integer degree n k ( m ); k is an integer index selected to arrange, in increasing order, the different roots n for a given m in the boundary conditions. The number of coefficients depends on the maximum spatial and temporal indices of the expansion, Kint and Q, respectively.

The geomagnetic components X, Y and Z are obtained as appropriate spatial derivatives of Eq. (1) in spherical coordinates, since the potential is non-observable. To overcome the non linearity problem that arises when combining vector measurements with total field measurements, we used a first order Taylor expansion of the total magnetic field intensity, as a square root function of the X, Y and Z components (Haines and Newitt, 1997).

After many tests, the model parameters that best represent the input data in the period between 1990 and 2010 were found to be Kint = 2 and Q = 3. The and coefficients were obtained through a least squares regression procedure.

The model so defined was chosen with basis functions defined over a cap with a semi-angle of 16°, in order to represent the main field and its secular variation, including the most significant harmonics of the regional geomagnetic field (the minimum and maximum degrees are approximately equal to 6.1 and 13.2, respectively). In order to reduce the typical problems in SV modelling over small regions (e.g., Torta et al., 2006), we chose to use a larger cap than that actually containing the data. Of course this choice is made at expenses of losing the orthogonality of the basis functions over the data interval, since the real data are mostly concentrated within the central 4° half-angle. However we still have the advantage of reducing significantly the number of model coefficients and improving the quality of the model. The coefficients of the Albanian-Italian Geomagnetic Reference Model developed by SCHA for this region are shown in the Table 1. The final model values are obtained with the addition of IGRF-11 values at 2000.0.

Table 1. Model coefficients and in nT/yearq of the Albanian-Italian Geomagnetic Reference Model developed by SCHA. Final model values are obtained adding the IGRF-11 field values at 2000.0.

After this SCHA model, here we present also another regional model that we consider a good compromise between the SHA and SCHA techniques. It takes advantage of both the original concepts of the SCHA and those typical of the most recent R-SCHA. We call this new simple method “Radially Simplified Spherical Cap Harmonic Analysis” (RS-SCHA). It consists in simplifying the expression in radial distance r of the geomagnetic potential V in Eq. (1) “emulating” the radial behaviour of the conventional spherical harmonics analysis (SHA) as follows:

((2))

Please note that although this form of geomagnetic field potential is not usual, nevertheless it satisfies Laplace’s equation separately for each spherical coordinate: indeed, for the radial variation it is the typical solution in the global case, while for the angular coordinates it is the typical SCHA solution. Given the differences in the radial function only with respect with the typical SCHA, we can affirm that our model satisfies the Laplace’s equation at Earth’s surface (r = a). When we consider satellite altitude we actually multiply all surface spherical cap functions by a given constant that it is taken into account in the final estimated model coefficients. Also here, the model is still solution of Laplace’s equation at the ground or satellite altitudes. We admit that there is a mixture of eigen-values m, n which is not correct. However this is analogous with what is normally done when removing a SHA model (typically IGRF) from observations and then applying SCHA inversion (and vice-versa when we synthesise the model values): the final potential is the sum of two different potentials, a global and a regional one, with two different n, m sets; the former is characterised by integer n values (SHA) while the latter is expressed with non-integer n k values (SCHA or even R-SCHA). Our proposal is such to say that SHA radial functions works better than SCHA radial functions, at least for the area and the period of interest.

To check the validity of this technique, we applied the RS-SCHA to the region under investigation (Albania, Southern part of Italy and surrounding seas) using the same parameters as above (SCHA model). The final RS-SCHA coefficients of the Albanian-Italian Geomagnetic Reference Model are shown in Table 2. As in the previous case, the final model values are obtained with the addition of IGRF-11 values at 2000.0.

Table 2. Model coefficients and in nT/yearq of the Albanian Geomagnetic Reference Model developed by RS-SCHA. Final model values are obtained adding the IGRF-11 field values at 2000.0.

4. Results

Regional models usually represent the spatial and temporal variations in a particular region in more detail and accuracy compared to global geomagnetic models (e.g. IGRF-10 and IGRF-11). This is confirmed also by our analysis. Table 3 shows the root mean square (RMS) fits of SCHA, RS-SCHA, IGRF-10 and IGRF-11 models to the analysed ground and satellite data for the magnetic field (in nT) and its secular variation (in nT/yr). From a statistical analysis in terms of RMS we see that both SCHA and RS-SCHA techniques provide better results for the field in comparison with both IGRFs. In particular the RS-SCHA model improves the fit of all ground and satellite components better then IGRF-s and better than SCHA. Moreover, we see that IGRF most recent version, IGRF-11, has better performances than the previous version, especially for F satellite data.

Table 3. Root mean square fits of SCHA, RS-SCHA, IGRF-10 and IGRF-11 models to the analysed ground data and satellite data for magnetic field (nT) (for the first four lines) and the secular variation (nT/yr) (for the second four lines). Underlined bold values are the best RMS values among the models.

Figure 2 shows the regional X, Y, Z charts for the analysed area for the epoch 2012.0 at sea level obtained from SCHA, RS-SCHA and IGRF-11. The evaluation in time from 1995 to 2010 of the geomagnetic field, especially of the Y component (Fig. 3), shows a clear fast westward drift of the geomagnetic field (e.g. Barraclough and Malin, 1999) of 0.4°–0.5°/year.

Fig. 2.
figure 2

Maps (in nT) for X (top left), Y (top right), Z (bottom left), and F (bottom right) elements for epoch 2012.0 at sea level obtained from SCHA (a), RS-SCHA (b) and IGRF-11 (c) models developed on an 16° half-angle cap.

Fig. 3.
figure 3

Maps (in nT) for Y component for epochs 1995.0 (top left), 2000.0 (top right), 2005.0 (bottom left) and 2010.0 (bottom right) at sea level obtained from SCHA model. It is clear a fast westward drift of 0.4°–0.5°/year.

To estimate the SV at some epoch we used the differences between two close epochs, one year before and one year after; however the numerical results would have been the same by performing numerical derivation of the polynomial part of the regional models. All SV models (SCHA, RS-SCHA and IGRF-11) show practically the same RMS for all the components (Table 3) and confirm that the region under study presents low SV values for this period (see, Figs. 4 and 5), with less than 50 nT/year for all components (e.g. Gubbins, 1990). We notice that the larger differences in both the field and SV between IGRF-11 and our regional models are in Albania. To better visualise this fact, in Fig. 6 we show the differences between observed data, SCHA, RS-SCHA and IGRF-11 values in X, Y and Z components for the sites of Tirana (Albania) and Masseria Maserino (Italy). As it can see, the IGRF-11 does not fit the observations in the Albanian station as well as in the Italian station. This could be probably ascribed to some crustal bias present in the Albanian station but absent in Italian stations, details that only a regional model is able to grasp in comparison with those provided by a global model. This is confirmed by both SCHA and RS-SCHA models, especially in Y- and Z-components when compared with IGRF-11.

Fig. 4.
figure 4

Maps (in nT/year) for the secular variation of X (top left), Y (top right), Z (bottom left), and F (bottom right) magnetic elements centred at 2005 (deduced from differences from 2004 to 2006) at sea level obtained from the SCHA (a), RS-SCHA (b) and IGRF-11 (c) models.

Fig. 5.
figure 5

Maps (in nT/year) for the secular variation of X (top left), Y (top right), Z (bottom left), and F (bottom right) magnetic elements centred at 2010 (deduced from differences from 2009 to 2011) at sea level obtained from the SCHA (a), RS-SCHA (b) and IGRF-11 (c) models.

Fig. 6.
figure 6

Comparison between observed (stars), SCHA, RS-SCHA and IGRF-11 X, Y and Z components and total intensity F for sites of Tirana (Albania) and Masseria Maserino (Italy).

The SV differences between observed data, SCHA, RS-SCHA and IGRF-11 values in X, Y and Z and F for the sites of Tirana (Albania) and Masseria Maserino (Italy) are shown in Fig. 7. We notice that both regional models provide a better polynomial interpolation among the observed SV components in the two repeat stations than the abrupt IGRF-11 SV changes at each 5-year interval.

Fig. 7.
figure 7

Comparison between observed (stars), SCHA, RS-SCHA and IGRF-11 values of the secular variation of the X, Y and Z components and total intensity F for sites of Tirana (Albania) and Masseria Maserino (Italy).

5. Conclusions

In this paper, we present two regional models for SouthEast Italy, Albania and surrounding sea area, both based on a spherical cap harmonic expansion of the potential but with different radial functions. The introduction of a radially simplified version of SCHA, i.e. RS-SCHA, characterised by a SHA-like radial variation of the field, allows us to model both ground and satellite data better than SCHA. The validity in time of both regional models is in the period 1990–2010 with prediction to 2012. One of our objectives was also to show how IGRF-11 works in this area of the world. The comparison between RMS fits of the regional models to real data and those of IGRF allows us to confirm that both SCHA and RS-SCHA models represent an improvement with respect to the global model when representing the field and its secular variation, probably because some of the most recent data were not used to construct IGRF-11 model. Thus in conclusion, the regional models can be used to estimate the values of the geomagnetic field (and its secular variation) all over the region considered (South-East Italy and Albania, seas included). They can be utilized as well for reducing magnetic survey data taken in the area of interest in the period of validity of the model.