Etn@ref: a geodetic reference frame for Mt. Etna GPS networks

In volcanology, geodetic data are one of the most important instruments for the scientific community interested in modeling physical processes related to magma movements in the shallow crust. Since the end of the 1980s, GPS surveys and continuous GPS stations have greatly improved the possibility of measuring such movements with high time and space resolution. However, physical modeling requires that any external influence on the data that is not directly related to the investigated quantity must be filtered. One major tricky factor in determining a deformation field using GPS displacement vectors and velocities is the correct choice of a stable reference frame. In this study, we defined a local reference frame using more than a decade of GPS measurements, to refer the Mt. Etna ground deformation pattern to a rigid block. In particular, we used a weighted least-squares inversion to estimate the Euler pole for the rigid block by minimizing the adjustments to two horizontal components of GPS velocity at 13 «fiducial» sites located within a 350-km radius of Mt. Etna. The inversion inferred a Euler pole located at 38.450 ̊N and −107.702 ̊ E, and a rotation rate of 0.263 deg/Myr.


Introduction
Routine use of GPS for monitoring ground deformation started at Mt. Etna in 1988, when a network of 18 benchmarks was surveyed by both GPS and EDM (Electronic Distance Measurement) techniques [Briole et al. 1992].Since that first survey, the network has been improved and measured by GPS in survey mode almost yearly.In late 2000, the installation of a continuous GPS permanent network on Mt.Etna began.The network geometry was gradually upgraded in the following years to reach the current configuration of 30 permanent stations that provide a dense cover of most areas of the volcano edifice (Figure 1).The Mt. Etna volcano GPS network (Etn@net) allows monitoring of the volcanic deformation at different time scales.We can measure slow (days or months) or fast (minutes or hours) changes in site positions using different techniques in relation to data acquisition and processing.In the first case, we process GPS data collected at 30-s sampling rates on a daily basis [Patanè et al. 2005], while in the second we perform high-frequency (1 Hz) instantaneous GPS positioning [Mattia et al. 2004].Global and regional solutions for GPS networks are typically analyzed in a global reference frame which is defined at a particular time by choosing: (i) an origin representing the point where three axes intersect; (ii) a scale representing the unit of measurement; and (iii) directions for three orthogonal Cartesian axes.The origin or geocenter is often chosen to coincide with the center of mass of the Earth.The typical unit of measurement is the meter.The axes are usually chosen to provide alignment with historical Earth orientation estimates.Each of the seven reference-frame parameters has a rate representing drift in the frame over time.
The best currently available geocentric reference frame is the International Terrestrial Reference Frame (ITRF).The ITRF is maintained by the International Earth Rotation Service, which monitors the Earth orientation parameters through a global network of observing stations.This is carried out using space geodesy techniques, such as: very long baseline interferometry, lunar laser ranging, satellite laser ranging, Doppler orbitography and radio positioning integrated by satellite, and GPS observations (see Altamimi et al. [2007] for an overview).These regional solutions in the ITRF can then be transformed into a «rigid plate» or other physically meaningful reference frames, either by using published Euler vectors that describe plate motions in the ITRF, or by including stations in an analysis that can be used to define a rigid plate or block.In this study, we defined a local reference frame from more than a decade of GPS measurements, to refer the Mt.Etna ground deformation pattern to a rigid block.In particular, by using site velocities and associated errors, we estimated the Euler pole of rotation by minimizing the adjustments to velocities of a dense set of well-distributed stations located around the volcano edifice.

Why a geodetic reference frame for Mt. Etna?
On Mt.Etna, the use of a consistent and stable reference frame to correctly isolate the ground deformation has become critically important since the first GPS measurements were carried out on the volcano edifice and surrounding

Subject classification:
GPS monitoring, Euler pole, GPS velocity, Mt.Etna, Rigid block.areas.To isolate volcano deformation, the researchers working on Mt.Etna have referred GPS measurements to: i) an external network located outside the volcanic area (Figure 1, black circles; Palano et al. [2008] and reference therein); ii) a specific station located inside the volcanic area, close to: 1. the summit area (an older benchmark located very close to the EPLU station; Murray [1994]), 2. the middle southern flank of the volcano (the ENIC station; Bonaccorso et al. [2002]), 3. the basal area (the EIIV station; Aloisi et al. [2003]; Patanè et al. [2005]); iii) the NOTO/NOT1 station located about 90 km south of Mt.Etna (Figure 1, inset;Houlie et al. [2006]).
The use of these different reference frames over time makes the direct combination and comparison of results very difficult, especially for computing valid velocity estimates.This problem can be overcome by using a consistent and stable reference frame over the time by adopting the same rigorous approach used to describe plate motions.As Mt.Etna is located close to the Eurasia-Nubia plate boundary, geodynamic studies of this plate boundary usually rotate GPS velocity fields to a Eurasian and/or Nubian reference frame.In a Eurasian frame, the GPS sites located in Sicily largely show a northward motion with an average velocity of about 5 mm/yr, in agreement with the Eurasia-Nubia convergence process [Hollenstein et al. 2003, D'Agostino and Selvaggi 2004, Mattia et al. 2009].In a Nubian frame, the GPS sites show a systematic eastward motion with an average velocity of about 2 mm/yr, revealing a different motion with respect to the Nubia plate [D'Agostino andSelvaggi 2004, Serpelloni et al. 2005].In light of this, it is clear that both the Eurasia and Nubia reference frames are not suited to correctly isolate the «volcanic» ground deformation from the background «geodynamic» pattern.Only the establishing of a local reference frame that can be used to: (i) isolate in a rigorous mode the volcanic deformation (or other local deformation); and (ii) directly allow comparisons of ground deformation over time, will prove useful to researchers wishing to analyze ground deformation on Mt.Etna, both for research and monitoring activities.For these reasons, in the following we computed a local reference frame for Mt.Etna.

GPS data and processing
To define a local reference frame to which the solutions of the Etn@net network can be referred, and to test this, we analyzed continuous GPS data collected by the following GPS networks in southern Italy between 1996.00 and 2009.75: • The Italian National GPS Network RING (http:// ring.gm.ingv.it/), a permanent, integrated and real-timemonitoring GPS network, which was implemented by the Istituto Nazionale di Geofisica e Vulcanologia (INGV) starting from 2004.The RING consists of about 120 stations distributed in Italy, with a greater concentration in southern Italy.
• The Sicili@net network, a permanent GPS network implemented by the INGV, sezione di Catania (INGV-CT) since 2005.At present this comprises 20 stations covering Sicily, with a larger concentration in its eastern part.
• The Etn@net network, a permanent monitoring network implemented by the INGV-CT since 2000 on Mt.Etna volcano.The network geometry was gradually upgraded in the following years, reaching the current configuration of 30 permanent stations that provide a dense cover of most areas of the volcano edifice.
To provide a consistent and stable reference frame over the past years, we also processed GPS survey-mode data from some benchmarks located around Mt. Etna (Figure 1) that have been used as reference frames in the past (see Palano et al. [2008] and references therein).
We analyzed all of the available GPS data using the GAMIT/GLOBK software [Herring et al. 2006] in a two-step approach [Dong et al. 1998].In the first step, we used GPS phase observations to estimate loosely constrained station coordinates, while in the second step we used these estimates from each day and from each network or survey as quasiobservations in a Kalman filter, to estimate a consistent set of coordinates and velocities.The quasi-observations are estimated parameters and associated covariance matrices for the station coordinates, and Earth-rotation parameters and source positions generated from analyses of the daily GAMIT solution [Herring et al. 2006].
In particular, the data were processed with the International Global Navigation Satellite Systems (GNNS) Service (IGS; http://igscb.jpl.nasa.gov)precise ephemerides and Earth orientation parameters from the International Earth Rotation Service (http://www.iers.org)Bulletin B. We tied the regional measurements to an external global reference frame by including data from nine continuously operating IGS stations (Figure 1; AJAC, CAGL, GRAS, GRAZ, LAMP, MATE, MEDI, NOT1 and NOTO) in the regional analysis.The regional quasi-observations were then combined with quasi-observations from global solutions provided by the Scripps Orbit and Permanent Array Center (ftp://garner.ucsd.edu/pub/hfiles)for the IGS1, IGS2, IGS3, IGS4, EURA and EUREF networks, to create a daily combined network solution aligned to the ITRF2005 reference frame [Altamimi et al. 2007].Before estimating THE Etn@ref GPS REFERENCE FRAME Table 1.Details of the sites used to estimate the Euler pole.GPS velocities (Vel.) in the ITRF2005 reference frame (as determined in this study), 1-sigma uncertainties (v), correlations between the east and north components of velocity (RHO), and residual velocity (res.) with respect to the Euler pole, as defined in Table 2, are also reported.Long., longitude; La., latitude.

Site
Long velocities in the second step of our analysis, we examined all of the position-time series for outliers or jumps.Only time series with an interval longer than 2.5 years were considered for velocity estimation, to mitigate the effects of bias introduced by a periodic annual signal [Blewitt and Lavallèe 2002].The reference frame for our velocity was defined in the second step, in which we applied generalized constraints [Dong et al. 1998] while estimating a seven-parameter transformation (three network rotations, three network translations, and one scaling parameter).In particular, we defined the reference frame by minimizing the horizontal velocities of 20 global IGS stations with respect to the ITRF2005 reference frame [Altamimi et al. 2007].Velocities in this reference frame are shown in Figure 2a.

The Etn@ref geodetic reference frame
The use of a consistent and rigid reference frame is critically important for terrestrial measurements.Two observers using different frames cannot directly combine or compare their results.One observer using different frames over time cannot compute valid velocity estimates.Scalar quantities, such as baseline length, have a minimal dependence on the reference frame and are affected by the choice of unit, and not by the origin or orientation of the frame [Heflin et al. 2002].Relative vector quantities, such as north, east, and vertical baseline components, or the relative sea level depend on both the unit and orientation, and not on the choice of origin.Absolute vector quantities, such as position, velocity, or absolute sea level, depend upon all of the reference frame parameters.Taking into account the ITRF2005 site velocities and associated errors, we estimated the three Euler vector parameters (latitude and longitude of pole, rotation rate) for a block, by minimizing the adjustments to two horizontal components of GPS velocity at each intra-block site with a weighted least-squares inversion.The Euler pole can be visualized as the fixed point on the surface of the Earth (not generally within the block itself ) around which the block rotates.This rotation model essentially constrains the block to move rigidly on the Earth surface (no radial motion).

Euler pole estimation
In a mathematical formulation, let a point P on the Earth surface be defined by its Cartesian coordinates X, Y and Z.We have chosen to work in Earth-centred-Earth-fixed (ECEF) Cartesian coordinates because the plate motion calculations are much easier in the ECEF frame using rotation vectors.In this ECEF frame, at a given location P, for a plate rotation defined by the rotation vector X (~x, ~y, ~z), the velocity vector in geocentric coordinates V (v x , v y , v z ) is given by the cross-product: (1) where t is the mean Earth radius (~6378 km).
Site coordinates and the resulting velocities are then in the ECEF frame, so they need to be converted to the local NEU frame.A valid velocity conversion at the corresponding geodetic location (m, {, h) is a combination of the three rotations needed to align the ECEF frame with the NEU frame: (2 where R is the rotation matrix: (3) and L (N, E, U) are the components of the velocity in the local frame.
Regarding the inverse problem, to calculate the pole parameters by knowing the velocities of some reference sites, the forward problem shown in Equation ( 1) can be written in a matrix form.In particular, for a given site, and taking into account the measurement errors, the equation is: (4) If there are at least two sites with known velocities, the problem is over-determined and can be solved by minimizing the residue model errors through the use of the least-squares approach [Tarantola and Vallette 1982]: (5 where W is the data weight matrix that is typically chosen as the inverse of the measurement error matrix.The model covariance matrix is then: The robustness of the inversion solution and the validity of the results are assured by the hypothesis of a Gaussian probability density for the parameters X and by introducing a Monte Carlo resampling step to estimate the mean and covariance of the parameters considered.Operationally, the Euler vector parameters were estimated by taking into account the ITRF2005 velocities of 13 sites located within a 350-km radius of Mt.Etna (Table 1).Defining a block-fixed reference frame requires careful selection of the sites that best represent the overall motion of the block.Formally, small standard deviations of the velocities (errors <0.2 mm/yr), observation histories longer than 2.5 years, no systematic effects in the coordinate time series (e.g.jumps), and geographic cover with respect to Mt. Etna were the criteria adopted for selecting the 13 «fiducial» sites for the referenceframe definition.To estimate the consistency of the fitting procedure results, the correlation coefficient r was used to evaluate the statistical significance of the correlation: (7) where O i and P i represent the observed and predicted values, respectively, and N represents the number of datasets [Palano et al. 2009].The r values can take on any value between 0 and 1, with a value closer to 1 indicating a better fit.The meaning of r can be explained by the coefficient of determination, r 2 , which indicates the percentage of data O i that, in turn, can be explained by data P i .According to statistical decision theory, r 2 >70% (r = 0.837) is an acceptable   test for prediction [Mulargia and Gasperini 1995].Furthermore, we are aware that the selected stations are located across the Eurasia-Nubia plate boundary; however, the RMS of residuals (about 0.7 mm/yr and 1.1 mm/yr for the east and north components, respectively) and the fitting procedure result in r 2 = 96.6%(r = 0.983), which indicates that a rigid plate is an adequate approximation given the location of the GPS stations.To verify that the predicted Euler pole is not critically affected by the selected sites, we performed some inversions that excluded for each computation three GPS sites, and repeating the estimation of the Euler vectors.The resulting poles did not significantly differ with respect to the full dataset estimates; in particular, we obtained differences less than 0.001 deg/Myr and 3.5˚ for the rotation rate and the location, respectively.
The observed and predicted velocities, misfit vectors, and best-fit Euler pole are shown in Figure 2. Table 2 provides a listings of the Euler parameters, their uncertainties, and the data misfit statistics.
Finally, the MatlabTM software package for the pole-rotation computation was implemented.The software package was developed to work readily with GAMIT/GLOBK velocity field files (org files) and velocity «apriori» files (apr file).This allows to calculate the expected velocity value for any point located on the Earth for a given Euler pole (direct problem), or to infer the Euler pole parameters by inverting the observed velocity at a set of sites located on a rigid block (inverse problem).This software package, called PlatEmotion (Figure 3), is available for the scientific community, and can be freely obtained by simply contacting the corresponding author.

Some examples of time series and velocity fields referred to the Etn@ref frame
In the following, we report two examples referred to the Etn@ref reference frame.As a first step, the daily positions for each site are combined to form time series of geodetic positions using the GLORG module of GLOBK.In particular, the unconstrained daily solution has been aligned to the Etn@ref reference frame by means of a weighted least-squares Helmert transformation.Figure 4 shows the NEU time series of the ENIC station.The large variations and the drift in the horizontal components are related to the volcanic activity, and in particular to the cycle of inflationdeflation due to magma uprising before, and discharge after, an eruption.
As a second step, the daily positions for each site are combined to estimate a long-term average velocity in an Etn@ref reference frame, using GLORG.Figure 5 shows an example of a velocity field observed on Mt.Etna and referred to the Etn@ref frame.This velocity field spans the 2009.00-2009.70time interval.

Discussion and conclusions
Although the GPS satellites provide a natural dynamic framework for ground-based geodesy, the doubly differenced phase observations (or as an equivalent, the undifferenced phase with clocks estimated or provided) do not tie a ground station to the orbital constellation at the millimeter level we require for scientific studies.Rather, researchers define and realize a precise terrestrial frame by applying constraints to one or more sites in their network.The simplest approach, called «finite constraints», is to keep the coordinates of one or more sites fixed over the time.A more rigorous approach to frame realization is through «generalized constraints», in which the adjustments of the coordinates of the framedefining sites are minimized while estimating the Helmert parameters [Herring et al. 2006].With this approach, all of the reference sites are free to adjustment.
In this study, more than a decade of GPS measurements have been used to define a local reference frame for a Mt.Etna ground-deformation pattern definition with respect to a rigid block.According to the «generalized constraints approach», we estimated the Euler pole for the rigid block by using a weighted least-squares inversion to minimize the adjustments to two horizontal components of the GPS velocity at each intra-block site.To this end, we selected 13 «fiducial» sites located within a 350-km radius of Mt.Etna.
In our case, the accuracy of the estimation of the Euler pole might prove problematic, mainly because the rigid block we considered is too small, and also because Mt.Etna is located in eastern Sicily, which is very close to the Eurasia-Nubia plate boundary, such that plate boundary deformation effects and/or differences in site motions in one or the other plates are possible.Roughly speaking, the angular-velocity parameters for any given plate would be best constrained by a network that maximally spans the rigid plate interior.For very small plates, as in our case, the motion can be characterized by a horizontal translation within the sensitivity of geodetic measurements.In this case, there is a high correlation between the rate of rotation and the location of the Euler pole normal to the direction of plate motion, and so the concepts of rate of rotation and Euler pole essentially lose their meaning.Nevertheless, what is important to geophysical processes is not the precision of the Euler pole and the rate of rotation, but rather the precision to which the relative motion is known across plate boundaries.We have constrained this parameter well, considering a set of «fiducial» sites that have good geometry with respect to the boundary of the rigid block.Finally, the block rigidity is assessed by the low misfit of the model to the data (coefficient of determination = 96.6%),which indicates that the Euler pole estimation is highly accurate.In this way, we can affirm that the objective of a stable frame for the Mt.Etna GPS network has been reached.This result is of direct application to the Etn@net network, and strengthens the role of the Etn@net network in both research and monitoring activities.It can also be considered an important starting point for every researcher working with GPS data on Mt.Etna.
Finally, the RING and other important projects have planned for the implementation of a denser and larger network of permanent GPS stations to improve deformation process studies over the whole of Sicily.This network improvement will allow us to estimate a new Etn@ref reference frame within the next three years.

Figure 1 .
Figure 1.Sketch map of Mt.Etna.The permanent GPS network stations are indicated as white circles.The older survey-mode external reference network stations are shown as black circles.The elevation contour interval is 0.2 km.Inset, the IGS stations used in this study.

Figure 2 .
Figure 2. (a) GPS velocities in ITRF2005 reference frame.The uncertainties have been excluded from the figure because their values are small with respect to the scale of the displacements (see Table 1 for details).Sites used to estimate the Euler pole are shown in yellow.The gray lines show the trace of the Plio-Pleistocene subduction front between the Eurasia and Nubia plates.(b, c) Predicted velocities (b, white arrows) and misfit vectors (c, black arrows) with 95% confidence ellipses.The estimated errors also take into account the Euler pole error, computed through a Monte Carlo method.(d) Euler pole derived in this study.

Figure 4 .
Figure 4. NEU time series of ENIC series referred to an Etn@ref reference frame.

Figure 3 .
Figure 3. Overview of the PlatEmotion Matlab TM software package.

Figure 5 .
Figure 5. Velocity field (black arrows) at Mt. Etna during the 2009.00-2009.70time interval referred to an Etn@ref reference frame.The permanent GPS network stations are indicated as white circles.

Table 2 .
Parameters of the Etn@ref Euler pole.The associated errors are also given.