###### Abstract  
This study evaluates the capability of an unconstrained tomographic algorithm to capture 3D water vapor density variability throughout 2017 in Houston, U.S. The algorithm relies solely on Global Navigation Satellite System (GNSS) observations and does not require an initial guess or other specific constraints regarding water vapor density variability within the tomographic domain. The test domain, featuring 9 km horizontal, 500 m vertical, and 30 min temporal resolutions, yielded remarkable results when compared to data retrieved from the ECMWF Reanalysis v5 (ERA5), regional Weather Research and Forecasting Model (WRF) data, and GNSS-Radio Occultation (RO). For the first time, a time series of Preciptable Water Vapor maps derived from the Interferometric Synthetic Aperture Radar (InSAR) technique was used to validate the spatially integrated water vapor computed by GNSS tomography. Tomographic results clearly indicate the passage of Hurricane Harvey, with integrated water vapor peaking at 60 kg/m\({}^{2}\) and increased humidity at altitudes up to 7.5 km. Our findings suggest that GNSS tomography holds promise as a reliable source of atmospheric water vapor data for various applications. Future enhancements may arise from denser and multi-constellation networks. GNSS tomography; GNSS meteorology; water vapor density; precipitable water vapor; tomographic algorithm; atmospheric water vapor +
Footnote †: journal: Remote Sensing  
0000-0001-5000-0002-9999-9999-9999-9999-99999-9999-9999-99999-99999-99999-99999-99999-99999-99999-99999-999999-999999-999999-999999-999999-999999-999999-99999  
with incomplete data due to the heterogenous distribution of observations through the 3D tomography grid. One way to mitigate this limitation is to densify the GNSS network; however, that solution is costly and may not solve the problem depending on the grid horizontal resolution and the GNSS elevation cut-off mask used, typically 7-15\({}^{\circ}\) to avoid the delay error induced by incorporating low elevation observations. Another way consists of adding some inter-voxel physical constraints to add information to voxels not crossed by any ray (e.g., [2; 3; 4; 5]) or adding a priori information from ground-based, air-borne, or space-borne instruments (e.g., [6; 7; 8; 9]) from Numerical Weather Prediction models (NWP) data (e.g., [8; 10; 11; 12; 13]) or even creating virtual data (e.g., [14; 15; 16; 17; 18; 19]). Constraints and external information become necessary to overcome the rank deficiency problem introduced by empty voxels. However, the uncertainty introduced by that information can affect the tomography solution, restricting the mathematical model's ability to fully simulate the dynamics of water vapor in the atmosphere and compromising the model's reliability. Another tomography challenge is the inversion algorithm. Some studies apply the generalized inverse method or a derivative [2; 17; 20; 21; 22] but these approaches can reveal numerical problems in the case of huge matrices and require large computational resources [23]. An alternative that avoids matrix inversion and even the creation of large matrices is the Algebraic Reconstruction Technique (ART) family. ART is an iterative algorithm that solves a system of linear equations by iteratively refining an initial guess [23]. A drawback of this technique is the requirement of solid initialization data, a general relaxation factor, and a criterion for stopping the iteration, which can depend on a specific configuration, leading to a degrading or oscillating solution [24]. Another technique applied in several studies is the Kalman Filter (KF), which combines the last tomographic solution and the observations using the corresponding error covariance matrices, providing an update of the error covariance for the tomographic solution, but realistic error estimates required in the KF are challenging to obtain [24; 25; 26; 27]. An intercomparison between algorithms can be seen in [10; 23; 28; 29; 30].