New GNSS tomography of the atmosphere method – proposal and testing Michal Kačmařík1, Lukáš Rapant2 1Institute of Geoinformatics Faculty of Mining and Geology VŠB–Technical University of Ostrava 17. listopadu 15, 708 00, Ostrava-Poruba, Czech Republic michal.kacmarik@vsb.cz 2Department of Applied Mathematics Faculty of Electrical Engineering and Computer Science VŠB–Technical University of Ostrava 17. listopadu 15, 708 00, Ostrava-Poruba, Czech Republic lukas.rapant@vsb.cz Abstract Paper is focused on GNSS meteorology which is generally used for the determi- nation of water vapour distribution in the atmosphere from GNSS measurements. Water vapour in the atmosphere is an important parameter which influences the state and development of the weather. At first, the paper presents basics of the GNSS meteorology and tomography of the atmosphere and subsequently introduces a new GNSS tomography method which doesn't require an extensive network of GNSS receivers, but uses only a few receivers situated in a line. After a theoreti- cal concept describing this method and used mathematical background, the results from a real experiment are shown and discussed. Unfortunately the results indicate that presented method is not able to provide credible outputs. Possibly the main problem lies in an insufficient number of available signals from current global nav- igation satellite systems (GPS and GLONASS) where the improvement could be expected after the start of Galileo and Compass. Potential ways how to improve the results without increasing the number of satellites are outlined in the last section. Keywords: GNSS tomography, GNSS meteorology, slant wet delay 1. Introduction Signal from the GNSS satellite is during its path through the atmosphere influenced by troposphere and ionosphere. The influence of the ionosphere is possible to eliminate from the GNSS measurements using a suitable method of dual-frequency data processing. The troposphere cannot be simply eliminated and, due to its dependence on the variable amount of water vapour in the atmosphere, it cannot be modelled with sufficient accuracy based on external data. However, it is possible to determine the troposphere and thus also the water vapour from GNSS measurements. Currently used methods for atmospheric water vapour determination (e.g. meteorological radiosondes, remote sensing satellites, water vapour radiometer) are connected with some limitations and disadvantages and these could be completed by processing GNSS data. On Geoinformatics FCE CTU 9, 2012 63 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method the other side, even this method has its own problems which cause a limitation when used in numerical weather prediction (NWP) models. The main disadvantage is that it can produce only an absolute value of atmospheric water vapour in a zenith direction above the GNSS receiver in any time (Zenith Total Delay, ZTD). Thus, the vertical distribution of water vapour above the receiver, which is the information very important for data assimilation in NWP models, is not known. From this reason the method of GNSS tomography, which would allow a three-dimensional study of atmospheric water vapour, has been developed since 2001 (Flores et al., 2001, Noguchi et al., 2001). At the present day there exists a few approaches how to solve this task, but all of them require a dense network of GNSS reference stations to achieve a reasonable spatial and temporal resolution. Those networks had to often been built artificially in the past for the needs of tomographic research projects. However, nowadays a possibility of using national reference station networks for creating tomographic systems over area of a whole country becomes an interesting topic. Both local and national tomographic systems require accessibility to the data from a large number of GNSS reference stations and necessary technical background for its processing. Because all contemporary tomographic projects are connected with the necessity of dense network of GNSS receivers, authors attempted to propose a new easier method. Its core concept is to use only one cross-section throw the atmosphere instead of solving a whole three- dimensional network of voxels. Thanks to that, the necessity of large number of receivers in a network would be reduced to a few receivers situated in a single line in a terrain. The output would be then limited to the 2D tomographic grid instead of the three-dimensional water vapour reconstruction. 2. GNSS meteorology It has been shown many times that the GPS system is useful for tropospheric parameters estimation - the approach called GPS meteorology (Bevis et al., 1992, Duan et al., 1996). In the past only the signals from system GPS were used for the determination of tropospheric parameters. But nowadays a combination of signals from more global navigation satellite systems (GPS + GLONASS) can be used for such a processing. This approach is analogically called GNSS meteorology. The standard value of the GNSS signal delay due to the atmosphere in the zenith direction is at zero sea level 2.3 m (or 8 ns). Troposphere influences GNSS signals in two different ways. Firstly, signal bends in the response to gradients of the index of atmospheric refraction, travelling a curved path instead of a straight line. Secondly, the signal travels throw an environment with some density slower than would travel in a vacuum. The total delay is the sum of those two effects. For the purposes of GNSS meteorology, however, ZTD can be divided in a different way – to the part caused by the hydrostatic components of the atmosphere (Zenith Hydrostatic Delay, ZHD) and part caused by the non-hydrostatic components of the atmosphere (Zenith Wet Delay, ZWD). While ZHD depends mainly on the atmospheric pressure, ZWD depends mainly on the water vapour (Zhengdong, 2004). During the GNSS data processing both components are modelled separately and the relation is expressed in (1). From a general point of view the hydrostatic characteristics of the atmosphere are rather stable in a time and it is relatively easy to model them. On the other side, the atmospheric Geoinformatics FCE CTU 9, 2012 64 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method water vapour amount is temporally highly instable and thus much more difficult for modelling. ZTD = ZHD + ZWD (1) On the basis of Zenith Total Delays (ZTD) and meteorological parameters observed in situ of GNSS measurements it is possible to precisely determinate the atmospheric water vapour (Precipitable Water Vapour, PWV). Significant advantage of the GPS/GNSS meteorology is its possible high spatial and temporal resolution. By processing data from twenty receivers and deriving PWV values in an hour interval we get 480 values for studied territory and a day. The contemporary used meteoro- logical radiosondes for the whole Czech Republic provide only 6 PWV measurements a day at two distinct places only. The stable and reliable results even under very severe atmospheric conditions including extreme storms, gales or torrential rains is another advantage of GNSS meteorology as confirmed by Liou and Huang (2000), Song and Grejner-Brzezinska (2009), Kačmařík (2010). First attempts to assimilate ZTD or PWV values from GPS measurements to the NWP models have been already done in the last years of the 20th century. Many researches shown that this step has a positive impact on a short-term prediction of rain localization and its power (Guerova, 2003, Koizumi and Sako, 2003, Nakamura et al., 2004, Peng a Zou 2004, Vedel a Huang, 2004, Shoji et al., 2009). In Europe, ZTD products from GNSS data processing are assimilated in NWP models in Great Britain (Bennitt, 2008) and in France (Moll et al., 2008). 3. Tomography of the Atmosphere For the modern generation of NWP models with a horizontal resolution of a few kilometres the existing network of sites launching meteorological radiosondes is not able to provide corresponding information about the vertical distribution of atmospheric water vapour. Since the classical GNSS meteorology cannot provide vertical profile, but GNSS tomography of the atmosphere can what has been already proven by many projects (e.g. Gradinarsky a Jarlemark, 2004, Champollion et al., 2004, Troller, 2004, Bender et al., 2011, Notarpietro et al., 2011). The tomography method is generally known and with a success it is used in many directions of human activity as medicine, rock and material studies, archaeology etc. The principle is a spa- tial discretization of the space above the network of GNSS receivers into a three-dimensional network of cells called voxels. Information about the water vapour is reconstructed for each voxel from all accessible slant wet tropospheric delays (SWD) observed between GNSS re- ceivers and satellites. Usually, it is assumed that within one voxel the water vapour is constant. Slant total delay of the GNSS signal caused by atmosphere can be written (Bevis et al., 1992): STD = 10−6 ∫ S Nds (2) where Geoinformatics FCE CTU 9, 2012 65 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method STD slant total delay N atmospheric refraction s length of the signal path In GNSS tomography slant total delays are provided by independent observations to all visible GNSS satellites and the refraction parameters are thus reconstructed from a large number of such delays. This defines an inverse problem, which could results in more than one solution. The primary problem is a limited number of signals given with the number of ground receivers and with a time variable number of distant satellites. Using a discretization of space above the ground network of GNSS receivers into voxels, the relation (2) can be written Ax = m (3) where m vector of i observations, i.e. slant path delays (SWD) x vector of j refractivities (N) in j voxels A kernel matrix, the mapping of the state x on the observations m. The matrix elements aij are the subsections of the i-th slant path in the j-th grid cell. (Bender et al., 2011) For the tomography purposes the path of signal is assumed to be a straight line and an influence of signal bending is neglected therefore. Thanks to this assumption a whole problem of inversion becomes linear and leads to a system of linear equations. 3.1. Input data As above mentioned, the slant tropospheric delays in various azimuthal and elevation angles are used for the purpose of GNSS tomography and not just ZTDs which are estimated in classical GNSS processing model. More concretely, the SWDs caused by the non-hydrostatic characteristics of the atmosphere are used instead of STD. The relationship between ZWD and SWD is given as follows: SWD(ε,θ) = M(e) × (ZWD + cot ε× (GN cos θ + GE sin θ)) + R (4) where ε satellite elevation angle θ satellite azimuthal angle M mapping function GN,GE horizontal gradients for north (N) and east (E) direction R post-fit residual Mapping function defines a relationship between tropospheric delay in a zenith angle and those observed in all elevation angles. The basis of any mapping function is 1/ sin ε, which corresponds to the simplest case. Currently, the most used mapping function is defined by Niell (1996). Horizontal gradients are considered to bring in the model the information about an azimuthal anisotropy in the atmosphere and in most accurate application these are usually estimated too. Post-fit residual represents a difference between a measured observation and an observation adjusted during the data processing. It is expected to contain the anisotropic part of the Geoinformatics FCE CTU 9, 2012 66 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method wet delay, but also various other unmodeled errors which influence the GNSS measurements: multipath, errors of antenna phase centre variations or signal noise. Some tomography projects use the tropospheric horizontal gradients and/or post-fit residuals as additional information of water vapour when reconstructing the anisotropic part of the atmosphere (Flores et al., 2001, Noguchi et al., 2001, Gradinarsky a Jarlemark, 2004, Bender et al., 2011), but other authors (Gradinarsky and Jarlemark, 2004, Nilsson et al., 2005) are sceptical of such an approach. Kačmařík et al. (2012) used various approaches for deriving SWD values from GNSS reference station GOPE and confronted them with in situ water vapour radiometer measurements. The comparison shown that using the tropospheric horizontal gradients or post-fit residuals (with or without correcting for multipath) did not give improved estimated SWDs. Adding the post-fit residuals to the GPS SWD had only a very small positive impact on bias, while it negatively influenced the standard deviation. The situation slightly improved after using the stacking maps for removing the multipath effect from the residuals. The results show that the post-fit residuals, even after multipath elimination, do not represent only the anisotropic part of the water vapour distribution but represent further insufficiencies in the model. Generally, a process of slant wet delay determination can be written as a sequence of these steps (Skone et al. 2003, Miidla et al. 2008): • ZTD computation with optional determination of horizontal gradients or post-fit resid- uals • ZHD computation using an appropriate hydrostatic model and values of meteorological quantities measured at the GNSS receiver’s location, derivation of ZWD from ZTD and ZHD • Computation of SWD for particular observations using appropriate mapping function and optionally horizontal gradients or post-fit residuals with or without multipath cor- rection 3.2. General Quality of Tomographic Reconstructions Quality evaluation of GNSS tomographic results is possible to perform on the basis of com- parison with vertical profiles from meteorological radiosondes or NWP models. Regardless of used ground network of GNSS receivers or applied tomographic solution so far, results are confronted with two common trends: • during standard atmospheric conditions, when the water vapour descends approximately exponentially with the elevation, the GNSS tomography achieves very good results comparable with radiosondes measurements, • in cases of atmospheric inversion or with distinctively different vertical development of water vapour from the standard atmospheric conditions, the tomography can represent this trend only with strong difficulties. This situation is partially given by the necessity of adding constraints in the tomographic system defining mutual relations between neighbouring voxels. If the values between voxels are more smoothed, it is impossible to catch such an inverse trend. On the other way, if the Geoinformatics FCE CTU 9, 2012 67 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method weights in the tomographic system are set more freely, the system becomes less stable and is more prone to provide unreliable results. 4. Theoretical concept of new tomographic method As already mentioned in the introduction, the main idea of new tomography method is to use only a limited number of GNSS receivers situated in a line and reconstruct vertical 2D water vapour distribution. The most signals are crossing throw the 2D voxels situated above the centre receiver and therefore the quality of tomographic reconstruction in these voxels should be the best. Therefore the primary output from the method would be a vertical profile for a receiver situated in the middle of the line. Vertical profile is an output produced also by meteorological radiosondes, which is perfectly suitable for NWP assimilation requirements. The orientation of used line of receivers in a terrain is very important for the new tomographic method. It directly influences the number of visible satellites in a thin atmospheric strip demarcated parallely with the line of receivers. It is possible to use standard applications for planning GNSS measurements to propose a suitable orientation of the line and demonstrate an affordable state. Fig. 1 shows a GPS satellites constellation for an interval 14.30 – 16.00 UTC for selected day at the reference sta- tion VSBO. This graphical output was acquired from Topcon Occupation Planning tool. The orientation of notional line was chosen in the east-west direction and simulation of obstruc- tions was provided. Thus, only parts of the sky with azimuthal angles between 88◦ – 92◦ and 268◦ - 272◦ remained visible. Six various GPS satellites were accessible in the narrow strip and selected interval. Similar, but mostly worse situations held for different parts of a day or different line orientations. Therefore proposed method would be generally able to produce results only during a few intervals in a day depending on the actual satellite constellation. Also radiosondes measurements are limited with such a phenomenon while being launched each six or twelve hours. Besides the maximum number of visible satellites, very important is their elevation. The lowest atmospheric layer above the ground contains the most of water vapour and therefore the signals at low elevations are of the highest interest for GNSS tomography. While choosing an ideal line orientation and time interval for the reconstruction, it is important to try to achieve a state with the highest number of visible satellites distributed on various elevations including the lowest ones. 5. Terrain measurements, data selection On August 2nd 2012 an 8-hour long testing measurement was done at six GNSS receivers situated almost in a line with a northwest-southeast direction. Five of those receivers were GNSS reference stations (TVID, BISK, TKRN, VSBO, CFRM) while the sixth receiver (Top- con Hiper GD) was stabilized for the measurement’s time only. All receivers worked in a GPS and GLONASS mode. GLONASS measurements were included to raise a number of satellites and increase an accessible reconstruction quality. Total length of the line was 113 km and distances between two adjacent stations ranged from 21 to 25 km. The line orientation wasn´t ideal but came out from the accessibility of GNSS measurements from selected reference stations and overall feasibility to realize such Geoinformatics FCE CTU 9, 2012 68 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method Figure 1: GPS satellites constellation, GNSS reference station VSBO, time interval April 11th 2012 14.30 - 16.00 UTC, dark grey colour represents notional obstructions an experiment. Generally, an ideal orientation shouldn´t deviate a lot from the east-west direction. For processing GNSS measurements Bernese GPS Software (Dach et al., 2007) was used. Accurate coordinates of temporarily stabilized receiver were determined using a Precise Point Positioning technique (Zumberge, 1997). Consequently a double differencing technique was used for ZTD processing using settings summarized in a Table 1 and a Europe-wide network containing approximately 30 reference stations. For ZHD determination a Saastamoinen model and measured values of meteorological pa- rameters were used. Slant wet delays were computed from ZWD estimated in 30-second intervals from GPS and GLONASS observations. Horizontal gradients or post-fit residuals weren´t used for SWD computation. Values of atmospheric air pressure and temperature were measured only at the station VSBO and in a place of the temporarily stabilized receiver. For locations of rest reference stations the values of atmospheric air pressure were re-counted using the Babinet's formula given in (5) and values of pressure and temperature measured at the location of the nearest receiver or at a close synoptic meteorological station operated by the Czech Hydrometeorological Institute. Geoinformatics FCE CTU 9, 2012 69 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method z = 16 000 × (1 + t/273) × ((p1 −p2)/(p1 + p2)) (5) where z difference between elevations of point 1 and 2 [m] t atmospheric air temperature measured at point 1 or 2 [◦C] p1 atmospheric air pressure at point 1 [hpa] p2 atmospheric air pressure at point 2 [hpa] After SWD computations while horizontal gradients or post-fit residuals had not been added the observations from a narrow strip in a vertical direction above the line of receivers were selected. It was done in a geometrical way. A plane was held through a centre of the Earth and two ending points of the line. Plane’s equation was computed from the coordinates of these three points in a geocentric coordinate system ITRF. A test whether a distance between a particular satellite and the given plane is in a concrete time less than 1 160 km for GPS satellite and less than 1 113 km for the GLONASS satellite was held. Values defining those distances were obtained by delimiting an angle of 2.5◦ on a distance of 26 560 km which represents a distance between the centre of the Earth and a GPS satellite on its orbit. For GLONASS satellites the distance of 25 510 km was used. List of satellites satisfying given condition and therefore located in a defined strip was saved into a text file. Satellites from this file were consequently confronted with the real data from the RINEX observation files to exclude satellites satisfying given condition but being located on the opposite site of the Earth during the testing terrain campaign. Positions of all satellites in ITRF coordinate system in 15-minute interval were acquired from files with precise ephemerides used previously for ZTD processing. Finally, an interval between 18.45 and 19.45 UTC was chosen to test the tomographic recon- struction. During this selected period one GPS satellite and five various GLONASS satellites were observed at elevation angles between 7◦ and 168◦. Generally it would be better to select a constellation with more GPS satellites because there are still some difficulties with process- ing GLONASS measurements (mainly it is not possible to fix the ambiguities of GLONASS Ephemerides, Satellite clocks CODE Rapid Sampling rate 30 s Elevation cut-off angle 3◦ Mapping function Niell Phase Centre Correction IGS model applied Ocean Loading Applied Observables Double differences ZTD, gradient estimation interval 30 min Pole information CODE rapid Differential Code Bias information CODE 30-day solution Table 1: Basic Characteristics of ZTD Processing in Bernese GPS SW Geoinformatics FCE CTU 9, 2012 70 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method measurements and the quality of some input products like satellites ephemerides and clocks is still lower for GLONASS than for GPS). SWD values for GPS satellites observations should be therefore on higher quality level. Unfortunately according to measured data it was not possible to select any hour interval where more than one GPS satellite had been observed. 6. Mathematical background of the reconstruction Tomographic reconstruction is primarily a mathematical task. In our case certain numerical tools were used to compute 2D tomographic grid. Beside slant wet delays the ZWD values for particular receivers were used as input. First of all the discretize of tomographic system was done by dividing atmosphere above the stations into rectangular network where vertical lines are placed between the stations and horizontal lines are placed at desired height levels. Coordinates of the receivers for this purpose were used and also the vertical restriction in the height of 10 km where amount of water vapour was considered to be null. Signals that left the tomographic system below this height of 10 km were excluded from the reconstruction. Generally such signals can be included in the reconstruction only if it is possible to determine impact of water vapour on the part of the signal beyond the tomographic system. This operation requires an external source of information about the vertical water vapour distribution. Because of relatively high distance between two adjacent stations in the presented solution there were only a few signals which had to be excluded due to this problem. All of them were measured on low elevation angles on receivers situated by the line ends. From the discretization and slant SWD measurements, a set of linear equations defined by equation (3) can be created. This gives a matrix A =   a11 · · · a1ls ... ... ... ak1 · · · akls   , (6) where k is number of slant SWD measurements, l is number of horizontal layers and s is number of measuring stations. Value of each element of the matrix is determined by length of trajectory of slant signal passing through corresponding rectangle. Rectangles are num- bered from the bottom and the left (ak1 correspond to bottom leftmost rectangle, ak(s+1) corresponds to the leftmost rectangle in the second row and so on). Right side of this system is created from the SWD of the corresponding measurement m =   SWD1 ... SWDk   . (7) This gives us a system of equations from (3) from which we can compute amount of water vapours in each rectangle. Graphical illustration about the tomographic discretization and input data presents Fig. 2. Geoinformatics FCE CTU 9, 2012 71 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method However, in the computations there are more slant measurements than the number of rect- angles in the discretization. This creates a problem of over determined system of equations. Such a system has no exact solution. Usual method for solving over stated systems of linear equations is the method of least squares (Björck, 1996). The method of least squares solves the problem Ax ≈ b, A ∈ Rm×n, x ∈ Rm, b ∈ Rn, (8) where b is not in the field matrix values of A, i.e. there is no such a vector x that solves the equation. The method of least squares tries to find an alternative solution xLS, which minimizes min x∈Rm ||Ax − b||2. (9) It can be rewritten equivalently with vector e, which modifies the right side, as min e∈Rn ||e||2 (10) and where (b + e) is in the field matrix values of A. This transforms system of linear equation into task of finding elements of vector e that are minimal in the Euclidean norm and satisfy previous condition. It can be done by various means such as QR or SVD decompositions and by iterative methods. In our solutions we used LSQR iterative method. After using the method of least squares we obtain approximate solution of linear equations system (3). Figure 2: Basic graphical illustration about the tomographic system 7. Results Using computed SWD values and specified system the tomographic reconstruction was done. However its results stated in Table 2 can´t represent a real water vapour distribution. The Geoinformatics FCE CTU 9, 2012 72 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method table shows ZWD values reconstructed for particular cells. Their horizontal centres were situated to locations of GNSS receivers and used tomographic system was divided only into three layers vertically (thickness of two lower vertical layers was 2.5 km and 5 km for the highest layer). Generally tomographic projects use about ten vertical layers to represent the vertical distribution of water vapour in the troposphere. The classification into only three vertical layers would be insufficient for real meteorological applications and was chosen just to test the potential of the new tomographic method. If the results were promising a grid with much higher vertical resolution would be tested consequently. Because during standard atmospheric conditions approximately 50 % of water vapour is sit- uated in the surface layer of the atmosphere below 1.5 km from the ground and above 5 km there is only 5 % of its amount (Seidel, 2002), reconstructed ZWD with the highest values in the top layer can’t represent a real state of the atmosphere prevailing during the test measurements. 10.5 10.5 11 11.5 11 9.5 3.25 3.25 1.25 0.75 2.25 5.5 3.75 -0.75 5 5.25 2.75 0 POINT 1 POINT 2 POINT 3 POINT 4 POINT 5 POINT 6 Table 2: Tomographic reconstruction results, ZWD values in [cm] for particular cells. The bottom layer is in the lowest row of the table. But the conditions of convergence for successful run of the least squares method were not fulfilled what denotes an insufficient amount of input SWD data. Whole situation indicates that two contemporary accessible Global navigation satellite systems GPS and GLONASS are not able to provide a sufficient number of signals for a possibly meaningful usage of proposed tomographic method. 8. Future insights Potential improvement in described problem with low number of GNSS signals can be ex- pected after start of Galileo and Compass systems when the number of satellites in space will be practically doubled. Except the increase in an absolute number of GNSS signals it should also theoretically bring more time intervals during the day, when higher number of satellites will be situated in a narrow strip above the line of receivers. More detailed work at optimization and broaden of the mathematical apparatus used for the first testing purposes would also help. The usage of iterative component trying to find the most ideal reconstruction result based on the input data and a given initial state of the atmosphere would possibly help and should be implemented in a near future. The initial state will be derived from the standard atmospheric conditions when the water vapour decreases exponentially with the increasing height. Some tomographic projects use external source of data about the water vapour (e.g. Numerical Weather Prediction models, radiosonde profiles) to support the GNSS tomographic reconstruction but we do not consider doing that in future. It would bring new problems with the necessity of assimilation such a data in the solution and also made the suggested method dependent on external sources. Geoinformatics FCE CTU 9, 2012 73 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method 9. Conclusion GNSS tomography of the atmosphere is a modern method allowing spatial and temporal study of the water vapour distribution in the atmosphere. So far presented projects were focused on three-dimensional reconstructions above large networks of GNSS receivers. A new tomographic method based on measurements from a line of a few receivers and reconstructing a 2D vertical grid above the line was proposed and described in the paper. Unfortunately the results based on described terrain measurements and mathematic solution are not too much promising. Situation indicates that GPS and GLONASS systems can´t provide sufficient number of observations for this method. Authors expect a possible improvement after an optimization of mathematic apparatus, considering further constraints or applying a priori data and using Galileo and Compass systems in future. If it is proven that new method can provide credible information about vertical distribution of water vapour it could became an effective source of data for NWP models and possibly reduce the necessity of meteorological radiosondes. Acknowledgment This work was supported by the European Regional Development Fund in the IT4Innovations Centre of Excellence project (CZ.1.05/1.1.00/02.0070). Data acquisition of the GNSS continuous operating reference station VSBO is supported by the project of large research infrastructure CzechGeo/EPOS, Grant No.LM2010008. References [1] Bender, M., Dick, G., Ge, M., Deng, Z., Wickert, J., Kahle, H.-G., Raabe, A., Tetzlaff, G. Development of a GNSS water vapour tomography system using algebraic reconstruction techniques, Advances in Space Research, Vol. 47, Issue 10, pp. 1704-1720, 2011 [2] Bennitt, G. Use of Ground based GNSS data in NWP at UK Met Office, E-GVAP workshop, Copenhagen, Denmark, 6th November 2008 [3] Bevis, M., Businger, S., Herring, T.A., Rocken, C., Anthes, R.A., Ware, R.H. GPS meteorology – remote-sensing of atmospheric water-vapor using the global positioning system, Journal of Geohysical Research - Atmospheres, Vol. 97, Issue D14, pp. 15787- 15801, 1992 [4] Björck, Å. Numerical methods for least squares problems. Society for Industrial and Applied Mathematics, USA, ISBN 0-89871-360-9, 1996 [5] Champollion, C., Masson, F., Bouin, M.-N., Walpersdorf, A., Doerflinger, E., Bock, O., Van Baelen, J. GPS water vapour tomography: preliminary results from the ESCOMPTE field experiment, Atmospheric research, vol. 74, p. 253-274, 2004 [6] Dach, R., Hugentobler, U., Fridez, P., Meindl, M. GPS Bernese Software, Version 5.0, Astronomical Institute, University of Berne, Berne, 2007 [7] Duan, J., Bevis, M., Fang, P., Bock, Y., Chiswell, S., Businger, S., Rocken, C., Solheim, F., van Hove, T., Ware, R., McClusky, S., Herring, T., King, R. GPS meteorology: Direct Geoinformatics FCE CTU 9, 2012 74 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method estimation of the absolute value of precipitable water. Journal Appl. M., 24(24), 830–838, 1996 [8] Flores, A., Rius, A., Vilá-Guearou, J., Escudero, A. Spatio-Temporal Tomography of the Lower Troposphere Using GPS Signals, Phys. Chem. Earth (A), Vol. 26, No. 6-8, pp. 405-411, 2001 [9] Gradinarsky, L. P., Jarlemark, P. Ground-Based GPS Tomography of Water Vapor: Analysis of Simulated and Real Data, Journal of the Meteorological Society of Japan, Vol. 82, No. 1B, pp. 551-560, 2004 [10] Guerova, G. Application of GPS derived water vapour for Numerical Weather Prediction in Switzerland, PhD thesis, University of Bern, 2003 [11] Kačmařík, M., Douša, J., Zapletal, J. Comparison of GPS slant wet delays acquired by different techniques, Acta Geodynamica et Geomaterialia, v. 9, No. 4(168), 2012, in print [12] Kačmařík, M. Monitoring of precipitable water vapour by GPS under extreme weather conditions, Advances in Geoinformation Technologies 2010, VŠB – Technical University of Ostrava, pp. 151-161. ISBN 978-80-248-2357-7, 2010 [13] Koizumi, K., Sato, Y. Impact of GPS and TMI Precipitable Water Data on Mesoscale Numerical Weather Prediction Model Forecasts, Journal of the Meteorological Society of Japan, Vol. 82, No. 1B, pp. 453-457, 2004 [14] Liou Y.-A., Huang Ch.-Y. GPS observations of PW during the passage of a typhoon, Earth, planets and space, vol. 52, pp. 709-712, 2000 [15] Moll, P., Poli, P., Ducrocq, V. Use of ground based GNSS data in NWP at Météo-France, E-GVAP workshop, Copenhagen, Denmark, 6th November 2008 [16] Miidla, P., Rannat, K., Uba, P. Simulated studies of water vapour tomography, WSEAS Transactions on Environment and Development, Issue 3, vol. 4, March 2008 [17] Nakamura, H., Koizumi, K., Mannoji, N. Data Assimilation of GPS Precipitable Water Vapor into the JMA Mesoscale Numerical Weather Prediction Model and its Impact on Rainfall Forecasts, Journal of the Meteorological Society of Japan, Vol. 82, No. 1B, pp. 441-452, 2004 [18] Niell, A. E. Global Mapping Functions for the Atmosphere Delay at Radio Wavelengths. Journal of Geophysical Research, Vol. 101, B2, pp. 3227-3246, 1996 [19] Nilsson, T., Gradinarsky, L., Elgered, G. Assessment of tomographic methods for esti- mation of atmospheric water vapour using ground-based GPS. Chalmers University of Technology, Göteborg, Sweeden, 2005 [20] Noguchi, W., Yoshihara, T., Tsuda, T., Hirahara, K. Time-Height Distribution of Water Vapor Derived by Moving Cell Tomography During Tsukuba GPS Campaigns, Journal of the Meteorological Society of Japan, Vol. 82, No. 1B, pp. 561-568, 2004 [21] Notarpietro, R., Cucca, M. Gabella, M., Venuti, G., Perona, G. Tomographic reconstruc- tion of wet and total refractivity fields from GNSS receiver networks, Advances in Space Research, Vol. 47, pp. 898–912, 2011 Geoinformatics FCE CTU 9, 2012 75 Kačmařík, M., Rapant, L.: New GNSS tomography of the atmosphere method [22] Seidel, D. J. Water Vapor: Distribution and Trends. The Earth System: Physical and Chemical Dimensions of Global Environmental Change, John Wiley & Sons, Ltd, 2002 [23] Shoji, Y., Kunii, M., Saito, K. Assimilation of Nationwide and Global GPS PWV Data for a Heavy Rain Event on 28 July 2008 in Hokuriku and Kinki, Japan, Scientific Online Letters on the Atmosphere, Vol. 5, pp. 45-48, 2009 [24] Skone, S.H., Shrestha, S.M.: Strategies for 4-D Regional Modeling of Water Vapour Using GPS, Wuhan University Journal of Natural Sciences, Vol. 8, No. 2B, 627-635, 2003 [25] Song, D.-S., Grejner-Brzeninska, D. Remote sensing of atmospheric water vapor varia- tion from GPS measurements during a severe weather event, Earth, planets and space, vol. 61, n. 10, pp. 1117-1125, 2009 [26] Troller, M. GPS based Determination of the Integrated and Spatially Distributed Water Vapor in the Troposphere, PhD thesis, Swiss Federal Institute of Technology, Zurich, Switzerland, 2004 [27] Vedel, H., Huang, X. Impact of Ground Based GPS Data on Numerical Weather Pre- diction, Journal of the Meteorological Society of Japan, Vol. 82, No. 1B, pp. 459-472, 2004 [28] Zhengdong, B. Near-Real Time GPS Sensing of Atmospheric Water Vapour, PhD thesis, Queensland University of Technology, Australia, 2004 [29] Zumberge, J. F., Heflin, M. B., Jefferson, D. C., Watkins, M. M., Webb, F. H. Pre- cise Point Positioning for the Efficient and Robust Analysis of GPS Data from Large Networks. Journal of Geophysical Research, 102(B3):5005–5017, 1997 Geoinformatics FCE CTU 9, 2012 76