GNSS Signals: A Powerful Source for Atmosphere and Earth’s Surface Monitoring

It is well known that Global Navigation Satellite Systems signals (which include for example the U.S. GPS and its modernization, the Russian GLONASS, the future European Galileo, the Chinese COMPASS), commonly processed for navigation purposes, can also be used to characterize media where they propagate in. In the last decade, GNSS atmospheric and Earth’s surface remote sensing become more and more important, thanks to technical improvements applied to the processing of such “free-of-charge”, everywhere available and weather insensitive signals.


Introduction
It is well known that Global Navigation Satellite Systems signals (which include for example the U.S. GPS and its modernization, the Russian GLONASS, the future European Galileo, the Chinese COMPASS), commonly processed for navigation purposes, can also be used to characterize media where they propagate in.In the last decade, GNSS atmospheric and Earth's surface remote sensing become more and more important, thanks to technical improvements applied to the processing of such "free-of-charge", everywhere available and weather insensitive signals.
For example, remote sensing of wet part of troposphere is possible "extracting" the atmospheric delays from GNSS observations.These delays are associated to water vapour and are accumulated by the signal along its propagation path.In the double difference phase observation adjustment (a standard GNSS signal pre-processing) it is possible and quite easy to estimate the wet contribution to atmospheric total delay mapped into the zenith direction, the so-called Zenith Wet Delay.From one side the estimate of propagation delays is essential to improve the accuracy of the height determination in the geodetic positioning framework (Kleijer, 2004).From the remote sensing point of view, Zenith Wet Delay may be then transformed into the so-called Integrated Precipitable Water Vapour (IPWV).Therefore, the knowledge of the temporal behaviour of IPWV above a GPS receiver network allows meteorologists to know the evolution of total water vapour content in atmosphere, which is one of the variable operatively used in Numerical Weather Prediction Models.These aspects are described in section 2.
A second important application allows to add vertical variability information to the atmospheric parameter distribution with respect to the previous one, which represents an "integrated" quantity.The amplitude and phase variations experienced by GNSS signal crossing the atmospheric "limb" and received on-board a Low Earth Orbit satellite, can be used to infer temperature and water vapor profiles, thanks to the GNSS Radio Occultation technique (Melbourne et al., 1994;Ware et al., 1996;Kursinski et al., 1997;Hajj, 2002).Even if aspects related to such very important Remote Sensing technique are not treated in the present chapter (a comprehensive tutorial can be found in Liou et al (2010), while review of results Interferometric Synthetic Aperture Radar (InSAR).Several techniques are well established to derive the vertically Integrated Precipitable Water Vapor (IPWV)1 , in particular using ground-based and spaced-based radiometers, radiosonde observations and GNSS receivers.
Radiosonde observations produce an accurate measurement of the water vapour profile, but the temporal and spatial resolution is rather poor.Radiosondes are typically launched every 6 to 12 hours, which may cause significant variations in water vapour to go undetected.
Ground-based microwave radiometers show problems during periods of rain fall and spacebased radiometer observations can be degraded in the presence of clouds.This prevents reliable measurements during periods where changes in water vapour could be quite great.Besides these limitations, all systems involve considerable costs.
The technique to estimate IPWV by means of GNSS receivers is based on measurements of the tropospheric delay time of navigation signals.Therefore the delay, regarded as a nuisance parameter by geodesists, can be directly related to the amount of water vapour in the atmosphere, and hence is a product of considerable value for meteorologists.Furthermore, water vapour estimation with ground-based GNSS receivers is not affected by rain fall and clouds, and can therefore be considered an all-weather system.So, GNSS is a valuable complement to radiosondes and radiometers, taking into account that GNSS IPWV estimates come from an existing GNSS infrastructure and frequently from quite dense receiver networks.

Description of observables, theoretical basis and retrieval technique
The use of GNSS receivers to estimate IPWV is based on measurements of the delay affecting the navigation signals during their propagation in troposphere (neutral atmosphere) from the GNSS satellites to the receivers on ground.The dispersive ionospheric effect can be removed with a good level of accuracy by a linear combination of dual frequency data.
Such a technique is founded on the non-dispersive refractive characteristics of the neutral atmosphere, governed by its composition.The water vapour molecules in atmosphere are polar in nature possessing a permanent dipole moment.All the other gases are non-polar molecules and a dipole moment is induced among these gases when microwave propagates through atmosphere.These molecules reorient themselves according to the polarity of propagating wave.In the retrieval technique to be described the atmosphere is considered as the sum of a dry component (mainly due to O 2 ) and a wet component.
Consequently, the neutral delay due to the troposphere can be decomposed into the hydrostatic delay associated with the induced dipole moment of the atmosphere constituents and the wet delay associated with the permanent dipole moment of water vapour (Askne & Nordius, 1987;Brunner & Welsch, 1993;Treuhaft & Lanyi, 1987).The zenith hydrostatic delay (ZHD) has a typical magnitude of about 2.4 m at sea level, and it grows with increasing zenith angle reaching about 9.3 m for elevation angle of 15°.With simple models and accurate surface pressure measurements, it is usually possible to predict accurately the ZHD.The zenith wet delay (ZWD) can vary from a few millimeters in very arid condition to more than 350 mm in very humid condition, and it is not reliable to predict the wet delay with an useful degree of accuracy from surface measurements of pressure, temperature and humidity.
Therefore, from GNSS radio signals the total tropospheric delay is provided and, measuring the ZHD, it is possible to retrieve the remaining ZWD, incorporating mapping functions which describe the dependence on path orientation.The ZWD time series are then directly transformed into an estimate of IPWV: GNSS receivers can estimate IPWV with a temporal resolution of 30 min or better and with an accuracy better than 0.15 cm.

Retrieval algorithm
In this section the retrieval algorithm used for the estimation of IPWV from GNSS observations is presented.
Using GNSS methods of path delay correction, developed for geodetic applications, it is possible to estimate time-varying atmospheric zenith neutral delay ZTD (excess path length due to signal travel in the troposphere at zenith) defined as: 6 10 ( ) where ds has units of length in the zenith, H is the surface height and N(s), usually expressed in parts per million (ppm), is the refractivity of air given by (Thayer, 1974): where P d is the dry air pressure (hPa), T is the air temperature (K), e is the partial pressure of water vapour (hPa), Z d and Z w are the dry air and water vapour compressibility factors, that consider the departure of air from an ideal gas.Values for inverse dry and wet compressibility factors differ from unity of about one part per thousand, and are given by: ( ) where T c is temperature in Celsius.
In eq. 2, the first two terms of N are due to the induced dipole effect of the neutral atmospheric molecules (dry gases and water vapour), and the third term is caused by the permanent dipole moment of the water vapour molecule.Therefore, the hydrostatic part is described by: 61 1 10 and the wet part is:


(5) In the estimation algorithm of IPWV we can identify four principal steps: 1. We start with the estimation of the neutral zenith path delay from GNSS observations (Bevis et al., 1992), which are elaborated using a specific GNSS software (e.g.Bernese GPS software or others).The neutral radio path delay has to be estimated using precise orbit ephemerides, choosing a proper cut-off angle (e.g. 15 degrees), resolution time (e.g. 30 minutes), and a suitable law as dry and wet mapping functions.Different kinds of mapping functions exist and they are different in number of meteorological parameters involved (Herring, 1992;Ifadis, 1986;Niell, 1996).2. Computation of the ZHD component of the atmosphere, that is the greater component in magnitude of ZTD but it is less variable with respect to ZWD.
If atmospheric profiles of temperature and dry pressure are available near the GPS station, the ZHD can be computed using eq. 4. Since such an availability is difficult in time and space, alternative and more simple procedures can be adopted for a reliable estimation of the hydrostatic delay.
If surface pressure is known with an accuracy of 0.3 hPa or better, ZHD can be estimated through simple models to better than 1 mm (Elgered et al., 1991), e.g. using the Saastamoinen model (Saastamoinen, 1972): where ZHD depends on actual surface pressure P s (hPa), on latitude λ (rad) and on the surface height H (km).The error introduced by the assumption of hydrostatic equilibrium in the model formulation is typically of the order of 0.01%, corresponding to 0.2 mm in the zenith delay.
3. Then ZWD is computed by subtracting ZHD from ZTD. 4. Finally, it is possible to retrieve IPWV using the relationship: Typical values for the parameter П are approximately 0.16, so 6 mm of ZWD is equivalent to about 1 mm of IPWV.
The parameter П is a function of various physical constants and of the weighted mean temperature T m of the atmosphere (Askne & Nordius, 1987;Davis et al., 1985): where ρ is the density of liquid water, R ν is the specific gas constant for water vapour, m is the ratio of molar masses of water vapour and dry air, and k 1 , k 2 , k 3 are the constants defined previously.
The transformation described in eq. 8 assumes that the wet path delay is entirely due to water vapour and that liquid water and ice do not contribute significantly to the wet delay (Duan et al., 1996).
Although the IPWV retrieval algorithm from ZTD measurements is well-established, different strategies were adopted for the time-varying parameter П. Anyway, П can be estimated with such an accuracy that very little uncertainty is introduced during the computation of eq. 8. Bevis et al. (1994) provided an error budget for П and showed that in most practical conditions the uncertainty for this parameter is essentially due to the uncertainty for T m (usually predicted from the surface temperature T s on the basis of regressions), leading to a relative error in П of the order of 2%.In fact, exact calculations of T m require profiles of atmospheric temperature and water vapor, as from radiosoundings or analysis from Numerical Weather Prediction Models (e.g the global European model, ECMWF).Since those data are not easily available, T m is commonly estimated using station data of surface air temperature with empirical linear or more complicated relationship (the so-called T m -T s relationship) that can be site-dependent and may vary seasonally and diurnally (Bevis et al., 1994).
A simple and alternative approach can be considered for П estimation: the use of a linear regression (ZWD and IPWV as predictors and predictands, respectively) from historical data base of radiosoundings or ECMWF available near the site of interest for the water vapour estimation, leading again to a relative error in П just above 2%.Considering monthly averages of П the uncertainty is around 1.5% (Basili et al., 2001).This approach does not need measurements of surface temperature for each computation of П.
Estimation of water vapour features by GNSS is valuable from the point of view of climate monitoring, atmospheric research, and other applications such as ground-based and satellite-based sensor calibration and validation.GNSS tropospheric delays are also useful for operational weather prediction models (Gutman & Benjamin, 2001;Macpherson et al., 2008;Smith et al., 2000).
The IPWV retrieval by means of a GNSS ground-based receiver can be used to monitor in situ water vapour time series, or to compare the IPWV values estimated by co-located ground-based sensors (e.g.microwave radiometer, photometer).Networks of GNSS receivers can be used to monitor the water vapour field, mapping its horizontal distribution.
The possibility of mapping IPWV measured by GNSS networks has been explored (de Haan et al., 2009;Morland & Matzler, 2007), also combining IPWV data retrieved from GNSS receivers and from satellite-based radiometers to produce IPWV maps over extended areas (Basili et al., 2004;Lindenbergh et al, 2008).

Results
The degree of accuracy in IPWV estimation by GNSS receivers exploiting the tropospheric propagation delay at L-band is usually around 0.10-0.20 cm.The horizontal resolution of zenith columnar water vapour associated to a single receiver using standard methods (azimuthally symmetric weighting functions) is in the order of tens of kilometers, roughly corresponding to the aperture of the cone which includes all the lines of sight of the various GNSS satellites observed at different elevation angles.
Besides GNSS, several techniques are well established to derive the vertically IPWV, such as ground-based microwave radiometers (MWR), radiosonde observations (RAOBs), analysis data from Numerical Weather Prediction Models (e.g ECMWF).Some examples of IPWV comparisons among different techniques during experimental campaigns are reported in this sub-section.
For instance, during an experimental campaign in Rome, Italy (20 September -3 October, 2008), different instruments managed by the Sapienza University of Rome were operative at the same site: a GPS receiver (included in the Euref Permanent Network) a MWR (a dualchannel type, 23.8 and 31.4GHz, model WVR-1100, Radiometrics) and six RAOBs (Pierdicca et al., 2009).Also, analysis data from ECMWF nearest the site were considered.The IPWV time series for the entire campaign are plotted in Fig. 1.The IPWV root mean square (rms) difference of GPS compared with MWR is 0.10 cm, with RAOBs and ECMWF is around 0.15 cm.
With reference to an Italian ground-based network of GPS receivers, managed by the Italian Space Agency (ASI), another experimental campaign was conducted in Cagliari (Italy), during the whole 1999 (Basili et al., 2001).The experimental site was selected at the Cagliari GPS station where a ground-based dual-channel microwave radiometer (WVR-1100) was operated for the whole campaign of measurements.Also, data from RAOBs released at Cagliari every six hours were available.Results of the experiment for the whole 1999 are shown in Fig. 2, gathered in non-precipitating conditions to avoid problems with the radiometer measurements.The comparison is performed considering a sampling time of 6 hours, in coincidence with RAOB releases.
This long-term comparison has shown a fairly good agreement among the two remote sensors and the RAOBs, with an error standard deviation similar to other experiments reported in literature.

Wet atmospheric refractivity maps through tomography
As it has already been shown in Section 2, the remote sensing of "wet" troposphere is possible by estimating the wet contribution to atmospheric total delay mapped into the zenith direction, the ZWD, in the general adjustment of double difference phase observations.Following a step ahead, it is also possible to try to extract some information on the three dimensional distribution of atmospheric parameters, from total delay observations taken by different line of sights.Tomography deals with the inversion of integral measurements collected from a great variety of directions, for the extraction of nonhomogeneous signatures inside the analyzed volume.Requirements necessary to make tomographic inversion procedures effective are well known.The geometry of the signal paths is crucial for the stability of the inversion procedure.All the voxels (volume pixels) have to be crossed by a lot of rays coming from different directions.Horizontal resolution can be improved only considering quite dense GNSS networks.Vertical resolution can be improved if receivers are deployed on a sloped area.This section presents results already published by Notarpietro et al. (2011), results obtained applying a tomographic inversion to real observations taken on October 2010 in Italy, by a dense network of GNSS receivers.

State of the art
Several activities were carried out in the past in the field of neutral atmospheric tomography based on observations performed on GNSS signals.Starting from one of the first concept description given by Elosegui et al. (1999), the effectiveness of a 4-dimensional (4D) water vapour field tomographic reconstruction was assessed by Flores et al. (2000) on a 20x20x15 km atmospheric domain against ECMWF (European Centre for Medium-Range Weather Forecast) data.After that, several methods were applied to different kind of real or simulated GPS observables (obtained by more or less dense receiver networks), demonstrating the effectiveness of water vapour field reconstructions on different atmospheric volume sizes, with different resolutions, against radiosonde data, Numerical Weather Prediction models or other independent water vapour dataset.Some reference papers (the list is not exhaustive) are that of Hirahara (2000), Gradinarsky and Jarlemark (2004), Champollion (2005Champollion ( , 2009)), Bi et al. (2006), Troller et al. (2006), Nilsson and Gradinarsky (2006).
In the framework of the European Space Agency project METAWAVE (Mitigation of Electromagnetic Transmission errors induced by Atmospheric Water Vapour Effects), we applied a new approach to the ZWDs estimated from the observations collected by a local network of GPS geodetic receivers deployed over a small area around the city of Como.Such new approach is based on an algorithm previously developed, on which we shown, from a simulative point of view only, the possibility to infer wet refractivity fields without using first guess atmospheric models and without adopting any a priori informations (Notarpietro et al., 2008).Such an algorithm has been applied to real measurements collected by a local network of GPS receivers.In what follows we will summarize the results we obtained.

Theoretical basis, retrieval technique, observables and validation approach
Basically, two different classes of algorithms can be applied to perform atmospheric tomography (and tomography in general).The first belongs to iterative reconstruction techniques (for example the Algebraic, the Multiplicative Algebraic or the Simultaneous Iterative Reconstruction Techniques, respectively called ART, MART and SIRT, see Herman 1980) which need a good first guess atmospheric model to converge at the "good" solution.
The second belongs to the Least Square Inversion (or Generalized Inversion) techniques, which are "one-step" algorithms and do not need a first guess.Notarpietro et al. (2008), shown the possibility to infer Wet Refractivity fields without using first guess atmospheric models.The algorithm accomplishes the reconstruction in two consecutive steps.The first step allows the retrieval of a "raw" three dimensional wet refractivity distribution directly from Slant Wet Delays (SWD) observables ΔΦ wet (defined as equivalent optical length), which in turn depend on the wet refractivity (N w ) distribution along the ray path, in the way defined by the following equation: Linearizing eq.11 and considering the entire observation dataset, the following matrix equation turns out: where L is the Data Kernel to be inverted to obtain the wet refractivity distribution, which is a matrix containing for each row, the lengths of each segment inside each voxel crossed by the generic rectilinear ray-path connecting the receiver and the satellite.This tomographic pre-processing step pertains to Least Square Inversion algorithms (Lawson & Hanson, 1974).It achieves the result through the constrained inversion (using Singular Value Decomposition) of the Tikonov-regularized Data Kernel matrix.Although the resolution obtainable with this pre-processing step is quite rough (the entire tropospheric volume has been divided into 2x2x20 voxels grid), this result is used as first guess for the algebraic technique used in the second phase of the proposed reconstruction algorithm.In particular we applied the SIRT technique to obtain the distribution of wet refractivity inside the tropospheric volume characterized by the final resolution (4x4x20 voxels grid).
With the aim of studying the potentialities of GNSS in the determination of local wet refractivity fields, needed for instance to correct InSAR derived landslide deformation maps, we used observations collected during a couple of weeks in 2008 by the MisT GPS network, defined by eight geodetic receivers that were deployed around the COMO Permanent Network station (which is placed in the North West part of Italy).This network was born for different purposes from the tomographic reconstruction of the wet refractivity field, and its design was not fully compliant with the requirements of this technique (details about each MisT station are reported in Table 1, while the MisT network topology is shown in Fig. 3).A daily multi-station adjustment of observations collected by the whole network was performed via the Bernese V0.5 software, to estimate jointly the station positions and the Hourly ZWDs parameters.These are basically averaged value of the tropospheric delay zenithal projection, affecting all the signals from the considered station to all the satellites in view, as they move along their orbits in 1 h time.Differences between the actual instantaneous slant delays and these averaged values projected back on the slant direction are to be found in the double difference adjustment residuals (this analysis is not described here).More precisely, carrier phase double differences were processed, all the single differences being formed with respect to the COMO reference station.The Bernese software models the tropospheric delay in each station-receiver phase measurement as the sum of a hydrostatic component and a wet one.The first can be modelled (and slanted toward the satellite position using the dry Niell's mapping function (Niell, 1996)) considering the Saastamoinen formulation (Davis et al., 1985) and interpolating surface pressure data (in time and space) obtained by 0.25°x0.25°ECMWF analysis.The second can be expressed as the product of an unknown parameter, the ZWD, by a known coefficient computed in our case from the wet Niell's mapping function.For each MisT station, input data were Hourly ZWDs, estimated during the week from October 12th to October 18th, 2008 and from November 13 th to November 19 th , 2008.Hourly ZWDs related to each MisT station, were then "geometrically" projected along the slant paths (using Niell's mapping functions) by upsampling at 1-min sample intervals the 15 min GPS satellites positions obtained from International GNSS Service (IGS) sp3 files and inverted using the developed tomographic procedure.
It has to be pointed out that the standard dataset adopted for tomographic reconstructions is built up by considering only 6 out of 9 MisT receivers.Firstly, COMO, ANZA and CAST are the three stations belonging to the so called MisT inner sub-network.We considered only ANZA among the three close stations of COMO, ANZA and CAST (deployed at distances less than 200 m from one-another), whose ZWDs are highly correlated (>95%).Moreover, ZWD data obtained processing NAND observations are used for self-consistency validation purposes ('leave-one-out' quality assessment) and are not included in the input dataset.
As we have previously stated, our tomographic approach is based on two consecutive reconstruction steps.The first one (data kernel generalized inversion) creates the first guess field for the second one (algebraic tomography), which doubles the horizontal resolution (from 2x2x20 to a 4x4x20 voxels grid, i.e. means 4.5x6.5x0.5 km 3 ).It has to be stressed that volume resolution is strictly related to the geometrical distribution of GNSS receivers and to the availability of observations.Higher resolutions would introduce an increasing number of voxels not crossed by any ray, thus worsening the final results.On the contrary, lower resolutions would imply a too coarse description of the field.
Considering the available observables we were able to obtain 168 or 144 Hourly wet refractivity maps (for the October or the November week respectively).Validation is carried out considering the difference between ZWD GNSS measurements taken over NAND receiver and corresponding ZWD estimates evaluated by vertically integrating the reconstructed wet refractivity maps.Considering the entire observing period, final statistics are thus based on 168 (144) ZWD differences (measured-estimated) distribution for the October (November) week and results are given in terms of their mean values and their rms values.

Results
In what follows, we will show results related to the so called baseline scenario and improvements obtained adding observations taken by mountainous receivers and from low elevation angles.Some hints about the impact of distance and height of the reconstruction error and about validation against independent data will be also given.

Baseline scenario results and effect of mountainous observation ingestion
The baseline scenario is that defined considering observations taken by the reduced MisT network formed by ANZA, BRUN, LAPR, PRCO, MGRA and DANI stations.For the October week, tomographic reconstructions were carried out considering ZWDs observed by the reduced MisT network observations taken during 12-18 October.The good agreement between measured and estimated ZWD time series evaluated above NAND during this period and for this scenario is shown in Fig. 4, while some statistics are given in column A of Table 2. Since data from the two mountainous receivers (BISB and BOLE in Fig. 3) were available only on 12 th October, 2008, between 9.00 am and 7.00 pm, comparisons of measured and estimated ZWDs above NAND receiver were performed also considering observations taken by the reduced MisT network in this smaller period (column B of Table 2).A bias decrease of 0.4 mm is observed adding BOLE (1199 m a.s.l.) observations (see column C, Table 2) and of 1 mm adding both BOLE and BISB (1373 m a.s.l.) data (see column D, Table 2) in the input dataset.This demonstrate the necessity of measurements collected at higher altitudes which allows a best reconstruction of vertical refractivity gradients characterizing the first three atmospheric layers.The high rms error with respect the one characterizing the baseline result given in column B, Table 2, is probably due to the more noisy data acquired by the two portable mountainous receivers (this is also evidenced by the decrease in correlation observed between NAND ZWDs measurements and estimates).Table 2. Statistics of the ZWD difference (measured-estimated after reconstruction) over the NAND reference station.

Ingestion of low elevation observations
Considering the baseline scenario described in paragraph 3.3.1, it is clear that the improvement in the reconstruction of lower layers is strictly related to the availability of trajectories crossing (and discriminating) the lower tropospheric layers.In our tomographic reconstruction, only rays exiting from the top boundary of the analyzed 18x26 km 2 x10 km volume were considered.In our case, the mean elevation angle was about 30°.Since the MisT network topography is fixed, to overcome this limit and therefore improving the retrieved field, we try to ingest also low elevation trajectories which enter from the lateral boundaries of the analyzed volume.Since SWDs associated to these rays contains both a contribution of the wet refractivity field inside the considered volume (namely, the inner volume) and outside the volume (the outer volume) up to 10 km height, we modelled and removed this last quantity from the SWDs associated to low elevation (< 30°) ray before entering the tomographic approach.The wet refractivity model considered in the outer volume was obtained considering three different approaches: a. from a very coarse tomographic reconstruction performed on a bigger volume using the same GNSS experimental data (considering as input data those observed by the entire MisT network except those taken by the NAND receiver); b. interpolating the CIRA-Q wet atmospheric climatologic model (Kirchengast et al., 1999) in the outer volume; c. considering data taken by ECMWF analysis (91 pressure levels, 0.25°x0.25°grid resolution), collocated in time and space with the centre of each voxel belongs to the outer volume (this was done by a bilinear space interpolation and a linear time interpolation of the meteorological data).
Results related to this analysis are summarized in Table 3.They confirm the importance of the availability of low elevation measurements issued from different altitudes to improve the estimation of vertical refractivity gradients in such a tomographic approach.It has to be noted that the availability of external independent information (atmospheric models or, better, meteorological data) for modelling the SWD component of low elevation observations in the outer volume seems to be necessary in this case.Because of the MisT network design (receivers not homogeneously distributed in the inner volume), the internal procedure based on the coarse tomographic reconstruction (case a)) is not very effective.
Table 3. Self-consistency results considering SWD derived by low elevation observations (taken during the October week) after the application of the outer volume wet refractivity modelling strategies a., b. and c.. Results are relative to the statistics of ZWD errors (measured-estimated after reconstruction) over the NAND reference station.Results related to the baseline scenario are reported in the first column as a reference.In the last column the evident outliers due to measurements (see blue dots in Fig. 4) were removed.

Impact of distance and height on reconstruction goodness
Results described previously are good, but are related to the baseline scenario.
Considering this scenario, the validation has been performed above NAND receiver, which is in a good position since its baseline from the COMO master station is between the nearest and the farthest stations.In this further analysis we have considered all the measurements (ZWDs or ZTDs) available from the MisT network (8 receivers) during the entire week (excluding only the COMO receiver, see Fig. 3).Then we have excluded data (ZWDs or ZTDs) observed by one receiver per time, keeping such data as reference for the self-consistency validation purpose for that receiver.For each case we have run our tomographic reconstruction considering all the 168 Hourly ZWDs (or ZTDs) available per each station for the October week, mapping them into the slant directions and including also low elevation observations (following the procedure described in paragraph 3.3.2).The obtained 168 Wet Refractivity Maps (considering ZWDs as input to the tomography) or Total Refractivity Maps (considering ZTDs as input) have then been used to evaluate the ZWD and ZTD estimates above the reference receiver, which are compared with the ZWD and ZTD observations above that receiver.This analysis has been repeated for each receiver of the MisT network.
Root Mean Squares of ZWDs and ZTDs differences (measured-estimated after reconstruction) are then reported in function of the distance of the station from the COMO master station or in function of the height of the station.Such results are plotted in Fig. 5.The same analysis has been performed considering data taken by the MisT network extended to the two mountainous receivers during 12 October from 9:00 AM to 7:00 PM (only 10 Hourly averaged ZWDs or ZTDs observations are contemporaneously available to any receivers of the "extended" network).In this case results are shown in Fig. 6.First of all this analysis confirms the impact of a good height displacement of receivers in the network.Even if MisT network topography has not been optimized for the geography of the analyzed area and for tomographic applications, if we consider the impact of height in the evaluation of propagation delays, we can say that the lack of receivers placed at higher altitudes will worsen final results.In particular, considering the original MisT network, where all the receivers are more or less placed in the same layer of the map (Fig. 5) we want to highlight that, if data observed at the highest receiver (namely BISB, which is placed in another vertical layer) are not given in input to the tomography, the rms of the difference between estimated and measured zenith delays (both Wet and Total) is generally doubled.
Things are better if we consider the MisT network plus the mountainous receivers.In this case the worsening is not so emphasized, since it is compensated by other receivers placed at similar altitudes (see Fig. 6).In both cases it seems that the results worsening follows a (more than) linear rule.It is absolutely not clear why the effects on the evaluation of Wet delays and Total delays are inverted, considering or not considering the mountainous receivers.It has to be noted that the network solution obtained for the mountainous receivers is not as accurate as that obtained for the other receivers, since the mountainous sensor positions have not been fixed.Moreover, results reflect 10 hours of observations instead of the entire week.
As far as the impact with distance is concerned, it is quite difficult to identify a clear relationship with results.Obviously if we exclude data observed by the nearest receivers (ANZA or CAST) to the reference one (COMO), results are better (rms is halved considering both the weekly data of the original MisT network and the 10 hours data of the MisT network plus mountainous receivers) than that we can obtain excluding one of the other (farther) receivers.But for all the other cases, it seems that final results are insensitive to distance.It is a surprising result since we expected a certain error correlation with distance.But the farthest receivers (MGRA and DANI) are placed in opposite positions with respect the map center and are the southest receivers (see Fig. 3).
If we take into account low elevation observations (even if such observations are averaged, since they are obtained simply mapping hourly averaged Zenith observations into slant directions), rays related to the northern receivers (all the others) anyway interest the atmospheric volume above the southest receivers (and not viceversa, given the orbital positions of GPS satellites).And this could probably compensate the "distance" effect.Anyway, also in this case, further analysis and measurements are necessary to better understand if there is a clear relationship.

Validation against independent data
In order to assess the goodness of inferred wet refractivity fields in different points of the grid considering independent data, we also did a comparison of ZWDs obtained vertically integrating wet refractivity fields derived after tomographic reconstruction along each column of retrieved maps with those derived by ECMWF analysis co-located in the same points (and times), even if the ECMWF horizontal resolution (0.25°x0.25°) and time resolution (6 h) are too coarse with respect those characterizing our final maps.
Statistical comparisons were performed considering the 168 wet refractivity maps obtained using data observed by the reduced MisT network (plus NAND receiver) collected during the October week and considering the 144 maps obtained for the November one.Results are shown in Fig. 7, where the time series of both ZWDs estimated after tomographic reconstruction (blue lines) and evaluated using ECMWF data (red lines) are plotted for each column of our volume discretization.We classified the areas accordingly to the corresponding rms values (computed for each ZWD difference time series, after the average bias removal) using green, yellow and red colors.As expected, the northern part is where the agreement is worse.In that area we had no receiver and less satellites were in view in the north direction.On the other hands, in the southern area, agreement is better even if no receivers were present, thanks to the availability of a higher number of rays.The best area is obviously the central one.meteorological events happened during the observing period.Anyway, an increase of wet refractivity (water vapour concentration) can be evidenced between the 100th (4 am, 16th October) and the 120th (midnight, 16th October) hours, with a peak around the 110th hour (2 pm, 16th October).The increase is well reproduced in terms of integrated wet refractivity along zenith (see ZWD evolution in Fig. 8(top)).Moreover, meteorological data (not shown here) confirmed an increase of cloud covering during that time interval.

GNSS reflectometry
Other than for atmosphere monitoring, GNSS signals may be used to characterize the Earth surface.In this section this kind of remote sensing technique is described, considering two scenarios of observation: ocean and land.
The exploitation of GNSS signals reflected off the oceans allows to obtain altimetry measurements (sea surface heights), surface roughness from which wind intensity and direction is determined, sea-ice topography and its stratification.Additionally, land observations are used to determine the soil moisture content and to monitor the surface snow cover.
The most of performed experiments are based on code measurements, since signal phase coherence after reflections is not many times maintained, because smooth surfaces are rarely found in reality.

Description of observables, theoretical basis and retrieval technique
For remote sensing purposes, the reflected and direct GNSS signals coming from the same satellite are collected on bistatic radar geometry; at least two antennas are required: the first RHCP (Right Hand Circularly Polarized) and zenith looking in charge of receiving the direct signal, the second LHCP (Left Hand Circularly Polarized) and nadir looking used to track the reflections.
In order to be more precise, the overall system could be considered as a multistatic observing system, since up to 6/7 GNSS transmitters are contemporary visible by the receiver antenna.
Each reflection is geo-referenced knowing the geometry of acquisition, looking at the point where the GNSS signal is reflected under specular condition; for doing this, the observer coordinates are necessary.Therefore the direct signal is used not only as a reference but also for computing the position of the receiver.
Three acquisition scenarios are possible: • Ground based: in this static configuration, the receiver is placed over mountains, towers and bridges and the collected measurements are used for testing the instrument functionalities and for monitoring small areas (i.e.coastal altimetry, local soil moisture content determination); • On aircraft: the sensor is placed on aircrafts or rarely on balloons to demonstrate its performances and to monitor small regions with higher spatial resolution than spacebased measurements.This dynamic configuration requires an evaluation of the Doppler shift due to the non-zero velocity of the aircraft; furthermore, this Doppler shift improves the resolution on the surface by means of iso-Doppler lines computation.

•
Space based: the sensor is placed on board a LEO satellite (400-800 km) with the aim of monitoring the entire Earth surface assuring a global coverage of the acquired reflections, which may be detected also very far from coastal zones (i.e. in the middle of the ocean); the Doppler shift experienced by the signal is the largest achievable among the three described scenarios.
The shape and extension of the footprint of the reflections depends on: the surface roughness, the sensor height above the Earth surface, the elevation of the reflected ray, the direction of the incidence plane respect to the receiver velocity.
The footprint must be considered lying on a plane tangent to the Earth surface in the specular reflection point.The distance of the specular reflection point from the receiver nadir increases when the elevation of the GNSS satellite decreases.
Inside the area interested by the reflection, the smallest resolution achievable from a geometrical point of view is determined by the cells generated by the intersections of the iso-delay and iso-Doppler lines.
Iso-delay lines are determined considering the points on the surface by which the reflected signal arrives at the receiver with the same delay.Generally speaking, these points are ellipses and are determined considering a single chip of the GNSS code as relative delay associated to each ellipse respect to the adjacent one (Martin-Neira, 1993).
Iso-Doppler lines are determined considering the hyperbolas on the surface where reflected signals come to the receiver with the same Doppler shift.The zero Doppler line is computed as the line passing through the receiver and orthogonal to its velocity direction (Martin-Neira, 1993).
Clearly, we cannot forget the antenna footprint, which acts as a filter in delay and Doppler on the surface looks.When the surface is smooth, the total power received is almost coming from the first Fresnel zone defined around the specular scattering point (Beckmann & Spizzichino, 1987).In this case, the computation of the cross-correlation between the reflected signal and the local GPS code replica gives a waveform simply delayed respect to the cross-correlation of the direct signal, but with the same triangle shape and a noise floor around.
When the surface is rough non-coherent reflections are expected and the use of the Fresnel zone to model the received power is ineffective.In this case, the glistening zone represents the source of scattered power (Beckmann & Spizzichino, 1987).
N scattering elements contained in the glistening zone are considered in determining the cross-correlation function where Λ is the triangle cross-correlation function and the m index indicates the quantities to the modelled signal generated with the local GPS code replica.Through this formulation Rp becomes a summation of triangle functions weighted with the amplitude of the n th element scattered field and delayed accordingly to the phase shift associated to each n th scattering element.The final correlation function shape in this case is shown in Fig. 9.
www.intechopen.comThe delay is used to determine the surface height, so is considered in case of GNSS signals reflected off water surfaces (Martin-Neira et al., 2001;Hajj & Zuffada, 2003).The height of the surface respect to the observer is retrieved in eq.14 through the relative delay Δτ, the speed of light c and the elevation angle of the reflection γ.
On the other hand, the reflected power is used to determine the surface reflectivity and the scattering cross section (Masters et al., 2004).
The surface reflectivity belongs to the coherent part of the scattered power that is measurable from the specular part of the received echo; it is used to determine the reflection coefficient that is related to the incident angle and dielectric constant.The dielectric constant is related to the soil composition and to its moisture content following empirical models or carefully calibrating the data (Masters et al., 2004).
The surface can be characterized looking at its roughness from the scattering cross section, since it contains the non-specular part of the reflected power.In this case we consider reflected power part calculated from the amplitude and the gradient of the correlation function on the right side of its maximum.
In order to retrieve surface winds over the sea, the shape of the non-specular echo is compared with a simulated one obtained using a sea surface model (Zavorotny & Voronovich, 2000;Elfouhaily et al., 2002).

State of the art
Winds retrieval and altimetry are the more consolidated applications, while soil moisture and ice monitoring are under-development.
Many instruments were developed up to now (Nogues-Correig et al.,2007) and the techniques of retrieval have been tested through many experimental activities.The early experiments deal basically with altimetry; measurements were collected either from a static position (Martin-Neira et al, 2001), from balloon (Cardellach et al., 2003) or from aircraft (Lowe et al, 2002).Other set of experiments were developed to retrieve the ocean surface state (Garrison et al., 2000), such as wind or sea roughness.Last but not least, the technique was demonstrated on board a small satellite, the UK-DMC (Gleason et al., 2005).
Nevertheless, nowadays no operative missions exist in this field.
From our point of view, during the SMAT-F1 project we developed a prototype based on a Software Defined Radio solution, using a navigation software receiver (Tsui, 2005).This is the NGene SW receiver, developed by NAVSAS group of Politecnico di Torino (Fantino et al., 2009).The instrument is highly reconfigurable, since collects raw I and Q IF samples of the incoming signals (direct and reflected).A sampling frequency of 8.1838 MHz is used, giving about 8 samples per C/A code chip.
Moreover, the small hardware architecture is made up of cheap COTS (Commercial Of The Shelf) components, with very low overall weight and power consumptions.These features make the system suitable to be easily placed on board aircrafts, also small U.A.V.s (Unmanned Aerial Vehicle) (Cucca et al., 2010).

Results
Using the described receiving system, we carried out two experiments.The first data collection has been made on a static position looking at the sea surface from a high cliff.The second was performed placing the receiver on an aircraft and acquiring GNSS signals reflected from rice fields.

Sea surface data collection
The first data collection was carried out on December 2010, from Sardinia Eastern coast at 157 m above the sea surface, near Cala Gonone.This region is characterized by high cliffs like those shown in Fig. 10.

RX H= 157
For satellite 16 (Fig. 12 (right)) two different echoes are visible.The first echo is coming from the coast (range closest to the observer and lower intensity due to scattering from the terrain); while the second is characterized by a greater delay, a different Doppler and successive returns lasting about 2 chips.In this case the correlation peak expires after 24 samples.This is a typical example of the capability to extract informations also from the uncoherent part of the signal.
Thus, our receiving system is able to track coherent and un-coherent reflections and to contemporary distinguish between echoes with different delays, Doppler shifts and intensity.

Rice fields data collection
During the second data collection of May 2011, an experiment performed flying over an area placed in the Piedmont region (north west part of Italy), the receiving system was placed on board a small aircraft in order to track reflections from rice fields.Since rice fields are flooded during this month, they are a perfect scenario to study reflection phenomena.
Like in the previous experiment, the geometry of reflections was analyzed and all the satellites with elevation lower than 33° were discarded, since below this elevation the specular reflections did not enter inside the -3 dB beam-width of the LHCP nadir looking antenna.
The signal to noise ratio detected from the reflected signal was normalized respect to the correspondent direct signal; moreover, we compute all the specular reflection points visible and to each point we associate the relative normalized signal to noise ratio.
On board the aircraft, a video camera was placed to see which fields were really flooded during the acquisition.The panoramic view extracted from the video was superimposed on ©Google Maps, together with the specular reflection points (Fig. 13).After the superposition, we have noticed a good agreement between the fields' state and the received power (see Fig. 13 and 14).The minimum received power correspondent to a low normalized signal to noise ratio is clearly associated to not flooded fields.

Conclusions and outlook
Scope of this chapter was to give an overview on some very powerful and quite recent Remote Sensing possibilities emerged exploiting GNSS observations, which complement the atmospheric and Earth's surface remote sensing traditionally performed by dedicated payloads and instrumentation.
GPS ground receivers can provide valuable and accurate information on integrated precipitable water vapor, considering that single receivers or fairly dense networks are available in many part of the world, providing a quite cheap and reliable source of information.As described previously, many investigations have been carried out in this respect to develop processing techniques, to validate the results through comparisons with independent sources and to exploit the final product.For instance, ZTD or IPWV data from a GPS ground based network can be assimilated into Numerical Weather Prediction models, or integrated with additional sources of IPWV to produce two-dimensional water vapour fields, leading to improved products.
As far as the tomographic approach for the retrieval of Neutral Atmospheric Refractivity maps is concerned, we demonstrated that it is possible (and with a good level of accuracy) as long as some tricks are taken into account.In particular it has to be outlined that, in order to make neutral atmospheric tomography more effective, the choice of the GNSS network topology is a key aspect.A good horizontal receiver's distribution guarantees a good retrieval of horizontal gradients.A good vertical receiver's distribution guarantees also a good retrieval of vertical gradients.Even if our network topology was not optimal for tomographic purposes, the inclusion of measurements (even if not very accurate) performed by two receivers placed at higher heights and of the low elevation observations, demonstrate this aspect.Since a suitable vertical receiver distribution is difficult to implement, the availability of quasi-horizontal observations is necessary.Then, limb sounding Radio Occultation observations are necessary in order to guarantee good observations coming also from low elevation angles (this aspect has already been demonstrated by Foelsche andKirchengast, 2001 andNotarpietro et al., 2008).
GNSS signals reflected off the Earth surface which represent an error source for navigation purposes, are instead useful for characterizing land and sea surfaces both from a monitoring and early-warning point of view.In particular the possibility of extracting information about the sea height and roughness, the soil moisture content, the snow and ice cover state have been successfully proven.Presently, no operative missions exist but many experimental activities have been carried out and the interest of national space agencies is constantly growing.From our point of view, we put some efforts in developing an instrument capable of collecting reflected GNSS signals, since we believe in the potentialities of this technique.
We definitely believe that the "expansion" of GNSS sources expected when also the European GALILEO, the Indian IRNSS and the Chinese BEIDOU navigation satellite systems will be deployed, together with the consequent availability of Radio Occultation observations, and the consequent availability of "vertical" and "horizontal" observations, will improve definitively all the techniques here presented.

Fig. 3 .
Fig. 3. Geographic distribution of the MisT network.The two "mountainous" GPS receivers are highlighted.The final volume discretization is also superimposed.

Fig. 4 .
Fig. 4. Time series of ZWDs measured (blue dots) and estimated (red dots) after reconstruction above NAND station, for the baseline experiment.

Fig. 5 .
Fig. 5. rms of the differences between ZWDs (blue dots) or ZTDs (red dots) observed and estimated above each reference receiver, excluding data of that receiver from the input dataset before the reconstruction.All data observed by the MisT network during the entire week are taken into account.(Left) rms are plotted against the distance of the reference receiver from COMO master station.(Right) rms are plotted against the height of the reference receiver above WGS84.The degraded results obtained excluding BRUN receiver (which is the highest one) are highlighted.

Fig. 6 .
Fig. 6.Like Fig. 5, but considering all data observed by the MisT network and by the two mountainous receivers during the 10 hours of 12 th October, 2008.

Fig. 7 .
Fig. 7. Time series of ZWD obtained integrating ECMWF (red) collocated and estimated with tomography (blue) wet refractivity maps.Left: October week data; right: November week observations.The black numbers shown the column "number" inside the map.Even if our main goal was to demonstrate the effectiveness in adopting tomographic reconstruction procedures for the evaluation of propagation delays inside water vapour fields, the real water vapour vertical variability and its time evolution is also well reproduced.Fig.8(bottom)  shows the time evolution of wet refractivity vertical profiles evaluated in the map centre (voxel 11 -see Fig.7) during the overall October week, considering data taken by all the available MisT receivers.Unfortunately, no meaningful

Fig. 9 .
Fig. 9. Shape of the correlation function for non-coherent reflections (black); each triangle refers to the signals received by an isorange, 8 samples equal to 1 C/A chip The basic observables are the delay of the reflected signal respect to the direct one, and the received power after reflection.Both observables are retrieved looking at the correlation function of the reflected signal and eventually comparing or normalizing it with the correspondent correlation of the direct signal.

Fig. 13 .
Fig. 13.Specular reflection points tracks for satellite 8 and 26 over Piedmont rice fields with relative normalized signal to noise ratio.See Fig. 13 for the red rectangle zoom.Furthermore, we compare the signals of two different satellites with similar elevation but different azimuth; we notice a high correlation between the two specular reflection point tracks both from the qualitative (Fig.13, Fig.14) and the quantitative (Fig.15) point of view.The quantitative comparison is performed considering the reflected power coming from the same longitude, considering a bean of 0.01°.Further investigations on this behavior are under development.