Node-Based Optimization of GNSS Tomography with a Minimum Bounding Box Algorithm

Abstract

Global Navigation Satellite Systems (GNSS) tomography plays an important role in the monitoring and tracking of the tropospheric water vapor. In this study, a new approach for improving the node-based GNSS tomography is proposed, which makes a trade-off between the real observed region and the complexity of the discretization of the tomographic region. To obtain dynamically the approximate observed region, the convex hull algorithm and minimum bounding box algorithm are used at each tomographic epoch. This new approach can dynamically define the tomographic model for all types of study areas based on the GNSS data. The performance of the new approach is tested by comparing it against the common node-based GNSS tomographic approach. Test data in May 2015 are obtained from the Hong Kong GNSS network to build the tomographic models and the radiosonde data as a reference are used for validating the quality of the new approach. The experimental results show that the root-mean-square errors of the new approach, in most cases, have a 38 percent improvement and the values of standard deviation reduce to over 43 percent compared with the common approach. The results indicate that the new approach is applicable to the node-based GNSS tomography.

1. Introduction

Water vapor plays an important role in the Earth’s weather and climate systems, the water cycle and the regulation of air temperature. Since the near-Earth surface includes nearly all of the atmosphere’s water vapor, remote sensing techniques focus on the area of the lowest layer of Earth’s atmosphere, which is the troposphere. Global Navigation Satellite Systems (GNSS) tomography is one of the primary techniques for determining the three-dimensional distribution of water vapor in the tomographic region. The GNSS tomography has proved its capacity to utilize satellite constellation and ground observation stations (i.e., GNSS networks) to build the water vapor field and the meteorological applications. The wet delay of the GNSS signals, caused by atmospheric water vapor, is used to estimate the amount of water vapor between a GNSS satellite and each GNSS station. The total integrated information, which includes spatial information about the water vapor, can be combined into a spatially resolved field through tomographic reconstruction techniques. Due to its ability to provide all-weather, high-efficiency, high-resolution vertical and horizontal measurements of water vapor, the GNSS tomography has potential value for severe weather study and short-term weather forecasting in the middle and small scale areas. The size of the areas used in this process depends on the distribution of GNSS signals between the ground-based GNSS network and the GNSS satellite constellation. Today, most tomographic regions used in scientific research are cuboid-shaped in their region. For example, Flores et al. [1] use the Kilauea GNSS network to cover a region of 20 × 20 km² in the horizontal direction and 15 km in altitude. Other examples of the approximate tomographic region dimensions from current research are shown in

Table 1. The dimensions of the Global Navigation Satellite Systems (GNSS) tomographic region.

Paper	Flores et al. [1]	Hirahara [10]	Seko et al. [11]	Champollion et al. [12]
Tomographic region (km3)	20 × 20 × 15	27 × 23 × 10	30 × 30 × 9	20 × 20 × 10
Paper	Troller et al. [13]	Notarpietro et al. [14]	Rohm et al. [15]	Perler [16]
Tomographic region (km3)	300 × 150 × 15	30 × 20 × 10	50 × 40 × 10	266 × 166 × 15
Paper	Xia et al. [17]; Ye et al. [18]; Ding et al. [19]; Yao [20]
Tomographic region (km3)	The Hong Kong Satellite Positioning Reference Station Network (SatRef) roughly covers an area of 75 × 60 × 10

. Although the cuboid is the most commonly used tomographic region, there are still some approaches to extend or refine the cuboid-shaped region based on the distribution of the GNSS signals. For example, Rohm et al. [2] propose the double models (i.e., inner model and outer model) for the GNSS tomography to extend the tomographic region and increase the information for estimating the results. Chen and Liu [3] compress the height of the tomographic region from 15 km to 8.5 km to introduce more GNSS signals and then move the tomographic region to the optimal location where the number of the voxels crossed by signals reaches its maximum. Ding et al. [4] propose an adaptive node parameterization approach to adjust the tomographic region at each tomographic epoch dynamically. The methods of using the signals penetrating from the side face of the tomography area [5,6,7] and the data of GNSS observations outside the study area are also proposed. Zhang et al. [8] optimize the dynamic scope of the tomographic area and propose a height factor model to utilize GNSS rays passing from the side boundary [9].

The tomographic region usually depends on the spatial distribution of the GNSS signal from the ground-based GNSS receivers. However, the most common methods used to determine the dimensions of the tomographic model only consider the GNSS signals’ distribution through the top boundary of the model. This is based on the assumption that, for the selected model, any GNSS signals that pierce the top of the model must also be distributed through the lower levels of the model. This assumption is based on the inverted cone-shaped region defined by the GNSS signal distribution from the Earth’s surface to the GNSS satellites. The selected rectangular tomographic region is then chosen to include the top-level GNSS signals approximately. This predetermined approximate tomographic region causes some problems. First, the box-shaped region includes several large empty spaces, which have no GNSS signals passing through them. Second, the voxels and their nodes (as the cell of parameterization) in the empty spaces are unnecessary for the tomographic model and become a computational burden. Third, the fixed region for all epochs has larger redundant space (i.e., the region must ensure that various distributions of signals at each epoch are covered), which does not fit the distribution of GNSS signals at a defined specific epoch.

To address the above problems, a node-based innovative approach to dynamically determine the dimensions of the tomographic region at each epoch by the minimum bounding box algorithm is proposed. In each epoch, the new approach defines the two-dimensional tomographic boundaries at different height of tomographic planes, which comprise the tomographic region in the three-dimensional space. Two steps determine the tomographic boundaries: (1) the first approximation model using the convex hull algorithm to define the boundary (i.e., convex polygon) of the collection of rays; (2) the secondary approximation model using the minimum bounding box algorithm to define the final tomographic boundary (i.e., minimum-area rectangle) that encases this convex polygon. By using the secondary approximation modeling method, the position of the nodes can be easily assigned in the tomographic region. The expected intervals, defined in Section 3.2, were used to adjust the nodes’ distribution at each tomographic plane. When the nodes have been determined, the discretization of the tomographic region is finished. Moreover, the new approach is not limited to the test in this study but can be utilized for all GNSS networks. In this paper, we focus on the improvement of the common node-based GNSS tomographic approach. The common approach was used for comparative analysis and validation of this new approach.

This paper mainly discusses the new approach that improves node-based GNSS tomography with a minimum bounding box algorithm. This paper is structured as follows. Section 2 introduces the GNSS tomography technique. This section details the observations used in the tomographic model and the common tomographic method. Section 3 proposes a new approach, which can dynamically define the tomographic model. Section 4 presents the test results and analysis. At the end, a few conclusions are drawn in Section 5.

2. Methodology

2.1. Observations Used in the Tomographic Model

For successful tomographic modeling, the slant wet delay (SWD) or slant water vapor (SWV) needs to be estimated in the tomographic equations. As part of the wet, the SWD is separated from the slant total delay (STD), which is produced due to the GNSS signals across the entire atmosphere. The SWD can be approximated by

$$\mathrm{S}\mathrm{W}\mathrm{D}(el,az)={\mathrm{m}}_{\mathrm{w}}(el)\cdot \mathrm{Z}\mathrm{W}\mathrm{D}+\left[{\mathrm{m}}_{\mathrm{G}}(el)\cdot \left({\mathrm{G}}_{\mathrm{W}}^{\mathrm{N}}\cdot \cos (az)+{\mathrm{G}}_{\mathrm{W}}^{\mathrm{E}}\cdot \sin (az)\right)+R_{\mathrm{c}}\right]$$

(1)

where SWD(el,az) is the function of the SWD, m_w (el) is the wet mapping function; ZWD is the zenith wet delay; m_G (el) is the gradient mapping function; ${\mathrm{G}}_{\mathrm{W}}^{\mathrm{N}}$ is the north–south gradient; ${\mathrm{G}}_{\mathrm{W}}^{\mathrm{E}}$ is the east–west gradient; R_c is the cleaning postfit residuals [21].

The SWD is used as the tomographic observation, or it can be further converted into SWV by

$$\mathrm{S}\mathrm{W}\mathrm{V}=Scal\cdot \mathrm{S}\mathrm{W}\mathrm{D}$$

(2)

where Scal is the conversion factor [22,23].

In GNSS tomographic solution process, the SWV is then used as the GNSS tomographic observation for the estimation of water vapor densities at each node in the tomographic model.

2.2. Classical Tomographic Model

The commonly used tomographic approaches are based on a predetermined cuboid-shaped region, which is a fixed region selected to approximately cover the complete distribution of GNSS signals for all epochs. Two spatial discretization models are adopted for the regional division. As shown in Figure 1, one is the voxel-based model (top panel) and another is the node-based model (bottom panel). The voxel-based model uses voxels to discretize the tomographic region and only within the yellow-shaded voxels do GNSS signals actually pass. Conversely, as mentioned before, the empty spaces in the same tomographic model will still contain many voxels (the non-yellow-shaded boxes), which should not be used in the tomographic solution. Moreover, SWV parameters must still be estimated in these voxels by the horizontal/smoothing constraint [24,25,26,27] for the tomographic modeling to be resolved at the expense of reduced accuracy.

The same problem exists when using a node-based model, which uses nodes at different tomographic planes (different heights) to parameterize the tomographic region. This approach also includes many problematic nodes (white circles in the bottom panel in Figure 1) without surrounding GNSS signals. However, in this instance, only the yellow nodes (yellow circles in the bottom panel) are needed in the tomographic solution. However, the layout of nodes is more flexible than that of voxels since the voxels in the tomographic region need to adjoin each other, but the nodes in the tomographic region do not have this constraint. Therefore, the node-based tomographic region can be improved easily to match the real observation region (i.e., the coverage areas of the GNSS signals).

3. New Node-Based Tomographic Approach

The present research presents an innovative tomographic modeling approach that takes advantage of the greater flexibility of the node-based GNSS tomographic approach. It also adapts the tomographic region and the position of nodes to make sure that the adjusted region fits the distribution of the GNSS signals at the specified epoch. This new approach also adjusts the nodes’ distribution as evenly distributed as a classical method. Further, the new tomographic model also significantly reduces the parameters to be solved by only solving nodes with nearby GNSS signals. This approach also increases not only the computing efficiency but also the accuracy of the tomographic results.

In the new approach, GNSS signal paths are assumed to be straight lines for simplification. A temporal resolution of 20-min GNSS observations is selected for the tomographic sampling epoch. Different vertical intervals between two consecutive tomographic planes at different altitudes are adopted due to the nonuniformness of water vapor distribution in the vertical dimension. The new node-based tomography modeling process is achieved through the following steps: (1) determine the tomographic region dimension at each epoch by the minimum bounding box algorithm; (2) determine the nodes position and density at each of the different tomographic planes; (3) construct and solve the tomographic equations. Each of the above steps is elaborated below.

3.1. Determination of a Tomographic Region

In this study, the dimensions of the tomographic region are determined by searching for the minimum bounding box enclosing the intersections between the paths of all GNSS signals and the plane at different height of the tomographic planes (the intersection points are named pierce points). Freeman and Shapira [28] demonstrate clearly that the minimum bounding box of a set of points is also the minimum area rectangle of the points’ convex hull polygon. Thus, the primary goal is to determine the convex hull of the pierce points on each tomographic plane.

3.1.1. Convex Hull of the Pierce Points

In each of the tomographic planes, the Graham scan [29] is used to find the so-called convex hull of the pierce points by efficiently utilizing a stack to detect and remove concavities in the boundary. As shown in Figure 2, the Graham scan procedure mainly includes four steps: (1) Determination of the start point with the lowest y-coordinate. As shown in the top-left panel, the P₀ is the start point. If there is more than one point that exists with the lowest y-coordinate, the point with lowest x-coordinate is the start point. (2) All points except P₀ should be sorted in increasing order of the polar angle. The P₀ is the origin of the polar coordinate, which is so-called pole and the other points have the subscripts of unique sequence number serving as the index (i.e., P₁, P₂, P₃ and P₄). (3) For each point, the approach needs to make sure that whether the two previously considered points (e.g., P₀ and P₁ in the top-left panel) to this point (e.g., P₂) is a “left turn” or a “right turn” and these three points constitute a system to detect the points of the convex hull. “Left turn” means that along the direction of the first two points (e.g., P₀ to P₁) of the system, the last point (e.g., P₂) on the left side of the direction. Similarly for “right turn” the last point of these three points on the right side of the remaining two points’ direction. If it is a left turn, the second-to-last point (e.g., P₁ in the top-left panel) is part of the convex hull. As shown in the top-right panel, the approach moves on to the next point in the sorted array and the detection system consists of P₁, P₂ and P₃. The same determination is then made for the system until a “right turn” is encountered. In the bottom-left panel, P₄ is a “right turn” and this means that the second-to-last point (P₃) is not at the convex hull and should be removed from the array (stack). (4) Again, determining whether three points constitute a “left turn” or a “right turn” and this process will eventually return to the start point P₀ and the stack now contains the points on the convex hull in counterclockwise order (i.e., P₀, P₁, P₂ and P₄ in the bottom-right panel).

Based on the above steps, the convex hull of the points can then be determined, which initially discretizes the tomographic region. First, a tomographic region for a particular epoch (top panel in Figure 3) is divided vertically into several horizontal planes with nonuniform vertical intervals from 300 to 3800 m. The blue points in the bottom panel represent the pierce points of GNSS signals (red lines). After the discretization of the tomographic region in the vertical direction, the convex hulls of the pierce points on all the tomographic planes (i.e., all the black loops in bottom panel) are then defined.

In classical node parameterization approaches, the rectangular shape for the tomographic region is adopted (shown in Figure 1), from which the nodes’ position can be determined easily. In contrast, the tomographic region constructed by the convex hull algorithm is irregular (shown in Figure 3). This irregularity leads to difficulty in uniformly presetting the distribution of the nodes. In order to determine the best balance between the proximity of the distribution of GNSS signals to the nodes and the complexity of the nodes’ layout, the minimum bounding box algorithm is used for building the tomographic region. The procedure for this approach is discussed in the next section.

3.1.2. Minimum Bounding Box Algorithm

After the above process, the tomographic region is initially discretized by convex hulls at different tomographic planes (see Figure 3). The next step is to build the secondary approximation region by a minimum bounding box algorithm. As shown in Figure 4, the minimum bounding boxes of these convex hull polygons are defined by the theorem defined by [28] that the rectangle of the minimum area enclosing a convex polygon has a side collinear with one of the edges of the polygon.

As an example, the convex hull polygon at the top of the tomographic region in Figure 4 is used to define the corresponding minimum bounding boxes. This process is illustrated in Figure 5 and includes three steps: (1) First, each edge of the convex hull is set as an axis to find the bounding box and the ordinate origin of the axis is the middle of the edge (e.g., top panel). This is due to the fact that each edge of the minimum bounding box must be collinear with the one of the convex hull’s edge. Second, the point with the longest vertical distance between the edge and the points of the convex hull (e.g., R1) and the two points with the absolute values of the maximum positive and negative coordinates (e.g., R2 and R3) are used to build the bounding box, which is also collinear with the one of convex hull’s edge. Third, all the areas of the bounding box are calculated and then the lowest area is selected as the minimum bounding box. In the bottom panel, the blue rectangle is the minimum bounding box comparing with the other bounding boxes (dotted rectangles).

3.2. Determination of the Position and Density of Nodes

For the tomographic region, the convex hulls in Figure 3b are considered to be the first approximation and the minimum bounding boxes in Figure 4 are the secondary approximation. The second approximation keeps the balance between the proximity of the nodes to the distribution of GNSS signals and the complexity of the node-set. In this section, the nodes with different water vapor parameters are assigned based on the shape of the minimum bounding box.

It should be noted that the areas of the minimum bounding boxes at different heights of the tomographic planes are different. For example, the rectangular area at the bottom of the tomographic region is 1450 km² and that of the top tomographic boundary is 4309 km². However, the number of GNSS signals that pierce through the different planes at each height of the tomographic region is the same. This implies that the density of the GNSS signals in the lower layers is higher than that in the higher layer. So, for different heights, the densities of the nodes should be adjusted to match the tightness of the GNSS signals. The density of nodes can be expressed by

$${\mathrm{F}}_{\mathrm{L}}^{\mathrm{i}}=\left[\begin{array}{c}{\mathrm{D}}_{\mathrm{L}}^{\mathrm{i}}/{\mathrm{f}}_{\mathrm{L}}^{\mathrm{i}}\end{array}\right]$$

(3)

$${\mathrm{F}}_{\mathrm{s}}^{\mathrm{i}}={\left[{\mathrm{D}}_{\mathrm{S}}^{\mathrm{i}}/{\mathrm{f}}_{\mathrm{S}}^{\mathrm{i}}\right]}_{}$$

(4)

where ${\mathrm{D}}_{\mathrm{L}}^{\mathrm{i}}$ and ${\mathrm{D}}_{\mathrm{S}}^{\mathrm{i}}$ are the length of the long side and short side of the minimum bounding box at i th tomographic plane; ${\mathrm{f}}_{\mathrm{L}}^{\mathrm{i}}$ and ${\mathrm{f}}_{\mathrm{S}}^{\mathrm{i}}$ are the expected intervals for the long side and short side; the brackets round ${\mathrm{D}}_{\mathrm{L}}^{\mathrm{i}}$/${\mathrm{f}}_{\mathrm{L}}^{\mathrm{i}}$ and ${\mathrm{D}}_{\mathrm{S}}^{\mathrm{i}}$/${\mathrm{f}}_{\mathrm{S}}^{\mathrm{i}}$ to the nearest integer; ${\mathrm{F}}_{\mathrm{L}}^{\mathrm{i}}$ and ${\mathrm{F}}_{\mathrm{S}}^{\mathrm{i}}$ are the segments for the long side and short side, respectively. As shown in Figure 6, when the ${\mathrm{F}}_{\mathrm{L}}^{\mathrm{i}}$ and ${\mathrm{F}}_{\mathrm{S}}^{\mathrm{i}}$ are defined, the nodes (black square) can be determined easily at the tomographic plane.

As mentioned before, the densities of the nodes must be appropriate for the size of the minimum bounding boxes at a different height. Therefore, the variables ${\mathrm{f}}_{\mathrm{L}}$ and ${\mathrm{f}}_{\mathrm{S}}$ vary at different heights of the tomographic region (details seen in

Table 2. The size of ${\mathrm{f}}_{\mathrm{L}}^{\mathrm{i}}$ and ${\mathrm{f}}_{\mathrm{S}}^{\mathrm{i}}$ at different heights of the minimum bounding boxes.

Height (m)	Below 2400	2800–4500	5000–5500	6500–8500	11000
${\mathrm{f}}_{\mathrm{L}}^{\mathrm{i}}$ and ${\mathrm{f}}_{\mathrm{S}}^{\mathrm{i}}$ (m)	7000	8000	9000	9800	15000

). When the nodes at each of the tomographic planes have been determined, the discretization of the region is finished. These nodes are herein named tomographic nodes.

3.3. Construction and Solution of Tomographic Equations

Next, the tomographic equations are introduced to establish the relationships between the tomographic observations and tomographic nodes. These tomographic equations are solved by the algebraic reconstruction technique (ART). ART is used in computed tomography and it is also successfully applied to GNSS tomography [30]. ART is considered an excellent iterative solver of a system of linear equations since its high numerical stability, even under adverse conditions.

The procedure for the construction of the tomographic equations includes four main steps: (1) The SWV quantity in the tomographic observation can be expressed by section integration. The integration paths (parts of the ray path) of which are divided by the adjacent tomographic planes. (2) Each of these section integrations is first estimated approximately by the Newton–Cotes formulas (i.e., formulas for quadrature based on evaluating the integrand at equally spaced points). (3) The mathematical relationship between the equally spaced points and tomographic nodes is established by the distance inverse weight method and vertical interpolation function in both the horizontal and vertical direction, respectively. Finally, the tomographic equations are constructed to link the tomographic observations and the tomographic nodes. These steps are elaborated in the Ding et al. [31] and the tomographic equations are defined by

NX = B

(5)

where the N is the coefficient matrix of the tomographic model; X is the vector of the unknown parameters of water vapor densities; B is the vector of the SWVs. The access order scheme based on prime number decomposition (PND) and ART [19] is used to define the iteration sequences of (5) and then to estimate the tomographic results.

4. Test Results and Analysis

Data from the Hong Kong Satellite Positioning Reference Station Network (SatRef) were processed by the new node-based tomographic approach and the common node-based tomographic approach, for verification. These results were also compared against water vapor data, from King’s Park Meteorological Station (HKKP) to assess their quality. The horizontal and vertical distributions of these stations are presented in Figure 7. The initial data contain information on 17 reference stations from the SatRef, which are evenly distributed in Hong Kong (see Figure 7). In this study, a window, on the top of the tomographic region, is used as a sampling method for data acquisition. The window size is shown in Figure 7, which is an attempt to utilize the GNSS signals with large satellite elevation angles in the denser regions. Due to these GNSS data, the most important datasets generally have relatively small errors (i.e., relatively large satellite elevation angles). The new approach can build a high-quality tomographic model. The tomographic data use the Day-of-Year (DOY) numbering system. This is a common format used in research data. The DOY is the sequential day number starting with day 1 on 1 January. The data throughout May 2015 (i.e., DOY 121-151) with a fixed sampling interval (the 30s). For simplicity, the results from the common node-based tomographic approach and the new node-based tomographic approach are called CNT and NNT, respectively. The water vapor data from HKKP are called WVD for short. Some statistics, e.g., the root-mean-square error (RMSE) and bias, are applied to quantify the results of CNT and NNT. The values of the tomographic results at all WVD sampling points were calculated first using the interpolation method, then the statistics between the interpolated values and WVD observations at all the sampling points of the WVD profile from the ground surface to 11000m at each epoch was calculated.

For the new and the common node-based tomographic modeling, the SWV observations were collected with a sampling rate of 5min. The tomographic results at Coordinated Universal Time (UTC) 12 throughout the whole month are compared with the water vapor data acquired by a sounding balloon at HKKP. Since the reference data from HKKP only provides one moisture profile at each epoch, the interpolated tomographic results at the position of the HKKP are used for validation of the new node-based tomographic model.

Figure 8 shows the RMSE estimates for the CNT and NNT results at UTC 12 every day in May. This is calculated by subtracting the WVD from the tomographic results (CNT or NNT) at all sampling points of the sounding balloon below 10,800 m. The figure reveals that the NNT, for most of the time, outperformed the CNT technique except for DOY’s 127, 131 and 141, which are marked by yellow circles in Figure 8. Overall, the mean percentage improvement of RMSE between the CNT to NNT is 38 percent. Although CNT fared better than NNT in three days of this month, the mean percentage decline of RMSE from NNT to CNT is only 9 percent.

Focusing more clearly on these DOYs (see Figure 9), there is no apparent difference between the NNT and CNT derived water vapor profiles for these three days. DOY 139 (bottom right panel), shows the most significant differences in RMSEs values between NNT and CNT, which reveals that the results from CNT are poorer on this day. This may also suggest that in some instances, the GNSS tomographic solutions from CNT are largely beyond the predicted results and that the tomographic model of the CNT may be unstable.

Table 3. Statistics of two approaches throughout the entire month of May.

	Statistics	RMSE (g/m3)	Bias (g/m3)	PCC (%)
Approach
CNT		1.11	−0.145	96.4
NNT		0.981	−0.061	98.7
Improvement (%)		11.6	57.9	2.39

contains the RMSE, bias and the Pearson correlation coefficient (PCC) values for the CNT and NNT results at UTC 12 for the whole month of May. It can be seen that NNT outperforms CNT in all statistics. The improvement of RMSE is 11.6 percent from CNT to NNT, which indicates that the accuracy of NNT is higher than that of the CNT. The statistical bias result, which is a frequently used measure of the systematic error, also shows that the model from the NNT is more stable than the CNT model. The Pearson correlation coefficient of the NNT shows a weaker correlation than the CNT. This may be because the PCC has a great cardinal number. Even if the CNT technique produces a poor PCC result (e.g., 90 percent), the improvement from the CNT to NNT technique cannot exceed 10 percent.

Figure 10 presents the scatterplots of the calculated WVD value from the CNT and NNT tomographic techniques tomographic results (CNT/NNT), respectively. The standard deviation (SD) also present in Figure 10. The performances of the CNT and NNT technique can be evaluated directly by considering the correlation of the values with the baseline for the accuracy of the tomographic results (blue line). Values that are closer to this represent better tomographic results. The figure implies that the NNT technique (bottom panel) has higher accuracy in comparison with the CNT (top panel). Alongside, the statistics presented in Table 3, these figures represent a high degree of consistency in the CNT results.

5. Conclusions

A new node-based tomographic (NNT) approach is proposed to define the tomographic model at each epoch dynamically. This contains two parts: (1) determination of the tomographic region, and (2) determination of the position of the tomographic node at each tomographic epoch. First, the convex hull algorithm and the minimum bounding box algorithm are used to dynamically adjust the tomographic region, thus ensuring the region consistent with the distribution of the GNSS signals. Second, the tomographic node position is defined based on the expected intervals at each tomographic plane, which are used to match the density of the GNSS signals at different tomographic planes.

NNT is compared against the common node-based tomographic approach (i.e., CNT) for accuracy verification. Radiosonde data from King’s Park Meteorological Station were introduced for assessing the quality of NNT and CNT technique. Overall, daily and monthly statistics (at UTC 12 every day) of the NNT and CNT results throughout the entire month of May 2015 indicated that the NNT technique outperformed the CNT technique, except for three days (i.e., DOY 127, 131 and 141). However, the difference (9 percent in RMSE) between the NNT and CNT results can be neglected compared with that of the other days (above the average of 38 percent in RMSE). NNT also makes about a 43 percent improvement in SD. NNT was shown to have a smaller systematic error, evident by the bias statistic (in Table 3) and the plots in Figure 10b.

In the future, the work will be focusing on improving the mesh generation of the GNSS tomographic region, which is extremely important for GNSS tomography.