WO2021094753A1 - Configuration of global satellite navigation system infrastructure - Google Patents

Configuration of global satellite navigation system infrastructure Download PDF

Info

Publication number
WO2021094753A1
WO2021094753A1 PCT/GB2020/052874 GB2020052874W WO2021094753A1 WO 2021094753 A1 WO2021094753 A1 WO 2021094753A1 GB 2020052874 W GB2020052874 W GB 2020052874W WO 2021094753 A1 WO2021094753 A1 WO 2021094753A1
Authority
WO
WIPO (PCT)
Prior art keywords
receivers
water vapour
sub
pwv
region
Prior art date
Application number
PCT/GB2020/052874
Other languages
French (fr)
Inventor
Chen Yu
Nigel PENNA
Zhenhong Li
Original Assignee
University Of Newcastle Upon Tyne
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by University Of Newcastle Upon Tyne filed Critical University Of Newcastle Upon Tyne
Publication of WO2021094753A1 publication Critical patent/WO2021094753A1/en

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/01—Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/03—Cooperating elements; Interaction or communication between different cooperating elements or between cooperating elements and receivers
    • G01S19/07—Cooperating elements; Interaction or communication between different cooperating elements or between cooperating elements and receivers providing data for correcting measured positioning data, e.g. DGPS [differential GPS] or ionosphere corrections
    • G01S19/073—Cooperating elements; Interaction or communication between different cooperating elements or between cooperating elements and receivers providing data for correcting measured positioning data, e.g. DGPS [differential GPS] or ionosphere corrections involving a network of fixed stations

Definitions

  • Certain examples of the present disclosure provide a technique for configuring (e.g. refining and/or optimizing) Global Satellite Navigation System (GNSS) infrastructure, and in particular for configuring a network of GNSS receivers (or stations) with regards to tropospheric effects.
  • GNSS Global Satellite Navigation System
  • the present disclosure provides a technique for configuring GNSS infrastructure for use in geophysical, meteorological and/or engineering applications.
  • GNSS Global Navigation Satellite Systems
  • river channel monitoring e.g., Brasington, J., Rumsby, B. T. & McVey, R. A. Monitoring and modelling morphological change in a braided gravel-bed river using high resolution GPS- based survey. Earth Surface Processes and Landforms 25, 973-990 (2000)),
  • Network RTK may be applied worldwide, used for all the abovementioned RTK applications, and also for developing technologies which require seamless, reliable, instantaneous positioning with centimetre-level precison and accuracy, such as intelligent transport systems (e.g., Meng, X., Roberts, S., Cui, Y, Gao, Y, Chen, Q., Xu, C., He, Q., Sharpies, S. & Bhatia, P. Required navigation performance for connected and autonomous vehicles: Where are we now and where are we going? Transportation Planning and Technology 41, 104-118 (2016)), autonomous vehicles and the Internet of Things.
  • intelligent transport systems e.g., Meng, X., Roberts, S., Cui, Y, Gao, Y, Chen, Q., Xu, C., He, Q., Sharpies, S. & Bhatia, P. Required navigation performance for connected and autonomous vehicles: Where are we now and where are we going? Transportation Planning and Technology 41, 104-118 (2016)), autonomous vehicles
  • Network RTK works on the premise that the application of the transmitted, interpolated errors from the CORS enable the same centimetre-level user positional quality to be obtained as from local reference station RTK positioning.
  • the established CORS are optimally spaced with regards to the ensuing positional accuracy at the user’s location, and there is a dearth of information on the optimal CORS configuration needed and whether dominant error sources such as atmospheric delays may be sufficiently interpolated and mitigated with existing CORS infrastructure.
  • a method for configuring a Global Satellite Navigation System (GNSS) infrastructure comprising a plurality of GNSS receivers.
  • the method comprises: obtaining a water vapour map comprising water vapour values at different points; and performing a configuration (e.g. refinement) process comprising: defining one or more sub-regions based on an initial set of receivers, such that one or more receivers are located on (e.g.
  • GNSS Global Satellite Navigation System
  • each sub-region obtaining, from the water vapour map, a first water vapour value for each of a plurality of points on the sub-region; determining a second water vapour value for each of the plurality of points by using interpolation of water vapour values, obtained from the water vapour map, corresponding to positions of two or more receivers; determining a statistical value based on the first and second water vapour values for each point; and if the statistical value satisfies a predetermined criterion, increasing the number of receivers on the sub-region to obtain a modified set of receivers.
  • an apparatus for configuring a Global Satellite Navigation System (GNSS) infrastructure comprising a plurality of GNSS receivers, the apparatus comprising a processor configured to: obtain a water vapour map comprising water vapour values at different points; and perform a configuration (e.g. refinement) process comprising: defining one or more sub-regions based on an initial set of receivers, such that one or more receivers are located on (e.g.
  • GNSS Global Satellite Navigation System
  • each sub-region obtaining, from the water vapour map, a first water vapour value for each of a plurality of points on the sub-region; determining a second water vapour value for each of the plurality of points by using interpolation of water vapour values, obtained from the water vapour map, corresponding to positions of two or more receivers; determining a statistical value based on the first and second water vapour values for each point; and if the statistical value satisfies a predetermined criterion, increasing the number of receivers on the sub-region to obtain a modified set of receivers.
  • the GNSS infrastructure may be configured with regards to the mitigation of tropospheric effects.
  • the two or more receivers used for the interpolation comprise one or more receivers on the sub-region and/or one or more receivers on one or more neighbouring sub- regions.
  • the method further comprises, or the processor is further configured for, repeating the configuration process such that the modified set of receivers obtained at the end of one iteration is used as the initial set of receivers in the next iteration.
  • the configuration process is repeated until a termination criterion is satisfied.
  • the termination criterion is based on one or more of: the statistical value does not satisfy the predetermined criterion for every sub-region; the receiver spacing is less than a certain threshold; and the number of iterations has reached a certain threshold.
  • the water vapour map is obtained from Moderate Resolution Imaging Spectroradiometer (MODIS) data.
  • MODIS Moderate Resolution Imaging Spectroradiometer
  • the interpolation is performed based on an Iterative Tropospheric Decomposition (ITD) model.
  • ITD Iterative Tropospheric Decomposition
  • the initial set of receivers comprises a new or pre-existing regular or irregular array of receivers.
  • a sub-region comprises a polygon, and a receiver is located on one or more vertices of the polygon, on one or more edges of the polygon, and/or at the centre of the polygon.
  • the polygon comprises one or more of: square; rectangle; triangle; and hexagon.
  • the sub-regions comprise a regular array of rectangular or square sub- regions, or a Triangulated Irregular Network (TIN) obtained using Delaunay triangulation.
  • TIN Triangulated Irregular Network
  • the plurality of points on the sub-region have the same spatial resolution as the water vapour map.
  • the statistical value is the root mean square (RMS) of the differences between the first and second water vapour values.
  • the predetermined criterion is that the statistical value is greater than a predetermined threshold.
  • increasing the number of receivers on the sub-region comprises adding a receiver to the geometric centre of the sub-region.
  • the method further comprises, or the processor is further configured for, generating a report including position, spacing, and/or density information of the modified set of receivers.
  • the method further comprises, or the processor is further configured for, configuring the GNSS infrastructure according to the modified set of receivers.
  • GNSS Global Satellite Navigation System
  • Figure 1 illustrates the simulated global GNSS network distribution optimization for boreal summer (June - Aug 2016) and boreal winter (Dec 2015 - Feb 2016);
  • Figure 2 illustrates the relationship between the precipitable water vapour (PWV) and station spacing (blue line), and between the PWV variation and station spacing (red line), for local summer (a) and local winter (b), respectively.
  • the PWV or PWV variation is a global average of all pixels at the corresponding spacing, and is normalized linearly between [0, 1];
  • Figure 3 illustrates original and densified (optimally refined) regional networks during local summer of 2016: July 2016 for the British Isles and Italy, January 2016 for New Zealand and Java;
  • Figure 4 illustrates time series of optimal station spacing from January 2012 to December 2016 for each of the British Isles, Italy, New Zealand and Java networks
  • Figure 5 is a flow diagram of a global optimization process according to an example of the present disclosure
  • Figure 6 is a flow diagram of a local optimization process according to an example of the present disclosure.
  • FIG. 7 is a flow diagram of an Iterative Tropospheric Decomposition (ITD) model according to an example of the present disclosure
  • Figure 8 illustrates the evolution of a GNSS network during operation of the techniques disclosed herein in the case of global optimization and local optimization.
  • Figures 9-18 disclose Figures 1-10 of the Yu et al. (2017) paper referred to below.
  • examples of the present disclosure can be realized in the form of hardware, software or any combination of hardware and software. Any such software may be stored in any suitable form of volatile or non-volatile storage device or medium, for example a ROM, RAM, memory chip, integrated circuit, or an optically or magnetically readable medium (e.g. CD, DVD, magnetic disk or magnetic tape).
  • Certain embodiments of the present disclosure may provide a computer program comprising instructions or code which, when executed, implement a method, system and/or apparatus in accordance with any aspect, claim, example and/or embodiment disclosed herein. Certain embodiments of the present disclosure provide a machine-readable storage storing such a program.
  • Such an apparatus and/or system may be configured to perform a method according to any aspect, embodiment, example or claim disclosed herein.
  • Such an apparatus may comprise one or more elements, for example one or more of receivers, transmitters, transceivers, processors, controllers, modules, units, and the like, each element configured to perform one or more corresponding processes, operations and/or method steps for implementing the techniques described herein.
  • an operation/function of X may be performed by a module configured to perform X (or an X-module).
  • the one or more elements may be implemented in the form of hardware, software, or any combination of hardware and software.
  • GNSS Global Navigation Satellite Systems
  • CORS network real-time kinematic
  • existing CORS networks have often been established simply by attempting regular spacing or in clusters around cities, with little consideration of weather, climate and topography effects, which represent a substantial GNSS positional error and prevent homogeneous Network RTK accuracy attainment.
  • Atmospheric water vapour may also be expressed as precipitable water vapour (PWV) and it is a highly variable quantity both spatially and temporally (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 Geophysical Research 97, 15787-15801 (1992)).
  • PWV precipitable water vapour
  • the optimal CORS spacing for PWV interpolation may be considered a proxy for Network RTK CORS infrastructure design with regards to the mitigation of GNSS tropospheric effects.
  • Certain examples of the present disclosure provide centimetre- 1 eve I homogeneous positioning in real-time with GNSS network infrastructure.
  • Certain examples of the present disclosure provide a technique for optimizing the design of CORS ground infrastructure (e.g. optimizing CORS spacing for PWV interpolation corresponding to optimizing Network RTK CORS infrastructure design) to enable consistent user centimetre-level positioning globally with regards to tropospheric (PWV) effects.
  • PWV tropospheric
  • it is shown that denser networks are typically required in the tropics and in mountainous areas.
  • Certain examples of the present disclosure also provide a technique for improving existing national GNSS network infrastructure.
  • configuring a network of GNSS receivers may be understood as meaning determining, defining, designing and/or creating an arrangement of receivers.
  • the arrangement of receivers may be defined, for example, by specifying positions of individual receivers, by specifying the spacing between receivers and/or by specifying the density of receivers in a certain area. If there is an existing network, then configuring the network may be understood to include modifying, refining and/or optimizing the existing network. On the other hand, if no network exists, then configuring a network may be understood to include creating or defining a new network that may be modified, refined and/or optimized. The skilled person will understand that references herein to “optimizing” or “refining” may be taken to also include “configuring” according to the above definition.
  • the water vapour data is obtained from Moderate Resolution Imaging Spectroradiometer (MODIS) data, although the skilled person will appreciate that the water vapour data may be obtained from any additional or alternative source.
  • MODIS comprises two satellites providing high radiometric sensitivity (12 bit) in 36 spectral bands ranging in wavelength from 0.4 pm to 14.4 pm.
  • the present example uses the near-infrared band water vapour product which applies a near-infrared algorithm over clear land areas of the globe and above clouds over both land and ocean.
  • any other suitable band may be used.
  • Global water vapour reference grid is obtained from Moderate Resolution Imaging Spectroradiometer
  • the present example method uses the near-infrared band water vapour MODIS product over cloud-free land areas of the globe (Gao, B.-C. & Kaufman, Y. J. Water vapor retrievals using Moderate Resolution Imaging Spectroradiometer (MODIS) near-infrared channels. Journal of Geophysical Research 108 (D13), 4389 (2003), hereinafter Gao and Kaufman (2003)).
  • MODIS Moderate Resolution Imaging Spectroradiometer
  • this MODIS product provides a PWV measurement accuracy and precision of about 1-2 mm, based on comparisons with GPS and radiosondes (e.g. Li, Z., Muller, J. P. & Cross, P.
  • the PWV data may be generated at about 1-km spatial resolution, with 1 to 2 days required for satellite measurements of the whole Earth’s surface. If the ⁇ 1-km data are not uniformly spaced, the present method may sample them to a regular 0.05 by 0.05-degree grid ( ⁇ 5-km spacing), taking as the grid node value the nearest original (1-km) data point. To reduce the impact of meridional convergence on the uniformity of the grid point spacing, the present method may only use in tests grid points from the equator to latitude 60 degrees north/south, at which 0.05-degree grid spacing corresponds to ⁇ 5.6 km and ⁇ 2.8 km, respectively.
  • Figure 1 The resulting latitude 60 degrees north to south PWV distributions over land are shown in Figure 1 for 3-month averages for boreal summer (1 June 2016 to 31 August 2016) and for boreal winter (1 December 2015 to 29 February 2016).
  • Figure 1(a) illustrates average PWV content from MODIS (0.05 degree resolution);
  • Figure 1(b) illustrates PWV variations, defined as the PWV standard deviation of all points within a 0.8x0.8 degree window centered about each 0.05 degree MODIS pixel (0.05 degree resolution);
  • Figure 1(c) illustrates optimal station spacing (0.8 degree resolution), including percentage of grid cells with different spacing;
  • Figure 1(d) illustrates the correlation map between the optimal station spacing and PWV variation (0.8 degree resolution). Topography variations are denoted by the shaded background.
  • This averaging per grid point may be done to reduce the effect of weather variations during each season, but also to minimize the effect of cloud-induced data availability gaps.
  • PWV measurements over water bodies may be excluded using the Shuttle Radar Topography Mission 90 m digital elevation model (DEM) (Farr, T. G. et al. The shuttle radar topography mission. Reviews of Geophysics 45, RG2004 (2007)) and a DEM corresponding to the 0.05 degree PWV grid may be generated to enable the assessment of topographic effects.
  • DEM digital elevation model
  • the PWV values may be interpolated, for example using the ITD model developed by Yu et al. (2017).
  • the contents of Yu et al. (2017) are disclosed in the annex to the description of the present application, in combination with Figures 9-18 of the present application which correspond to Figures 1-10 of Yu et al. (2017).
  • the disclosure of Yu et al. (2017) in the annex to the description and in Figures 9-18 forms part of the disclosure of the present application.
  • the PWV is assumed to comprise a stratified (topography/elevation dependent) component, and a turbulent component, representing topography-independent PWV signals (Equation 1):
  • PWV i S(H i ) + T( i ) + e i (1)
  • S(H j ) L 0 e ⁇ ⁇ Hi which represents the stratified component at position vector x ⁇ ; i is the location considered; H, is the height above sea level; U is the stratified PWV at sea level; b is an exponential coefficient; T is the turbulent PWV component; e, is the residual PWV.
  • the residuals may then be interpolated, for example using inverse distance weighting to obtain an initial estimate of the turbulent component at the location of interest, and repeated to obtain initial values for the turbulent component for all grid points, and subtracted from the initial values of the stratified delay per grid point.
  • the stratified delay may then be re-estimated to obtain updated values for the exponential coefficients and subsequently updated residuals, which may then be interpolated again to obtain updated values for the turbulent component.
  • This iterative process may be repeated until the estimates for the stratified and turbulent PWV components converge, and the final value of the interpolated PWV may be obtained by summing the stratified and turbulent components, plus the residual.
  • six or seven iterations may be needed, for example as demonstrated in Yu et al. (2017), where further details of an exemplary method may be obtained.
  • an interpolated PWV value (PWVo) at a specific pixel (Po) may be determined based on PWV values (PWVi, PWW, PWV3, ...) at a plurality of other pixels (Pi, P2, P3, ).
  • the other pixels (Pi, P 2 , P 3 , ...) may correspond to positions of receivers, for example receivers within a defined range (e.g. 150km) from Po. Any receiver at Po may be referred to as a “roving” receiver (or station) or “rover” receiver (or station).
  • the PVW values (PVWi, PVW 2 , PVW 3 , ...) at the plurality of other pixels (Pi, P 2 , P 3 , ...) may be obtained from the global water vapour reference grid generated from MODIS data.
  • the interpolation may be performed using the ITD model of Yu et al. (2017). Accordingly, in certain examples, once the MODIS data is obtained, the algorithm for configuring a network of receivers may be performed on a computer without requiring any further measurements taken by receivers or other apparatus.
  • the PVW values (PVWi, PVW 2 , PVW 3 , %) at the plurality of other pixels (Pi, P 2 , P 3 , ...) may be obtained from the global water vapour reference grid generated from MODIS data.
  • one or more, or all, of the PVW values (PVWi, PVW 2 , PVW 3 , ...) at the plurality of other pixels (Pi, P 2 , P 3 , ...) may be obtained through measurements made by the receivers at those pixels and/or any other suitable apparatus (e.g. GPS or meteorological apparatus).
  • the optimal network distribution across the land areas of the globe may be evaluated by computing the station spacing needed in each cell of a 0.8x0.8 degree grid. This may be chosen because 0.8 degrees represents about 80 km, and may be assumed to be the coarsest CORS spacing used for existing RTK networks.
  • Each 0.8x0.8 degree cell may be sub-divided into 256 (16x16) pixels which match the 0.05 degree MODIS PVW reference grid, and then PVW values per pixel may be obtained by interpolation from a particular defined station spacing.
  • the RMS of the differences between the interpolated and reference value for each of the 256 pixels per 0.8x0.8 degree cell may be computed. This leads to a set of RMS PVW interpolation errors per cell of a 0.8 degree grid globally.
  • the station spacing used may be considered insufficient for that particular cell and densified until the RMS error for each and every cell is less than 1.5 mm. This may be realized by first densifying such grid cells to use 0.4 degree reference PVWs for interpolation, or a station spacing of ⁇ 40 km. The interpolation to the 256 pixels in the 0.8 degree cell may then be repeated, but using the 0.4 degree grid points within the 0.8 degree cell, plus the 0.8 degree grid points outside it. If the RMS error for the same 256 pixels is still greater than 1.5 mm, the 0.8 degree grid cell may be further densified to a simulated station spacing of 0.2 degrees (20 km), and the interpolation repeated.
  • the station spacing needed per 0.8 degree grid cell may then be deemed the coarsest grid spacing needed for the RMS error to be 1.5 mm or lower.
  • the resulting station spacing map shown in Figure 1 consequently comprises a station spacing per 0.8 degree cell ranging from 80, 40, 20, 10 to 5 km.
  • the standard deviation of the 256 0.05 degree grid PVW values within each 0.8 degree cell may be used to represent the PVW variation, and correlation coefficients between station spacing and PVW variation are shown in Figure 1. These are computed for each 0.8 degree grid cell, using the station spacing and PWV variation for the cell itself and the eight 0.8 degree cells immediately surrounding it.
  • a Triangulated Irregular Network may be formed, for example based on Delaunay triangulation.
  • Another difference may be that 1 -month averages of the MODIS PWV daily values per point of the global 0.05 degree grid may be used, not 3-month averages.
  • the procedure may comprise:
  • Figure 5 is a flow diagram of a global optimization process according to an example of the present disclosure. The skilled person will appreciate that features of the embodiment of Figure 5 may be incorporated into the methods described above, and vice versa.
  • the optimal spacings may be obtained by first thinning the PWV values from the MODIS 5 km global resolution reference data set to 80 km (0.8 degrees) resolution.
  • Each 0.8x0.8 degree cell i.e. sub region
  • Each 0.8x0.8 degree cell may be sub-divided into 256 (16x16) pixels which matches the 0.05 degree MODIS PWV reference grid, and then PWV values per pixel may be obtained by interpolation from a particular defined station spacing, for example using the ITD model of Yu et al. (2017) which considers the DEM and applies the water body mask.
  • the RMS of the differences between the interpolated and reference value for each sub region may be computed. This leads to a set of RMS PWV interpolation errors per sub region globally (i.e.
  • the station spacing used may be considered insufficient for that particular sub-region and will be halved to form a refined global non-uniform network. If the refined network is different from the previous one, the procedure may iterate until the RMS error for each sub-region is lower than 1.5 mm, or the 0.05 degree spacing is reached.
  • a DEM may be a 90 metre resolution global digital elevation model provided by the Shuttle Radar Topography Mission.
  • a water body mask may be a 90 metre resolution grid data provided along with the DEM which masks out the water body, for example including oceans, lakes and inland rivers.
  • Figure 6 is a flow diagram of a local optimization process according to an example of the present disclosure.
  • the local optimization may start from an initial network infrastructure which may be, for example, either an existing GNSS network with an irregular distribution of historic stations, or a newly defined network on some key reference locations.
  • the MODIS 5 km global resolution water vapour map may be used for reference.
  • the procedure may comprise the following operations:
  • Figure 7 is a flow diagram of an ITD model according to an example of the present disclosure.
  • ITD Model The main procedure of ITD is to iteratively estimate the height scaling function and find the optimal scaling parameters. To start with, ITD supposes no turbulent delay exists so the initial scaling parameters can be obtained. The turbulent delay is then updated by the residuals ( T k ) using the inverse distance weighting (IDW) function f, and is removed from the total delay to generate an updated turbulent delay (T k+i ). Several iterations of these steps are made until the difference between the updated turbulent delay and the previous turbulent delay is smaller than the threshold. The final output is generated using the final turbulent signal and scaled by the final height scaling parameter.
  • IDW inverse distance weighting
  • height scaling may refer to the use of an exponential function to simulate the ZTD:
  • IDW Inverse Distance Weighting
  • wi j is determined by the distance between reference points / and j.
  • the interpolation for each station may be performed using available stations that fall within a certain region (e.g. having a certain size, shape and/or position) defined with respect to each station.
  • the interpolation for a station may be performed using all the available stations that fall within a certain shape (e.g. a circle) having a certain size (e.g. 150 km radius) centered at that station.
  • at least three stations may be required for the interpolation.
  • the interpolation may use a relatively large number of stations surrounding the point that is to be interpolated, not just the three vertices.
  • the interpolation may be performed based on stations at the vertices only, for example if the vertices are located far enough away from each other and/or from the station being interpolated, for example over a certain threshold (e.g. over 150 km).
  • a certain threshold e.g. over 150 km.
  • the interpolation process may search for all the surrounding grids/stations within a certain shape having a certain size (e.g. a 150 km radius circle), not just the stations in its own grid (e.g. 80-km grid).
  • the same interpolation procedure may be used in both cases.
  • Figure 8 illustrates the evolution of a GNSS network during operation of the techniques disclosed herein in the case of global optimization and local optimization. As shown, stations are iteratively added to grid elements of the network until the RMS value for a grid element is lower than a certain threshold.
  • stations may be positioned on vertices of the grid cells and/or within (e.g. at the centre of) grid cells.
  • a station may be located at the centre of each square (or cell or region) (i.e. each square represents one station).
  • each cell e.g. 0.8x0.8 degree cell
  • 2x2 sub-cells e.g. 0.4 degree spacing.
  • This may be continued, if necessary, using successively smaller spacing (e.g. 0.2 degree, then 0.1 degree, then 0.05 degree etc.) until the RMS error for the pixels (e.g. 16x16 pixels) is lower than a certain threshold (e.g. 1.5 mm), or until a certain termination criteria is satisfied (e.g. a certain degree spacing, such as 0.05 degree spacing, is reached).
  • a certain degree spacing such as 0.05 degree spacing
  • cells may have any suitable regular or irregular shape (e.g. other than square or triangle, such as a hexagon); the cells may have different shapes and/or sizes; the cell array may have any suitable regular or irregular form (e.g. other than a regular square array or TIN); any suitable statistical value (e.g. other than an RMS value) and any suitable threshold value may be used to determine whether or not to densify a cell; receivers may be located on one or more vertices of a cell, on one or more edges of a cell, and/or in the interior (e.g. at the centre) of a cell; to densify a cell any suitable number of receivers (e.g.
  • one receiver or more than one receiver may be added to a cell in any suitable locations (e.g. within the cell (e.g. at the cell centre), on one or more cell boundaries (e.g. in the middle of a cell edge) and/or on one or more cell vertices) and in a uniform or non-uniform manner;
  • the receivers used for interpolation may comprise one or more receivers of a cell and/or one or more receivers of one or more neighbouring cells;
  • the interpolation may be perfomed in any suitable manner (e.g. other than using the ITD model);
  • the iterative process may be terminated according to any suitable termination criteria (e.g. once a maximum number of iterations has been reached and/or once a certain minimum receive spacing is reached); the water vapour map may be obtained from any suitable source (e.g. other than MODIS data).
  • MODIS data any suitable source
  • the optimal CORS station spacing globally with regards to tropospheric delay mitigation is shown in Figure 1 for both boreal summer and winter, and ranges from a defined, fixed upper limit of 80 km (the typical sparsest spacing between Network RTK CORS in many countries) and then where needed, increasingly dense halvings of 40, 20, 10 to 5 km.
  • the optimal spacings were obtained by first thinning the PWV values from the MODIS 5 km global resolution reference data set to 80 km (0.8 degrees) resolution then, as described in the Methods section above, interpolating using the ITD model of Yu et al. (2017) and cross-validating against the 5 km reference values.
  • This PWV error cut-off was chosen as it is equivalent to a zenith tropospheric delay error of about 1 cm (Bevis M., Businger, S., Chiswell, S., Herring, T, Anthes, R., Rocken, C. & Ware, R. H. GPS meteorology: Mapping zenith wet delays onto precipitable water. Journal of Applied Meteorology 33, 379-386 (1994)), which maps to a GNSS height error of about 2-3 cm (Brunner, F. K. & Welsch, W. M. Effect of the troposphere on GPS measurements.
  • FIG. 1 To understand the cause of the variation in the optimal station spacing, included in Figure 1 are maps of 3-monthly average PWV per 0.05 degree MODIS pixel. It can be seen that many of the areas with the densest station spacing occur where the PWV content is very high (approaching 10 cm), notably in the tropics and across India and south-east Asia in local summer. However, there are also regions which have a high PWV content, such as central Africa and the Amazon basin, yet the station spacing needed to accurately interpolate PWV is very sparse.
  • Figure 1 indicates that a dense station spacing is needed in north-east Canada in local summer, the large PWV variations there are not considered reliable because the MODIS PWV values are degraded by the large number of lakes, whose average size is several kilometres and close to the MODIS spatial resolution (Gao and Kaufman, 2003). In local winter however, both the average PWV content (-0.4 cm) and its variation (-0.08 cm) are low, meaning that this effect can be negligible and the maximum 80 km spacing is always sufficient.
  • the optimal station spacing for example with regards to tropospheric delay errors, globally has been considered by assuming a regular station grid and without considering existing geodetic GNSS infrastructure.
  • CORS stations already existing and, while in some cases they have been designed to form as regular a station spacing as possible (e.g., about 60-80 km in Germany), invariably they are irregularly spaced, and designers do not seem to have considered topographic or climatic variations. Consequently, existing networks do not necessarily provide optimal PWV mitigation performance.
  • Figure 3 illustrates original and densified (optimally refined) regional networks during local summer of 2016: July 2016 for the British Isles and Italy, January 2016 for New Zealand and Java.
  • the first and second columns represent the original and the densified networks, respectively, with the shaded triangles indicating areas where the interpolation RMS exceeds 1.5 mm.
  • the background underneath the TIN of the first two columns shows the interpolation differences between the reference MODIS water vapour map and the interpolated map using the reference stations.
  • the third column shows the topography and the additional stations inserted to obtain the optimal spacing.
  • Figure 4 illustrates time series of optimal station spacing from January 2012 to December 2016 for each of the British Isles, Italy, New Zealand and Java networks. Shown are the resultant monthly average station spacing and the number of stations added to the original network to be optimal. Dots denote the actual station numbers per month per network. Solid lines show the best fitting annual sinusoid per time series.
  • the first example considered is the British Isles (www.bigf.ac.uk), as the bulk of the network comprises national CORS geodetic infrastructure which is used for Network RTK service provision and has a near regular station spacing of about 60-80 km (Figure 3). It exhibits moderate PWV variations and fairly flat topography in England, but with more mountainous terrain in Scotland and Wales, albeit only rising to a maximum elevation of about 1300 m.
  • the third example is New Zealand’s GeoNet (https://www.geonet.orQ.nz) ⁇ which incorporates dense station spacing ( ⁇ 25 km) in the east of the north island. However, spacing is much sparser ( ⁇ 90 km) across the Southern Alps in the south island which exhibit a large range in topography from about 3500 m elevation mountains to plains close to sea level.
  • the fourth example is for western Java, which has very large PWV variations, large topographic variations from sea level to -2000 m elevation and, from Figure 1 , very dense suggested station spacings of towards 5 km being needed. It is also chosen as an area with few existing CORS and as Network RTK services are not apparently operating here, instead of showing how an existing network could best be densified, an initial hypothetical network with a regular 80 km spacing is assumed, and then densified.
  • the MODIS 0.05 degree grid point PWV values were averaged per month and CORS station PWV values were taken as those of the nearest grid point.
  • the PWV values from the CORS locations were then interpolated on to the 5 km PWV MODIS values and the RMS difference (error) with respect to the reference values computed per triangle. If the RMS error in any triangle exceeded 1.5 mm, extra stations were added (see Methods section above) until all triangles across the network had RMS errors of under 1.5 mm.
  • the number of additional stations needing to be inserted per network per month are shown in Figure 4, together with the resultant average station spacing across the network. Clear periodic, seasonal variations in the number of added stations and resultant average station spacing needed for optimal PWV interpolation are seen.
  • the amplitude of the seasonal variation is dependent on the type of climate, with mountainous climates such as Italy exhibiting an amplitude of about 140 stations, while the tropical climate of Java has an amplitude of only 30 stations.
  • the amplitude also provides an indication of how close the original network is to being optimal for Network RTK provision, as a small amplitude denotes that few stations need to be added, although it could also indicate that in local winter the original network is denser than required in places and could potentially be thinned.
  • a framework has been disclosed to optimize GNSS CORS network distributions with regard to the distribution of the variations in PWV and topography, for example in order to help facilitate centimetre-level homogeneous accuracy for Network RTK positioning.
  • the optimal station spacings were found globally using MODIS PWV values as reference, and cross-validating after interpolation from station spacings ranging from 80 km to 5 km. More detailed local examples have also been presented, demonstrating how assumed original existing networks could be densified through a Delaunay triangulation approach in order to optimize the interpolation of PWV errors, with a clear seasonal variation demonstrated in the optimal network.
  • the station spacing has been shown to depend on variations in PWV across the network considered, rather than the total amount of PWV.
  • the largest variations, and hence the densest station spacings tend to be in the tropical regions which have large PWV content, although topographic effects also have a major influence.
  • topographic effects For example, where large variations in topography arise such as in the Alps, a dense station spacing is needed.
  • Local extreme weather also contributes, such as the monsoon climate of India which requires a very dense network in the summer despite being fairly flat, but only a fairly sparse network is needed in winter.
  • the framework has been disclosed for the optimization of Network RTK GNSS positioning with regards to mitigating tropospheric delay errors.
  • the resulting optimal CORS network designs may also be important for repeat-pass Interferometic SAR (InSAR) data processing and applications, providing indications of whether existing GNSS stations will be able to provide sufficient quality corrections to reduce atmospheric effects on InSAR observations and hence improve estimates of post- and inter-seismic crustal movements (e.g. Yu etal., 2018a).
  • the ability to generate precise maps of PWV from GNSS CORS data may also be useful in meteorology in the nowcasting of rainfall (e.g., Bonafoni, S. & Biondi, R.
  • GNSS Global Navigation Satellite Systems
  • Atmospheric Research 167, 15-23 (2016) and assimilation into numerical weather models for precipitation forecasting (e.g., Rohm, W., Yuan, Y, Biadeglgne, B., Zhang, K. & Le Marshall, J. Ground-based GNSS ZTD/IWV estimation system for numerical weather prediction in challenging weather conditions.
  • the framework concept disclosed herein may be employed to optimize other monitoring networks, e.g. those for soil moisture, air pollution, and tree species monitoring, although their own dominating error source should be identified and the corresponding spatial and temporal variations should be considered.
  • the framework concept should be further developed to deal with multiple variables.
  • MODIS data used in the examples disclosed herein were obtained from the Land Processes Distributed Active Archive Centre (LP DAAC), managed by the NASA Earth Science Data and Information System (ESDIS) project.
  • LP DAAC Land Processes Distributed Active Archive Centre
  • ESDIS NASA Earth Science Data and Information System
  • MODIS near-IR water vapour products may be used to generate a global background dataset.
  • the skilled person will appreciate that the present disclosure is not limited to this example.
  • High-resolution European Centre for Medium-Range Weather Forecasts (ECMWF), or any other suitable high-resolution global numerical weather models, can also be used to generate such a background dataset.
  • any suitable high- resolution regional or local numerical weather models e.g. Weather Research and Forecasting model; Unified Model
  • Weather Research and Forecasting model can be used to generate the background dataset.
  • Unified Model any suitable high- resolution regional or local numerical weather models
  • one or more boundaries may be defined.
  • specific areas e.g. countries, provinces, counties, etc.
  • optimization of a GNSS network may be performed on one or more specific defined areas.
  • specific areas e.g. water bodies and special infrastructure
  • specific areas may be masked out while optimizing a GNSS network. For example, in the algorithms described above, calculations corresponding to positions falling within a masked area may be skipped, and insertion of new reference stations inside a masked area may be prohibited.
  • a GNSS network may be optimized for monitoring the stability of linear infrastructure (e.g. a high-speed railway line or long sea-crossing bridge).
  • the stations may have fixed positions. However, in other examples, one or more or all of the stations may be moveable. In addition, in certain examples of the present invention, the stations may be positioned onshore. However, in other examples, one or more or all of the stations may be positioned offshore. For example, the techniques described above may be used to optimize a GNSS buoy network, or to optimize the distribution of GNSS stations on oil platforms. Certain examples may use a combination of onshore and offshore stations, which may help reduce edge effects (interpolation versus extrapolation of errors) that can result from the lack of stations on one side of a boundary or area.
  • stations may be deployed for a relatively long period of time (e.g. they may be permanent stations). However, in other examples, one or more or all of the stations may be deployed for a relatively short period of time (e.g. they may be temporary stations).
  • Temporary stations may be used for campaign design, for example for monitoring of a particular process such as plate tectonics and certain engineering projects.
  • the stations may all be the same. However, in other examples the stations may vary in their quality or performance (e.g. accuracy, sensitivity, etc.), and the construction and/or optimization of a GNSS network may be performed in consideration of the quality and/or performance of the stations.
  • the stations may vary in their quality or performance (e.g. accuracy, sensitivity, etc.), and the construction and/or optimization of a GNSS network may be performed in consideration of the quality and/or performance of the stations.
  • Network performance in some regions may be subject to temporal (e.g. seasonal) variations.
  • temporal e.g. seasonal
  • the optimization of a GNSS network may be adjusted over time to account for such variations in performance. For example, local densification of the network may be performed during certain times of the year.
  • a first level report may show only relatively general statistics
  • a second level report may provide information on the approximate positions or distribution of stations
  • a third level report may provide detailed information giving accurate positions of the stations. For example, rather than specify the exact position of a station, an approximate position of the station may be given.
  • information specifying the station density e.g. the number of stations to be positioned within a certain area
  • Pointwise GPS measurements of tropospheric zenith total delay can be interpolated to provide high resolution water vapor maps which may be used for correcting SAR images, for numeral weather prediction and for correcting Network Real-time Kinematic GPS observations.
  • Several previous studies have addressed the importance of the elevation dependency of water vapor, but it is often a challenge to separate elevation-dependent tropospheric delays from turbulent components.
  • Time-varying two-dimensional precipitable water vapor (PWV) fields are used in meteorological nowcasting, including the identification of events dominated by horizontal advection [Benevides et al., 2015]; for assessing moisture transport in the lower troposphere [e.g., Mengistu Tsidu et al., 2015]; for relating humidity fields to precipitation events [e.g.
  • Such corrections are essential for centimeter level positioning, particularly heights, and enable (subject to sufficient GNSS base station coverage) Network RTK to be used for geophysical and engineering applications that have normally only used local base station RTK, such as river channel mapping [e.g., Brasington et ai, 2003], glacier flow and debris mapping [e.g., De Paoli and Flowers, 2009], coastal erosion [e.g., Lee et ai, 2013], crustal deformation [e.g., Genrich and Bock., 2006] and structural [e.g., Im et ai, 2013] monitoring, precision farming [e.g., Perez-Ruiz et ai, 2011], embankment instability and landslide monitoring [e.g., Gili et ai, 2000; Zabuski et ai., 2015]
  • river channel mapping e.g., Brasington et ai, 2003
  • glacier flow and debris mapping e.g., De Paol
  • Two-dimensional PWV maps/fields may be obtained directly from the National Aeronautics and Space Administration (NASA) Moderate Resolution Imaging Spectroradiometer (MODIS) on board each of the Terra and Aqua satellites (launched in 1999 and 2002, respectively).
  • NASA National Aeronautics and Space Administration
  • MODIS Moderate Resolution Imaging Spectroradiometer
  • MODIS swaths do not coincide in time with images from SAR satellites such as ERS-1/2 (nominally 30 min of time difference), Envisat (nominally 30 min of time difference), or Sentinel-1 (nominally 4.5 hours of time difference), and therefore for SAR atmospheric corrections temporal interpolation errors will always prevail.
  • SAR satellites such as ERS-1/2 (nominally 30 min of time difference), Envisat (nominally 30 min of time difference), or Sentinel-1 (nominally 4.5 hours of time difference), and therefore for SAR atmospheric corrections temporal interpolation errors will always prevail.
  • the Envisat SAR satellite included the Medium Resolution Imaging Spectrometer (MERIS), a passive push-broom imaging instrument measuring the solar radiation reflected from the Earth’s surface and clouds in the visible and near-IR spectral range during the daytime.
  • MERIS Medium Resolution Imaging Spectrometer
  • MERIS provided near-IR water vapor products above land or ocean surfaces under cloud-free conditions [Bennartz and Fischer, 2001] or above the highest cloud level under cloudy conditions [Albert et ai, 2001] at two nominal spatial resolutions: 0.3 km for full-resolution mode and 1.2 km for reduced-resolution mode.
  • Atmospheric Infrared Sounder (AIRS) instrument also on-board the Aqua satellite, provides global coverage amounting to 324,000 PWV profiles with a vertical resolution of 2 k and a spatial resolution at nadir of 45 km, but also only operates in clear skies and in partly cloudy conditions [Aumann et al., 2003] Comparison between AIRS products and GPS PWV showed RMS differences of about 3-4 mm and 15% differences in the 2 km layers between AIRS and radiosondes [Divakarla et al., 2006; Raja et al., 2008] PWV fields may also be obtained from global atmospheric models (both reanalysis and rapid update cycle predictions), but suffer from their coarse spacing, low temporal resolution (e.g., 3 or 6 hours for the European Centre for Medium- Range Weather Forecasts (ECMWF) model) and failures during periods of atmospheric turbulence [Jolivet et al., 2014]
  • EMWF European Centre for Medium- Range Weather Forecasts
  • GPS provides a means of obtaining PWV with continuous temporal resolution without any cloud or weather dependence, albeit at discrete points where the GPS receivers are located.
  • ZTD time- varying tropospheric zenith total delay
  • ZWD zenith wet delay
  • the ZWD may then be readily converted to PWV using estimates of the mean temperature of the atmosphere, based on surface temperature measurements [Sews et al., 1992]
  • the pointwise GPS PWV measurements agree to those from radiosondes and microwave radiometers with standard deviations of about 1-2 mm [e.g., Ohtani and Naito., 2000; Liou et al., 2001 ; Li et al., 2003], and may then be interpolated to provide PWV fields, which has been attempted using different models.
  • Jarlemark and Emardson [1998] evaluated a gradient model, a linear regression in time model that ignored observational directions, and a turbulence model that yielded at least 10% improved RMS error than the other two models.
  • the aim of this study is therefore to improve the accuracy of GPS interpolated tropospheric water vapor maps by accounting for the coupling effect of both the terrain elevation dependency and tropospheric turbulence, demonstrate this over varying terrain, and compare with the SKIm+Onn model.
  • PPP point positioning processing
  • This region was selected as there is a dense spacing of continuous GPS stations whose data have been extensively analyzed previously (including post- processed ZTDs readily available for quality control), it is comparable to the size of a SAR image swath (e.g., 100 km for Envisat Stripmap mode and 250 km for Sentinel-1 Interferometric Wide Swath mode), there are large variations in topography, and it experiences relatively cloud-free conditions as needed for validation with MODIS.
  • SAR image swath e.g., 100 km for Envisat Stripmap mode and 250 km for Sentinel-1 Interferometric Wide Swath mode
  • Tropospheric delays especially the part due to atmospheric water vapor, vary both vertically and laterally over short distances, and are often considered as the sum of (i) a stratified component highly correlated with topography which therefore delineates the vertical tropospheric profile, and (ii) a turbulent component resulting from disturbance processes (e.g., severe weather) in the troposphere which trigger uncertain patterns in space and time [Williams et al., 1998; Cho et al., 2003; Li et al., 2006a] We present an Iterative Tropospheric Decomposition (ITD) model to effectively separate the turbulent and elevation dependent ZTD components.
  • ITD Iterative Tropospheric Decomposition
  • the turbulent component usually consists of medium-to-long wavelength signals that can be interpolated by an IDW method. If n GPS stations are used in the region considered, the IDW model reads as: where w u denotes the interpolation coefficient; u and i are indices for the user and reference stations, respectively; d Ui represents the horizontal distance from the user to reference station. Reference stations at distances more than 100 km from the user station show limited correlation [Emardson and Johansson, 1998] so are not used.
  • the residuals s k which are the summation of the unmodeled errors and the turbulent component, are computed by subtracting per station the stratified delay (as modeled by the estimated exponential coefficients) from the ZTD.
  • Steps (ii) to (iv) are repeated until the exponential cofficients b and L 0 converge.
  • the final outputs are the exponential coefficients for the decorrelated range limit considered, plus the turbulent delay component and residuals per reference station.
  • the ZTD at the location of interest is obtained by interpolation of the reference station turbulent component and residuals, and added to the stratified delay computed using the final values of the exponential coefficients.
  • the ITD model by ZTDs in zero difference mode, it is also suitable for interpolating differenced ZTD or PWV/ZWD.
  • the input to the ITD model may be, as well as GPS-ZTDs, other water vapor products such as those from the ECMWF numerical weather model, MODIS, or their combination.
  • the pointwise ZTD values for the 41 study region GPS stations used for the generation of the interpolated maps were estimated using real-time mode GPS PPP processing.
  • PPP relies on highly accurate satellite orbits and clocks [Zumberge et al., 1997], which are usually held fixed in post-processed PPP.
  • IGS International GNSS Service
  • the satellite clocks have unpredictable behavior which makes their real-time prediction challenging [Li etai., 2014], so we did not fix these to the real time product values but estimated corrections to them using satellite clock parameters with the Gurium and Koch [2002] robust estimation method.
  • res(dt k ) dt k - dt k,RTS (5)
  • dt k RTS is the initial value of the satellite clock given by the real-time product and acts as a pseudo-observation
  • res(dt k ) and dt k are the satellite clock residual and value, respectively.
  • the satellite clock parameters were estimated as white noise parameters with a sigma of 0.001 ns, and the error messages contained in the real-time satellite clock product were used to determine the weights of the pseudo-observations in equation (4).
  • An iterative process was used to identify some clock outliers which were hence ignored or assigned less weight in subsequent iterations [Gurium and Koch, 2002]
  • Figure 1 Mean (a) and RMS (b) differences per station between RTPPP and GIPSY ZTDs, estimated every 5 minutes on 2 September 2015 for all 41 GPS stations in the California study region. For all stations, the mean difference and the mean RMS difference are 1.9 mm and 10.1 mm respectively, with this day being indicative of the median differences for all days of 2015. The white area represents the Pacific Ocean and San Francisco Bay, and the background shows the elevation. Stations P177, P230 and S300 are labeled as they are considered in Figure 9.
  • the RMS difference is 12.5 mm, commensurate with the spatial RMS difference and further indicating an RTPPP ZTD precision of about 1 cm, which is commensurate with previous real-time studies [e.g. Li et al., 2014; Ahmed et al., 2016]
  • Validation was carried out by comparing the interpolated ZTDs with the RTPPP ZTD estimates themselves at 14:00 local time and computing the RMS difference and Mean Absolute Error (MAE) for all stations for each day. This was repeated using interpolation with the SKIm+Onn model.
  • MAE Mean Absolute Error
  • the ITD interpolation model When considered for all stations for the entire year, the ITD interpolation model has been shown to substantially improve on the SKIm+Onn model. However, in the cross validation, some large differences occurred (more than 2 cm magnitude for both models), which suggests that the interpolation result is influenced by the GPS reference station distribution. To investigate this, the annual MAE and RMS differences per station are plotted in Figure 5, for both the ITD and SKIm+Onn models. It can be seen that the smaller MAE and RMS differences occur where the station coverage is denser, but the SKIm+Onn MAE and RMS values show substantial degradation compared with ITD in the north-west of the region where there are fewer stations, e.g.
  • the RTPPP GPS pointwise ZTD values were converted to PWV, interpolated to 1 km pixels across the entire study region and compared with the MODIS near-IR PWV product.
  • Pressure and temperature data at 5 minute temporal resolution from co-located meteorological sensors were available at four of the 41 GPS stations and obtained from unavo.org. These were supplemented by 10 meteorological stations which were located up to 10 km outside the study region. The meteorological data were first interpolated to all 41 stations using the Li et al.
  • the Level 2 data were generated at the 1-km spatial resolution of the MODIS instrument using the near-IR algorithm [Gao and Kaufman, 2003] About 30% of days had severe cloud conditions so we excluded them as only a few grids can be obtained. Areas with cloud conditions or above water were also masked and only the cloud free land areas were used in the comparison, with the percentage of data points available and used in the comparisons shown in Figure 6(e), with the average cloud-free pixels over the entire year being 35%.
  • cross validation was carried out at each GPS reference station, using one GPS PWV per day per station, taken at the MODIS acquisition time.
  • MODIS PWV as ‘truth’ values to validate the interpolated PWV at each GPS station (cross test), and the MODIS PWV values were averaged over boxes of 3 c 3 pixels centered on the GPS station’s location if the centered pixel was missing.
  • the spatial resolution of the MODIS PWV maps was effectively reduced from 1 km x 1 km to 3 km x 3 km, this resolution is still much greater than the ⁇ 15 km GPS station spacing across the study region, so provided a compromise between resolution and maximizing the number of MODIS pixels for comparison.
  • the SKIm+Onn model also produced a greater mean difference than ITD during times of lower PWV content (i.e., 0-10 mm), with the more effective ITD modeling in the times with less water vapor (i.e., from mid-autumn to mid-spring) being consistent with the GPS cross validation results shown in Figure 3, i.e. the improved separation of medium to long wavelength turbulence signals and improved performance across mountainous areas.
  • any values over water areas have been removed since PWVs above water (bay, lake or ocean) share different characteristics from PWVs over land areas which cannot be well-described by the interpolation model [Sobrino et al., 2003]
  • the mountainous areas give more negative differences than flat terrain, showing that MODIS tends to overestimate PWV compared with GPS with increasing altitude (i.e., small PWV values), as previously found by Li et al. [2003]
  • the edge areas with fewer GPS stations produce larger differences than central areas, confirming as discussed in section 4 that improved GPS station coverage will improve the quality of interpolated PWV maps.
  • the ITD model produces smoother difference maps than SKIm+Onn and has a lower percentage of large differences. ITD also performs much better than SKIm+Onn in coastal areas where the PWV is more changeable and gives more complicated turbulent components.
  • FIG 8 we plot ITD-based PWVs over mountainous areas, shown as (a) from 37° 09’ to 38° 30’ in latitude and -122° 12’ to -121° 00’ in longitude, and (b) from 37° 09’ to 37° 50’ in latitude and -122° 30’ to -121° 00’ in longitude, for 3 September 2015.
  • FIG 8 we include PWV profiles (smoothed using a tenth averaging window) for both ITD RTPPP and MODIS along lines of constant latitude and longitude, over both mountainous and flat areas, which enable detailed comparisons of the ITD PWV gradients with respect to topography.
  • the ITD PWV profiles change in a similar tendency with MODIS and share similar gradients.
  • the overall RMS difference between MODIS and ITD PWV for the eight profiles considered is 1.51 mm and the RMS differences for mountain (gray polygon in Figure 8) and flat areas are 1.57 mm and 1.47 mm, respectively.
  • FIG. 6 Cross validation of (a) ITD and (b) SKIm+Onn interpolated RTPPP GPS PWV with MODIS, using daily values at MODIS acquisition time on all 41 GPS stations for all of 2015. (c) and (d) display all the available pixels between MODIS PWV and ITD/SKIm+Onn PWV maps for year 2015.
  • the color scale represents the density of occurrence and the red boxes highlight instances of particular improvement (see section 5).
  • the daily cloud free MODIS PWV pixel density is displayed by (e) in which the vertical bar represents the total available pixel numbers divided by the maximum amount.
  • the color scale represents the average daily PWV of all available pixels.
  • (a,b,c,d) are for 3 September 2015 and (e,f,g,h) for 19 November 2015.
  • MODIS maps are MODIS maps
  • (b) and (f) are ITD maps
  • (c) and (g) are ITD difference maps
  • (d) and (h) are SKIm+Onn difference maps.
  • the large differences (red pixels) in the north east of (c) and (d) are likely due to the presence of thin clouds which are not labeled in the cloud mask product.
  • the PWV profile series are in the same order as the line segments in the PWV map, and are averaged by a tenth average window.
  • the gray polygon areas represent the mountain area.
  • the overall RMS difference between MODIS and ITD PWV along the eight profiles is 1.51 mm and the RMS difference for the mountainous (gray polygon) and flat areas are 1.57 mm and 1.47 mm, respectively.
  • the principal aim of the ITD model is to separate the elevation dependent ZTD/PWV component from the turbulent component, which is the most variable and uncertain part, and can easily bias the vertical ZTD scaling, making the separation of the two components challenging. Due to the constraints of the density of GPS stations, only medium-to-long wavelength turbulent signals are expected to be successfully modeled using ITD. To illustrate the size and variation of the turbulent component, and the importance of iterating the solution until convergence arises, a sample three GPS stations were considered: P177, P230 and S300 ( Figure 1), which are in different parts of the study region, are at different elevations, and are at varying distances from the nearest other GPS reference stations.
  • the convergence tendency in Figure 9 reveals that the turbulence effect can be reduced by separating the stratified and turbulent components through iteration.
  • the decoupling procedure acts very similarly to robust estimation, in which the ZTDs from stations exhibiting strong turbulence will contribute less in height scaling but account for more in the turbulent delay interpolation.
  • the iteration also allows for the detection of ZTD outliers, which is a not uncommon occurrence in real-time PPP due to the unpredictable behavior of the satellite clocks. Consequently, the ITD model enables both fitting of the tropospheric vertical profiles and also models the turbulence processes.
  • the fourth column in Figure 9 suggests, successfully accounting for this results in the overall RMS difference of the cross validation test reducing and converging.
  • FIG. 9 RTPPP ZTD turbulent component separated by ITD at each iteration step.
  • the first, second and third columns represent stations P177, P230 and S300, respectively, and the fourth column represents the ZTD cross validation RMS difference for all 41 stations on the corresponding day. Shown for sample days in each of spring, summer and autumn.
  • an iterative tropospheric decomposition model has been developed for the generation of high resolution regional tropospheric PWV fields (maps) from the interpolation of pointwise GPS ZTDs, without any data differencing.
  • a California study region of around 150 km x 150 km the approach of decoupling the terrain elevation dependency and the tropospheric turbulence contribution to ZTD in an iterative procedure (typically 4-6 iterations were required) led to improved accuracy interpolated tropospheric water vapor maps over those based on previous studies, such as the tropospheric turbulence and elevation dependent model SKIm+Onn of Xu et al.
  • the cross validation ITD model RMS and mean differences are similar for both mountainous and flatter terrain, and also similar for both coastal and inland areas.
  • the cross validation improvements using ITD are smallest in the summer months. Spatially, we generated PWV values for 1 km pixels for all land-covered parts of the region and compared with daily MODIS PWV near-IR product values, with the RMS difference for the year being improved from 1.96 mm using the SKIm+Onn model to 1.71 mm using ITD. Furthermore, the spatial PWV gradients using the ITD model and MODIS across a variety of topography were nearly identical to each other.
  • the overall RMS difference between MODIS and ITD PWV profiles is 1.51 mm and the RMS differences for mountain and flat areas are 1.57 mm and 1.47 mm, respectively.
  • the ITD PWV fields are also able to reveal detailed water vapor information over varying terrains.
  • FIG. 10 RTPPP PWV fields across the California study region every 2 hours on 2 November 2015 during a rainfall process from (a) 10:00 to (f) 20:00 UTC. Arrows represent PWV increasing (upwards) or decreasing (downwards) during the preceding 2 hours.
  • the generated spatially-dense PWV fields with continuous, high (5 minute) temporal resolution are not only suitable for correcting atmospheric effects in SAR images at the instant of acquisition, but they also will ensure the identification of water vapor variation from ground motion between image acquisition times [Foster et al., 2006]
  • the high performance of the dense PWV maps using the ITD model is especially useful for mitigating the effects of water vapor for SAR measurements in mountainous areas, which usually suffer from vertical stratification and turbulent mixing due to the orography [Wadge et al., 2002]
  • GMF Global Mapping Function
  • Emardson, T. R., and J. M. Johansson (1998), Spatial interpolation of the atmospheric water vapor content between sites in a ground-based GPS Network, Geophys. Res. Lett., 25(17), 3347-3350.
  • Li, Z., E. J. Fielding, and P. Cross (2009a) Integration of InSAR time series analysis and water vapour correction for mapping postseismic deformation after the 2003 Bam, Iran Earthquake, IEEE T. Geosci. Remote., 47, 3220-3230.

Landscapes

  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Position Fixing By Use Of Radio Waves (AREA)

Abstract

A method for configuring a Global Satellite Navigation System (GNSS) infrastructure comprising a plurality of GNSS receivers with regards to mitigating tropospheric delay errors, is disclosed. The method comprises: obtaining a water vapour map comprising water vapour values at different points; and performing a configuration process comprising: defining one or more sub-regions based on an initial set of receivers, such that one or more receivers are located on (e.g. on the boundary of, on one or more vertices of, and/or within) each sub-region; and for each sub-region: obtaining, from the water vapour map, a first water vapour value for each of a plurality of points on the sub-region; determining a second water vapour value for each of the plurality of points by using interpolation of water vapour values, obtained from the water vapour map, corresponding to positions of two or more receivers; determining a statistical value based on the first and second water vapour values for each point; and if the statistical value satisfies a predetermined criteria, increasing the number of receivers on the sub-region to obtain a modified set of receivers.

Description

CONFIGURATION OF GLOBAL SATELLITE NAVIGATION SYSTEM INFRASTRUCTURE
Field
Certain examples of the present disclosure provide a technique for configuring (e.g. refining and/or optimizing) Global Satellite Navigation System (GNSS) infrastructure, and in particular for configuring a network of GNSS receivers (or stations) with regards to tropospheric effects. For example, the present disclosure provides a technique for configuring GNSS infrastructure for use in geophysical, meteorological and/or engineering applications.
Background
Instantaneous Global Navigation Satellite Systems (GNSS) positioning with centimetre-level precision and accuracy is commonly used and needed in:
• utilities mapping (e.g., Li, S., Cai, H., Abraham, D. M. & Mao, P. Estimating features of underground utilities: Hybrid GPR/GPS approach. Journal of Computing in Civil Engineering 30, 04014108 (2016)),
• landslide monitoring (e.g., Jiao, Y, Wang, Z., Wang, X., Adoko, A. C. & Yang, Z. Stability assessment of an ancient landslide crossed by two coal mine tunnels. Engineering Geology 159, 36-44 (2013)),
• coastal erosion monitoring (e.g., Harley, M. D., Turner, I. L., Short, A. D. & Ranasinghe, R. Assessment and integration of conventional, RTK-GPS and image-derived beach survey methods for daily to decadal coastal monitoring. Coastal Engineering 58, 194-205 (2011)),
• glacier and permafrost monitoring (e.g., Pan, B., Cao, B., Wang, J., Zhang, G., Zhang, C., Hu, Z. & Huang, B. Glacier variations in response to climate change from 1972 to 2007 in the western Lenglongling mountains, northeastern Tibetan Plateau. Journal of Glaciology 58, 879-888 (2012); Goyanes, G., Vieira, G., Caselli, A., Cardoso, M., Marmy, A., Santos, R, Bernando, I. & Hauck, C. Local influences of geothermal anomalies on permafrost distribution in an active volcanic island (Deception Island, Antarctica). Geomorphology 225, 57-68 (2014)),
• river channel monitoring (e.g., Brasington, J., Rumsby, B. T. & McVey, R. A. Monitoring and modelling morphological change in a braided gravel-bed river using high resolution GPS- based survey. Earth Surface Processes and Landforms 25, 973-990 (2000)),
• precision agriculture (e.g., Stoll, A. & Kutzbach, H. D. Guidance of a forage harvester with GPS. Precision Agriculture 2, 281-291 (2000); Perez-Ruiz, M., Slaughter, D. C., Gliever, C. & Upadhyaya, S. K. Tractor-based Real-time Kinematic-Global Positioning System (RTK-GPS) guidance system for geospatial mapping of row crop transplant. Biosystems Engineering 111, 64-71 (2012), hereinafter Perez-Ruiz et al., 2012),
• highway construction (e.g., Pirti, A. Performance analysis of the real time kinematic GPS (RTK GPS) technique in a highway project (stake-out). Survey Review 39, 43-53 (2007)) structural monitoring (e.g., Nickitopoulou, A., Protopsalti, K. & Stiros, S. Monitoring dynamic and quasi-static deformations of large flexible engineering structures with GPS: Accuracy, limitations and promises. Engineering Structures 28, 1471-1482 (2006), hereinafter Nickitopoulou et al. , 2006),
• flood modelling (e.g., Narama, C., Duishonakunov, M., Kaab, A., Daiyrov, M. &Abdrakhmatov,
K. The 24 July 2008 outburst flood at the western Zyndan glacier lake and recent regional changes in glacier lakes of the Teskey Ala-Too range, Tien Shan, Kyrgyzstan. Natural Hazards and Earth System Sciences 10, 647-659 (2010)), and
• photogrammetric mapping and LiDAR surveys from aeroplanes/helicopters (e.g., Javernick,
L, Brasington, J. & Caruso, B. Modeling the topography of shallow braided rivers using
Structure-from-Motion photogrammetry. Geomorphology 213, 166-182 (2014)) and unmanned aerial vehicles (e.g., Rabah, M., Basiouny, M., Ghanem, E. & Elhadary, A. Using RTK and VRS in direct geo-referencing of the UAV imagery. NRIAG Journal of Astronomy and Geophysics 7, 220-220 (2018)).
Most of these applications and studies have realized such centimetre-level positioning using local, single baseline real-time kinematic (RTK) GNSS. However, this has the limitation that the user must establish their own GNSS reference station within about 10 km of the user’s roving receiver (Wanninger, L. Virtual reference stations (VRS). GPS Solutions 7, 143-144 (2003)), so that common user and reference station GNSS errors such as atmospheric delays can be eliminated by data differencing.
In parallel with the development of new GNSS that complement and enhance the original Global Positioning System (GPS), notably GLONASS, Galileo and Beidou, networks of Continuously Operating GNSS Reference Stations (CORS) have been increasingly established in many parts of the world. This has led to the development of the Network Real-time Kinematic (Network RTK) GNSS positioning technique, which adopts the working principle that the CORS at known ground locations continuously compute GNSS observational error corrections in real-time. These corrections are then transmitted to and applied by a single user receiver in order to improve the accuracy of its position to centimetre-level (e.g., Hu, G., R., Khoo, H. S., Goh, P. C. & Law, C. L. Development and assessment of GPS virtual reference stations for RTK positioning. Journal of Geodesy 77, 292-302 (2003)). Therefore with a sufficient global coverage of CORS, Network RTK may be applied worldwide, used for all the abovementioned RTK applications, and also for developing technologies which require seamless, reliable, instantaneous positioning with centimetre-level precison and accuracy, such as intelligent transport systems (e.g., Meng, X., Roberts, S., Cui, Y, Gao, Y, Chen, Q., Xu, C., He, Q., Sharpies, S. & Bhatia, P. Required navigation performance for connected and autonomous vehicles: Where are we now and where are we going? Transportation Planning and Technology 41, 104-118 (2018)), autonomous vehicles and the Internet of Things.
Current CORS networks (potentially) used for Network RTK positioning have been established in a relatively ad-hoc manner. In some instances, such as in the British Isles, Sweden and Victoria (Australia), Network RTK service providers make use of governmental geodetic networks of CORS which have been set up with approximately even geographical spacing (e.g., 60-80 km in Britain). However, in many nations such CORS are clustered in urban areas or along highways (e.g., Australia and Canada), or to provide national geodetic infrastructure for goals such as determining plate motions, earthquake characteristics and volcano dynamics, and hence many such CORS are concentrated around faults or plate boundaries (e.g., the Plate Boundary Observatory: http://pboweb.unavco.org/). Many other Network RTK service providers have established their own specific commercial CORS networks, although these also tend to comprise either urban area clusters or an approximately even geographical spacing.
Regardless of the set up, Network RTK works on the premise that the application of the transmitted, interpolated errors from the CORS enable the same centimetre-level user positional quality to be obtained as from local reference station RTK positioning. However, little apparent consideration has been given to whether the established CORS are optimally spaced with regards to the ensuing positional accuracy at the user’s location, and there is a dearth of information on the optimal CORS configuration needed and whether dominant error sources such as atmospheric delays may be sufficiently interpolated and mitigated with existing CORS infrastructure.
The above information is presented as background information only to assist with an understanding of the present disclosure. No determination has been made, and no assertion is made, as to whether any of the above might be applicable as prior art with regard to the present disclosure.
Summary
It is an aim of certain examples of the present disclosure to address, solve, mitigate or obviate, at least partly, at least one of the problems and/or disadvantages associated with the related art, for example at least one of the problems and/or disadvantages mentioned herein. Certain examples of the present disclosure aim to provide at least one advantage over the related art, for example at least one of the advantages mentioned herein.
The present invention is defined in the independent claims. Advantageous features are defined in the dependent claims.
Embodiments or examples disclosed in the description and/or figures falling outside the scope of the claims are to be understood as examples useful for understanding the present invention.
In accordance with an aspect of the present disclosure, there is provided a method for configuring a Global Satellite Navigation System (GNSS) infrastructure comprising a plurality of GNSS receivers. The method comprises: obtaining a water vapour map comprising water vapour values at different points; and performing a configuration (e.g. refinement) process comprising: defining one or more sub-regions based on an initial set of receivers, such that one or more receivers are located on (e.g. on the boundary of, on one or more vertices of, and/or within) each sub-region; and for each sub-region: obtaining, from the water vapour map, a first water vapour value for each of a plurality of points on the sub-region; determining a second water vapour value for each of the plurality of points by using interpolation of water vapour values, obtained from the water vapour map, corresponding to positions of two or more receivers; determining a statistical value based on the first and second water vapour values for each point; and if the statistical value satisfies a predetermined criterion, increasing the number of receivers on the sub-region to obtain a modified set of receivers.
In accordance with another aspect of the present disclosure, there is provided an apparatus for configuring a Global Satellite Navigation System (GNSS) infrastructure comprising a plurality of GNSS receivers, the apparatus comprising a processor configured to: obtain a water vapour map comprising water vapour values at different points; and perform a configuration (e.g. refinement) process comprising: defining one or more sub-regions based on an initial set of receivers, such that one or more receivers are located on (e.g. on the boundary of, on one or more vertices of, and/or within) each sub-region; and for each sub-region: obtaining, from the water vapour map, a first water vapour value for each of a plurality of points on the sub-region; determining a second water vapour value for each of the plurality of points by using interpolation of water vapour values, obtained from the water vapour map, corresponding to positions of two or more receivers; determining a statistical value based on the first and second water vapour values for each point; and if the statistical value satisfies a predetermined criterion, increasing the number of receivers on the sub-region to obtain a modified set of receivers.
In certain examples, the GNSS infrastructure may be configured with regards to the mitigation of tropospheric effects.
In certain examples, the two or more receivers used for the interpolation comprise one or more receivers on the sub-region and/or one or more receivers on one or more neighbouring sub- regions.
In certain examples, the method further comprises, or the processor is further configured for, repeating the configuration process such that the modified set of receivers obtained at the end of one iteration is used as the initial set of receivers in the next iteration.
In certain examples, the configuration process is repeated until a termination criterion is satisfied.
In certain examples, the termination criterion is based on one or more of: the statistical value does not satisfy the predetermined criterion for every sub-region; the receiver spacing is less than a certain threshold; and the number of iterations has reached a certain threshold.
In certain examples, the water vapour map is obtained from Moderate Resolution Imaging Spectroradiometer (MODIS) data.
In certain examples, the interpolation is performed based on an Iterative Tropospheric Decomposition (ITD) model.
In certain examples, the initial set of receivers comprises a new or pre-existing regular or irregular array of receivers. In certain examples, a sub-region comprises a polygon, and a receiver is located on one or more vertices of the polygon, on one or more edges of the polygon, and/or at the centre of the polygon.
In certain examples, the polygon comprises one or more of: square; rectangle; triangle; and hexagon.
In certain examples, the sub-regions comprise a regular array of rectangular or square sub- regions, or a Triangulated Irregular Network (TIN) obtained using Delaunay triangulation.
In certain examples, the plurality of points on the sub-region have the same spatial resolution as the water vapour map.
In certain examples, the statistical value is the root mean square (RMS) of the differences between the first and second water vapour values.
In certain examples, the predetermined criterion is that the statistical value is greater than a predetermined threshold.
In certain examples, increasing the number of receivers on the sub-region comprises adding a receiver to the geometric centre of the sub-region.
In certain examples, the method further comprises, or the processor is further configured for, generating a report including position, spacing, and/or density information of the modified set of receivers.
In certain examples, the method further comprises, or the processor is further configured for, configuring the GNSS infrastructure according to the modified set of receivers.
In accordance with another aspect of the present disclosure, there is provided a Global Satellite Navigation System (GNSS) infrastructure configured according to the method of any of the preceding aspects and/or examples.
Other aspects, advantages, and salient features of the present disclosure will become apparent to those skilled in the art from the following detailed description, which, taken in conjunction with the accompanying drawings, disclose examples of the present disclosure.
Brief Description of the Figures
Figure 1 illustrates the simulated global GNSS network distribution optimization for boreal summer (June - Aug 2016) and boreal winter (Dec 2015 - Feb 2016);
Figure 2 illustrates the relationship between the precipitable water vapour (PWV) and station spacing (blue line), and between the PWV variation and station spacing (red line), for local summer (a) and local winter (b), respectively. The PWV or PWV variation is a global average of all pixels at the corresponding spacing, and is normalized linearly between [0, 1];
Figure 3 illustrates original and densified (optimally refined) regional networks during local summer of 2016: July 2016 for the British Isles and Italy, January 2016 for New Zealand and Java;
Figure 4 illustrates time series of optimal station spacing from January 2012 to December 2016 for each of the British Isles, Italy, New Zealand and Java networks; Figure 5 is a flow diagram of a global optimization process according to an example of the present disclosure;
Figure 6 is a flow diagram of a local optimization process according to an example of the present disclosure;
Figure 7 is a flow diagram of an Iterative Tropospheric Decomposition (ITD) model according to an example of the present disclosure;
Figure 8 illustrates the evolution of a GNSS network during operation of the techniques disclosed herein in the case of global optimization and local optimization; and
Figures 9-18 disclose Figures 1-10 of the Yu et al. (2017) paper referred to below.
Description
The following description of examples of the present disclosure, with reference to the accompanying drawings, is provided to assist in a comprehensive understanding of the present invention, as defined by the claims. The description includes various specific details to assist in that understanding but these are to be regarded as merely exemplary. Accordingly, those of ordinary skill in the art will recognize that various changes and modifications of the examples described herein can be made.
The terms and words used in this specification are not limited to the bibliographical meanings, but are merely used to enable a clear and consistent understanding of the present disclosure.
The same or similar components may be designated by the same or similar reference numerals, although they may be illustrated in different drawings.
Detailed descriptions of elements, features, components, structures, constructions, functions, operations, processes, characteristics, properties, integers and steps known in the art may be omitted for clarity and conciseness, and to avoid obscuring the subject matter of the present disclosure.
Throughout this specification, the words “comprises”, “includes”, “contains” and “has”, and variations of these words, for example “comprise” and “comprising”, means “including but not limited to”, and is not intended to (and does not) exclude other elements, features, components, structures, constructions, functions, operations, processes, characteristics, properties, integers, steps and/or groups thereof.
Throughout this specification, the singular forms “a”, “an” and “the” include plural referents unless the context dictates otherwise. For example, reference to “an object” includes reference to one or more of such objects.
By the term “substantially” it is meant that the recited characteristic, parameter or value need not be achieved exactly, but that deviations or variations, including for example, tolerances, measurement errors, measurement accuracy limitations and other factors known to those of skill in the art, may occur in amounts that do not preclude the effect the characteristic, parameter or value was intended to provide. Throughout this specification, language in the general form of “X for Y” (where Y is some action, process, function, activity, operation or step and X is some means for carrying out that action, process, function, activity, operation or step) encompasses means X adapted, configured or arranged specifically, but not exclusively, to do Y.
Elements, features, components, structures, constructions, functions, operations, processes, characteristics, properties, integers, steps and/or groups thereof described herein in conjunction with a particular aspect, embodiment, example or claim are to be understood to be applicable to any other aspect, embodiment, example or claim disclosed herein unless incompatible therewith.
It will be appreciated that examples of the present disclosure can be realized in the form of hardware, software or any combination of hardware and software. Any such software may be stored in any suitable form of volatile or non-volatile storage device or medium, for example a ROM, RAM, memory chip, integrated circuit, or an optically or magnetically readable medium (e.g. CD, DVD, magnetic disk or magnetic tape).
Certain embodiments of the present disclosure may provide a computer program comprising instructions or code which, when executed, implement a method, system and/or apparatus in accordance with any aspect, claim, example and/or embodiment disclosed herein. Certain embodiments of the present disclosure provide a machine-readable storage storing such a program.
The techniques described herein may be implemented using any suitably configured apparatus and/or system. Such an apparatus and/or system may be configured to perform a method according to any aspect, embodiment, example or claim disclosed herein. Such an apparatus may comprise one or more elements, for example one or more of receivers, transmitters, transceivers, processors, controllers, modules, units, and the like, each element configured to perform one or more corresponding processes, operations and/or method steps for implementing the techniques described herein. For example, an operation/function of X may be performed by a module configured to perform X (or an X-module). The one or more elements may be implemented in the form of hardware, software, or any combination of hardware and software.
Instantaneous Global Navigation Satellite Systems (GNSS) positioning with centimetre-level precision and accuracy is critical for activites ranging from flood modelling, construction, precision agriculture, utilities mapping, intelligent transport systems, autonomous vehicles, and the Internet of Things. This may be achieved via the network real-time kinematic (Network RTK) GNSS approach, which uses a single receiver and a network of continuously operating GNSS reference stations (CORS). However, existing CORS networks have often been established simply by attempting regular spacing or in clusters around cities, with little consideration of weather, climate and topography effects, which represent a substantial GNSS positional error and prevent homogeneous Network RTK accuracy attainment.
Little apparent consideration has been given to whether the established CORS are optimally spaced with regards to the ensuing positional accuracy at the user’s location, and there is a dearth of information on the optimal CORS configuration needed and whether substantial or dominant error sources such as atmospheric delays (e.g. tropospheric delay and/or ionospheric delay) may be sufficiently interpolated and mitigated with existing CORS infrastructure.
One of the two components of atmospheric delay arises from the troposphere, and this tropospheric delay is recognized as a substantial remaining Network RTK positional error source (e.g., Grejner-Brzezinska, D. A., Kashani, I. & Wielgosz, R On accuracy and reliability of instantaneous network RTK as a function of network geometry, station separation, and data processing strategy. GPS Solutions 9, 212-225 (2005)), especially the error induced by atmospheric water vapour (Ahn, Y. W, Lachapelle, G., Skone, S., Gutman, S. & Sahm, S. Analysis of GPS RTK performance using external NOAA tropospheric corrections integrated with a multiple reference station approach. GPS Solutions 10, 171-186 (2006)). Atmospheric water vapour may also be expressed as precipitable water vapour (PWV) and it is a highly variable quantity both spatially and temporally (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 Geophysical Research 97, 15787-15801 (1992)).
Therefore the optimal CORS spacing for PWV interpolation may be considered a proxy for Network RTK CORS infrastructure design with regards to the mitigation of GNSS tropospheric effects.
Certain examples of the present disclosure provide centimetre- 1 eve I homogeneous positioning in real-time with GNSS network infrastructure. Certain examples of the present disclosure provide a technique for optimizing the design of CORS ground infrastructure (e.g. optimizing CORS spacing for PWV interpolation corresponding to optimizing Network RTK CORS infrastructure design) to enable consistent user centimetre-level positioning globally with regards to tropospheric (PWV) effects. For example, as discussed further below, it is shown that denser networks are typically required in the tropics and in mountainous areas. Certain examples of the present disclosure also provide a technique for improving existing national GNSS network infrastructure.
This includes determining how many CORS stations are required over a given area, by computing the staton spacing density needed for the successful interpolation of PWV (and hence tropospheric delay) with homogeneous accuracy, including the dependence on geographic location, topography and climate. Ensuring such optimal CORS network design will better enable the successful uptake of Network RTK globally and subsequently its practical uptake for all of the above mentioned RTK applications and developing technologies.
Certain examples will now be described for configuring a network of GNSS receivers (or stations). However, the skilled person will appreciate that the present disclosure is not limited to these specific examples. For example, the skilled person will appreciate that the present disclosure is not limited to the particular values of variables and parameters referred to, and that any other suitable values may be used in alternative methods.
In the present disclosure, configuring a network of GNSS receivers (or stations) may be understood as meaning determining, defining, designing and/or creating an arrangement of receivers. The arrangement of receivers may be defined, for example, by specifying positions of individual receivers, by specifying the spacing between receivers and/or by specifying the density of receivers in a certain area. If there is an existing network, then configuring the network may be understood to include modifying, refining and/or optimizing the existing network. On the other hand, if no network exists, then configuring a network may be understood to include creating or defining a new network that may be modified, refined and/or optimized. The skilled person will understand that references herein to “optimizing” or “refining” may be taken to also include “configuring” according to the above definition.
In the following, consideration is given to the impact of PWV global spatial (climatic) and temporal (seasonal) distributions in formulating CORS network optimal design, as well as topographic variations. Scenarios for idealistic new CORS Network RTK design are considered, ranging from CORS spacings of 80 km to 5 km (assumed to be the maximum density of a CORS network anywhere), along with case studies from selected national CORS networks as to how existing infrastructure can be best densified. Data from the Moderate Resolution Imaging Spectrometer (MODIS) instrument on board NASA’s Terra and Aqua satellites are used to provide PWV values at approximately 5 km (0.05 degree grid) resolution globally, to provide a reference data set for testing the PWV interpolation success, which is carried out using an Iterative Tropospheric Decomposition (ITD) model (Yu, C., Penna, N. T. & Li, Z. Generation of real-time mode high-resolution water vapor fields from GPS observations. Journal of Geophysical Research: Atmospheres 122, 2008-2025 (2017), hereinafter Yu et al. (2017)). This model is able to decouple both elevation-dependent and lateral PWV variations and, importantly, from the analysis of a preliminary station spacing assessment over California, it is found with this model that a 5 km spacing resulted in a similar accuracy as using a 2 km spacing for the interpolation of PWV, with interpolation errors from both spacings at the noise level of the PWV measurements (Yu., C., Li, Z. & Penna, N. T. Interferometric synthetic aperture radar atmospheric correction using a GPS-based iterative tropospheric decomposition model. Remote Sensing of Environment 204, 109-121 (2018a), hereinafter Yu et al. (2018a)).
Methods
In the present example, the water vapour data is obtained from Moderate Resolution Imaging Spectroradiometer (MODIS) data, although the skilled person will appreciate that the water vapour data may be obtained from any additional or alternative source. MODIS comprises two satellites providing high radiometric sensitivity (12 bit) in 36 spectral bands ranging in wavelength from 0.4 pm to 14.4 pm. The present example uses the near-infrared band water vapour product which applies a near-infrared algorithm over clear land areas of the globe and above clouds over both land and ocean. However, the skilled person will appreciate that any other suitable band may be used. Global water vapour reference grid. To generate global PWV grids at 5 km spatial resolution, the present example method uses the near-infrared band water vapour MODIS product over cloud-free land areas of the globe (Gao, B.-C. & Kaufman, Y. J. Water vapor retrievals using Moderate Resolution Imaging Spectroradiometer (MODIS) near-infrared channels. Journal of Geophysical Research 108 (D13), 4389 (2003), hereinafter Gao and Kaufman (2003)). In cloud- free conditions, this MODIS product provides a PWV measurement accuracy and precision of about 1-2 mm, based on comparisons with GPS and radiosondes (e.g. Li, Z., Muller, J. P. & Cross, P. Comparison of precipitable water vapor derived from radiosonde, GPS, and Moderate-Resolution Imaging Spectroradiometer measurements. Journal of Geophysical Research 108 (D20), 4651 (2003); Liu, J., Liang, H., Sun, Z. & Zhou, X. Validation of the Moderate-Resolution Imaging Spectroradiometer precipitable water vapor product using measurements from GPS on the Tibetan Plateau. Journal of Geophysical Research: Atmospheres 111, D14103 (2006)). It provides higher accuracy than PWV from high-resolution global and regional numerical weather models (Li, H. & Chen, B. The evolution of precipitable water associated with the Asian and Australian monsoons as revealed from MODIS/SSMI, ECMWF and NCEP reanalysis data sets. Geophysical Research Letters 32, L10811 (2005); Prasad, A. K. & Singh, R. P. Validation of MODIS Terra, AIRS, NCEP/DOE AMIP-II Reanalysis-2, and AERONET Sun photometer derived integrated precipitable water vapor using ground-based GPS receivers over India. Journal of Geophysical Research 114, D05107 (2009)), particularly over regions where few PWV observations have been assimilated, such as Africa.
The PWV data may be generated at about 1-km spatial resolution, with 1 to 2 days required for satellite measurements of the whole Earth’s surface. If the ~1-km data are not uniformly spaced, the present method may sample them to a regular 0.05 by 0.05-degree grid (~5-km spacing), taking as the grid node value the nearest original (1-km) data point. To reduce the impact of meridional convergence on the uniformity of the grid point spacing, the present method may only use in tests grid points from the equator to latitude 60 degrees north/south, at which 0.05-degree grid spacing corresponds to ~5.6 km and ~2.8 km, respectively.
The resulting latitude 60 degrees north to south PWV distributions over land are shown in Figure 1 for 3-month averages for boreal summer (1 June 2016 to 31 August 2016) and for boreal winter (1 December 2015 to 29 February 2016). Figure 1(a) illustrates average PWV content from MODIS (0.05 degree resolution); Figure 1(b) illustrates PWV variations, defined as the PWV standard deviation of all points within a 0.8x0.8 degree window centered about each 0.05 degree MODIS pixel (0.05 degree resolution); Figure 1(c) illustrates optimal station spacing (0.8 degree resolution), including percentage of grid cells with different spacing; Figure 1(d) illustrates the correlation map between the optimal station spacing and PWV variation (0.8 degree resolution). Topography variations are denoted by the shaded background.
This averaging per grid point may be done to reduce the effect of weather variations during each season, but also to minimize the effect of cloud-induced data availability gaps. PWV measurements over water bodies may be excluded using the Shuttle Radar Topography Mission 90 m digital elevation model (DEM) (Farr, T. G. et al. The shuttle radar topography mission. Reviews of Geophysics 45, RG2004 (2007)) and a DEM corresponding to the 0.05 degree PWV grid may be generated to enable the assessment of topographic effects.
Interpolation of PWV. The PWV values may be interpolated, for example using the ITD model developed by Yu et al. (2017). The contents of Yu et al. (2017) are disclosed in the annex to the description of the present application, in combination with Figures 9-18 of the present application which correspond to Figures 1-10 of Yu et al. (2017). The disclosure of Yu et al. (2017) in the annex to the description and in Figures 9-18 forms part of the disclosure of the present application. In the present example, the PWV is assumed to comprise a stratified (topography/elevation dependent) component, and a turbulent component, representing topography-independent PWV signals (Equation 1):
PWVi = S(Hi) + T( i) + ei (1) where S(Hj) = L0e~^Hi which represents the stratified component at position vector x^; i is the location considered; H, is the height above sea level; U is the stratified PWV at sea level; b is an exponential coefficient; T is the turbulent PWV component; e, is the residual PWV.
All PWV values within a defined 150 km (or other suitable value) tropospheric decorrelation distance (as per Yu., C., Li, Z., Penna, N. T. & Crippa, P. Generic atmospheric correction model for interferometric synthetic aperture radar observations. Journal of Geophysical Research: Solid Earth 123, 9202-9222 (2018b)) of the location considered may be used to estimate an initial value for the stratified component, computing Lo and b and setting the turbulent component to zero. The residuals (differences between the PWV values per point and the initial value for the stratified component) may then be computed per grid point, which may be deemed to contain the turbulent component plus any unmodelled stratified PWV component. The residuals may then be interpolated, for example using inverse distance weighting to obtain an initial estimate of the turbulent component at the location of interest, and repeated to obtain initial values for the turbulent component for all grid points, and subtracted from the initial values of the stratified delay per grid point. The stratified delay may then be re-estimated to obtain updated values for the exponential coefficients and subsequently updated residuals, which may then be interpolated again to obtain updated values for the turbulent component. This iterative process may be repeated until the estimates for the stratified and turbulent PWV components converge, and the final value of the interpolated PWV may be obtained by summing the stratified and turbulent components, plus the residual. Typically six or seven iterations may be needed, for example as demonstrated in Yu et al. (2017), where further details of an exemplary method may be obtained.
In certain examples, an interpolated PWV value (PWVo) at a specific pixel (Po) may be determined based on PWV values (PWVi, PWW, PWV3, ...) at a plurality of other pixels (Pi, P2, P3, ...). The other pixels (Pi, P2, P3, ...) may correspond to positions of receivers, for example receivers within a defined range (e.g. 150km) from Po. Any receiver at Po may be referred to as a “roving” receiver (or station) or “rover” receiver (or station). The PVW values (PVWi, PVW2, PVW3, ...) at the plurality of other pixels (Pi, P2, P3, ...) may be obtained from the global water vapour reference grid generated from MODIS data. The interpolation may be performed using the ITD model of Yu et al. (2017). Accordingly, in certain examples, once the MODIS data is obtained, the algorithm for configuring a network of receivers may be performed on a computer without requiring any further measurements taken by receivers or other apparatus.
As mentioned above, in certain examples, the PVW values (PVWi, PVW2, PVW3, ...) at the plurality of other pixels (Pi, P2, P3, ...) may be obtained from the global water vapour reference grid generated from MODIS data. In certain alternative examples, one or more, or all, of the PVW values (PVWi, PVW2, PVW3, ...) at the plurality of other pixels (Pi, P2, P3, ...) may be obtained through measurements made by the receivers at those pixels and/or any other suitable apparatus (e.g. GPS or meteorological apparatus).
Global optimal station spacing. The optimal network distribution across the land areas of the globe may be evaluated by computing the station spacing needed in each cell of a 0.8x0.8 degree grid. This may be chosen because 0.8 degrees represents about 80 km, and may be assumed to be the coarsest CORS spacing used for existing RTK networks. Each 0.8x0.8 degree cell may be sub-divided into 256 (16x16) pixels which match the 0.05 degree MODIS PVW reference grid, and then PVW values per pixel may be obtained by interpolation from a particular defined station spacing. The RMS of the differences between the interpolated and reference value for each of the 256 pixels per 0.8x0.8 degree cell may be computed. This leads to a set of RMS PVW interpolation errors per cell of a 0.8 degree grid globally. If in any cell the RMS error is greater than 1.5 mm, then the station spacing used may be considered insufficient for that particular cell and densified until the RMS error for each and every cell is less than 1.5 mm. This may be realized by first densifying such grid cells to use 0.4 degree reference PVWs for interpolation, or a station spacing of ~40 km. The interpolation to the 256 pixels in the 0.8 degree cell may then be repeated, but using the 0.4 degree grid points within the 0.8 degree cell, plus the 0.8 degree grid points outside it. If the RMS error for the same 256 pixels is still greater than 1.5 mm, the 0.8 degree grid cell may be further densified to a simulated station spacing of 0.2 degrees (20 km), and the interpolation repeated. This may be continued if necessary using 0.1 then 0.05 degree spacing until the RMS error for all 256 pixels is lower than 1.5 mm, or the 0.05 degree spacing is reached. The station spacing needed per 0.8 degree grid cell may then be deemed the coarsest grid spacing needed for the RMS error to be 1.5 mm or lower. The resulting station spacing map shown in Figure 1 consequently comprises a station spacing per 0.8 degree cell ranging from 80, 40, 20, 10 to 5 km.
The standard deviation of the 256 0.05 degree grid PVW values within each 0.8 degree cell may be used to represent the PVW variation, and correlation coefficients between station spacing and PVW variation are shown in Figure 1. These are computed for each 0.8 degree grid cell, using the station spacing and PWV variation for the cell itself and the eight 0.8 degree cells immediately surrounding it.
Local GNSS network infrastructure optimization. To assess whether existing GNSS networks are optimally spaced, or what CORS densification should be carried out, a similar interpolation approach to the global uniform grid simulations may be used. However, as existing GNSS networks may have an irregular distribution of stations, instead of starting from a coarse uniform grid, a Triangulated Irregular Network (TIN) may be formed, for example based on Delaunay triangulation. Another difference may be that 1 -month averages of the MODIS PWV daily values per point of the global 0.05 degree grid may be used, not 3-month averages. The procedure may comprise:
(i) Construct a TIN using the original network, i.e. using all existing stations, for example as at December 2016, where the PWV value for the station may be taken as that of the nearest 0.05 degree grid point;
(ii) Interpolate the PWV from the nodes of the TIN on to the 0.05 degree grid reference points;
(iii) Compute the RMS error (difference between the interpolated and reference PWV values) of all the grid points contained within each triangle of the TIN and, if the RMS error exceeds 1.5 mm, insert a new station at its geometric centre;
(iv) Reconstruct the TIN based on the original and added stations, and repeat steps (i) to (iii) until the RMS error for each and every triangle is not greater than 1.5 mm.
Further Methods
Figure 5 is a flow diagram of a global optimization process according to an example of the present disclosure. The skilled person will appreciate that features of the embodiment of Figure 5 may be incorporated into the methods described above, and vice versa.
Global optimization: The optimal spacings may be obtained by first thinning the PWV values from the MODIS 5 km global resolution reference data set to 80 km (0.8 degrees) resolution. Each 0.8x0.8 degree cell (i.e. sub region) may be sub-divided into 256 (16x16) pixels which matches the 0.05 degree MODIS PWV reference grid, and then PWV values per pixel may be obtained by interpolation from a particular defined station spacing, for example using the ITD model of Yu et al. (2017) which considers the DEM and applies the water body mask. The RMS of the differences between the interpolated and reference value for each sub region may be computed. This leads to a set of RMS PWV interpolation errors per sub region globally (i.e. the 5 km resolution global water vapour bias map). If in any sub-region the RMS error is greater than 1.5 mm, then the station spacing used may be considered insufficient for that particular sub-region and will be halved to form a refined global non-uniform network. If the refined network is different from the previous one, the procedure may iterate until the RMS error for each sub-region is lower than 1.5 mm, or the 0.05 degree spacing is reached.
Herein, a DEM may be a 90 metre resolution global digital elevation model provided by the Shuttle Radar Topography Mission.
Herein, a water body mask may be a 90 metre resolution grid data provided along with the DEM which masks out the water body, for example including oceans, lakes and inland rivers.
Figure 6 is a flow diagram of a local optimization process according to an example of the present disclosure.
Local optimization: The local optimization may start from an initial network infrastructure which may be, for example, either an existing GNSS network with an irregular distribution of historic stations, or a newly defined network on some key reference locations. The MODIS 5 km global resolution water vapour map may be used for reference. The procedure may comprise the following operations:
(i) Construct a Triangulated Irregular Network (TIN) using the initial network, i.e. using all existing stations (e.g. as at December 2016), where the PWV value for the station may be taken as that of the nearest 0.05 degree grid point;
(ii) Interpolate the PWV from the nodes of the TIN on to the 0.05 degree grid reference points, for example using the ITD model of Yu et al. (2017), which considers the DEM and applies the water body mask, to produce a simulated high resolution water vapour map;
(iii) Compute the RMS error (difference between the interpolated and reference PWV values) of all the grid points contained within each triangle of the TIN (i.e. high resolution water vapour bias map) and, if the RMS error exceeds 1.5 mm, insert a new station at its geometric centre;
(iv) Reconstruct the TIN based on the initial and added stations, and if the refined network is different from the previous one, repeat steps (i) to (iii) until the RMS error for each and every triangle is not greater than 1.5 mm.
Figure 7 is a flow diagram of an ITD model according to an example of the present disclosure.
ITD Model: The main procedure of ITD is to iteratively estimate the height scaling function and find the optimal scaling parameters. To start with, ITD supposes no turbulent delay exists so the initial scaling parameters can be obtained. The turbulent delay is then updated by the residuals ( Tk ) using the inverse distance weighting (IDW) function f, and is removed from the total delay to generate an updated turbulent delay (Tk+i). Several iterations of these steps are made until the difference between the updated turbulent delay and the previous turbulent delay is smaller than the threshold. The final output is generated using the final turbulent signal and scaled by the final height scaling parameter.
Herein, height scaling may refer to the use of an exponential function to simulate the ZTD:
ZTD = L0exp (- bH )
Where and b are the height scaling parameters, and H is the elevation. Herein, an Inverse Distance Weighting (IDW) function may be defined according to:
Figure imgf000017_0001
Where wij is determined by the distance between reference points / and j.
In certain examples, the interpolation for each station may be performed using available stations that fall within a certain region (e.g. having a certain size, shape and/or position) defined with respect to each station. For example, the interpolation for a station may be performed using all the available stations that fall within a certain shape (e.g. a circle) having a certain size (e.g. 150 km radius) centered at that station. In certain examples, at least three stations may be required for the interpolation. In the case of a local network, in certain examples, the interpolation may use a relatively large number of stations surrounding the point that is to be interpolated, not just the three vertices. On the other hand, in some examples, the interpolation may be performed based on stations at the vertices only, for example if the vertices are located far enough away from each other and/or from the station being interpolated, for example over a certain threshold (e.g. over 150 km). In the case of a global network, in certain examples, for each grid, the interpolation process may search for all the surrounding grids/stations within a certain shape having a certain size (e.g. a 150 km radius circle), not just the stations in its own grid (e.g. 80-km grid). In certain examples, the same interpolation procedure may be used in both cases.
Figure 8 illustrates the evolution of a GNSS network during operation of the techniques disclosed herein in the case of global optimization and local optimization. As shown, stations are iteratively added to grid elements of the network until the RMS value for a grid element is lower than a certain threshold.
In various examples, stations may be positioned on vertices of the grid cells and/or within (e.g. at the centre of) grid cells.
For example, in the global optimization case illustrated in Figure 8, a station may be located at the centre of each square (or cell or region) (i.e. each square represents one station). When a cell is densified, each cell (e.g. 0.8x0.8 degree cell) may be sub-divided into 2x2 sub-cells (e.g. 0.4 degree spacing). This may be continued, if necessary, using successively smaller spacing (e.g. 0.2 degree, then 0.1 degree, then 0.05 degree etc.) until the RMS error for the pixels (e.g. 16x16 pixels) is lower than a certain threshold (e.g. 1.5 mm), or until a certain termination criteria is satisfied (e.g. a certain degree spacing, such as 0.05 degree spacing, is reached).
On the other hand, in the local optimization case illustrated in Figure 8, at the beginning of an iteration a station is located on each vertex of each triangle (or cell or region), and when a cell is densified, a new station is added to the centre of the triangle, which becomes a vertex at the beginning of the next iteration.
The skilled person will appreciate that variation and modifications of the above examples are possible. For example, the skilled person will appreciate that one or more of the generalizations defined in the above Summary section of the present disclosure may be applied to any of the specific examples disclosed in the Description section of the present disclosure.
For example, cells may have any suitable regular or irregular shape (e.g. other than square or triangle, such as a hexagon); the cells may have different shapes and/or sizes; the cell array may have any suitable regular or irregular form (e.g. other than a regular square array or TIN); any suitable statistical value (e.g. other than an RMS value) and any suitable threshold value may be used to determine whether or not to densify a cell; receivers may be located on one or more vertices of a cell, on one or more edges of a cell, and/or in the interior (e.g. at the centre) of a cell; to densify a cell any suitable number of receivers (e.g. one receiver or more than one receiver) may be added to a cell in any suitable locations (e.g. within the cell (e.g. at the cell centre), on one or more cell boundaries (e.g. in the middle of a cell edge) and/or on one or more cell vertices) and in a uniform or non-uniform manner; the receivers used for interpolation may comprise one or more receivers of a cell and/or one or more receivers of one or more neighbouring cells; the interpolation may be perfomed in any suitable manner (e.g. other than using the ITD model); the iterative process may be terminated according to any suitable termination criteria (e.g. once a maximum number of iterations has been reached and/or once a certain minimum receive spacing is reached); the water vapour map may be obtained from any suitable source (e.g. other than MODIS data). The skilled person will readily contemplate other variations and modifications.
Results
Some results obtained using the examples described above will now be described.
Global Optimal Network Distribution. The optimal CORS station spacing globally with regards to tropospheric delay mitigation is shown in Figure 1 for both boreal summer and winter, and ranges from a defined, fixed upper limit of 80 km (the typical sparsest spacing between Network RTK CORS in many nations) and then where needed, increasingly dense halvings of 40, 20, 10 to 5 km. The optimal spacings were obtained by first thinning the PWV values from the MODIS 5 km global resolution reference data set to 80 km (0.8 degrees) resolution then, as described in the Methods section above, interpolating using the ITD model of Yu et al. (2017) and cross-validating against the 5 km reference values. 3-month averages of MODIS near-infrared PWV were used to reduce biases caused by cloud-induced data gaps and daily and weekly weather variations, with 1 June 2016 to 31 August 2016 selected to represent boreal summer and 1 December 2015 to 29 February 2016 for boreal winter. The spacings needed were computed on an 80x80 km window- by-window basis across the global land masses. Within the window considered, if the RMS of the differences between the interpolated and reference PWV values of all 0.05 degree MODIS pixels exceeded 1.5 mm, the spacing was densified in that window until the RMS reduced to below 1.5 mm. This PWV error cut-off was chosen as it is equivalent to a zenith tropospheric delay error of about 1 cm (Bevis M., Businger, S., Chiswell, S., Herring, T, Anthes, R., Rocken, C. & Ware, R. H. GPS meteorology: Mapping zenith wet delays onto precipitable water. Journal of Applied Meteorology 33, 379-386 (1994)), which maps to a GNSS height error of about 2-3 cm (Brunner, F. K. & Welsch, W. M. Effect of the troposphere on GPS measurements. GPS World 4, 42-51 (1993)), the typical precision of RTK positioning (e.g., Nickitopoulou et al., 2006; Perez-Ruiz et al., 2012; He, H., Li, J., Yang, Y, Xu, J., Guo, H. & Wang, A. Performance assessment of single- and dual frequency BeiDou/GPS single-epoch kinematic positioning. GPS Solutions 18, 393-403 (2014)).
Clear seasonal variations are obtained in the optimal station spacing, with a denser network being required in summer than in winter. The densest station coverage of 5-10 km is needed in the Andes, Central America, the Sierra Pacaraima, West Africa and the East Indies throughout the year. Other regions with very dense spacings of 5-10 km in summer are India and south-east Asia (Burma, Thailand, Vietnam, Cambodia, Laos, southern China, and southern Japan), as well as parts of the Alps. Over much of Asia, Europe, North America, northern and southern Africa, southern South America and Australia, the sparsest station spacing of 80 km is sufficient. The average station spacing globally (defined as being the land masses from 60 degrees north to 60 degrees south) is 52 km in local summer and 70 km in local winter. The spacing, in general but not exclusively, tends to be denser in the tropics (average spacing of 46 km) than in cooler regions.
To understand the cause of the variation in the optimal station spacing, included in Figure 1 are maps of 3-monthly average PWV per 0.05 degree MODIS pixel. It can be seen that many of the areas with the densest station spacing occur where the PWV content is very high (approaching 10 cm), notably in the tropics and across India and south-east Asia in local summer. However, there are also regions which have a high PWV content, such as central Africa and the Amazon basin, yet the station spacing needed to accurately interpolate PWV is very sparse. By comparing the maps of the station spacing with maps of PWV variation (also shown in Figure 1), defined as the standard deviation of all 0.05 degree grid PWV values within an 80x80 km window centred on each 0.05 degree MODIS pixel considered, it can be seen that it is the PWV spatial variation that has more influence on the station spacing than the PWV magnitude itself. The influence of PWV variation on station spacing is confirmed from inspection of the correlation plots also shown in Figure 1. For example, high correlation values between station spacing and PWV variation are seen around the Andes, west Africa and south-east Asia, and the Alps in local summer. Meanwhile, for the Amazon basin there is little correlation because the PWV variation across the region is small, and the station spacing is sparse.
Furthermore, in Figure 2 plots of PWV and PWV variation (normalized linearly to [0, 1] against the respective global minimum and maximum) against station spacing are shown, collated for all global values for both local summer and winter. It can be seen that there is a direct relationship between the normalized PWV variation and station spacing for both summer and winter, but not so between the normalized PWV and station spacing. One of the principal causes of large spatial variations in PWV above the ground at a particular location, and hence the dense station spacing, is the topography. The Andes, located at a similar latitude to the Amazon basin, exhibit rapid and large changes in topography over distances of only a few kilometres, and therefore rapid PWV variations because PWV content scales with height. Similar effects are seen over the Alps in local summer (but not in adjacent flatter regions such as France and Germany), and in particular in Nepal and Bhutan where the Himalayas meet the Indian plain. The terrain is flatter across much of India, but the monsoon climate means extreme variations in PWV can rapidly occur, and hence the need for a dense station spacing. However, in local winter in the same parts of India, only the minimal sparse spacing of 80 km is needed, as the PWV variations are only moderate (0.14 cm). Although Figure 1 indicates that a dense station spacing is needed in north-east Canada in local summer, the large PWV variations there are not considered reliable because the MODIS PWV values are degraded by the large number of lakes, whose average size is several kilometres and close to the MODIS spatial resolution (Gao and Kaufman, 2003). In local winter however, both the average PWV content (-0.4 cm) and its variation (-0.08 cm) are low, meaning that this effect can be negligible and the maximum 80 km spacing is always sufficient.
Local GNSS Network Enhancement. The optimal station spacing, for example with regards to tropospheric delay errors, globally has been considered by assuming a regular station grid and without considering existing geodetic GNSS infrastructure. In practice, there are a large number of CORS stations already existing and, while in some cases they have been designed to form as regular a station spacing as possible (e.g., about 60-80 km in Britain), invariably they are irregularly spaced, and designers do not seem to have considered topographic or climatic variations. Consequently, existing networks do not necessarily provide optimal PWV mitigation performance. Now considered are some examples of how existing GNSS networks in different parts of the world could be densified to achieve optimal PWV interpolation performance with the fewest additional stations, without having to redesign the entire network. These are indicative examples based on existing CORS infrastructure: only in some cases do they currently underpin commercial and/or free Network RTK services. Nevertheless, they provide indications of how infrastructure may be densified if all stations were able to be used for Network RTK services.
Figure 3 illustrates original and densified (optimally refined) regional networks during local summer of 2016: July 2016 for the British Isles and Italy, January 2016 for New Zealand and Java. The first and second columns represent the original and the densified networks, respectively, with the shaded triangles indicating areas where the interpolation RMS exceeds 1.5 mm. The background underneath the TIN of the first two columns shows the interpolation differences between the reference MODIS water vapour map and the interpolated map using the reference stations. The third column shows the topography and the additional stations inserted to obtain the optimal spacing.
Figure 4 illustrates time series of optimal station spacing from January 2012 to December 2016 for each of the British Isles, Italy, New Zealand and Java networks. Shown are the resultant monthly average station spacing and the number of stations added to the original network to be optimal. Dots denote the actual station numbers per month per network. Solid lines show the best fitting annual sinusoid per time series.
The first example considered is the British Isles (www.bigf.ac.uk), as the bulk of the network comprises national CORS geodetic infrastructure which is used for Network RTK service provision and has a near regular station spacing of about 60-80 km (Figure 3). It exhibits moderate PWV variations and fairly flat topography in England, but with more mountainous terrain in Scotland and Wales, albeit only rising to a maximum elevation of about 1300 m.
The second example considers Italy’s RING network (http://ring.gm.ingv.it/). which has stations clustered as densely as every 25 km in the Appennines, but in the north the distribution is much sparser at around 80 km, especially in and close to the Alps where particularly large PWV variations occur in the summer.
The third example is New Zealand’s GeoNet (https://www.geonet.orQ.nz)· which incorporates dense station spacing (~25 km) in the east of the north island. However, spacing is much sparser (~90 km) across the Southern Alps in the south island which exhibit a large range in topography from about 3500 m elevation mountains to plains close to sea level.
The fourth example is for western Java, which has very large PWV variations, large topographic variations from sea level to -2000 m elevation and, from Figure 1 , very dense suggested station spacings of towards 5 km being needed. It is also chosen as an area with few existing CORS and as Network RTK services are not apparently operating here, instead of showing how an existing network could best be densified, an initial hypothetical network with a regular 80 km spacing is assumed, and then densified.
The MODIS 0.05 degree grid point PWV values were averaged per month and CORS station PWV values were taken as those of the nearest grid point. A Triangulated Irregular Network (TIN), based on Delaunay triangulation, was formed between all existing CORS. The PWV values from the CORS locations were then interpolated on to the 5 km PWV MODIS values and the RMS difference (error) with respect to the reference values computed per triangle. If the RMS error in any triangle exceeded 1.5 mm, extra stations were added (see Methods section above) until all triangles across the network had RMS errors of under 1.5 mm.
The original (existing) and densified networks are shown in Figure 3, considered for local summer (January or July 2016) for each network, when PWV variations are greatest.
For the British Isles, 144 stations were added to the original 108. However, across most of England, very few additional stations were needed with the original 60-80 km spacing almost sufficing. The additional stations are mainly in the mountainous areas of Wales and western Scotland, and also in Ireland where the assumed network was less dense originally. The new stations reduce the average inter-station spacing from 65 km to 48 km, with an overall PWV RMS error reduction of 2.4 mm to 1.4 mm. In Italy, a very dense network of about 34 km spacing is needed in the Alps compared with the existing 80 km spacing, but this has the effect of reducing the overall RMS error from 2.5 mm to 1.2 mm. Few additional stations are needed across the Appennines, as the PWV variations induced by the topographic variations are already well mitigated by the dense original network.
Similar to Italy, for New Zealand many (108) new stations (now with average spacing 40 km) are needed in the south island to obtain the optimal RMS error, while in the north island those added tend to be in the flatter areas of the west coast but where there is a dearth of existing stations.
For Java, a very dense network of 5 km station spacing is needed in many areas, but it can be seen from Figure 3 that there are local variations in the spacing needed, including across parts of the mountains and adjacent flatter areas where PWV spatial variations are lower.
A seasonal variation in the global optimal station spacing is suggested from the summer and winter PWV variation and station spacing maps shown in Figure 1. To investigate the extent of such seasonal variations in optimal station spacing for the local network examples considered, the Delaunay triangulation optimal densification of the assumed existing networks was repeated for each month for five years from 2012 through to 2016.
The number of additional stations needing to be inserted per network per month are shown in Figure 4, together with the resultant average station spacing across the network. Clear periodic, seasonal variations in the number of added stations and resultant average station spacing needed for optimal PWV interpolation are seen. The amplitude of the seasonal variation is dependent on the type of climate, with mountainous climates such as Italy exhibiting an amplitude of about 140 stations, while the tropical climate of Java has an amplitude of only 30 stations. The amplitude also provides an indication of how close the original network is to being optimal for Network RTK provision, as a small amplitude denotes that few stations need to be added, although it could also indicate that in local winter the original network is denser than required in places and could potentially be thinned.
Discussion and Conclusions
A framework has been disclosed to optimize GNSS CORS network distributions with regard to the distribution of the variations in PWV and topography, for example in order to help facilitate centimetre-level homogeneous accuracy for Network RTK positioning. The optimal station spacings were found globally using MODIS PWV values as reference, and cross-validating after interpolation from station spacings ranging from 80 km to 5 km. More detailed local examples have also been presented, demonstrating how assumed original existing networks could be densified through a Delaunay triangulation approach in order to optimize the interpolation of PWV errors, with a clear seasonal variation demonstrated in the optimal network.
The station spacing has been shown to depend on variations in PWV across the network considered, rather than the total amount of PWV. However, the largest variations, and hence the densest station spacings, tend to be in the tropical regions which have large PWV content, although topographic effects also have a major influence. For example, where large variations in topography arise such as in the Alps, a dense station spacing is needed. Local extreme weather also contributes, such as the monsoon climate of India which requires a very dense network in the summer despite being fairly flat, but only a fairly sparse network is needed in winter.
The framework has been disclosed for the optimization of Network RTK GNSS positioning with regards to mitigating tropospheric delay errors. However, the resulting optimal CORS network designs may also be important for repeat-pass Interferometic SAR (InSAR) data processing and applications, providing indications of whether existing GNSS stations will be able to provide sufficient quality corrections to reduce atmospheric effects on InSAR observations and hence improve estimates of post- and inter-seismic crustal movements (e.g. Yu etal., 2018a). Furthermore, the ability to generate precise maps of PWV from GNSS CORS data may also be useful in meteorology in the nowcasting of rainfall (e.g., Bonafoni, S. & Biondi, R. The usefulness of the Global Navigation Satellite Systems (GNSS) in the analysis of precipitation events. Atmospheric Research 167, 15-23 (2016)) and assimilation into numerical weather models for precipitation forecasting (e.g., Rohm, W., Yuan, Y, Biadeglgne, B., Zhang, K. & Le Marshall, J. Ground-based GNSS ZTD/IWV estimation system for numerical weather prediction in challenging weather conditions. Atmospheric Research 138, 414-426 (2014); Yan, X., Ducrocq, V., Poll, P, Hakam, M., Jaubert, G. & Walpersdorf, A. Impact of GPS zenith delay assimilation on convective-scale prediction of Mediterranean heavy rainfall. Journal of Geophysical Research: Atmospheres 114, D03104 (2009)). Also, the framework concept disclosed herein may be employed to optimize other monitoring networks, e.g. those for soil moisture, air pollution, and tree species monitoring, although their own dominating error source should be identified and the corresponding spatial and temporal variations should be considered. In addition, there might be multiple dominating error sources in some cases, whereby the framework concept should be further developed to deal with multiple variables.
Data availability
The MODIS data used in the examples disclosed herein were obtained from the Land Processes Distributed Active Archive Centre (LP DAAC), managed by the NASA Earth Science Data and Information System (ESDIS) project.
In certain examples of the present disclosure, in relation to the background water vapour map used in global optimization, MODIS near-IR water vapour products may be used to generate a global background dataset. The skilled person will appreciate that the present disclosure is not limited to this example. For example, High-resolution European Centre for Medium-Range Weather Forecasts (ECMWF), or any other suitable high-resolution global numerical weather models, can also be used to generate such a background dataset.
In certain examples of the present disclosure, in relation to the background water vapour map used in regional/local optimization, to optimize regional or local GNSS networks, any suitable high- resolution regional or local numerical weather models (e.g. Weather Research and Forecasting model; Unified Model) can be used to generate the background dataset.
In certain examples of the present disclosure, one or more boundaries may be defined. For example, specific areas (e.g. countries, provinces, counties, etc.) may be defined and optimization of a GNSS network may be performed on one or more specific defined areas. In other examples, specific areas (e.g. water bodies and special infrastructure) may be masked out while optimizing a GNSS network. For example, in the algorithms described above, calculations corresponding to positions falling within a masked area may be skipped, and insertion of new reference stations inside a masked area may be prohibited.
In other examples, a GNSS network may be optimized for monitoring the stability of linear infrastructure (e.g. a high-speed railway line or long sea-crossing bridge).
In certain examples of the present disclosure, the stations may have fixed positions. However, in other examples, one or more or all of the stations may be moveable. In addition, in certain examples of the present invention, the stations may be positioned onshore. However, in other examples, one or more or all of the stations may be positioned offshore. For example, the techniques described above may be used to optimize a GNSS buoy network, or to optimize the distribution of GNSS stations on oil platforms. Certain examples may use a combination of onshore and offshore stations, which may help reduce edge effects (interpolation versus extrapolation of errors) that can result from the lack of stations on one side of a boundary or area.
In certain examples of the present disclosure, stations may be deployed for a relatively long period of time (e.g. they may be permanent stations). However, in other examples, one or more or all of the stations may be deployed for a relatively short period of time (e.g. they may be temporary stations). Temporary stations may be used for campaign design, for example for monitoring of a particular process such as plate tectonics and certain engineering projects.
In certain examples of the present disclosure, the stations may all be the same. However, in other examples the stations may vary in their quality or performance (e.g. accuracy, sensitivity, etc.), and the construction and/or optimization of a GNSS network may be performed in consideration of the quality and/or performance of the stations.
Network performance in some regions (e.g. regions except the tropics) may be subject to temporal (e.g. seasonal) variations. Hence there may be a need for a different optimal network at different times (e.g. from season to season). In this case, the optimization of a GNSS network may be adjusted over time to account for such variations in performance. For example, local densification of the network may be performed during certain times of the year.
By applying the algorithms described above, the total number of stations and the precise positions of those stations may be obtained for an optimized network. This detailed information may be processed to obtain different types of output information (for example of varying levels of detail and/or accuracy). For example, a first level report may show only relatively general statistics, a second level report may provide information on the approximate positions or distribution of stations, and a third level report may provide detailed information giving accurate positions of the stations. For example, rather than specify the exact position of a station, an approximate position of the station may be given. As another example, rather than specify the positions of stations, information specifying the station density (e.g. the number of stations to be positioned within a certain area) may be given.
While the invention has been shown and described with reference to certain examples, it will be understood by those skilled in the art that various changes in form and detail may be made therein without departing from the scope of the invention, as defined by the appended claims.
Annex (Yu et al. (2017))
Generation of real-time mode high resolution water vapor fields from GPS observations
Chen Yu, Nigel T. Penna, and Zhenhong Li
School of Civil Engineering and Geosciences, Newcastle University, Newcastle upon Tyne, UK.
Corresponding author: Zhenhong Li, School of Civil Engineering and Geosciences, Newcastle University, Newcastle upon Tyne, NE1 7RU, UK. (Zhenhong.Li@newcastle.ac.uk)
Key points:
1. An Iterative Tropospheric Decomposition model for zero differenced ZTD interpolation;
2. Generation of high resolution real-time mode GPS-based ZTD and PWV maps;
3. Validation of ITD PWV maps with the MODIS near-IR water vapor product.
Summary
Pointwise GPS measurements of tropospheric zenith total delay can be interpolated to provide high resolution water vapor maps which may be used for correcting SAR images, for numeral weather prediction and for correcting Network Real-time Kinematic GPS observations. Several previous studies have addressed the importance of the elevation dependency of water vapor, but it is often a challenge to separate elevation-dependent tropospheric delays from turbulent components. In this paper, we present an iterative tropospheric decomposition interpolation model that decouples the elevation and turbulent tropospheric delay components. For a 150 km x 150 km California study region, we estimate real-time mode zenith total delays at 41 GPS stations over 1 year using the precise point positioning technique, and demonstrate that the decoupled interpolation model generates improved high resolution tropospheric delay maps compared with previous tropospheric turbulence and elevation dependent models. Cross validation of the GPS zenith total delays yields an RMS error of 4.6 mm with the decoupled interpolation model, compared with 8.4 mm with the previous model. On converting the GPS zenith wet delays to precipitable water vapor and interpolating to 1 km grid cells across the region, validations with the MODIS near-IR water vapor product show 1.7 mm RMS differences using the decoupled model, compared with 2.0 mm for the previous interpolation model. Such results are obtained without differencing the tropospheric delays or water vapor estimates in time or space, whilst the errors are similar over flat and mountainous terrain, as well as for both inland and coastal areas.
1. Introduction
Time-varying two-dimensional precipitable water vapor (PWV) fields (maps) are used in meteorological nowcasting, including the identification of events dominated by horizontal advection [Benevides et al., 2015]; for assessing moisture transport in the lower troposphere [e.g., Mengistu Tsidu et al., 2015]; for relating humidity fields to precipitation events [e.g. Boniface et al., 2009]; for assessing the severity of tropical cyclones [e.g., Shoji et ai, 2011]; for assessing the impact of new assimilated observations for forecasting precipitation [e.g., Yan et ai, 2009] Such maps are also essential for correcting SAR images for atmospheric effects to enable small (and long-wavelength) geophysical signals to be measured, including interseismic strain accumulation and postseismic motion, observations of which not only give insight into the mechanics of a fault, but also play key roles in estimating the likelihood of future earthquakes [Wright et ai, 2004; Gourmelen andAmelung, 2005; Fialko, 2006; Walters et ai, 2013; Wen et ai, 2012; Li et ai., 2009a] Dense PWV fields also enable GNSS Network Real-time Kinematic (RTK) observations to be corrected for signal delays due to water vapor on propagating from space through the Earth’s neutral atmosphere (‘troposphere’) to a ground-based receiver. Such corrections are essential for centimeter level positioning, particularly heights, and enable (subject to sufficient GNSS base station coverage) Network RTK to be used for geophysical and engineering applications that have normally only used local base station RTK, such as river channel mapping [e.g., Brasington et ai, 2003], glacier flow and debris mapping [e.g., De Paoli and Flowers, 2009], coastal erosion [e.g., Lee et ai, 2013], crustal deformation [e.g., Genrich and Bock., 2006] and structural [e.g., Im et ai, 2013] monitoring, precision farming [e.g., Perez-Ruiz et ai, 2011], embankment instability and landslide monitoring [e.g., Gili et ai, 2000; Zabuski et ai., 2015]
Two-dimensional PWV maps/fields may be obtained directly from the National Aeronautics and Space Administration (NASA) Moderate Resolution Imaging Spectroradiometer (MODIS) on board each of the Terra and Aqua satellites (launched in 1999 and 2002, respectively). They provide global coverage every 1-2 days with observations in 36 spectral bands at moderate resolution (250- 1000 m), five of which (near-infrared (IR)) are used for remote sensing of water vapor with 1-km spatial resolution over clear land areas, oceanic areas with Sun glint and/or above clouds over both land and ocean [Gao and Kaufman, 2003] MODIS swaths do not coincide in time with images from SAR satellites such as ERS-1/2 (nominally 30 min of time difference), Envisat (nominally 30 min of time difference), or Sentinel-1 (nominally 4.5 hours of time difference), and therefore for SAR atmospheric corrections temporal interpolation errors will always prevail. The Envisat SAR satellite included the Medium Resolution Imaging Spectrometer (MERIS), a passive push-broom imaging instrument measuring the solar radiation reflected from the Earth’s surface and clouds in the visible and near-IR spectral range during the daytime. MERIS provided near-IR water vapor products above land or ocean surfaces under cloud-free conditions [Bennartz and Fischer, 2001] or above the highest cloud level under cloudy conditions [Albert et ai, 2001] at two nominal spatial resolutions: 0.3 km for full-resolution mode and 1.2 km for reduced-resolution mode. Both MODIS and MERIS near-IR PWV measurements have been validated against radiosonde and GPS-based products at discrete points with 1 mm and 1.2 mm respective RMS agreements [Li et ai, 2009b; Li et ai, 2006b; Li et ai, 2005; Li et ai, 2003; Li et ai, 2006c] However, neither the MODIS nor MERIS instruments provide reliable PWV values above land or ocean surfaces when clouds are present, which restricts their use. Similarly, the Atmospheric Infrared Sounder (AIRS) instrument, also on-board the Aqua satellite, provides global coverage amounting to 324,000 PWV profiles with a vertical resolution of 2 k and a spatial resolution at nadir of 45 km, but also only operates in clear skies and in partly cloudy conditions [Aumann et al., 2003] Comparison between AIRS products and GPS PWV showed RMS differences of about 3-4 mm and 15% differences in the 2 km layers between AIRS and radiosondes [Divakarla et al., 2006; Raja et al., 2008] PWV fields may also be obtained from global atmospheric models (both reanalysis and rapid update cycle predictions), but suffer from their coarse spacing, low temporal resolution (e.g., 3 or 6 hours for the European Centre for Medium- Range Weather Forecasts (ECMWF) model) and failures during periods of atmospheric turbulence [Jolivet et al., 2014]
GPS provides a means of obtaining PWV with continuous temporal resolution without any cloud or weather dependence, albeit at discrete points where the GPS receivers are located. First, the time- varying tropospheric zenith total delay (ZTD) is estimated per receiver using the GPS data alone, and then the zenith hydrostatic ‘dry’ delay (which may be accurately modeled using surface pressure data) is subtracted to leave the zenith wet delay (ZWD), which is spatially and temporally much more variable. The ZWD may then be readily converted to PWV using estimates of the mean temperature of the atmosphere, based on surface temperature measurements [Sews et al., 1992] The pointwise GPS PWV measurements agree to those from radiosondes and microwave radiometers with standard deviations of about 1-2 mm [e.g., Ohtani and Naito., 2000; Liou et al., 2001 ; Li et al., 2003], and may then be interpolated to provide PWV fields, which has been attempted using different models. For example, Jarlemark and Emardson [1998] evaluated a gradient model, a linear regression in time model that ignored observational directions, and a turbulence model that yielded at least 10% improved RMS error than the other two models. Williams et al. [1998] used a structure function to model the water vapor variation in space, but with respect to a reference value, whilst Janssen et al. [2004] found that Inverse Distance Weighted (IDW) and Ordinary Kriging (OK) interpolation models perform better than spline interpolation, but also only considered double differenced ZTDs as they were considering InSAR atmospheric corrections only. A deficiency of all these models is that they did not consider the terrain elevation dependence of water vapor and hence the interpolated values may contain large errors in regions with highly varying topography [e.g., Walters et al., 2013] To deal with the water vapor stratification, Emardson and Johansson [1998] incorporated a height scaling function with a best linear unbiased estimator and suggested an interpolated ZWD accuracy of about 1 cm, but only one station in Sweden was considered and the height variation across the network was only about 200 m. Li et al. [2006a] proposed a GPS topography-dependent turbulence model (GTTM) but reported that interpolation models should be applied to ZTD/ZWD values differenced in time rather than the absolute ZTD/ZWD values themselves, as this can reduce the influence of topographic effects on the ZTD/ZWD variations. For interpolating undifferenced GPS ZWD point values, Onn and Zebker [2006] utilized a frozen flow hypothesis to model the water vapor variation in time. Then Xuetal. [2011] showed that incorporating this height scaling function approach with an interpolator model based on the estimator of simple Kriging with varying local means (we will refer to this model as SKIm+Onn) improved the ZWD interpolation RMS accuracy by 29% compared with using the Berrada Baby et al. [1988] semi- empirical height scaling function. In a different approach to account for variations with topography, Bekaert et al. [2015] employed an InSAR phase observation-based power law correction model which used a fixed reference at the relative top of the troposphere, and described how the phase delay varies with altitude. To separate deformation and tropospheric signals, a frequency band insensitive to deformation is required. Benevides et al. [2016] also attempted to constrain GPS PWV with InSAR-derived PWV maps containing the topography signal. However, these models did not take into account that the InSAR measurements themselves have uncertainties of up to several centimeters and are susceptible to not detecting geophysical signals such as volcano inflation/deflation and interseismic slip rate [Williams et al., 1998] Hence we consider the SKIm+Onn model to represent the current state-of-the-art for the generation of PWV maps.
Several previous studies have noticed the coupling effect of the tropospheric turbulence and terrain elevation dependency [e.g., Treuhaft and Lanyi, 1987; Li et al., 2006a; Xu et al., 2011 ; Benevides et al. 2016] However, in the presence of strong atmospheric turbulence, previous models are inadequate for correcting SAR images to be used for sub-centimeter level ground motion monitoring [e.g., Walters et al., 2013; Fattahi and Amelung, 2015], or for the highest Network RTK positional precisions when such variations are not eliminated by data differencing. The aim of this study is therefore to improve the accuracy of GPS interpolated tropospheric water vapor maps by accounting for the coupling effect of both the terrain elevation dependency and tropospheric turbulence, demonstrate this over varying terrain, and compare with the SKIm+Onn model. We generate temporally-continuous high spatial resolution maps using undifferenced (neither spatially nor temporally) ZTD values estimated pointwise at GPS stations using the precise point positioning processing (PPP) technique in real-time mode. Hence they are directly applicable for correcting SAR images, enabling rapid response to ground motion from large earthquakes and volcanoes; for removing tropospheric effects in Network RTK GPS processing; for meteorological nowcasting via potential assimilation in numerical weather models, noting that ZTD is the form of GPS-based water vapor observation used operationally by for example the Met Office [Bennitt and Jupp, 2012] and Meteo-France [PoU et al., 2007] We validate our interpolated absolute ZTD fields first by cross- validation with the GPS pointwise ZTDs, and then by converting to PWV and comparing with the MODIS near-IR water vapor product, which we choose due to its high spatial resolution and accuracy [Li et al., 2003] The MERIS near-IR water vapor product has similar quality with an even higher spatial resolution [Li et al., 2006c], but it is only available during the period from May 2002 and April 2012 before the International GNSS Real-time Service began (April 2013). We use data from all of 2015 (hence covering all four seasons) for a study region of California of about 150 km x 150 km, i.e. 37° 06' to 38° 30' N, 122° 40' to 121° 00' W, incorporating 41 Plate Boundary Observatory (PBO) GPS stations as shown in Figure 1. This region was selected as there is a dense spacing of continuous GPS stations whose data have been extensively analyzed previously (including post- processed ZTDs readily available for quality control), it is comparable to the size of a SAR image swath (e.g., 100 km for Envisat Stripmap mode and 250 km for Sentinel-1 Interferometric Wide Swath mode), there are large variations in topography, and it experiences relatively cloud-free conditions as needed for validation with MODIS.
2. Iterative Tropospheric Decomposition Model
Tropospheric delays, especially the part due to atmospheric water vapor, vary both vertically and laterally over short distances, and are often considered as the sum of (i) a stratified component highly correlated with topography which therefore delineates the vertical tropospheric profile, and (ii) a turbulent component resulting from disturbance processes (e.g., severe weather) in the troposphere which trigger uncertain patterns in space and time [Williams et al., 1998; Cho et al., 2003; Li et al., 2006a] We present an Iterative Tropospheric Decomposition (ITD) model to effectively separate the turbulent and elevation dependent ZTD components. The ITD model decouples the ZTD into a stratified delay and a turbulent delay, which enables the more accurate interpolation of dense ZTD fields from pointwise values from a set of GPS reference stations across a region. It is defined as:
Figure imgf000030_0001
where, for the ZTD at location k, T represents the turbulent component and xk is the station coordinate vector in the local topocentric coordinate system; the stratified component is represented with an exponential function with coefficient b in which L0 is, for the region considered, the stratified component delay at sea level, and hk = ( hk -
Figure imgf000030_0002
-h^) is the scaled height; sk represents the remaining unmodeled residual errors, including unmodeled stratified and turbulent signals. The turbulent component usually consists of medium-to-long wavelength signals that can be interpolated by an IDW method. If n GPS stations are used in the region considered, the IDW model reads as:
Figure imgf000030_0003
where wu denotes the interpolation coefficient; u and i are indices for the user and reference stations, respectively; dUi represents the horizontal distance from the user to reference station. Reference stations at distances more than 100 km from the user station show limited correlation [Emardson and Johansson, 1998] so are not used.
In order to decompose the ZTD into stratified and turbulent components, which can account for substantial amounts of the ZTD but behave very differently, an iterative separation procedure was used:
(i) The ZTDs from all GPS stations within the user to reference station decorrelated range limit are used to estimate initial values for the exponential coefficients b and L0 , assuming the turbulent component values in equation (1) are zero.
(ii) The residuals sk , which are the summation of the unmodeled errors and the turbulent component, are computed by subtracting per station the stratified delay (as modeled by the estimated exponential coefficients) from the ZTD.
(iii) The turbulent component, T in equation (1), is computed per reference station from the residuals sk using the I DW function wui given in equation (2):
Figure imgf000031_0001
(iv) The updated values for the turbulent component per reference station are subtracted from the ZTD per reference station in equation (1), and a new set of exponential coefficients are obtained.
(v) Steps (ii) to (iv) are repeated until the exponential cofficients b and L0 converge. The final outputs are the exponential coefficients for the decorrelated range limit considered, plus the turbulent delay component and residuals per reference station.
(vi) The ZTD at the location of interest is obtained by interpolation of the reference station turbulent component and residuals, and added to the stratified delay computed using the final values of the exponential coefficients.
It should be noted that the assumption of the ITD model is that the turbulent component obeys the IDW interpolation law and the stratified component obeys the exponential law and, importantly, that these two components are not tightly coupled together. We later show (section 6 and Figure 9) that the two components are indeed separable and the convergence state can be reached rapidly. Although here we present the ITD model by ZTDs in zero difference mode, it is also suitable for interpolating differenced ZTD or PWV/ZWD. Hence the input to the ITD model may be, as well as GPS-ZTDs, other water vapor products such as those from the ECMWF numerical weather model, MODIS, or their combination.
3. Real-time Mode GPS ZTD Estimation and Validation 3.1. Real-time Mode GPS Precise Point Positioning Method
The pointwise ZTD values for the 41 study region GPS stations used for the generation of the interpolated maps were estimated using real-time mode GPS PPP processing. We used a PPP software which is a highly self-modified version of RTKLIB (www.rtkslb.coni), employing an extended Kalman filter [Gelb, 1974] to estimate in the state vector the constant ambiguities and time varying receiver coordinates, receiver and satellite clocks (considered as white noise), whilst the ZWD was estimated as a random walk parameter as a correction to an a priori ZTD from the UNB3 global empirical model [Leandro etai, 2006], employing the Global Mapping Function [Boehm etai, 2006a], and east-west and north-south tropospheric gradients were estimated. We used the ionospherically- free pseudorange and carrier phase observables and applied absolute IGS satellite and receiver antenna phase center offset corrections. We also applied corrections for antenna phase wind up [Beyerie, 2009], relativistic effects [Kouba, 2009], pseudorange Differential Calibration Delays [Kouba, 2009], earth tide [McCarthy, 1996] and ocean tide loading effects using FES2004 coefficients obtained from http://holt.oso.chalmers.se/loadinq. Uncalibrated phase and pseudorange hardware delays were assumed to be absorbed by the (float) ambiguity parameters and estimated receiver clocks, respectively.
PPP relies on highly accurate satellite orbits and clocks [Zumberge et al., 1997], which are usually held fixed in post-processed PPP. For our real-time mode processing, we used fixed satellite orbits from the International GNSS Service (IGS) Real-time Service, which were generated by decoding the IGC01 solution streams (products.igs-ip.net), every 15 seconds to match the GPS observation sampling rate of the PBO stations used. However, the satellite clocks have unpredictable behavior which makes their real-time prediction challenging [Li etai., 2014], so we did not fix these to the real time product values but estimated corrections to them using satellite clock parameters with the Gundlich and Koch [2002] robust estimation method. Additional constraints were introduced to overcome the rank deficiency of the normal equations, namely:
Figure imgf000032_0001
res(dtk ) = dtk - dt k,RTS (5) where dtk RTS is the initial value of the satellite clock given by the real-time product and acts as a pseudo-observation; res(dtk ) and dtk are the satellite clock residual and value, respectively. The satellite clock parameters were estimated as white noise parameters with a sigma of 0.001 ns, and the error messages contained in the real-time satellite clock product were used to determine the weights of the pseudo-observations in equation (4). An iterative process was used to identify some clock outliers which were hence ignored or assigned less weight in subsequent iterations [Gundlich and Koch, 2002]
We processed the data from 1 January to 31 December 2015 per GPS station in daily, discrete 24 hour batches in real-time mode, with an elevation angle cut-off of 10°. The tropospheric delay was estimated every 5 minutes using a process noise of 5.0e-8 km/Vs. To enable the fastest PPP solution convergence and separation of the ambiguities from the other estimated parameters, which is particularly problematic when using real-time satellite orbits and clocks [e.g., Yao et al., 2014], we applied loose constraints to a priori receiver coordinate values obtained from the PBO GPS Station Position Time Series (http://pbo.unavco.org/data/gps). Nearly 70% of daily solutions converged within 30 minutes (convergence time here means, from the beginning epoch to an epoch whose horizontal component bias is less than 10 cm and height component bias is less than 15 cm, and the overall standard deviation of its next 20 consecutive epochs also satisfies this requirement), with 90% of daily solutions converging within 50 minutes. The results presented hereafter are based only on the ZTD values after convergence was attained, according to these criteria.
3.2. Validation of Real-time Mode GPS Pointwise ZTDs
To validate our real-time mode PPP (RTPPP) GPS ZTD estimates per station, we compared with post-processed 5 minute ‘truth’ values computed by the Geodesy Laboratory at Central Washington University using the NASA JPL / Caltech GIPSY software version 6.2 in PPP mode with fixed IGS final orbits and clocks. The values were obtained from ftp ^/data- out unayco.org/pub/producis/iroposphere/ and in their processing, the VMF1 gridded tropospheric mapping function was used together with ECMWF gridded a priori ZHDs and ZWDs [Boehm et al., 2006b], whilst estimating the ZWD and tropospheric gradients (east-west and north-south) with process noise values of 5.0e-8 km/Vs and 5.0e-9 km/Vs respectively. We computed the differences between our RTPPP and GIPSY ZTDs at the common 5 minute epochs, excluded all the outliers greater than three times the standard deviation, then for each station computed per day the mean difference and the Root Mean Square (RMS) of the differences to assess the quality of the real-time mode processing. These are shown for a sample day (2 September 2015) in Figure 1 , chosen as it is indicative of the median differences for all days of 2015. The mean of the per station differences across the whole network is 1.9 mm for the sample day, indicating that no significant systematic error exists between the RTPPP and GIPSY ZTDs, including stations in mountainous areas. The mean RMS is 10.1 mm and more than 80% of the stations have an RMS value smaller than 12 mm, which is deemed sufficient quality for assimilation into real-time meteorological models [Shoji et ai, 2011]
Figure 1. Mean (a) and RMS (b) differences per station between RTPPP and GIPSY ZTDs, estimated every 5 minutes on 2 September 2015 for all 41 GPS stations in the California study region. For all stations, the mean difference and the mean RMS difference are 1.9 mm and 10.1 mm respectively, with this day being indicative of the median differences for all days of 2015. The white area represents the Pacific Ocean and San Francisco Bay, and the background shows the elevation. Stations P177, P230 and S300 are labeled as they are considered in Figure 9.
Apart from considering the spatial distribution of the differences, it is also important to assess the RTPPP performance over time. We therefore in Figure 2 show the 5 minute RTPPP ZTDs plotted against the common epoch 5 minute GIPSY ZTDs, from all 41 stations for the entire year, and plot the differences as a histogram. A linear regression fit gave: GIPSY ZTD = 0.989(±0.002) c RTPPP ZTD + 0.024(±0.003), and the correlation coefficient between them was 0.99, demonstrating high consistency between the RTPPP and GIPSY ZTDs not just spatially but also temporally. About 82% of solutions show differences smaller than 15 mm with 73% below 12 mm. The RMS difference is 12.5 mm, commensurate with the spatial RMS difference and further indicating an RTPPP ZTD precision of about 1 cm, which is commensurate with previous real-time studies [e.g. Li et al., 2014; Ahmed et al., 2016]
Figure 2. Comparison between RTPPP and GIPSY ZTDs for all 41 GPS reference stations in the California study region, from 1 January to 31 December 2015 at an interval of 5 minutes (a) Correlation analysis with a linear model: GIPSY ZTD = Slope c RTPPP ZTD + Intercept, and (b) Histogram of the differences.
4. Cross Validation of Interpolated ZTD at GPS Stations
Cross validation was used to evaluate the performance of the ITD model for interpolating ZTDs, and compared with the SKIm+Onn model. Whilst Xu et al. [2011] applied the SKIm+Onn model to ZWD, as ZWD also dominates the spatiotemporal variations of ZTD (with the ZHD being readily determined with surface pressure) we may also apply it to ZTD interpolation. Since the dry and wet components are both important for some applications such as InSAR atmospheric corrections [Elliott et al., 2008; Jolivet et al., 2014], and the separation of ZWD from GPS derived ZTD will introduce additional uncertainties, it is better to use total delays rather than wet delays only. To provide an indication of the ZTD quality used for single epoch SAR corrections, and to be commensurate with subsequent MODIS validations (section 5), we adopted the approach of Xu etal. [2011] and used one ZTD value per day (that at 14:00 Pacific Standard Time, i.e. local time) per station for the whole of 2015. The RTPPP ZTDs were interpolated for each GPS station in turn, with the 40 other GPS station ZTD values providing the input. Hence for ITD, per epoch, U and b of equation (1) were estimated for the network and the turbulent component estimated per station. Validation was carried out by comparing the interpolated ZTDs with the RTPPP ZTD estimates themselves at 14:00 local time and computing the RMS difference and Mean Absolute Error (MAE) for all stations for each day. This was repeated using interpolation with the SKIm+Onn model.
Time series of the cross validated daily MAE and RMS differences from the ITD and SKIm+Onn ZTD interpolation models are shown in Figure 3. It is clear that the ITD model leads to both lower annual mean MAE and RMS difference values than SKIm+Onn, i.e. the MAE and RMS reduce from 6.2 mm to 3.2 mm and 7.4 mm to 4.1 mm, respectively. It can also be seen from Figure 3 that the improvement of ITD is greater than that of SKIm+Onn in colder seasons (e.g., between days 0 to 100 and 280 to 365), when the medium-to-long wavelength turbulent signals are greater than the short wavelength ones and can be effectively modeled by ITD. This in turn enables stratified delays to be separated from the turbulent component in the ITD iterative process and better determined using regression analysis. The performance of the two models is more similar in the summer (i.e., from around day of year 150 to 220), indicating that short wavelength water vapor effects are significant and variable and cannot be fully mitigated by either model. Figure 3. Time series of cross validated daily MAE and RMS differences, based on 14:00 Pacific Standard Time (i.e., local time) RTPPP interpolation at all 41 GPS stations using the ITD and SKIm+Onn models, and comparing with RTPPP ZTDs. The annual mean MAE values are 3.2 mm for ITD and 6.2 mm for SKIm+Onn, whilst the annual mean RMS differences are 4.1 mm for ITD and 7.4 mm for SKIm+Onn.
Cross comparisons of the daily interpolated ZTD values are shown in Figure 4 for both the ITD and SKIm+Onn interpolation models, for all 41 stations for all of 2015. As for the time series, substantial reductions in the scatter is observed for ITD compared with SKIm+Onn, i.e. the RMS difference decreases from 8.4 mm to 4.6 mm. Small correlation and linear fit improvements are also obtained with ITD compared with SKIm+Onn: the correlation coefficient increases from 0.98 to 0.99, the linear fit improves from a slope of 1.025 to 1.012 and the intercept from 0.026 to 0.004 m, as detailed in Table 1. Furthermore, the proportion of differences under (magnitude) 10 mm increases from 61% for SKIm+Onn to 89% for ITD, and increases from 32% to 69% for the proportion under 5 mm magnitude.
Figure 4. Cross validation of (a) ITD and (b) SKIm+Onn RTPPP ZTDs daily at 14:00 local time for 1 year for all the 41 GPS stations. Linear model (Initial RTPPP ZTD = Slope* Interpolated RTPPP ZTD + Intercept) was also applied.
When considered for all stations for the entire year, the ITD interpolation model has been shown to substantially improve on the SKIm+Onn model. However, in the cross validation, some large differences occurred (more than 2 cm magnitude for both models), which suggests that the interpolation result is influenced by the GPS reference station distribution. To investigate this, the annual MAE and RMS differences per station are plotted in Figure 5, for both the ITD and SKIm+Onn models. It can be seen that the smaller MAE and RMS differences occur where the station coverage is denser, but the SKIm+Onn MAE and RMS values show substantial degradation compared with ITD in the north-west of the region where there are fewer stations, e.g. ~12 mm MAE for SKIm+Onn compared with ~5 mm for ITD. Meanwhile, the largest RMS value for any station is only ~8 mm for ITD, improved from ~12 mm for SKIm+Onn. Such station distribution degradations are consistent with the interpolation law and can be mitigated as far as possible by using stations with as uniform a distribution as possible. In terms of terrain effects on the MAE and RMS, for the ITD model, stations in the mountainous areas show approximately comparable precision with those at lower altitudes, whereas with SKIm+Onn larger MAE and RMS values arise, and the same applies in coastal areas. We attribute this improvement to the height scaling model working under the assumption that the turbulent delays are insignificant and have short-wavelength compared with stratified delays, therefore the height scaling is easily biased when strong tropospheric turbulence with medium-to- long-wavelength signals occur. Figure 5. Spatial distribution of cross validation ZTD results, showing MAE and RMS of daily (14:00 local time) differences of interpolated versus RTPPP ZTDs, computed over all of 2015 per GPS station (a) MAE using ITD, with an overall mean of 3.6 mm; (b) MAE using SKIm+Onn, with overall mean of 6.1 mm; (c) RMS using ITD, with an overall RMS of 4.6 mm; (d) RMS using SKIm+Onn, with an overall RMS of 8.4 mm.
5. Validation of Interpolated PWV Maps With MODIS
To provide further validation of the ITD interpolation model and its improvement over the SKIm+Onn model, including a detailed spatial resolution assessment, and to provide an accuracy assessment with an independent data set, the RTPPP GPS pointwise ZTD values were converted to PWV, interpolated to 1 km pixels across the entire study region and compared with the MODIS near-IR PWV product. Pressure and temperature data at 5 minute temporal resolution from co-located meteorological sensors were available at four of the 41 GPS stations and obtained from unavo.org. These were supplemented by 10 meteorological stations which were located up to 10 km outside the study region. The meteorological data were first interpolated to all 41 stations using the Li et al. [2003] differential models and, according to their cross tests, the resulting pressure and temperature errors were less than 1 hPa and 2 K, respectively. The interpolated pressure measurements at each GPS station were used to directly compute ZHD using the Saastamoinen [1972] model and subtracted from the RTPPP ZTD estimates, with the resulting ZWD pointwise values converted to PWV using the Sews et al. [1992] model with the interpolated surface temperature measurements. To enable the comparisons, 1 year of MODIS Level 2 data from the Terra satellite were obtained across the study region, providing one PWV map at the Terra orbit track time of each day (during daytime, about 10:30 local time). The Level 2 data were generated at the 1-km spatial resolution of the MODIS instrument using the near-IR algorithm [Gao and Kaufman, 2003] About 30% of days had severe cloud conditions so we excluded them as only a few grids can be obtained. Areas with cloud conditions or above water were also masked and only the cloud free land areas were used in the comparison, with the percentage of data points available and used in the comparisons shown in Figure 6(e), with the average cloud-free pixels over the entire year being 35%.
As the first step, cross validation was carried out at each GPS reference station, using one GPS PWV per day per station, taken at the MODIS acquisition time. We used MODIS PWV as ‘truth’ values to validate the interpolated PWV at each GPS station (cross test), and the MODIS PWV values were averaged over boxes of 3 c 3 pixels centered on the GPS station’s location if the centered pixel was missing. Whilst in such cases the spatial resolution of the MODIS PWV maps was effectively reduced from 1 km x 1 km to 3 km x 3 km, this resolution is still much greater than the ~15 km GPS station spacing across the study region, so provided a compromise between resolution and maximizing the number of MODIS pixels for comparison. PWV cross comparisons for all daily values available for each GPS station for the whole of 2015 are shown in Figure 6 (a and b), with the RMS of the differences being 1.48 mm for the ITD model and 1.73 mm for SKIm+Onn, as also listed in Table 1. The ITD model also results in a better linear regression fit, with a slope of 0.97 and intercept of 0.33 m compared with respective values of 0.95 and 0.63 mm with the SKIm+Onn model, and the correlation coefficient increased from 0.97 to 0.98. The SKIm+Onn model also produced a greater mean difference than ITD during times of lower PWV content (i.e., 0-10 mm), with the more effective ITD modeling in the times with less water vapor (i.e., from mid-autumn to mid-spring) being consistent with the GPS cross validation results shown in Figure 3, i.e. the improved separation of medium to long wavelength turbulence signals and improved performance across mountainous areas.
The interpolated RTPPP PWV values and resulting maps were then compared spatially with the MODIS PWV maps. Results are shown for (RTPPP GPS minus MODIS) in Figure 7 for MODIS images acquired on 3 September and 19 November 2015, chosen as they are sample days which are virtually free from cloud conditions across the whole study region. The height for each grid cell was resampled by the 30 m Shuttle Radar Topography Mission DEM. Some large differences are visually apparent, mostly across the areas with frequent cloud masks and near San Francisco Bay, but MODIS PWVs above water areas also involve a different retrieval algorithm compared to those above the land, resulting in differences and discontinuities at the land edge. Furthermore, any values over water areas have been removed since PWVs above water (bay, lake or ocean) share different characteristics from PWVs over land areas which cannot be well-described by the interpolation model [Sobrino et al., 2003] On average, the mountainous areas give more negative differences than flat terrain, showing that MODIS tends to overestimate PWV compared with GPS with increasing altitude (i.e., small PWV values), as previously found by Li et al. [2003] It can also be seen in Figure 7 that the edge areas with fewer GPS stations produce larger differences than central areas, confirming as discussed in section 4 that improved GPS station coverage will improve the quality of interpolated PWV maps. The ITD model produces smoother difference maps than SKIm+Onn and has a lower percentage of large differences. ITD also performs much better than SKIm+Onn in coastal areas where the PWV is more changeable and gives more complicated turbulent components.
To further spatially assess the improvement of ITD RTPPP PWV values with MODIS over those from the SKIm+Onn model, differences between the RTPPP interpolated PWV and all MODIS pixels (cloud free) were computed across the study region bounded by the GPS stations for the entire year 2015 and for all weather conditions. The differences are shown in Figures 6c and 6d, for each of the ITD and SKIm+Onn models. In general, most PWV pairs are located along the 1 :1 line, implying a strong correlation between the RTPPP-based PWV and MODIS PWV maps. The SKIm+Onn PWV plot exhibits greater differences in comparison with the ITD PWV plot. The scatters are distributed rather unsymmetrically, especially when the PWV amounts fell between 5-15 mm (typically for mountainous areas) and 30-35 mm due to the lower scale factor (0.912 vs 0.934).
To illustrate the finer spatial detail of PWV and the performance of interpolated RTPPP PWV, in Figure 8 we plot ITD-based PWVs over mountainous areas, shown as (a) from 37° 09’ to 38° 30’ in latitude and -122° 12’ to -121° 00’ in longitude, and (b) from 37° 09’ to 37° 50’ in latitude and -122° 30’ to -121° 00’ in longitude, for 3 September 2015. In Figure 8 we include PWV profiles (smoothed using a tenth averaging window) for both ITD RTPPP and MODIS along lines of constant latitude and longitude, over both mountainous and flat areas, which enable detailed comparisons of the ITD PWV gradients with respect to topography. The ITD PWV profiles change in a similar tendency with MODIS and share similar gradients. The overall RMS difference between MODIS and ITD PWV for the eight profiles considered is 1.51 mm and the RMS differences for mountain (gray polygon in Figure 8) and flat areas are 1.57 mm and 1.47 mm, respectively. These agreements demonstrate that the ITD model is capable of retrieving detailed water vapor distributions over a wide region, thereby showing its potential application for monitoring local extreme weather events.
Figure 6. Cross validation of (a) ITD and (b) SKIm+Onn interpolated RTPPP GPS PWV with MODIS, using daily values at MODIS acquisition time on all 41 GPS stations for all of 2015. (c) and (d) display all the available pixels between MODIS PWV and ITD/SKIm+Onn PWV maps for year 2015. The color scale represents the density of occurrence and the red boxes highlight instances of particular improvement (see section 5). The daily cloud free MODIS PWV pixel density is displayed by (e) in which the vertical bar represents the total available pixel numbers divided by the maximum amount. The color scale represents the average daily PWV of all available pixels.
Figure 7. MODIS PWV and ITD RTPPP GPS PWV maps and RTPPP GPS minus MODIS PWV difference maps at 1 km spatial resolution, for both ITD and SKIm+Onn interpolations. (a,b,c,d) are for 3 September 2015 and (e,f,g,h) for 19 November 2015. (a) and (e) are MODIS maps, (b) and (f) are ITD maps, (c) and (g) are ITD difference maps, (d) and (h) are SKIm+Onn difference maps. The large differences (red pixels) in the north east of (c) and (d) are likely due to the presence of thin clouds which are not labeled in the cloud mask product.
Figure 8. ITD RTPPP map and MODIS PWV (red line) and ITD RTPPP PWV (blue line) profiles along certain latitudes and longitudes, after shifting a constant number, for 3 September 2015. The PWV profile series are in the same order as the line segments in the PWV map, and are averaged by a tenth average window. The gray polygon areas represent the mountain area. The overall RMS difference between MODIS and ITD PWV along the eight profiles is 1.51 mm and the RMS difference for the mountainous (gray polygon) and flat areas are 1.57 mm and 1.47 mm, respectively.
6. Discussion on Tropospheric Turbulence
The principal aim of the ITD model is to separate the elevation dependent ZTD/PWV component from the turbulent component, which is the most variable and uncertain part, and can easily bias the vertical ZTD scaling, making the separation of the two components challenging. Due to the constraints of the density of GPS stations, only medium-to-long wavelength turbulent signals are expected to be successfully modeled using ITD. To illustrate the size and variation of the turbulent component, and the importance of iterating the solution until convergence arises, a sample three GPS stations were considered: P177, P230 and S300 (Figure 1), which are in different parts of the study region, are at different elevations, and are at varying distances from the nearest other GPS reference stations. P177 is near the ocean, whilst S300 and P230 are in mountainous areas with elevations of ~500m and ~700m, respectively. Three epochs, from different seasons (spring, summer and autumn) were considered, and the variation of the turbulent component and its convergence with the number of iterations is illustrated in Figure 9.
It is clear from Figure 9 that the turbulent ZTD component can reach several centimeters and can be efficiently separated from the elevation dependent component. Although the first iteration enables the majority of the turbulent component to be determined, the subsequent iterations are needed for robust estimation. The far right hand column in Figure 9 further indicates the performance improvement with increasing number of iterations, with the RMS difference (computed through the RTPPP ZTD cross validation from all 41 stations at the corresponding epoch of each row) becoming smaller and tends towards convergence. Around six iterations are typically needed for convergence, after which sub-millimeter RMS changes arise, with improvements of typically less than 1 mm.
As the most important feature of the ITD model, the convergence tendency in Figure 9 reveals that the turbulence effect can be reduced by separating the stratified and turbulent components through iteration. The decoupling procedure acts very similarly to robust estimation, in which the ZTDs from stations exhibiting strong turbulence will contribute less in height scaling but account for more in the turbulent delay interpolation. In this way, the systematic patterns of turbulence resulting from local weather conditions or topography [Betts, 2001 ; Cho et al., 2003] can be better modeled. The iteration also allows for the detection of ZTD outliers, which is a not uncommon occurrence in real-time PPP due to the unpredictable behavior of the satellite clocks. Consequently, the ITD model enables both fitting of the tropospheric vertical profiles and also models the turbulence processes. As the fourth column in Figure 9 suggests, successfully accounting for this results in the overall RMS difference of the cross validation test reducing and converging.
Figure 9. RTPPP ZTD turbulent component separated by ITD at each iteration step. The first, second and third columns represent stations P177, P230 and S300, respectively, and the fourth column represents the ZTD cross validation RMS difference for all 41 stations on the corresponding day. Shown for sample days in each of spring, summer and autumn.
7. Conclusions and Outlook
In this study, an iterative tropospheric decomposition model has been developed for the generation of high resolution regional tropospheric PWV fields (maps) from the interpolation of pointwise GPS ZTDs, without any data differencing. For a California study region of around 150 km x 150 km, the approach of decoupling the terrain elevation dependency and the tropospheric turbulence contribution to ZTD in an iterative procedure (typically 4-6 iterations were required) led to improved accuracy interpolated tropospheric water vapor maps over those based on previous studies, such as the tropospheric turbulence and elevation dependent model SKIm+Onn of Xu et al. [2011] To be applicable to not only post-processed SAR atmospheric corrections, i.e. to also facilitate SAR for rapid response to monitoring earthquakes and volcanoes, Network RTK positioning and meteorological forecasting, we used real-time mode PPP GPS ZTD values estimated every 5 minutes (which were validated with post-processed GIPSY ZTDs with an overall RMS difference of 12.5 mm for all 41 stations for all of 2015) to generate the tropospheric maps. Cross validation of the GPS ZTD values resulted in 4.6 mm RMS differences using the ITD model compared with 8.4 mm using the SKIm+Onn model, using one value per station per day (14:00 local time) for all of 2015. Whereas the SKIm+Onn interpolation model has degraded performance over mountainous areas, the cross validation ITD model RMS and mean differences are similar for both mountainous and flatter terrain, and also similar for both coastal and inland areas. The cross validation improvements using ITD are smallest in the summer months. Spatially, we generated PWV values for 1 km pixels for all land-covered parts of the region and compared with daily MODIS PWV near-IR product values, with the RMS difference for the year being improved from 1.96 mm using the SKIm+Onn model to 1.71 mm using ITD. Furthermore, the spatial PWV gradients using the ITD model and MODIS across a variety of topography were nearly identical to each other. The overall RMS difference between MODIS and ITD PWV profiles is 1.51 mm and the RMS differences for mountain and flat areas are 1.57 mm and 1.47 mm, respectively. Hence the ITD PWV fields are also able to reveal detailed water vapor information over varying terrains.
To provide an indication of the potential of the real-time mode ITD model interpolated tropospheric maps for meteorological and geodetic applications, including revealing detailed information of local weather processes, we show in Figure 10 the detailed 2-hourly PWV information during a rainfall process over the study region on 2 November 2015 (10 am to 8 pm UTC). Arrows represent the PWV increasing (upwards) or decreasing (downwards) during each preceding two hours. One important fact is that the PWV over mountainous areas decreased during the whole rainfall process but other areas experienced increasing and decreasing PWV before and after the rainfall, respectively. We do not explain the patterns further as the focus of this paper is on showing the quality of RTPPP ZTD and PWV maps.
Figure 10. RTPPP PWV fields across the California study region every 2 hours on 2 November 2015 during a rainfall process from (a) 10:00 to (f) 20:00 UTC. Arrows represent PWV increasing (upwards) or decreasing (downwards) during the preceding 2 hours.
The generated spatially-dense PWV fields with continuous, high (5 minute) temporal resolution are not only suitable for correcting atmospheric effects in SAR images at the instant of acquisition, but they also will ensure the identification of water vapor variation from ground motion between image acquisition times [Foster et al., 2006] What is more, the high performance of the dense PWV maps using the ITD model is especially useful for mitigating the effects of water vapor for SAR measurements in mountainous areas, which usually suffer from vertical stratification and turbulent mixing due to the orography [Wadge et al., 2002]
In the near future, more research is planned to investigate the capacity of high resolution real-time tropospheric products using the ITD interpolation model. This includes the testing of the ITD model on ECMWF ZTD, MODIS and radiosonde PWV in addition to the GPS interpolations considered here. The abundant data sources will contribute to generate global coverage, all weather condition dense water vapor fields and to further investigate detailed tropospheric turbulence characteristics and processes.
Acknowledgements
The GPS data sets and orbit/clock products used in this paper were gratefully obtained from the PBO (http://pbo.unavco.org/data/qps) and the IGS Real-time Service, respectively. We also acknowledge UNAVCO for the operation of the PBO and for providing high accuracy station coordinate time series and tropospheric products, along with Central Washington University. Also thanks to NASA for MODIS water vapor data. This work was supported by a Chinese Scholarship Council studentship awarded to Chen Yu. Part of this work was also supported by the UK Natural Environmental Research Council (NERC) through the Centre for the Observation and Modelling of Earthquakes, Volcanoes and Tectonics (COMET, ref.: come30001) and the LICS and CEDRRIC projects (ref. NE/K010794/1 and NE/N012151/1, respectively), and the ESA-MOST DRAGON-4 projects (ref. 32244). The figures were generated using the Generic Mapping Tools software [Wessel et ai, 2013]
References
Albert, P., R. Bennartz, and J. Fischer (2001), Remote sensing of atmospheric water vapor from backscattered sunlight in cloudy atmospheres, J. Atmos. Ocean. Tech., 18, 865-874.
Ahmed, F., P. Vaclavovic, F. N. Teferle, J. Dousa, R. Bingley, and D. Laurichesse (2016), Comparative analysis of real-time precise point positioning zenith total delay estimates, GPS Solut. , 20(2), 187-199.
Aumann, H. H., M. T. Chahine, C. Gautier, M. D. Goldberg, E. Kalnay, L. M. McMillin, H. Revercomb, P. W. Rosenkranz, W. L. Smith, D. H. Staelin, L. L. Strow, and J. Susskind (2003), AIRS/AMSU/HSB on the Aqua mission: design, science objectives, data products, and processing systems, IEEE T. Geosci. Remote., 41, 253-264.
Bekaert, D. P. S., A. Hooper, and T.J. Wright (2015), A spatially variable power law tropospheric correction technique for InSAR data, J. Geophys. Res., 120(2), 1345-1356.
Benevides, P, J. Catalao, and P. M. A. Miranda (2015), On the inclusion of GPS precipitable water vapour in the nowcasting of rainfall, Nat. Hazards Earth Syst. Sci., 15, 2605-2616.
Benevides, R, G. Nico, J. Catalao, and P. M. A. Miranda (2016), Bridging InSAR and GPS tomography: A new differential geometrical constraint, IEEE T. Geosci. Remote., 54(2), 697-702.
Bennartz, R., and J. Fischer (2001), Retrieval of columnar water vapour over land from back- scattered solar radiation using the Medium Resolution Imaging Spectrometer (MERIS), Remote Sens. Environ., 78, 271-280.
Bennitt, G. V, and A. Jupp (2012), Operational assimilation of GPS zenith todal delay observations into the Met Office numerical weather prediction models, Mon. Wea. Rev., 140, 2706-2719.
Berrada Baby, H., P Gole, and J. Lavergnat (1988), A model for the tropospheric excess path length of radio waves from surface meteorological measurements, Radio Sci., 23(6), 1023-1038.
Betts, R. A. (2001), Biogeophysical impacts of land use on present-day climate: Near-surface temperature change and radiative forcing, Atmospheric Science Letters, 2(1-4), 39-51.
Bevis, M., S. Businger, T. Herring, C. Rocken, R. Anthes, and R. Ware (1992), GPS meteorology: Remote sensing of atmospheric water vapor using the Global Positioning System, J. Geophys. Res., 97(D14), 15787-15801.
Beyerle, G. (2009), Carrier phase wind-up in GPS reflectometry, GPS Solut., 13(3), 191-198.
Boehm, J., A. Niell, P Tregoning, and H. Schuh (2006a), Global Mapping Function (GMF): A new empirical mapping function based on numerical weather model data, Geophys. Res. Lett., 33, L07304.
Boehm, J., B. Werl, and H. Schuh (2006b), Troposphere mapping functions for GPS and very long baseline interferometry from European Centre for Medium-Range Weather Forecasts operational analysis data, J. Geophys. Res., 111, B02406.
Boniface, K., V. Ducrocq, G. Jaubert, X. Yan, P Brousseau, F. Masson, C. Champollion, J. Ch'ery, and E. Doerflinger (2009), Impact of high-resolution data assimilation of GPS zenith delay on Mediterranean heavy rainfall forecasting, Ann. Geophys., 27, 2739-2753.
Brasington, J., J. Langham, and B. Rumsby (2003), Methodological sensitivity of morphometric estimates of coarse fluvial sediment transport, Geomorphology, 53, 299-316.
Cho, J. Y. N., R. E. Newell, B. E. Anderson, J. D. W. Barrick, and K. L. Thornhill (2003), Characterizations of tropospheric turbulence and stability layers from aircraft observations, J. Geophys. Res.: Atmos., 108(D20), 8784, .
De Paoli, L., and G. E. Flowers (2009), Dynamics of a small surge-type glacier using one dimensional geophysical inversion. J. Glaciol., 55(194), 1101-1112.
Divakarla, M. G., C. D. Barnet, M. D. Goldberg, L. M. McMillin, E. Maddy, W. Wolf, L. Zhou, and X. Liu (2006), Validation of Atmospheric Infrared Sounder temperature and water vapor retrievals with matched radiosonde measurements and forecasts, J. Geophys. Res.: Atmos., 111 , D09515.
Elliott, J. R., J. Biggs, B. Parsons, and T. J. Wright (2008), InSAR slip rate determination on the Altyn Tagh Fault, northern Tibet, in the presence of topographically correlated atmospheric delays, Geophys. Res. Lett., 35(12).
Emardson, T. R., and J. M. Johansson (1998), Spatial interpolation of the atmospheric water vapor content between sites in a ground-based GPS Network, Geophys. Res. Lett., 25(17), 3347-3350.
Fattahi, H., and F. Amelung (2015), InSAR bias and uncertainty due to the systematic and stochastic tropospheric delay, J. Geophys. Res.: Sol. Ea., 120(12), 8758-8773.
Fialko, Y. (2006), Interseismic strain accumulation and the earthquake potential on the southern San Andreas fault system, Nature, 441, 968-971, doi:910.1038/nature04797.
Foster, J., B. Brooks, T. Cherubini, C. Shacat, S. Businger, and C. L. Werner (2006), Mitigating atmospheric noise for InSAR using a high resolution weather model, Geophys. Res. Lett., 33(16).
Gao, B. C., and Y. J. Kaufman (2003), Water vapor retrievals using Moderate Resolution Imaging Spectroradiometer (MODIS) near-infrared channels, J. Geophys. Res., 108, 4389, doi:4310.1029/2002 JD003023.
Gelb, A. (1974), Applied optimal estimation. MIT press.
Genrich, J. F., and Y. Bock (2006), Instantaneous geodetic positioning with 10-50 Hz GPS measurements: Noise characteristics and implications for monitoring networks, J. Geophys. Res.: Sol. Ea., 111(B3).
Gili, J. A., J. Corominas, and J. Rius (2000), Using Global Positioning System techniques in landslide monitoring, Eng. Geol., 55(3), 167-192.
Gourmelen, N., and F. Amelung (2005), Postseismic mantle relaxation in the Central Nevada Seismic Belt, Science, 310, 1473 - 1476, DOI: 1410.1126/science.1119798.
Gundlich, B., and K. R. Koch (2002), Confidence regions for GPS baselines by Bayesian statistics, J. Geod., 76(1), 55-62.
Im, S., S. Hurlebaus, and Y. Kang (2013), Summary review of GPS technology for structural health monitoring, J. Struct. Eng., 139(10), 1653-1664.
Janssen, V., L. Ge, and C. Rizos (2004), Tropospheric corrections to SAR interferometry from GPS observations, GPS Solut., 8, 140-151.
Jarlemark, P. O., and T. R. Emardson (1998), Strategies for spatial and temporal extrapolation and interpolation of wet delay, J. Geod., 72(6), 350-355.
Jolivet, R., P. S. Agram, N. Y. Lin, M. Simons, M.-P. Doin, G. Peltzer, and Z. Li (2014), Improving InSAR geodesy using Global Atmospheric Models, J. Geophys. Res.: Sol. Ea., 119, 2324-2341, doi: 10.1002/2013J B010588.
Kouba, J. (2009), A guide to using International GNSS Service (IGS) products, http://iqscb.jpj nasa qov/saecb/resource/pube/UsinalGSProductsVer21.pdf.
Leandro, R., M. C. Santos, and R. B. Langley (2006), UNB neutral atmosphere models: Development and performance, Proceedings of ION NTM, 52(1), 564-73.
Lee, J.-M., J.-Y. Park, and J.-Y. Choi (2011), Evaluation of sub-aerial topographic surveying techniques using total station and RTK-GPS for applications in macrotidal sand beach environment, J. Coastal Res., 65, 535-540.
Li, Z., E. J. Fielding, P. Cross, and J.-P. Muller (2006a), Interferometric synthetic aperture radar atmospheric correction: GPS topography-dependent turbulence model, J. Geophys. Res., 111 , B02404, doi: 10.1029/2005J B003711. Li, Z., E. J. Fielding, and P. Cross (2009a), Integration of InSAR time series analysis and water vapour correction for mapping postseismic deformation after the 2003 Bam, Iran Earthquake, IEEE T. Geosci. Remote., 47, 3220-3230.
Li, Z., E. J. Fielding, P. Cross, and J.-P. Muller (2006b), Interferometric synthetic aperture radar atmospheric correction: Medium Resolution Imaging Spectrometer and Advanced Synthetic Aperture Radar integration, Geophys. Res. Lett., 33, L06816, doi:06810.01029/02005GL025299.
Li, Z., E. J. Fielding, P. Cross, and R. Preusker (2009b), Advanced InSAR atmospheric correction: MERIS/MODIS combination and stacked water vapour models, Int. J. Remote Sens., 30, 3343- 3363.
Li, Z., J. P. Muller, and P. Cross (2003), Comparison of precipitable water vapor derived from radiosonde, GPS, and Moderate-Resolution Imaging Spectroradiometer measurements, J. Geophys. Res.: Atmos., 108(D20).
Li, Z., J.-P. Muller, P. Cross, P. Albert, J. Fischer, and R. Bennartz (2006c), Assessment of the potential of MERIS near-infrared water vapour products to correct ASAR interferometric measurements, Int. J. Remote Sens., 27, 349-365, doi: 310.1080/01431160500307342.
Li, Z., J. P. Muller, P. Cross, and E. J. Fielding (2005), Interferometric synthetic aperture radar (InSAR) atmospheric correction: GPS, Moderate Resolution Imaging Spectroradiometer (MODIS), and InSAR integration, J. Geophys. Res.: Sol. Ea., 110(B3).
Li, X., G. Dick, M. Ge, S. Heise, J. Wickert, and M. Bender (2014), Real-time GPS sensing of atmospheric water vapor: Precise point positioning with orbit, clock, and phase delay corrections, Geophys. Res. Lett., 41, 3615-3621.
Liou, Y. A., Y. T. Teng, T. Van Hove, and J. C. Liljegren (2001), Comparison of precipitable water observations in the near tropics by GPS, microwave radiometer, and radiosondes, J. Appl. Meteorol. , 40(1), 5-15.
McCarthy, D. D. (1996). IERS Conventions (1996), IERS Technical Note 21. US Naval Observatory.
Mengistu Tsidu, G., T. Blumenstock, and F. Hase (2015), Observations of precipitable water vapour over complex topography of Ethiopia from ground-based GPS, FTIR, radiosonde and ERA-lnterim reanalysis, Atmos. Meas. Tech., 8, 3277-3295.
Ohtani, R., and I Naito (2000), Comparisons of GPS-derived precipitable water vapors with radiosonde observations in Japan, J. Geophys. Res.: Atmos., 105(D22), 26917-26929.
Onn, E, and H. A. Zebker (2006), Correction for interferometric synthetic aperture radar atmospheric phase artifacts using time series of zenith wet delay observations from a GPS network, J. Geophys. Res.: Sol. Ea., 111 (B9).
Perez-Ruiz, M., J. Carballido, J. Aguera, and J. A. Gil (2011), Assessing GNSS correction signals for assisted guidance systems in agricultural vehicles, Precis. Agric., 12(5), 639-652.
Poli, P., P. Moll, F. Rabier, G. Desroziers, B. Chapnik, L. Berre, S. B. Healy, E. Andersson, and F.- Z. El Guelai (2007), Forecast impact studies of zenith total delay data from European near real time GPS stations in Meteo France 4DVAR, J. Geophys. Res.: Atmos., 112, D06114.
Raja, M. K. R. V, S.l. Gutman, J. G. Yoe, L. M. McMillin, and J. Zhao (2008), The validation of AIRS retrievals of integrated precipitable water vapor using measurements from a network of ground- based GPS receivers over the contiguous United States, J. Atmos. Ocean. Tech., 25, 416-428.
Saastamoinen, J. (1972), Atmospheric correction for the troposphere and stratosphere in radio ranging of satellites. The Use of Artificial Satellites for Geodesy, Geophys. Monogr., 15, Amer. Geophys. Union, 247-251.
Shoji, Y, M. Kunii, and K. Saito (2011), Mesoscale data assimilation of Myanmar cyclone Nargis Part II: Assimilation of GPS-derived precipitable water vapor, J. Meteorol. Soc. Jpn, 89(1), 67-88.
Sobrino, J. A., J. El Kharraz, and Z. L. Li (2003), Surface temperature and water vapour retrieval from MODIS data,. Int. J. Remote Sens., 24(24), 5161-5182.
Treuhaft, R. N., and G. E. Lanyi (1987), The effect of the dynamic wet troposphere on radio interferometric measurements, Radio Sci., 22(2), 251-265.
Wadge, G., P. W. Webley, I. N. James, R. Bingley, A. Dodson, S. Waugh, T. Veneboer, G. Puglisi, M. Mattia, D. Baker, S. C. Edwards, S. J. Edwards, and P. J. Clarke (2002), Atmospheric models, GPS and InSAR measurements of the tropospheric water vapour field over Mount Etna, Geophys. Res. Lett., 29(19), 1905.
Walters, R. J., J. R. Elliott, Z. Li, and B. Parsons (2013), Rapid strain accumulation on theAshkabad fault (Turkmenistan) from atmosphere-corrected InSAR, J. Geophys. Res., 118, 1-17.
Wen, Y, Z. Li, C. Xu, I. Ryder, and R. Burgmann (2012), Postseismic motion after the 2001 Mw 7.8 Kokoxili earthquake in Tibet observed by InSAR time series, J. Geophys. Res., 117, B08405.
Wessel, P, W. H. F. Smith, R. Scharroo, J. Luis, and F. Wobbe (2013), Generic Mapping Tools: Improved version released, Eos Trans. AGU, 94(45), 409-210.
Williams, S., Y. Bock, and P. Fang (1998), Integrated satellite interferometry: Tropospheric noise, GPS estimates and implications for interferometric synthetic aperture radar products, J. Geophys. Res., 103(B11), 27051-27067.
Wright, T. J., B. Parsons, P. C. England, and E. J. Fielding (2004), InSAR observations of low slip rates on the major faults of western Tibet, Science, 305, 236-239.
Xu, W. B., Z. W. Li, X. L. Ding, and J. J. Zhu (2011), Interpolating atmospheric water vapor delay by incorporating terrain elevation information, J. Geod., 85(9), 555-564.
Yan, X., V. Ducrocq, G. Jaubert, P. Brousseau, P. Poli, C. Champollion, C. Flamant, and K. Boniface (2009), The benefit of GPS zenith delay assimilation to high-resolution quantitative precipitation forecasts: A case-study from COPS IOP 9, Q. J. Roy. Meteor. Soc., 135(644), 1788-1800.
Yao, Y, C. Yu, and Y. Hu (2014), A new method to accelerate PPP convergence time by using a global zenith troposphere delay estimate model, J. Navigation, 67(5), 899-910.
Zabuski, L., W. Swidzihski, M. Kulczykowski, T. Mrozek, and I. Laskowicz (2015), Monitoring of landslides in the Brda river valley in Koronowo (Polish Lowlands), Environmental Earth Sciences, 73(12), 8609-8619.
Zumberge, J. E, M. B. Heflin, D. C. Jefferson, M. M. Watkins, and F. H. Webb (1997), Precise point positioning for the efficient and robust analysis of GPS data from large networks, J. Geophys. Res.: Sol. Ea., 102(B3), 5005-5017. 44 able 1. Summary of ITD and SKIm+Onn Model Interpolation Performance From Cross Validation of all ZTDs and all Common Epoch RTPPP an ODIS PWVs for all of 2015.
>0 y1 No. Data Points Slope2 Intercept2 Correl. Coeff. RMS (mm)
ITD ZTD RTPPP ZTD 14883 1.012 ± 0.002 0.004 ± 0.003 m 0.99 4.6 SKIm+Onn ZTD RTPPP ZTD 14883 1.025 ± 0.002 0.026 ± 0.003 m 0.98 8.4 ITD PWV MODIS PWV3 8523 0.971 ± 0.002 0.329 ± 0.004 mm 0.98 1.48 SKIm+Onn PWV MODIS PWV3 8523 0.948 ± 0.002 0.634 ± 0.004 mm 0.97 1.73 ITD Map MODIS Map4 288 million 0.934 ± 0.003 1.223 ± 0.004 mm 0.97 1.71 SKIm+Onn Map MODIS Map4 288 million 0.912 ± 0.002 2.101 ± 0.004 mm 0.96 1.96 The linear model is Y = Slope * X + Intercept ncertainties are 95% confidence MODIS PWV pixels with co-located GPS stations AII available MODIS pixels

Claims

Claims
1. A method for configuring a Global Satellite Navigation System (GNSS) infrastructure comprising a plurality of GNSS receivers, the method comprising: obtaining a water vapour map comprising water vapour values at different points; and performing a configuration process comprising: defining one or more sub-regions based on an initial set of receivers, such that one or more receivers are located on (e.g. on the boundary of, on one or more vertices of, and/or within) each sub-region; and for each sub-region: obtaining, from the water vapour map, a first water vapour value for each of a plurality of points on the sub-region; determining a second water vapour value for each of the plurality of points by using interpolation of water vapour values, obtained from the water vapour map, corresponding to positions of two or more receivers; determining a statistical value based on the first and second water vapour values for each point; and if the statistical value satisfies a predetermined criteria, increasing the number of receivers on the sub-region to obtain a modified set of receivers.
2. A method according to claim 1, wherein the two or more receivers used for the interpolation comprise one or more receivers on the sub-region and/or one or more receivers on one or more neighbouring sub-regions.
3. A method according to claim 1 or 2, further comprising repeating the configuration process such that the modified set of receivers obtained at the end of one iteration is used as the initial set of receivers in the next iteration.
4. A method according to claim 3, wherein the configuration process is repeated until a termination criterion is satisfied.
5. A method according to claim 4, wherein the termination criterion is based on one or more of: the statistical value does not satisfy the predetermined criterion for every sub-region; the receiver spacing is less than a certain threshold; and the number of iterations has reached a certain threshold.
6. A method according to any preceding claims, wherein the water vapour map is obtained from Moderate Resolution Imaging Spectroradiometer (MODIS) data.
7. A method according to any preceding claim, wherein the interpolation is performed based on an Iterative Tropospheric Decomposition (ITD) model.
8. A method according to any preceding claim, wherein the initial set of receivers comprises a new or pre-existing regular or irregular array of receivers.
9. A method according to any preceding claim, wherein a sub-region comprises a polygon, and wherein a receiver is located on one or more vertices of the polygon, on one or more edges of the polygon, and/or at the centre of the polygon.
10. A method according to claim 9, wherein the polygon comprises one or more of: square; rectangle; triangle; and hexagon.
11. A method according to claim 9, wherein the sub-regions comprise a regular array of rectangular or square sub-regions, or a Triangulated Irregular Network (TIN) obtained using Delaunay triangulation.
12. A method according to any preceding claim, wherein the plurality of points on the sub-region have the same spatial resolution as the water vapour map.
13. A method according to any preceding claim, wherein the statistical value is the root mean squared (RMS) of the differences between the first and second water vapour values.
14. A method according to any preceding claim, wherein the predetermined criterion is that the statistical value is greater than a predetermined threshold.
15. A method according to any preceding claim, wherein increasing the number of receivers on the sub-region comprises adding a receiver to the geometric centre of the sub-region.
16. A method according to any preceding claim, further comprising generating a report including position, spacing, and/or density information of the modified set of receivers.
17. A method according to any preceding claim, further comprising configuring the GNSS infrastructure according to the modified set of receivers.
18. An apparatus for configuring a Global Satellite Navigation System (GNSS) infrastructure comprising a plurality of GNSS receivers, the apparatus comprising a processor configured to: obtain a water vapour map comprising water vapour values at different points; and perform a configuration process comprising: defining one or more sub-regions based on an initial set of receivers, such that one or more receivers are located on (e.g. on the boundary of, on one or more vertices of, and/or within) each sub-region; and for each sub-region: obtaining, from the water vapour map, a first water vapour value for each of a plurality of points on the sub-region; determining a second water vapour value for each of the plurality of points by using interpolation of water vapour values, obtained from the water vapour map, corresponding to positions of two or more receivers; determining a statistical value based on the first and second water vapour values for each point; and if the statistical value satisfies a predetermined criteria, increasing the number of receivers on the sub-region to obtain a modified set of receivers.
19. A Global Satellite Navigation System (GNSS) infrastructure configured according to the method of any of claims 1 to 17.
PCT/GB2020/052874 2019-11-11 2020-11-11 Configuration of global satellite navigation system infrastructure WO2021094753A1 (en)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
GBGB1916381.5A GB201916381D0 (en) 2019-11-11 2019-11-11 Configuration of global satellite navigation system infrastructure
GB1916381.5 2019-11-11

Publications (1)

Publication Number Publication Date
WO2021094753A1 true WO2021094753A1 (en) 2021-05-20

Family

ID=69062179

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/GB2020/052874 WO2021094753A1 (en) 2019-11-11 2020-11-11 Configuration of global satellite navigation system infrastructure

Country Status (2)

Country Link
GB (1) GB201916381D0 (en)
WO (1) WO2021094753A1 (en)

Cited By (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113269805A (en) * 2021-06-11 2021-08-17 北京大学 Remote sensing rainfall inversion training sample self-adaptive selection method guided by rainfall event
CN114019584A (en) * 2021-10-11 2022-02-08 武汉大学 VRS resolving method for high-precision CORS network in large-altitude-difference area
CN115015970A (en) * 2022-05-17 2022-09-06 中铁第一勘察设计院集团有限公司 Method for constructing strain field for Beidou continuous operation reference station coverage area
CN115099159A (en) * 2022-07-20 2022-09-23 武汉大学 MODIS water vapor inversion method based on neural network and considering earth surface difference
CN115343734A (en) * 2022-10-13 2022-11-15 武汉地铁集团有限公司 GNSS deformation monitoring method based on bilinear interpolation hemisphere model
CN117148398A (en) * 2023-10-31 2023-12-01 中国测绘科学研究院 Two-network-integrated station distribution geometric configuration assessment method, system and equipment
US11880014B2 (en) * 2022-06-09 2024-01-23 Shandong University Prediction method and system for changes of physical parameters after earthquakes

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20190166542A1 (en) * 2016-08-08 2019-05-30 Huawei Technologies Co., Ltd. Method and apparatus for updating network rtk reference station network

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20190166542A1 (en) * 2016-08-08 2019-05-30 Huawei Technologies Co., Ltd. Method and apparatus for updating network rtk reference station network

Non-Patent Citations (105)

* Cited by examiner, † Cited by third party
Title
AHMED, F.P. VACLAVOVICF. N. TEFERLEJ. DOUSAR. BINGLEYD. LAURICHESSE: "Comparative analysis of real-time precise point positioning zenith total delay estimates", GPS SOLUT., vol. 20, no. 2, 2016, pages 187 - 199, XP035895534, DOI: 10.1007/s10291-014-0427-z
AHN, Y. W.LACHAPELLE, G.SKONE, S.GUTMAN, SSAHM, S.: "Analysis of GPS RTK performance using external NOAA tropospheric corrections integrated with a multiple reference station approach", GPS SOLUTIONS, vol. 10, 2006, pages 171 - 186, XP019388114, DOI: 10.1007/s10291-005-0017-1
ALBERT, P.R. BENNARTZJ. FISCHER: "Remote sensing of atmospheric water vapor from backscattered sunlight in cloudy atmospheres", J. ATMOS. OCEAN. TECH., vol. 18, 2001, pages 865 - 874
AUMANN, H. H.M. T. CHAHINEC. GAUTIERM. D. GOLDBERGE. KALNAYL. M. MCMILLINH. REVERCOMBP. W. ROSENKRANZW. L. SMITHD. H. STAELIN: "AIRS/AMSU/HSB on the Aqua mission: design, science objectives, data products, and processing systems", IEEE T. GEOSCI. REMOTE., vol. 41, 2003, pages 253 - 264, XP011073315
BEKAERT, D. P. S.A. HOOPERT.J. WRIGHT: "A spatially variable power law tropospheric correction technique for InSAR data", J. GEOPHYS. RES., vol. 120, no. 2, 2015, pages 1345 - 1356
BENEVIDES, P.G. NICOJ. CATALAOP. M. A. MIRANDA: "Bridging InSAR and GPS tomography: A new differential geometrical constraint", IEEE T. GEOSCI. REMOTE., vol. 54, no. 2, 2016, pages 697 - 702, XP011595768, DOI: 10.1109/TGRS.2015.2463263
BENEVIDES, P.J. CATALAOP. M. A. MIRANDA: "On the inclusion of GPS precipitable water vapour in the nowcasting of rainfall", NAT. HAZARDS EARTH SYST. SCI., vol. 15, 2015, pages 2605 - 2616
BENNARTZ, R.J. FISCHER: "Retrieval of columnar water vapour over land from backscattered solar radiation using the Medium Resolution Imaging Spectrometer (MERIS", REMOTE SENS. ENVIRON., vol. 78, 2001, pages 271 - 280
BENNITT, G. V.A. JUPP: "Operational assimilation of GPS zenith todal delay observations into the Met Office numerical weather prediction models", MON. WEA. REV., vol. 140, 2012, pages 2706 - 2719
BERRADA BABY, H.P. GOLEJ. LAVERGNAT: "A model for the tropospheric excess path length of radio waves from surface meteorological measurements", RADIO SCI., vol. 23, no. 6, 1988, pages 1023 - 1038
BETTS, R. A.: "Biogeophysical impacts of land use on present-day climate: Near-surface temperature change and radiative forcing", ATMOSPHERIC SCIENCE LETTERS, vol. 2, no. 1-4, 2001, pages 39 - 51
BEVIS M.BUSINGER, S.CHISWELL, S.HERRING, T.ANTHES, R.ROCKEN, C.WARE, R. H.: "GPS meteorology: Mapping zenith wet delays onto precipitable water", JOURNAL OF APPLIED METEOROLOGY, vol. 33, 1994, pages 379 - 386
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 GEOPHYSICAL RESEARCH, vol. 97, 1992, pages 15787 - 15801
BEVIS, M.S. BUSINGERT. HERRINGC. ROCKENR. ANTHESR. WARE: "GPS meteorology: Remote sensing of atmospheric water vapor using the Global Positioning System", J. GEOPHYS. RES., vol. 97, no. D14, 1992, pages 15787 - 15801
BEYERLE, G: "Carrier phase wind-up in GPS reflectometry", GPS SOLUT., vol. 13, no. 3, 2009, pages 191 - 198, XP019665541
BOEHM, J.A. NIELLP. TREGONINGH. SCHUH: "Global Mapping Function (GMF): A new empirical mapping function based on numerical weather model data", GEOPHYS. RES. LETT., vol. 33, 2006, pages L07304
BOEHM, J.B. WERLH. SCHUH: "Troposphere mapping functions for GPS and very long baseline interferometry from European Centre for Medium-Range Weather Forecasts operational analysis data", J. GEOPHYS. RES., vol. 111, 2006, pages B02406
BONAFONI, S.BIONDI, R.: "The usefulness of the Global Navigation Satellite Systems (GNSS) in the analysis of precipitation events", ATMOSPHERIC RESEARCH, vol. 167, 2016, pages 15 - 23
BONIFACE, K.V. DUCROCQG. JAUBERTX. YANP. BROUSSEAUF. MASSONC. CHAMPOLLIONJ. CH'ERYE. DOERFLINGER: "Impact of high-resolution data assimilation of GPS zenith delay on Mediterranean heavy rainfall forecasting", ANN. GEOPHYS., vol. 27, 2009, pages 2739 - 2753
BRASINGTON, J.J. LANGHAMB. RUMSBY: "Methodological sensitivity of morphometric estimates of coarse fluvial sediment transport", GEOMORPHOLOGY, vol. 53, 2003, pages 299 - 316
BRASINGTON, J.RUMSBY, B. T.MCVEY, R. A.: "Monitoring and modelling morphological change in a braided gravel-bed river using high resolution GPS-based survey", EARTH SURFACE PROCESSES AND LANDFORMS, vol. 25, 2000, pages 973 - 990
BRUNNER, F. K.WELSCH, W. M.: "Effect of the troposphere on GPS measurements", GPS WORLD, vol. 4, 1993, pages 42 - 51
CHO, J. Y. N.R. E. NEWELLB. E. ANDERSONJ. D. W. BARRICKK. L. THORNHILL: "Characterizations of tropospheric turbulence and stability layers from aircraft observations", J. GEOPHYS. RES.: ATMOS., vol. 108, no. D20, 2003, pages 8784
DE PAOLI, L.G. E. FLOWERS: "Dynamics of a small surge-type glacier using one-dimensional geophysical inversion", J. GLACIOL., vol. 55, no. 194, 2009, pages 1101 - 1112
DIVAKARLA, M. G.C. D. BARNETM. D. GOLDBERGL. M. MCMILLINE. MADDYW. WOLFL. ZHOUX. LIU: "Validation of Atmospheric Infrared Sounder temperature and water vapor retrievals with matched radiosonde measurements and forecasts", J. GEOPHYS. RES.: ATMOS., vol. 111, 2006, pages D09515
ELLIOTT, J. R.J. BIGGSB. PARSONST. J. WRIGHT: "InSAR slip rate determination on the Altyn Tagh Fault, northern Tibet, in the presence of topographically correlated atmospheric delays", GEOPHYS. RES. LETT., vol. 35, no. 12, 2008
EMARDSON, T. R.J. M. JOHANSSON: "Spatial interpolation of the atmospheric water vapor content between sites in a ground-based GPS Network", GEOPHYS. RES. LETT., vol. 25, no. 17, 1998, pages 3347 - 3350
FARR, T. G. ET AL.: "The shuttle radar topography mission", REVIEWS OF GEOPHYSICS, vol. 45, 2007, pages RG2004, XP055729368, DOI: 10.1029/2005RG000183
FATTAHI, HF. AMELUNG: "InSAR bias and uncertainty due to the systematic and stochastic tropospheric delay", J. GEOPHYS. RES.: SOL. EA., vol. 120, no. 12, 2015, pages 8758 - 8773
FIALKO, Y.: "Interseismic strain accumulation and the earthquake potential on the southern San Andreas fault system", NATURE, vol. 441, 2006, pages 968 - 971
FOSTER, J.B. BROOKST. CHERUBINIC. SHACATS. BUSINGERC. L. WERNER: "Mitigating atmospheric noise for InSAR using a high resolution weather model", GEOPHYS. RES. LETT., vol. 33, no. 16, 2006
GAO, B. C.Y. J. KAUFMAN: "Water vapor retrievals using Moderate Resolution Imaging Spectroradiometer (MODIS) near-infrared channels", J. GEOPHYS. RES., vol. 108, 2003, pages 4389
GAO, B.-C.KAUFMAN, Y. J.: "Water vapor retrievals using Moderate Resolution Imaging Spectroradiometer (MODIS) near-infrared channels", JOURNAL OF GEOPHYSICAL RESEARCH, vol. 108, no. D13, 2003, pages 4389
GELB, A.: "Applied optimal estimation", 1974, MIT PRESS
GENRICH, J. F.Y. BOCK: "Instantaneous geodetic positioning with 10-50 Hz GPS measurements: Noise characteristics and implications for monitoring networks", J. GEOPHYS. RES.: SOL. EA., vol. 111, no. B3, 2006, XP009143930, DOI: 10.1029/2005JB003617
GILI, J. A.J. COROMINASJ. RIUS: "Using Global Positioning System techniques in landslide monitoring", ENG. GEOL., vol. 55, no. 3, 2000, pages 167 - 192, XP004185115, DOI: 10.1016/S0013-7952(99)00127-1
GOURMELEN, N.F. AMELUNG: "Postseismic mantle relaxation in the Central Nevada Seismic Belt", SCIENCE, vol. 310, 2005, pages 1473 - 1476
GOYANES, G.VIEIRA, G.CASELLI, A.CARDOSO, M.MARMY, A.SANTOS, F.BERNANDO, I.HAUCK, C: "Local influences of geothermal anomalies on permafrost distribution in an active volcanic island (Deception Island, Antarctica", GEOMORPHOLOGY, vol. 225, 2014, pages 57 - 68
GREJNER-BRZEZINSKA, D. A.KASHANI, I.WIELGOSZ, P.: "On accuracy and reliability of instantaneous network RTK as a function of network geometry, station separation, and data processing strategy", GPS SOLUTIONS, vol. 9, 2005, pages 212 - 225, XP019357642
GUNDLICH, B.K. R. KOCH: "Confidence regions for GPS baselines by Bayesian statistics", J. GEOD., vol. 76, no. 1, 2002, pages 55 - 62
HARLEY, M. D.TURNER, I. L.SHORT, A. D.RANASINGHE, R.: "Assessment and integration of conventional, RTK-GPS and image-derived beach survey methods for daily to decadal coastal monitoring", COASTAL ENGINEERING, vol. 58, 2011, pages 194 - 205, XP027568119
HE, H.LI, J.YANG, Y.XU, J.GUO, H.WANG, A.: "Performance assessment of single- and dual-frequency BeiDou/GPS single-epoch kinematic positioning", GPS SOLUTIONS, vol. 18, 2014, pages 393 - 403
HU, G., R.KHOO, H. S.GOH, P. C.LAW, C. L.: "Development and assessment of GPS virtual reference stations for RTK positioning", JOURNAL OF GEODESY, vol. 77, 2003, pages 292 - 302
IM, S.S. HURLEBAUSY. KANG: "Summary review of GPS technology for structural health monitoring", J. STRUCT. ENG., vol. 139, no. 10, 2013, pages 1653 - 1664
JADE SRIDEVI ET AL: "Water vapor study using MODIS and GPS data at 64 continuous GPS stations (2002-2017) in indian subcontinent", JOURNAL OF ATMOSPHERIC AND SOLAR-TERRESTRIAL PHYSICS, PERGAMON, AMSTERDAM, NL, vol. 196, 1 October 2019 (2019-10-01), XP085919710, ISSN: 1364-6826, [retrieved on 20191001], DOI: 10.1016/J.JASTP.2019.105138 *
JANSSEN, V.L. GEC. RIZOS: "Tropospheric corrections to SAR interferometry from GPS observations", GPS SOLUT., vol. 8, 2004, pages 140 - 151
JARLEMARK, P. O.T. R. EMARDSON: "Strategies for spatial and temporal extrapolation and interpolation of wet delay", J. GEOD., vol. 72, no. 6, 1998, pages 350 - 355
JAVERNICK, L.BRASINGTON, J.CARUSO, B.: "Modeling the topography of shallow braided rivers using Structure-from-Motion photogrammetry", GEOMORPHOLOGY, vol. 213, 2014, pages 166 - 182, XP028842202, DOI: 10.1016/j.geomorph.2014.01.006
JIAO, Y.WANG, Z.WANG, X.ADOKO, A. C.YANG, Z.: "Stability assessment of an ancient landslide crossed by two coal mine tunnels", ENGINEERING GEOLOGY, vol. 159, 2013, pages 36 - 44, XP028531696, DOI: 10.1016/j.enggeo.2013.03.021
JOLIVET, R.P. S. AGRAMN. Y. LINM. SIMONSM.-P. DOING. PELTZERZ. LI: "Improving InSAR geodesy using Global Atmospheric Models", J. GEOPHYS. RES.: SOL. EA., vol. 119, 2014, pages 2324 - 2341
KOUBA, J., A GUIDE TO USING INTERNATIONAL GNSS SERVICE (IGS) PRODUCTS, 2009, Retrieved from the Internet <URL:http://igscb.jpl.nasa.gov/igscb/resource/pubs/UsingIGSProductsVer21.pdf>
LEANDRO, R.M. C. SANTOSR. B. LANGLEY: "UNB neutral atmosphere models: Development and performance", PROCEEDINGS OF ION NTM, vol. 52, no. 1, 2006, pages 564 - 73, XP056004938
LEE, J.-M.J.-Y. PARKJ.-Y. CHOI: "Evaluation of sub-aerial topographic surveying techniques using total station and RTK-GPS for applications in macrotidal sand beach environment", J. COASTAL RES., vol. 65, 2011, pages 535 - 540
LI, H.CHEN, B: "The evolution of precipitable water associated with the Asian and Australian monsoons as revealed from MODIS/SSMI, ECMWF and NCEP reanalysis data sets", GEOPHYSICAL RESEARCH LETTERS, vol. 32, 2005, pages L10811
LI, S.CAI, H.ABRAHAM, D. M.MAO, P: "Estimating features of underground utilities: Hybrid GPR/GPS approach", JOURNAL OF COMPUTING IN CIVIL ENGINEERING, vol. 30, 2016, pages 04014108
LI, XG. DICKM. GES. HEISEJ. WICKERTM. BENDER: "Real-time GPS sensing of atmospheric water vapor: Precise point positioning with orbit, clock, and phase delay corrections", GEOPHYS. RES. LETT., vol. 41, 2014, pages 3615 - 3621
LI, Z.E. J. FIELDINGP. CROSS: "Integration of InSAR time series analysis and water vapour correction for mapping postseismic deformation after the 2003 Bam, Iran Earthquake", IEEE T. GEOSCI. REMOTE., vol. 47, 2009, pages 3220 - 3230
LI, Z.E. J. FIELDINGP. CROSSJ.-P. MULLER: "Interferometric synthetic aperture radar atmospheric correction: GPS topography-dependent turbulence model", J. GEOPHYS. RES.,, vol. 111, 2006, pages B02404
LI, Z.E. J. FIELDINGP. CROSSJ.-P. MULLER: "Interferometric synthetic aperture radar atmospheric correction: Medium Resolution Imaging Spectrometer and Advanced Synthetic Aperture Radar integration", GEOPHYS. RES. LETT., vol. 33, 2006, pages L06816
LI, Z.E. J. FIELDINGP. CROSSR. PREUSKER: "Advanced InSAR atmospheric correction: MERIS/MODIS combination and stacked water vapour models", INT. J. REMOTE SENS., vol. 30, 2009, pages 3343 - 3363
LI, Z.J. P. MULLERP. CROSS: "Comparison of precipitable water vapor derived from radiosonde, GPS, and Moderate-Resolution Imaging Spectroradiometer measurements", J. GEOPHYS. RES.: ATMOS., vol. 108, no. D20, 2003
LI, Z.J. P. MULLERP. CROSSE. J. FIELDING: "Interferometric synthetic aperture radar (InSAR) atmospheric correction: GPS, Moderate Resolution Imaging Spectroradiometer (MODIS), and InSAR integration", J. GEOPHYS. RES.: SOL. EA., vol. 110, no. B3, 2005
LI, Z.J.-P. MULLERP. CROSSP. ALBERTJ. FISCHERR. BENNARTZ: "Assessment of the potential of MERIS near-infrared water vapour products to correct ASAR interferometric measurements", INT. J. REMOTE SENS., vol. 27, 2006, pages 349 - 365
LI, Z.MULLER, J. PCROSS, P: "Comparison of precipitable water vapor derived from radiosonde, GPS, and Moderate-Resolution Imaging Spectroradiometer measurements", JOURNAL OF GEOPHYSICAL RESEARCH, vol. 108, no. D20, 2003, pages 4651
LIOU, Y. A.Y. T. TENGT. VAN HOVEJ. C. LILJEGREN: "Comparison of precipitable water observations in the near tropics by GPS, microwave radiometer, and radiosondes", J. APPL. METEOROL., vol. 40, no. 1, 2001, pages 5 - 15
LIU, J.LIANG, H.SUN, Z.ZHOU, X: "Validation of the Moderate-Resolution Imaging Spectroradiometer precipitable water vapor product using measurements from GPS on the Tibetan Plateau", JOURNAL OF GEOPHYSICAL RESEARCH: ATMOSPHERES, vol. 111, 2006, pages D14103
MAOLIN TANG ET AL: "Area-Oriented Reference Station Placement for Network RTK", COMPUTER SCIENCE AND SOFTWARE ENGINEERING, 2008 INTERNATIONAL CONFERENCE ON, IEEE, PISCATAWAY, NJ, USA, 12 December 2008 (2008-12-12), pages 919 - 922, XP031377881, ISBN: 978-0-7695-3336-0 *
MCCARTHY, D. D.: "IERS Conventions", IERS TECHNICAL NOTE 21. US NAVAL OBSERVATORY, 1996
MENG, X.ROBERTS, S.CUI, Y.GAO, Y.CHEN, Q.XU, C.HE, Q.SHARPLES, S.BHATIA, P.: "Required navigation performance for connected and autonomous vehicles: Where are we now and where are we going?", TRANSPORTATION PLANNING AND TECHNOLOGY, vol. 41, 2018, pages 104 - 118
MENGISTU TSIDU, G.T. BLUMENSTOCKF. HASE: "Observations of precipitable water vapour over complex topography of Ethiopia from ground-based GPS, FTIR, radiosonde and ERA-Interim reanalysis", ATMOS. MEAS. TECH., vol. 8, 2015, pages 3277 - 3295
NARAMA, C.DUISHONAKUNOV, M.KAAB, A.DAIYROV, M.ABDRAKHMATOV, K.: "Natural Hazards and Earth System Sciences", vol. 10, 24 July 2008, TIEN SHAN, article "outburst flood at the western Zyndan glacier lake and recent regional changes in glacier lakes of the Teskey Ala-Too range", pages: 647 - 659
NICKITOPOULOU, A.PROTOPSALTI, K.STIROS, S: "Monitoring dynamic and quasi-static deformations of large flexible engineering structures with GPS: Accuracy, limitations and promises", ENGINEERING STRUCTURES, vol. 28, 2006, pages 1471 - 1482, XP025095318, DOI: 10.1016/j.engstruct.2006.02.001
OHTANI, R.I NAITO: "Comparisons of GPS-derived precipitable water vapors with radiosonde observations in Japan", J. GEOPHYS. RES.: ATMOS.,, vol. 105, no. D22, 2000, pages 26917 - 26929
ONN, F.H. A. ZEBKER: "Correction for interferometric synthetic aperture radar atmospheric phase artifacts using time series of zenith wet delay observations from a GPS network", J. GEOPHYS. RES.: SOL. EA.,, vol. 111, no. B9, 2006
PAN, B.CAO, B.WANG, J.ZHANG, G.ZHANG, C.HU, Z.HUANG, B.: "Glacier variations in response to climate change from 1972 to 2007 in the western Lenglongling mountains, northeastern Tibetan Plateau", JOURNAL OF GLACIOLOGY, vol. 58, 2012, pages 879 - 888
PEREZ-RUIZ, M.J. CARBALLIDOJ. AGUERAJ. A. GIL: "Assessing GNSS correction signals for assisted guidance systems in agricultural vehicles", PRECIS. AGRIC., vol. 12, no. 5, 2011, pages 639 - 652, XP019947379, DOI: 10.1007/s11119-010-9211-4
PEREZ-RUIZ, M.SLAUGHTER, D. C.GLIEVER, C.UPADHYAYA, S. K.: "Tractor-based Real-time Kinematic-Global Positioning System (RTK-GPS) guidance system for geospatial mapping of row crop transplant", BIOSYSTEMS ENGINEERING, vol. 111, 2012, pages 64 - 71, XP028343097, DOI: 10.1016/j.biosystemseng.2011.10.009
PIRTI, A: "Performance analysis of the real time kinematic GPS (RTK GPS) technique in a highway project (stake-out", SURVEY REVIEW, vol. 39, 2007, pages 43 - 53
POLI, P.P. MOLLF. RABIERG. DESROZIERSB. CHAPNIKL. BERRES. B. HEALYE. ANDERSSONF.-Z. EL GUELAI: "Forecast impact studies of zenith total delay data from European near real-time GPS stations in Meteo France 4DVAR", J. GEOPHYS. RES.: ATMOS., vol. 112, 2007, pages D06114
PRASAD, A. K.SINGH, R. P: "Validation of MODIS Terra, AIRS, NCEP/DOE AMIP-II Reanalysis-2, and AERONET Sun photometer derived integrated precipitable water vapor using ground-based GPS receivers over India", JOURNAL OF GEOPHYSICAL RESEARCH, vol. 114, 2009, pages D05107
RABAH, M.BASIOUNY, M.GHANEM, E.ELHADARY, A.: "Using RTK and VRS in direct geo-referencing of the UAV imagery", NRIAG JOURNAL OF ASTRONOMY AND GEOPHYSICS, vol. 7, 2018, pages 220 - 220
RAJA, M. K. R. V.S.I. GUTMANJ. G. YOEL. M. MCMILLINJ. ZHAO: "The validation of AIRS retrievals of integrated precipitable water vapor using measurements from a network of ground-based GPS receivers over the contiguous United States", J. ATMOS. OCEAN. TECH., vol. 25, 2008, pages 416 - 428
ROHM, W.YUAN, Y.BIADEGLGNE, B.ZHANG, K.LE MARSHALL, J.: "Ground-based GNSS ZTD/IWV estimation system for numerical weather prediction in challenging weather conditions", ATMOSPHERIC RESEARCH, vol. 138, 2014, pages 414 - 426
SAASTAMOINEN, J.: "The Use of Artificial Satellites for Geodesy, Geophys. Monogr.", vol. 15, 1972, AMER. GEOPHYS. UNION, article "Atmospheric correction for the troposphere and stratosphere in radio ranging of satellites", pages: 247 - 251
SHOJI, Y.M. KUNIIK. SAITO: "Mesoscale data assimilation of Myanmar cyclone Nargis Part II: Assimilation of GPS-derived precipitable water vapor", J. METEOROL. SOC. JPN, vol. 89, no. 1, 2011, pages 67 - 88
SOBRINO, J. A.J. EL KHARRAZZ. L. LI: "Surface temperature and water vapour retrieval from MODIS data", INT. J. REMOTE SENS., vol. 24, no. 24, 2003, pages 5161 - 5182
STOLL, A.KUTZBACH, H. D.: "Guidance of a forage harvester with GPS", PRECISION AGRICULTURE, vol. 2, 2000, pages 281 - 291, XP019214609, DOI: 10.1023/A:1011842907397
TREUHAFT, R. N.G. E. LANYI: "The effect of the dynamic wet troposphere on radio interferometric measurements", RADIO SCI., vol. 22, no. 2, 1987, pages 251 - 265
WADGE, G.P. W. WEBLEYI. N. JAMESR. BINGLEYA. DODSONS. WAUGHT. VENEBOERG. PUGLISIM. MATTIAD. BAKER: "Atmospheric models, GPS and InSAR measurements of the tropospheric water vapour field over Mount Etna", GEOPHYS. RES. LETT., vol. 29, no. 19, 2002, pages 1905
WALTERS, R. J.J. R. ELLIOTTZ. LIB. PARSONS: "Rapid strain accumulation on theAshkabad fault (Turkmenistan) from atmosphere-corrected InSAR", J. GEOPHYS. RES., vol. 118, 2013, pages 1 - 17
WANNINGER, L.: "Virtual reference stations (VRS", GPS SOLUTIONS, vol. 7, 2003, pages 143 - 144
WEN, Y.Z. LIC. XUI. RYDERR. BURGMANN: "Postseismic motion after the 2001 Mw 7.8 Kokoxili earthquake in Tibet observed by InSAR time series", J. GEOPHYS. RES., vol. 117, 2012, pages B08405
WESSEL, P.W. H. F. SMITHR. SCHARROOJ. LUISF. WOBBE: "Generic Mapping Tools: Improved version released", EOS TRANS. AGU, vol. 94, no. 45, 2013, pages 409 - 210
WILLIAMS, S.Y. BOCKP. FANG: "Integrated satellite interferometry: Tropospheric noise, GPS estimates and implications for interferometric synthetic aperture radar products", J. GEOPHYS. RES., vol. 103, no. B11, 1998, pages 27051 - 27067
WRIGHT, T. J.B. PARSONSP. C. ENGLANDE. J. FIELDING: "InSAR observations of low slip rates on the major faults of western Tibet", SCIENCE, vol. 305, 2004, pages 236 - 239
XU, W. B.Z. W. LIX. L. DINGJ. J. ZHU: "Interpolating atmospheric water vapor delay by incorporating terrain elevation information", J. GEOD., vol. 85, no. 9, 2011, pages 555 - 564, XP019942007, DOI: 10.1007/s00190-011-0456-0
YAN, X.DUCROCQ, V.POLI, P.HAKAM, M.JAUBERT, G.WALPERSDORF, A.: "Impact of GPS zenith delay assimilation on convective-scale prediction of Mediterranean heavy rainfall", JOURNAL OF GEOPHYSICAL RESEARCH: ATMOSPHERES, vol. 114, 2009, pages D03104
YAN, X.V. DUCROCQG. JAUBERTP. BROUSSEAUP. POLIC. CHAMPOLLIONC. FLAMANTK. BONIFACE: "The benefit of GPS zenith delay assimilation to high-resolution quantitative precipitation forecasts: A case-study from COPS IOP 9", Q. J. ROY. METEOR. SOC., vol. 135, no. 644, 2009, pages 1788 - 1800
YAO, Y.C. YUY. HU: "A new method to accelerate PPP convergence time by using a global zenith troposphere delay estimate model", J. NAVIGATION, vol. 67, no. 5, 2014, pages 899 - 910
YONG WON AHN ET AL: "Analysis of GPS RTK performance using external NOAA tropospheric corrections integrated with a multiple reference station approach", GPS SOLUTIONS, SPRINGER, BERLIN, DE, vol. 10, no. 3, 5 January 2006 (2006-01-05), pages 171 - 186, XP019388114, ISSN: 1521-1886, DOI: 10.1007/S10291-005-0017-1 *
YU, C.PENNA, N. T.LI, Z.: "Generation of real-time mode high-resolution water vapor fields from GPS observations", JOURNAL OF GEOPHYSICAL RESEARCH: ATMOSPHERES, vol. 122, 2017, pages 2008 - 2025
YU.C.LI, Z.PENNA, N. T.: "Interferometric synthetic aperture radar atmospheric correction using a GPS-based iterative tropospheric decomposition model", REMOTE SENSING OF ENVIRONMENT, vol. 204, 2018, pages 109 - 121, XP085288248, DOI: 10.1016/j.rse.2017.10.038
YU.C.LI, Z.PENNA, N. T.CRIPPA, P: "Generic atmospheric correction model for interferometric synthetic aperture radar observations", JOURNAL OF GEOPHYSICAL RESEARCH: SOLID EARTH, vol. 123, 2018, pages 9202 - 9222
ZABUSKI, L.W. SWIDZIRISKIM. KULCZYKOWSKIT. MROZEKI. LASKOWICZ: "Monitoring of landslides in the Brda river valley in Koronowo (Polish Lowlands", ENVIRONMENTAL EARTH SCIENCES, vol. 73, no. 12, 2015, pages 8609 - 8619
ZUMBERGE, J. F.M. B. HEFLIND. C. JEFFERSONM. M. WATKINSF. H. WEBB: "Precise point positioning for the efficient and robust analysis of GPS data from large networks", J. GEOPHYS. RES.: SOL. EA., vol. 102, no. B3, 1997, pages 5005 - 5017, XP055482091, DOI: 10.1029/96JB03860

Cited By (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113269805A (en) * 2021-06-11 2021-08-17 北京大学 Remote sensing rainfall inversion training sample self-adaptive selection method guided by rainfall event
CN113269805B (en) * 2021-06-11 2022-06-28 北京大学 Rainfall event guided remote sensing rainfall inversion training sample self-adaptive selection method
CN114019584A (en) * 2021-10-11 2022-02-08 武汉大学 VRS resolving method for high-precision CORS network in large-altitude-difference area
CN115015970A (en) * 2022-05-17 2022-09-06 中铁第一勘察设计院集团有限公司 Method for constructing strain field for Beidou continuous operation reference station coverage area
US11880014B2 (en) * 2022-06-09 2024-01-23 Shandong University Prediction method and system for changes of physical parameters after earthquakes
CN115099159A (en) * 2022-07-20 2022-09-23 武汉大学 MODIS water vapor inversion method based on neural network and considering earth surface difference
CN115099159B (en) * 2022-07-20 2023-03-07 武汉大学 MODIS water vapor inversion method based on neural network and considering earth surface difference
CN115343734A (en) * 2022-10-13 2022-11-15 武汉地铁集团有限公司 GNSS deformation monitoring method based on bilinear interpolation hemisphere model
CN115343734B (en) * 2022-10-13 2023-01-17 武汉地铁集团有限公司 GNSS deformation monitoring method based on bilinear interpolation hemisphere model
CN117148398A (en) * 2023-10-31 2023-12-01 中国测绘科学研究院 Two-network-integrated station distribution geometric configuration assessment method, system and equipment
CN117148398B (en) * 2023-10-31 2023-12-29 中国测绘科学研究院 Two-network-integrated station distribution geometric configuration assessment method, system and equipment

Also Published As

Publication number Publication date
GB201916381D0 (en) 2019-12-25

Similar Documents

Publication Publication Date Title
Yu et al. Generation of real‐time mode high‐resolution water vapor fields from GPS observations
WO2021094753A1 (en) Configuration of global satellite navigation system infrastructure
Yu et al. Interferometric synthetic aperture radar atmospheric correction using a GPS-based iterative tropospheric decomposition model
Benevides et al. On the inclusion of GPS precipitable water vapour in the nowcasting of rainfall
Chae et al. Landslide prediction, monitoring and early warning: a concise review of state-of-the-art
Frappart et al. Radar altimetry backscattering signatures at Ka, Ku, C, and S bands over West Africa
Baade et al. TanDEM-X IDEM precision and accuracy assessment based on a large assembly of differential GNSS measurements in Kruger National Park, South Africa
Li et al. Water level changes of Hulun Lake in Inner Mongolia derived from Jason satellite data
Sabel et al. Development of a Global Backscatter Model in support to the Sentinel-1 mission design
Morel et al. Validity and behaviour of tropospheric gradients estimated by GPS in Corsica
Gonçalves et al. Monitoring Local Shoreline Changes by Integrating UASs, Airborne LiDAR, Historical Images and Orthophotos.
Frappart et al. The 2013 Ibiza calibration campaign of Jason-2 and SARAL altimeters
Zhao et al. An optimal tropospheric tomography approach with the support of an auxiliary area
Darvishi et al. Performance evaluation of phase and weather-based models in atmospheric correction with Sentinel-1data: Corvara landslide in the Alps
Sutterley et al. Evaluation of reconstructions of snow/ice melt in Greenland by regional atmospheric climate models using laser altimetry data
Klemas et al. Remote sensing of soil moisture: An overview in relation to coastal soils
Orlandi et al. Vertical accuracy assessment of the processed SRTM data for the Brazilian territory
Phoeurn et al. Assessment of satellite rainfall estimates as a pre-analysis for water environment analytical tools: A case study for Tonle Sap Lake in Cambodia
Pirmoradian et al. Performance evaluation of IMERG and TMPA daily precipitation products over CONUS (2000–2019)
Zhang et al. Reduction of atmospheric effects on InSAR observations through incorporation of GACOS and PCA into small baseline subset InSAR
Zhang et al. Assessment of ERA-Interim and ERA5 reanalysis data on atmospheric corrections for InSAR
Rawat et al. Effect of DEM data resolution on low relief region sub-watershed boundaries delineating using of SWAT model and DEM derived from CARTOSAT-1 (IRS-P5), SRTM and ASTER
Mohamad et al. Monitoring groundwater depletion due to drought using satellite gravimetry: A review
Yue et al. Accuracy assessment of SRTM V4. 1 and ASTER GDEM V2 in high-altitude mountainous areas: A case study in Yulong Snow Mountain, China
Worqlul et al. Comparison of TRMM, MPEG and CFSR rainfall estimation with the ground observed data for the Lake Tana Basin, Ethiopia

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 20811050

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 20811050

Country of ref document: EP

Kind code of ref document: A1