ITalian Geomagnetic Reference Field (ITGRF): update for 2000 and secular variation model up to 2005 by autoregressive forecasting

The updated version of the ITalian Geomagnetic Reference Field (ITGRF) for 2000.0 and its secular variation model up to 2005.0 are presented in this paper. The main field model is based on a simple polynomial approximation in latitude and longitude of the geomagnetic field elements computed from IGRF on a 12° ¥ 11° grid centred over Italy. The annual means from L’Aquila observatory were used to determine the baseline level, imposing a constant observatory anomaly bias. This procedure gives a set of 6 coefficients every 5 years from 1960 to 2005 for the horizontal H, total field F, vertical Z and declination D elements of the geomagnetic field. The extrapolation of ITGRF to 2005 is based on an autoregressive forecasting of the L’Aquila observatory annual means. Comparison of the field values computed from the model with those recorded at the other Italian observatory (Castello Tesino) shows that the ITGRF improves the fit of the secular variation pattern with respect to the global IGRF model by a factor of 3. The ITGRF represents a reliable alternative to global models when reducing magnetic surveys to a common reference epoch over the Italian region. Mailing address: Dr. Angelo De Santis, Istituto Nazionale di Geofisica e Vulcanologia, Via di Vigna Murata 605, 00143 Roma; e-mail: desantisag@ingv.it


Introduction
It is well known that regional models are more accurate than global models for describing the main features of the geomagnetic field over a limited region of the Earth.Global models, like the International Geomagnetic Reference Field (IGRF), are developed to represent magnetic features with wavelengths larger than regional dimensions, and they are not able to represent high spatial frequency magnetic features.Regional models, on the other hand, are obtained from data measured within the region of interest and they are constrained to represent the behaviour of the field only over this area without considering what is happening outside.They are therefore able to model shorter wavelengths.
Many countries have developed their own regional magnetic reference model, and several studies have been carried out to represent the behaviour of the geomagnetic field over limited regions, with areas covering from a few squared geographical degrees to continental scales.Different methods are used to generate the sets of coefficients that constitute such reference fields (Haines, 1990), ranging from simple procedures such as second degree polynomial fitting to more sophisticated techniques like the relatively new method of Spherical Cap Harmonic Analysis (SCHA; Haines, 1985) and its derivatives (De Santis, 1991a, 1992).
In the case of Italy, both simple and complex methods have been applied during recent years to create a regional geomagnetic reference field.Molina and De Santis (1987) proposed a reference model based on a polynomial fitting.This model, valid in its first version for the period 1965-1984, has been called the ITalian Geomagnetic Reference Field (ITGRF).In its last revision (De Santis, 1991b), the model was valid up to 1990 with a secular variation model until 1995.The polynomial approximation was chosen because of its simplicity, and a second-degree polynomial was considered adequate to represent the field, following the empirical Bullard's rule (Bullard, 1967;see also De Santis and Meloni, 1991), and taking into account that the mean errors of the approximation were very low with respect to the corresponding IGRF model.
Another model is that created from the measurements taken at the different points of the Italian repeat station network.In Coticchia et al. (2001), a description of the field measurements and procedures, maps and data for the latest 2000.0Italian survey, can be found.By using the SCHA technique, De Santis et al. (1990) obtained a set of 100 coefficients able to represent the main features of the field with wavelengths down to 330 km.The main disadvantage of this method is the need for a computer with the appropriate software to obtain the values of the field at a given location and epoch, and this is sometimes is not available when carrying out magnetic surveys, where a simple reference model is more useful.
In this paper, we present an up-to-date version of the ITGRF which is valid for epochs from 1960 to 2000 with a secular variation model valid until 2005.The model, characterised by a set of six coefficients for each element of the field for each five-year interval was produced by taking as the base level the L'Aquila observatory annual means and assuming a spatial gradient equal to that given by the IGRF.The ITGRF's validity was checked by comparing its generated field values with real data not used during its development, such as the annual means of Castello Tesino observatory and the values of the Italian network of repeat stations.

Data used and model generation
To a latitude-longitude grid of the IGRF values for the horizontal H and vertical Z components, declination D and total intensity F of the magnetic field for epochs 1960-2005, we fit a second-degree polynomial of the form where E is the value of the considered element at latitude f and longitude l, and a 0 ' , a 1 ,..., a 5 are the coefficients to be determined by the least squares technique.Latitude and longitude, expressed in degrees, are referred to a central point (12°E, 42°N) of the area concerned.The coefficients a 1 , a 2 (in nT/degree) and a 3 , a 4 , a 5 (in nT/degree 2 ) determine the spatial gradient for each element, and they are obtained assuming that the variation in the field values is that predicted by the IGRF; a 0 ' (nT) represents the base level of the model, and it will be modified to a 0 (see below) using the annual mean values for each element at L'Aquila observatory (AQU; 13.32°E, 42.38°N), in order to have the same difference between the model and AQU.
The spatial gradients of the ITGRF model are taken to be the same as those given by the IGRF because the latter is a global model based on a global data set and is therefore likely to represent these gradients well, avoiding boundary effects.This is not the case with regard to the absolute values of the field; so, for instance, a comparison of the values predicted by the IGRF with those measured at L'Aquila shows some clear differences (De Santis and Molina, 1988).
To correct the values at L'Aquila for the typical local anomaly of the observatory, we define the coefficient a 0 = a 0 ' + e in order to satisfy, for each year, the following equality: (2.3) Coefficients k 0 , k 1 ,..., k m have been estimated by least squares.Regressions of different order m were made, and m = 4 was finally chosen because the corresponding autoregressive model showed a better Root Mean Square (RMS) deviation with respect to the actual data.In the case of epoch 1960.0,we employed the so-called retrodiction, in the sense that we inverted time as we were observing the phenomenon in the opposite direction (i.e. from 2000.5 to 1960.5) and applied (2.3) again.In the same way, a synthetic grid of values from 1960.0 to 2005.0 at five-year intervals constructed from the eighth generation IGRF (IAGA, 2000) for each element has also been used as input data.This grid extends from 6.5° to 18.5° for east longitudes and between 36.5° and 47.5° for north latitudes at 1° intervals, covering the whole of Italy and adjacent areas.
The model coefficients are referred to sea level (altitude h = 0).In this sense, the mean values of each element at L'Aquila observatory (situated at h = 0.682 km) have been reduced to sea level by applying a correction that takes into account only the dipolar term of the geomagnetic potential (i.e.contributions with degree n equal to 1; e.g., Langel, 1987).For any intensive field element with value E at a height h, the correction DE for computing its value at h = 0 is expressed as where R = 6371.2km is the Earth's mean radius.
In conclusion, the ITGRF model provides the values of each E element as follows: (2.5) where k = ij, and the coefficients are defined as c 00 = a 0 , c 10 =a 1 , c 01 = a 2 , c 20 = a 3 , c 02 = a 4 , and c 11 = a 5 .

Model analyses and results
Values of the ITGRF coefficients a 1 ,..., a 5 computed using the least squares method at five-year intervals with the associated misfit error with respect to the corresponding IGRF are shown in table I.The a 0 coefficients (table II) are computed for every year from 1960 to 2005 applying eq.(2.1) and taking into account the local anomaly introduced in eq.(2.2).The value of the field for an intermediate epoch is computed by interpolation between sets of neighbouring epochs.
Since the ITGRF model was developed using independent polynomials for each element, the self-consistency of the model must be checked.For instance, we could check if the total intensity  0 ; this is called the geometrical condition.Another method is to impose the absence of vertical electric currents in the atmosphere (Chapman, 1942), which in terms of the gradients of the horizontal components, has the form .
(3.1)This is called the physical condition.Figure 1a,b shows both conditions for the ITGRF at epoch 2000.As we can see, DF values over the Italian territory are small, being lower than 8 nT, and lower than 4 nT per degree for the condition (3.1) on the horizontal gradients; so we can conclude that the ITGRF is geometrically and physically self-consistent.
With the obtained coefficients, the values of the field elements can be computed for a given location in Italy and compared with those obtained using the IGRF. Figure 2 shows the behaviour of the total field intensity at L'Aquila from 1960 to 2005 (annual mean values until 2000, values after 2000 are obtained by the autoregressive method) and those predicted by the ITGRF and the IGRF.It appears that, during the complete period, the ITGRF follows the observed F values better than the IGRF.Quantifying this difference in terms of RMS values, we found RMS Aquila-ITGRF = 13.2 nT (this difference being the same for every epoch because of the way the model has been designed) whereas RMS Aquila-IGRF = 15.4 nT.
This better behaviour of the ITGRF was an expected result, because of the determination of the a 0 coefficient from the annual values measured at L'Aquila.To test the validity of the regional model, it is necessary to compare its values with data not used in the development of the model, for example with the annual means of Castello Tesino observatory (CTS; 11.65°E, 46.05°N).  2 but for Castello Tesino F values.Apparently, in this case, the global model behaves than the regional one when comparing the two models with the observed values.This is due to the fact that data from Castello Tesino are directly used in the development of the IGRF, whereas they are not for the ITGRF.Nevertheless, the regional model follows the temporal variation of the different elements better than the IGRF does.Figure 4a presents the differences between the F annual mean values at Castello Tesino and the synthetic values obtained from the ITGRF and the IGRF.While the difference is practically constant with respect to the ITGRF, it changes sharply for the IGRF because the latter model has constant secular variation over given fiveyear intervals.The ITGRF takes into account the changes in secular variation by the inclusion of the a 0 coefficient, allowing the regional model to follow accurately the secular variation as illustrated in fig.4b.
This fact is confirmed if we compare the variation of the annual means of the magnetic elements measured at CTS with respect to the respective mean values over the interval 1965 to 1995, and the same for the values obtained from the ITGRF and the IGRF models.For example, if we consider the H values (fig.5), we obtain an RMS of 3.0 nT for the ITGRF variations with regard to those of CTS, while an RMS of 6.4 nT is reached for the IGRF.This difference increases when considering the total field intensity, with an RMS of 3.7 nT for the ITGRF and 12.3 nT for the IGRF.
We can compare the values obtained from the ITGRF with those computed from the INGRF2000 model, a reference field created from the measurements taken at the different locations of the Italian repeat stations network (Coticchia et al., 2001).This model has the same form as a second-degree polynomial in latitude and longitude, valid at sea level for the whole Italian territory except marine areas.Differences between a 0.5° ¥ 0.5° grid of values generated from the INGRF2000 and grids with same locations created using the ITGRF and the IGRF are similar, and in the case of F for epoch 2000 they are less than 20 nT over more than 90% of the Italian territory (fig.6).
Once the self-consistency and validity of the ITGRF have been demonstrated, and it has been shown that the ITGRF secular variation follows the values measured at L'Aquila and Castello

Conclusions
In this paper we have presented the so-called ITGRF2000, an up-to-date version of the Italian Geomagnetic Reference Field model that allows the determination of the different elements of the geomagnetic field over Italy from 1960 to 2005.The model is a second-degree polynomial in longitude and latitude and it consists of sets of six coefficients for each magnetic element at five-year intervals.
The design of the model, based on the spatial pattern of the IGRF global model but with a base level computed from annual mean values measured at L'Aquila observatory, makes it an accurate model for representing the secular variation of the field element values in Italy.This fact, together with its simplicity, makes the ITGRF a useful model for reducing magnetic surveys carried out at different epochs over the Italian territory, including future surveys, up to 2005, to a common epoch.
The intention of the Geomagnetism Unit of the Istituto Nazionale di Geofisica e Vulcanologia, the body responsible for geomagnetic field monitoring in Italy, is to revise the model at five-year intervals and to allow access to the ITGRF2000 through the homepage of the Institute (http:// www.ingv.it)for the computation of field values at any location and for a given epoch.
the field value E ITGRF (AQU) generated by the ITGRF model at the location of L'Aquila observatory is equal to the observatory value E(AQU) corrected by an anomaly e.This anomaly is computed as the difference between the values of each element predicted by IGRF and those measured at L'Aquila for epoch 1979.0, when the global model was partly based on data measured by the Magsat satellite and is therefore expected to be one of the most accurate IGRFs(Peddie, 1982).Actually, Magsat was in orbit from November 1979 till June 1980, so the chosen epoch was only at the very beginning of the satellite data coverage.However, for the sake of continuity, followingMolina and De Santis (1987), we have retained this epoch as the reference for the anomaly bias determination.The resulting differences between the values computed from the IGRF model and L'Aquila observatory annual means for 1979.0 were 13.2 nT for F, 32.0 nT for H, -4.0 nT for Z, and -4.4 arcmin for D. The input data for producing the model are the annual mean values for F, H, Z, and D measured at L'Aquila observatory from 1960.5 to 2000.5, obtained from the observatory yearbooks.Annual mean values from 2001.0 to 2005.0 have been computed by extrapolation assuming an autoregressive secular variation to define the coefficient a 0 up to epoch 2005, i.e. we imposed the following condition relating a L'Aquila annual mean E t (AQU), taken at time t (in years), with the m previous annual means: Fig. 1a,b.a) Difference (nT) between total field intensity and the square root of the sum of squares of the horizontal and vertical components at epoch 2000; b) difference (in nT/degree) between the gradients of the horizontal components.The star marks the location of L'Aquila observatory.

Figure 3
Figure 3 is analogous to fig. 2 but for Castello Tesino F values.Apparently, in this case, the global model behaves than the regional one when comparing the two models with the

Fig. 2 .
Fig. 2. Annual mean values of F measured at L'Aquila (circles) together with the values predicted by the ITGRF (solid line) and the IGRF (dashed line).

Fig. 3 .
Fig. 3. Annual mean values of F measured at Castello Tesino (circles) together with the values predicted by the ITGRF (solid line) and the IGRF (dashed line).

FigFig. 5 .
Fig. 4a,b.a) Difference between annual mean values for the F component at Castello Tesino and those obtained from the ITGRF (solid line) and IGRF (dashed line) models, respectively; b) behaviour of the secular variation for the F component at Castello Tesino (circles), for the IGRF (dashed line) and the ITGRF (solid line) models.

Fig. 6 .
Fig. 6.Differences between the INGRF2000 and the ITGRF for F at epoch 2000.

Fig
Fig. 7a-d.Maps for: a) F; b) Z; c) H, and d) D for epoch 2005.Units are nanotesla for F, Z and H, and arcmin for D.