Introduction

Hurricanes are some of the most potent hydrometeorological hazards and can cause severe damage to the coastal regions they strike. Hurricanes are thermally powered, large-scale storms with wind speeds exceeding 33 m/s that circulate around a well-defined core and are characterized by strong near-surface torrential rainfall (Emanuel 2017a, 1989). Over the years, the dissipation of power in hurricanes has increased as their strengths have intensified (Elsner et al. 2008; Klotzbach 2006; Bhatia et al. 2019). More moisture evaporates into the atmosphere over a warmer ocean, which invigorates and enlarges hurricanes (Trenberth et al. 2018). Some studies have suggested that rising temperatures are to blame for the increased intensity of hurricanes that have been observed in recent years (USGCRP 2017). However, predicting the intensity of a hurricane remains challenging (Emanuel 2017b; Elsberry 2014). The already complex underlying physics of hurricanes are further compounded by the combination of warm ocean water, moist air, and winds (Emanuel et al. 2006). Meteorologists attempt to forecast the intensity and locations of tropical cyclones days in advance. Nevertheless, the reliability of a prediction is complicated when the physical mechanisms that control the creation of a hurricane are not well understood. For several decades, atmospheric researchers have sought to unravel the secrets behind the nature of hurricanes (Rios-Berrios et al. 2014). Various thermodynamic parameters are necessary for the formation of hurricanes, but the conditions sufficient for the development of hurricanes remain poorly understood (Emanuel 1989). Nevertheless, new technological land-based observations and recent advances have increased the accuracy of tropical cyclone (TC) forecasting assessments (Stith et al. 2019). Further, advances in satellite remote sensing have enabled additional inroads in TC tracking (Ackerman et al. 2019). Infrared measurements obtained from space platforms to study TC characteristics and intensity have also been improved (Olander and Velden 2007). Moreover, for higher spatial and temporal resolutions, further research on TCs may be conducted as a result of the recent Cyclone Global Navigation Satellite System (CYGNSS) satellite missions (Ruf et al. 2016). These satellites receive forward-reflected GPS signals and measure ocean surface roughness, from which wind speed is recovered without appreciable attenuation by rain. However, intensity forecasting continues to be challenged by an inadequate understanding of the boundary layer and air–ocean exchange physics and by the deficiencies of real-time measurements of ocean mixed layer properties affecting the thermal signature of the subsurface ocean heat content under high winds (Emanuel 2018).

It is well established that TCs, i.e., hurricanes, draw heat energy from warm ocean areas through sensible and latent heat fluxes (Emanuel 1991). This energy transfer occurs mainly through the evaporation of water into the atmosphere, consequently boosting radial circulation within the TC by conserving the rotational momentum of the rotating air surface (Smith 2000). The band of air and wind that surrounds the eye of a hurricane acts as a conduit for vertical wind, which drags water vapor up from the warm ocean surface (Emanuel 1989). Hence, warm ocean surface temperatures are necessary for the development of a hurricane, and hurricane prediction schemes could be improved by monitoring atmospheric water vapor (Rios-Berrios et al. 2014). Atmospheric water vapor has been identified as an essential indicator for deep convection, leading to the formation of a hurricane once a certain threshold is met (Wang et al. 2010a, b). In the western Northern Hemisphere, hurricanes form mainly off the Gulf and East Coasts of the USA; these hurricanes are instigated by easterly winds that form over western Africa (Mo et al. 2001). Hurricane activity commonly moves from these formation locations upward and eastward on route to the mainland. The energy of hurricanes invigorates severe weather events when they discharge their energy to the mainland. The 2017 Atlantic hurricane season was one of the most active on record, resulting in the third-highest number of major hurricanes in a single year over the past century, surpassed only by the 1961 and 2005 seasons (Lim et al. 2018). Among the 17 TCs that occurred during the 2017 hurricane season, 10 developed rapidly to a hurricane level and had an unusually long lifespan, and six became major hurricanes (Gert, Harvey, Irma, Jose, Lee, and Maria). Hurricanes Harvey and Irma were category 4 hurricanes that produced more rainfall along the Gulf and East Coasts of the USA than any other storm that season, and thus, these hurricanes provide the conditions for an outstanding case study. These hurricanes proved to be natural catastrophes for the nation, resulting in numerous deaths, displacing several thousands of people, and causing widespread socioeconomic loss (Blake and Zelinsky 2018).

On August 25, 2017, Hurricane Harvey devastated large parts of Texas with a severity that was second only to that of Hurricane Katrina in 2005. The effect of Hurricane Harvey’s landfall was exacerbated by recording high ocean temperatures that fuelled and sustained the hurricane’s intensity; the land surface experienced subsidence due to the weight of the storm water (Milliner et al. 2018). Immediately following the cessation of Hurricane Harvey, Hurricane Irma followed and quickly grew into the strongest hurricane since Hurricane Katrina. On September 6, 2017, Hurricane Irma crossed several islands and made three landfalls in a single day in Barbuda, St. Martin and the British Virgin Islands. Hurricane Irma continued making landfalls in the Bahamas, Cuba, the Florida Keys and southwest Florida between September 9 and 10, 2017 (Rahmstorf 2017). Extreme meteorological events such as Hurricanes Harvey and Irma emphasize the need to track TCs using further observational tools to mitigate their risks.

The Global Positioning System (GPS) accurately captures water vapor distributions in the atmosphere under all weather conditions (Bevis et al. 1992). As the theoretical basis for these measurements, the propagation delay that affects GPS signals depends on the amount of water vapor in the atmosphere. This water vapor sensitivity of GPS delay measurements has been used to derive integrated water vapor (IWV) in the atmosphere (Bevis et al. 1992). In regions affected by dense clouds and heavy precipitation, where visible, infrared, and microwave-based satellite measurements are extensively contaminated, GPS sensors can provide independent measurements of atmospheric processes (Vergados et al. 2013). Furthermore, meteorological instruments such as radiosondes and spectroradiometers cannot measure the water vapor field with the fine resolution necessary to characterize its variability (Rocken et al. 1995; Gao and Kaufman 2003).

Currently, GPS tropospheric products are employed mainly for the investigation of mesoscale weather systems and near-real-time applications (Falvey and Beavan 2002; Nakamura et al. 2004; Kawabata et al. 2013). These tropospheric products are gradually being assimilated into numerical weather prediction (NWP) models for forecasting purposes (Wilgan et al. 2015), which has had a positive impact on short-range moisture field forecasts, with the ultimate goal of improving precipitation predictions for severe rainfall events (Bennitt and Jupp 2012). Developments in real-time GPS and in the now-casting NWP model make these types of applications possible by providing IWV measurements that complement meteorological predictions at some level (Li et al. 2015). Liou and Huang (2000) and Sapucci et al. (2016) reported a sharp increase in GPS-IWV prior to a rainfall event. In addition, recent studies have investigated the accuracy of GPS tropospheric products by identifying factors that affect the quality of GPS processing (Kačmařík et al. 2017; Klos et al. 2018, Ejigu et al. 2019). Furthermore, in areas where a dense network of GPS stations is available, it is possible to produce optimal GPS-IWV distribution maps to investigate spatiotemporal changes in water vapor. Thus, the use of GPS tropospheric products for severe storm prediction is a popular research focus with applications that the geodetic community wishes to promote among meteorological research groups. At present, the Gulf and East Coasts of the USA are two of the few regions covered by a dense array of GPS stations; accordingly, this study took advantage of the freely available high-quality GPS data from this array to estimate the zenith total delay (ZTD). Milliner et al. (2018) were the first to monitor the amount of water drained by the ground from hurricanes across the Gulf and East Coasts of the USA. More recently, Graffigna et al. (2019) showed a significant shift in ZTD gradients prior to the arrival of a hurricane. Ejigu et al. (2020) assessed the connection between GPS-derived IWV and its influence on TC structure during Hurricane Florence in 2018.

The IWV was employed for the first time in this study to track the major hurricanes of the 2017 season. We compared the temporal distribution of the GPS-IWV content to the amount of precipitation for days when a hurricane appeared in the area under hurricane influence. The primary motivation for this work was to track and characterize hurricane movements using GPS-IWV maxima patterns. We created spaghetti plots using GPS-derived atmospheric water vapor fields and determined how these fields could help hurricane track locations. We constructed a GPS-IWV distribution map within a 6-h interval. From maps of these GPS-IWV fields, we constructed pairs of geodetic coordinates representing the GPS-IWV maxima for August 24–31, 2017 for Hurricane Harvey and for September 05–12, 2017 for Hurricane Irma at 6-h intervals. Then, we employed the introduced spaghetti plot to provide forecasting information on where a hurricane might occur. This study will contribute to the quantitative assessment of complementary GPS-IWV measures used in hurricane path prediction.

Data and methods

This section describes the methods used to derive GPS-IWV levels reached during major hurricane events that hit the Gulf and East Coasts of the USA in 2017. We provide an analysis of the formal errors approach used for our GPS-derived IWV measurements. We introduce spaghetti plots based on the maximum GPS-derived IWV field for tracking the path moment of hurricanes.

Retrieval of IWV fields from the GPS ZTD

To assess the impacts of major hurricanes in 2017, we selected a dense array of GPS stations (Fig. 1a). We partitioned the GPS stations into two subregional networks based on the areas that the hurricanes impacted. The first subregional network affected by Hurricane Harvey includes 360 GPS stations with data spanning over one month from August 14 to September 14, 2017. These stations are distributed nearly 20 km across southern Texas and Louisiana. The second subregional network is located in the area along the path of Hurricane Irma with data spanning over one month from August 28 to September 28, 2017 and encompasses 562 GPS stations.

Fig. 1
figure 1

GPS-derived IWV measured during Hurricanes Harvey and Irma. a shows the distribution of GPS stations shown as green circles. Red and orange lines denote the actual paths of Hurricanes Harvey and Irma, respectively. The numbers shown on these lines indicate the category of each hurricane. b, c show six-hourly stacked GPS-IWV time series drawn from individual (black point sequence) and averaged (red lines) GPS stations superimposed onto daily GPM/IMERG precipitation (blue lines) time series for the two hurricanes: Hurricanes Harvey (within the red rectangle in a) and Irma (within the orange rectangle in a). Magenta regions denote the periods in which each hurricane shows maximum variations in GPS-IWV and precipitation. d, e show scatter plots presenting the regression between TRMM-derived precipitation and GPS-IWV. The red line denotes the estimated linear regression, and the light green shaded region denotes the 95% confidence interval. The regression was carried out with a three-hour resolution over different periods (August 25–31, 2017 for Hurricane Harvey over Houston and September 6–12, 2017 for Hurricane Irma over Florida). f shows the radiosonde station at Corpus Christi (CRP) and GPS station TXCC. g shows the radiosonde station at Lake Charles (LCH) in Houston and GPS station LAC1. The red line denotes the estimated linear regression. The number of data samples (N), Pearson correlation, slope, intercept, and root-mean-square error is given in the legends. Note that a description of the precipitation, grid satellite, and radiosonde datasets is given in the supporting information (Texts S1-S4)

In this study, we used ZTD Tropo SINEX products obtained from the Nevada Geodetic Laboratory (NGL) using GIPSY/OASIS-II software with a precise point positioning (PPP) strategy (Blewitt et al. 2018). An elevation cutoff angle of 7° was applied for 30-s raw observation. Jet Propulsion Laboratory (JPL) final orbit and clock products are applied. NGL prior wet and dry tropospheric delays were interpolated from the gridded reanalysis of European Center for Medium‐Range Weather Forecast (ECMWF) reanalysis data included through the Vienna Mapping Function 1 (VMF1) along with their respective VMF1 mapping functions (Landskron and Böhm 2018; Böhm et al. 2006). The final NGL ZTD and horizontal gradients were estimated as random walks at a 5-min temporal resolution with uncertainty levels of typically 1.5 to 4.5 mm (Bock et al. 2016). The ionospheric first-order effect was removed by carrier phase linear combination (LC) and code pseudorange linear combination (PC), and second- and third-order effects were modelled using IONEX data with IGRF12 (ftp://cddis.gsfc.nasa.gov/pub/gps/products/ionex). Further details on this GPS observation processing strategy are available at http://geodesy.unr.edu/gps/ngl.acn.

The ZTD can be further split into the zenith hydrostatic delay (ZHD) and zenith wet delay (ZWD). The ZHD can be determined precisely by measuring the surface pressure and surface air temperature at a station. We calculated ZHDs using pressure data available from VMF gridded files (Landskron and Böhm 2018; Böhm et al. 2006) as follows:

$${\text{ZHD}} = \left( {2.2768 \pm 0.0015} \right)\frac{{P_{s} }}{{1 - 2.66 \cdot 10^{ - 3} \cos \left( {2\phi } \right) - 2.8 \cdot 10^{ - 7} H}}$$
(1)

where ZHD is measured in units of mm, \(\phi\) is the station latitude in degrees and H is the station height in meters above the ellipsoid. According to Davis et al. (1985), the uncertainties constant in Eq. (1) (\(\pm 0.0015\)) is calculated by considering all uncertainties of the \({P}_{s}\) and H parameters.

In GPS meteorology, it is common practice to extract the ZWD directly by subtracting the ZHD from the estimated ZTD. Furthermore, we calculated the IWV using the following formula (Bevis et al. 1992):

$${\text{IWV }} = { }\frac{{{\text{ZWD}}}}{\Pi }$$
(2)
$$\Pi = 10^{ - 6} \left( {K_{2 }^{{\prime }} + \frac{{K_{3} }}{{T_{m} }}} \right)R_{v} \rho$$

where \(K_{2 }^{^{\prime}} = 22.1{ } \pm 2.2\,K/{\text{hPa}}\) and \(K_{3} = 3.739 \pm 0.012 \cdot 105\,K^{2} /{\text{hPa}}\) are physical atmospheric reflectivity constants, \(R_{v} = 461\,\left( {{\text{J}}\,{\text{kg/K}}} \right)\) represents the ideal gas constant for water vapor, \(\rho\) is the density of water vapor, and Tm is the atmospheric column mean weighted temperature. Tm is an essential parameter for retrieving IWV estimates from the ZWD of propagating GPS signals. The VMF-observed surface pressure and Tm were derived from the ECMWF Reanalysis-Interim. Then, we applied a bilinear interpolation technique to extract Tm at the station height. Furthermore, we estimated the uncertainty of GPS-IWV before conducting the analysis based on an approach similar to that used by Ning et al. (2013).

$$\sigma_{{{\text{PWV}}}} = \sqrt {\left( {\frac{{\sigma_{{{\text{ZTD}}}} }}{{f\left( {\phi ,h} \right)\Pi }}} \right)^{2} + \left( {\frac{{\sigma_{{{\text{ZHD}}}} }}{{f\left( {\phi ,h} \right)\Pi }}} \right)^{2} + \left( {\frac{{\sigma_{\Pi } }}{{f\left( {\phi ,h} \right)\Pi }}} \right)^{2} }$$
(3)

where \(\sigma_{{{\text{ZTD}}}}\) is the uncertainty of the ZTD, \(\sigma_{{{\text{ZHD}}}}\) is the uncertainty of the ZHD and \(\sigma_{\Pi }\) is the uncertainty of conversion factor \(\Pi\) defined in Eq. (1). \(\sigma_{{{\text{ZTD}}}}\) is derived during GPS processing, which implicitly covers satellite orbit errors, unmodeled ionospheric delays, signal multipath, antenna-related errors, e.g., antenna phase center variations (Ejigu et al. 2019), and mapping function errors. \(\sigma_{{{\text{ZHD}}}}\) is determined from uncertainties in surface pressure and the constant (2.2767 mm/hPa) as follows:

$$\sigma_{{{\text{ZHD}}}} = \sqrt {\left( {\frac{{2.2767\sigma_{Po} }}{{f\left( {\phi ,h} \right)}}} \right)^{2} + \left( {\frac{{0.0015P_{o} }}{{f\left( {\phi ,h} \right)}}} \right)^{2} }$$
(4)

where \(\sigma_{Po}\) is uncertainty in surface pressure. An error of 1 hPa in surface pressure produces an error of approximately 2.3 mm (Davis 1986). Taking into account the uncertainty of surface pressure \(\sigma_{ZHD}\) is very important, especially in extreme weather conditions. Under severe weather conditions, vertical acceleration can reach 1% of gravity acceleration, which can cause errors in the ZHD of approximately 20 mm (Davis et al. 1985). From comparisons between surface synoptic observations and local ground measurements made at 48 globally distributed GPS stations covering an 8-year time period, an averaged root-mean-square difference of 1.7 hPa was found by Wang et al. (2007). Additionally, Heise et al. (2009) and Ning et al. (2013) compared interpolated ground pressure levels obtained from the ECMWF analysis to local ground measurements. Heise et al. (2009) found a standard deviation of 0.9 hPa from more than 60 globally distributed GPS stations using one year of data while Ning et al. (2013) found 0.6 hPa for the single GPS station at the Onsala Space Observatory using over 10 years of data. We adopt the mean \(P_{s}\) uncertainty of \(\sigma_{{P_{s} }} = 0.6\,{\text{hPa}}\). \(\sigma_{\Pi }\) is determined by the uncertainties of \(\rho_{w} ,R_{w} ,T_{m} ,K_{3}\), and \(K_{2}^{{\prime }}\). However, the contributions of uncertainty from \(\rho_{w} ,\) and \(R_{w}\) represent \(<\) 0.1% of the total uncertainty and were considered to be insignificant in this study, depending mainly on \(T_{m} ,K_{3}\), and \(K_{2}^{{\prime }}\) as shown below:

$$\sigma_{\Pi } = 10^{ - 6} R_{v} \rho \sqrt {\left( {\frac{{\sigma_{{K_{3} }} }}{{T_{m} }}} \right)^{2} + \sigma_{{K_{2}^{{\prime }} }}^{2} + \left( {k_{3} \frac{{{\upsigma }_{{T_{m} }} }}{{T_{m}^{2} }}} \right)^{2} }$$
(5)

Wang et al. (2005) found a root-mean-square bias of 1.3 K in Tm based on a comparison between radio sound measurement data and the reanalysis numerical weather predicting (NWP) model NCEP/NCAR (National Center for Environmental Prediction/National Center for Atmospheric Research of USA weather service). Additionally, we adopted the mean uncertainty of Tm, \({\upsigma }_{{T_{m} }} = 1.3\,{\text{K}}.\)

Spaghetti plot

Spaghetti plots are developed to reflect potential hurricane paths. A spaghetti plot represents an “ensemble” forecast where NWP models with slightly different initial boundary conditions are applied to respond to model uncertainty. Hurricane track prediction is an ensemble forecast application that captures the intrinsic spread in hurricane tracks. Knipp (2016) showed that ensemble forecasting has a long history in meteorology and that the space weather community has increasingly adopted the spaghetti approach. To contribute to the advancement of current operational hurricane forecasting efforts, we examined water vapor fields estimated from GPS data. First, we constructed a GPS-IWV map every 6 h. Our GPS-IWV distribution map was constructed using the Generic Mapping Tools (GMT version 6) software “surface” function that uses spline interpolation with an adjustable tension value (Smith and Wessel 1990). To describe the sensitivity of the smoothens of IWV maps, we compared different tension factors of the GMT “surface” tool. Here, we considered tension values of zero with a root-mean-square difference with an actual IWV of 0.29 mm. Then, we extracted the GPS-IWV maximum \(\left( {{\text{IWV}}_{\max } } \right)\) from each map by estimating the value from a location given by the coordinates \(\left( {{\text{lon}}_{0} ,{\text{ lato}}_{0} } \right)\) of a zero crossing of the first derivative of the IWV field, (\({\text{IWV}} \propto f({\text{lon}}, {\text{lat}}, t = 6\,{\text{h}})\)) in six-hour forward intervals if:

$$\left. {\frac{{\partial {\text{IWV}}\left( {{\text{lon}}_{0} ,{\text{lat}}_{0} } \right)}}{{\partial {\text{lon}}}}} \right|_{{t = 6\,{\text{h}}, {\text{priori}}}} = \left. {\frac{{\partial {\text{IWV}}\left( {{\text{lon}}_{0} ,{\text{lat}}_{0} } \right)}}{{\partial {\text{lat}}}}} \right|_{{t = 6\,{\text{h}},{\text{priori}}}} = 0$$
(6)

where \({\text{IWV}}_{\max } = {\text{ IWV}} \left( {{\text{lon}}_{0} ,{\text{lat}}_{0} } \right).\)

Further, we employed the second derivative when we found more than one local maximum at and/or very close to the coordinate of local maxima to determine a single IWV maximum and its coordinates. Based on the maxima of the IWV fields, we plot the spaghetti line (Padilla et al. 2017) to link all of the maxima over several days. In a spaghetti plot, each maximum is represented as a geodetic coordinate pair. To introduce sufficient variation, we add white noise to perturb the positions of the IWV maxima:

$$f \sim f \left( {{\text{lon}}, {\text{lat}}} \right) = r\left( {{\text{lon}}_{0} ,{\text{lat}}_{0} } \right) + \xi$$
(7)

where r represents the geodetic coordinates of \({\text{IWV}}_{\max }\) and ξ is assumed to be normally distributed Gaussian noise with a specified mean and variance.

By adding various instances of normally distributed noise to the geodetic coordinates representing the IWV maxima, we plot different lines linking the IWV maxima at different points in time, and in the next section, we compare this line with the actual cyclone path to assess the tropical cyclone tracking ability of this GPS-based method. Finally, additional lines are obtained by perturbing the positions of the IWV maxima and linking the perturbed positions together. Each line in this spaghetti plot can be visualized as a single strand of spaghetti that is formed by connecting each maximum.

Results

In this section, we examine in detail the performance of GPS-derived IWV fields in monitoring and tracking hurricane paths. First, we provide an analysis of GPS-derived IWV and its agreement with other external precipitation and infrared (IR) brightness temperature data. Second, the tracking of hurricane paths using GPS-derived IWV is characterized, and in this way, we show how IWV may be used to predict the paths of future hurricanes.

Monitoring hurricanes using GPS-derived integrated water vapor fields

Figure 1b, c presents the stacked GPS-IWV time series for the two hurricanes together with the Global Precipitation Measurement (GPM/IMERG) satellite rainfall product time series associated with the two events. Stations located within the perimeters of the corresponding colored rectangles shown in Fig. 1a were used for stacking to understand the IWV characteristics present as each hurricane passed. The effects of these hurricanes were maximal in these areas. We stacked the GPS-IWV time series in each hurricane’s vicinity to obtain a clear understanding of the IWV characteristics therein. In the 2017 hurricane season, Hurricane Harvey greatly affected southern Texas and Louisiana. Hence, we stacked all available GPS station IWV time series around the Houston area that caused extensive damage, particularly in the northeastern Caribbean and Florida Keys. For Hurricane Irma, we similarly stacked all available GPS station IWV time series around Florida for two months. The GPS-IWV time series demonstrate maximum values reached during the passage of these hurricanes over their respective areas. The results demonstrate that the IWV time series are consistent among all of the stations. GPS-IWV exceeded 74 mm on August 27, 2017 when Hurricane Harvey crossed into Houston and on September 11, 2017 when Hurricane Irma crossed into the Florida Keys. We note that a separate IWV increase for Hurricane Harvey before August 12, 2017 corresponds to a lesser-known storm occurring around the Houston area.

We also show the daily GPM/IMERG precipitation attributable to each hurricane from August 19 to September 29, 2017 for comparison (Fig. 1b, c, blue line). Most of the precipitation occurred from August 24 to 29, 2017 for Hurricane Harvey and from September 5 to 12, 2017 for Hurricane Irma. Daily GPM/IMERG precipitation reached 32 cm over Houston as Hurricane Harvey passed over the city. Figure 1d, e shows scatter plots of GPS-IWV against Tropical Rainfall Measuring Mission (TRMM) precipitation for peak periods in the regions affected by the two hurricanes. A correlation analysis was performed on GPS-IWV and precipitation values to determine the strength of their relationship. These parameters exhibited strong coupling and followed a similar footprint associated with the hurricanes, with a Pearson correlation reaching 66%. Further, comparisons drawn between the radiosonde and GPS-IWV solutions for nearby sites show very good agreement. These values are in agreement with other comparisons drawn with similar distances (Torres et al. 2010). Nevertheless, we found RMS errors of between 2.9 and 3.5 mm and a bias of 2.5 mm between the radiosonde and GPS techniques, consistent with previous studies (Deblonde et al. 2005). This is not unexpected due to the spatial distances separating radiosonde and GPS stations and the uncertainties of radiosonde data. The nighttime GPS-IWV observations show better agreement with the radiosonde-derived IWV estimates than the daytime observations with a less than 1 mm bias. A large bias in daytime radiosonde observations is likely caused by solar heating of the radiosonde sensor (Miloshevich et al. 2009).

The construction of water vapor maps using a well-distributed GPS station network helps us produce IWV maps of the Gulf and East Coasts of the USA; these maps vary significantly in time. To quantify how well GPS-IWV data capture a hurricane, we compared the temporal GPS-IWV distribution with the precipitation for days for which a hurricane appeared in the area of influence. Figure 2 depicts daily average GPS-IWV evolution and daily cumulative GPM/IMERG precipitation for Hurricane Harvey from August 25 to September 1, 2017. Between August 17 and 23, 2017, Hurricane Harvey crossed the Caribbean while traveling northwest and then began to intensify and develop into a hurricane on August 24, 2017.

Fig. 2
figure 2

Comparison of GPS-derived IWV and GPM precipitation. Distribution maps of daily average GPS-IWV (first row ad and third row il) and daily cumulative GPM/IMERG precipitation (second row eh and fourth row mp) during Hurricane Harvey between August 25 and September 1, 2017. The IWV and precipitation maps were created using Generic Mapping Tools (GMT) by spherical surface spline gridding. The contour interval of the GPS-IWV data is set to 5 mm. Precipitation values are drawn from the GPM/IMERG mission. Hurricane Harvey’s path is plotted as the cyan line, and the hurricane symbol is plotted in brown

On August 25, the IWV values started to rise at a fairly rapid rate. GPS-IWV considerably increased when Hurricane Harvey reached the Houston area (approximately 26–August 27, 2017). Significant and intense precipitation of between 20 and 30 cm was produced over Houston to the north and east of Hurricane Harvey’s center. In particular, on August 27, 2017, significant precipitation (30–32 cm) was observed over and east of Houston, corresponding to Hurricane Harvey’s path. Later, on August 28, 2017, the maximum IWV intensity rose to 71 mm and shifted to the heart of the hurricane as it migrated back toward the middle of the Texas coastline. GPS-IWV rose even further while centered on Hurricane Harvey’s location and then moved to and from the center of the hurricane. As soon as the center passed over warm water again, Hurricane Harvey was reinvigorated and then moved northeast. Although Hurricane Harvey became slightly reintensified as it returned via western Louisiana, IWV increased on August 29, 2017. We found IWV values exceeding 70 mm (Fig. 2d) in western Louisiana, which is the area where Hurricane Harvey subsequently made landfall again on August 30, 2017. Following this landfall and after the passage of Hurricane Harvey, the IWV magnitude decreased, particularly over the western part of the region, and daily precipitation in Houston fell dramatically (to less than 1 mm). Two weeks later, early in the day on September 6, 2017, another hurricane of similar magnitude, Hurricane Irma, passed over Antigua and Barbuda. Distribution maps of Hurricane Irma’s GPS-IWV evolution and corresponding daily GPM/IMERG precipitation are shown in Figure S1. IWV increased significantly by September 7, 2017 when Hurricane Irma crossed over the Caribbean islands (throughout which some important GPS stations are available) and over the next few days. On September 10–11, 2017, Hurricane Irma moved rapidly northward, reaching just inland along the western coast of Florida. IWV intensified considerably by the time Hurricane Irma reached mainland Florida, particularly near Marco Island in southwestern Florida and surged to 74 mm in the southern part of the state. The daily commutative precipitation values were also greatly elevated (up to 20 cm, see Fig. S2). Further, GPS-IWV and satellite images of infrared brightness temperature from Hurricane Harvey are shown in Fig. 3. Here, one can see that GPS-IWV and IR brightness temperatures capture similar center locations and shapes for Hurricane Harvey that are both consistent. The IR brightness temperatures clearly indicate the region of intense convection and the possible region of the eye wall, indicating relatively strong low-level convection over both the land area and ocean. In general, lower brightness temperatures have moderate to strong updrafts in GPS-IWV maxima. These findings indicate that complementary GPS-IWV products can provide valuable information for tracking extreme meteorological events such as hurricanes.

Fig. 3
figure 3

Temporal evolution of IWV and infrared brightness as influenced by hurricanes. Hourly distribution maps of GPS-IWV (wysiwyg-coloured scale) and infrared (IR) brightness temperature (polar-coloured scale) centered on Hurricane Harvey. Dates are provided at the top, and hours are shown to the left of each panel

Additionally, the 6-hourly estimated spatiotemporal evolution of the IWV for various hurricanes (Hurricanes Harvey and Irma) is depicted in Figures S2–S4. GPS-IWV and IR brightness temperatures for Hurricane Irma are shown in Figure S6. Furthermore, animation videos of 6-hourly GPS-IWV and hurricane tracks (in the supporting materials Animation S1-S2) show that high IWV is associated with Hurricanes Harvey and Irma because every six-hour interval in one week is given as dynamic content in the supporting materials.

Tracking hurricane paths using GPS

With hurricane forecasts, scientists are responsible for providing officials with factual and visual data that will inform their decisions regarding whether (and when) they should call for compulsory evacuations and how best to allocate disaster preparedness resources. Hurricane track forecasting efforts continue to improve. Since 2005, National Hurricane Center (NHC) three-day track forecasts have been 33% more accurate than they were in the late 1980s (McAdie and Lawrence 2000). Nonetheless, forecasts of TC intensity remain difficult to successfully achieve (Broad et al. 2007) due to a lack of frequent sampling within the inner cores of such TCs (Rogers et al. 2013). In addition, many satellite missions are unable to examine heavy rainfall events that accompany the passage of hurricanes’ inner core regions.

For both hurricanes studied here, we found GPS-IWV maxima exceeding 57 mm. The GPS-IWV maxima were found to reach their highest magnitudes over our track predictions, which are similar to the paths of the final tracks from the NHC model. The NHC model uses a cone of uncertainty to determine the general region that may be influenced by a hurricane’s path. All of our estimated GPS-IWV maxima were found to fall within this cone of influence as provided by the NHC model. Based on the maxima of the IWV fields, we plot spaghetti lines for the path forecasts of Hurricanes Harvey and Irma. In this case, our input values for the spaghetti plot are geodetic coordinates representing the GPS-IWV maxima. Implicit to our geodetic coordinates are the time, pressure, and temperature. To produce different initial inputs for the spaghetti plot, we add Gaussian noise comprising 200 initial seeds with a standard error of 0.1° and a zero mean to our geodetic coordinates representing the maxima. This procedure created 200 different geodetic coordinates that reflect the virtual maximum IWV values along the hurricane trajectory, thereby creating 200 different spaghetti lines. More specifically, the spaghetti plots are obtained from the location of the IWV maximum at each step (every 6 h), and then a single spaghetti line is drawn to link all of these maxima over several days, after which additional lines are obtained by perturbing the positions of the maxima and linking the perturbed position together. Figure 4a displays spaghetti lines representing 200 different perturbed coordinates where the average isolines envelop the maximum and are plausible.

Fig. 4
figure 4

Spaghetti line plots derived from our GPS-estimated IWV levels. The top panels show spaghetti line plots (blue lines) transecting the GPS-IWV maxima for Hurricanes Harvey (left) and Irma (right). Coloured circles show the GPS-IWV maxima estimated for every six-hour interval. Light magenta polygons show the best tracks from the NHC model based on a poststorm analysis of all available data presented at six-hourly intervals. Middle panels show GPS-IWV crossing profiles along a straight line and sampled every 1 km over a grid for August 27, 2017 at 00 UTC (bottom left) and for September 10, 2017 at 18 UTC (bottom right). Each colour corresponding to the colour shown inside the circles in the bottom panels and the GPS-IWV profile is calculated within 45° of the path (colour in the circles)

For Hurricane Harvey, we extracted the different members of the GPS-IWV maxima from gridded GPS-IWV distribution maps starting at 06 UTC on August 25, 2017 (i.e., an intensified GPS-IWV location approximately 200 km away from the location of Hurricane Harvey’s development) until 18 UTC on August 30, 2017. Similarly, we extracted the different members of the GPS-IWV maxima for Hurricane Irma starting at 00 UTC on September 6, 2017 and continued this process every 6 h until 18 UTC on September 11, 2017. Early on September 6, 2017, before Hurricane Irma reached the Caribbean Islands, the GPS-IWV maxima showed a clear pattern of moving westward. Over the next few days, as Hurricane Irma continued moving westward, skirting the Caribbean Islands, we observed an IWV exceeding 65 mm along its path within a six-hour interval prior to the hurricane arriving in the area; this observed value is much higher than the typical IWV level of 40 mm. Late on September 9, 2017, Irma made its anticipated turn northwest into the Straits of Florida. We found a maximum GPS-IWV value in excess of 70 mm in the six-hour interval prior to the hurricane crossing the Florida Straits. We then found a significant maximum IWV of greater than 75 mm along Florida’s west coast late on September 10, 2017. According to the NHC’s TC report, maximum sustained winds were measured in western Florida from September 10–11, 2017.

After the initial sources of the GPS-IWV maxima (the colored circles in Fig. 4a) were identified, we added Gaussian noise and plotted each of the trajectories (differently colored lines in Fig. 4a). The results show that each spaghetti line was projected near the actual path. More importantly, the projected GPS-IWV maxima along the spaghetti plots indicate the most likely path along which the center of a storm will travel next.

To further elucidate the potential to predict hurricane tracks solely from GPS-IWV estimates, we selected two locations for both hurricanes subject to very complex conditions (see Fig. 4b, c). The potential tracks from the NHC for the prestorm period are not consistent. For instance, the NHC predicted that Hurricane Harvey would move eastward toward Louisiana at 10 UTC on August 25, 2017; however, by 04 UTC on August 26, 2017, Hurricane Harvey had moved further west across Texas Hill County. After various rotations, the update on August 27, 2017 moved the model back to Louisiana, predicting that the storm would move offshore. Here, we found the median GPS-IWV content for both hurricanes to be significantly higher in the storm center. Additionally, both hurricanes showed maximum GPS-IWV levels to be roughly parallel to the storm track axis and similar to the final track (Fig. 4b).

Discussion

A detailed analysis of the GPS-IWV time series was carried out for Hurricanes Harvey and Irma, which devastated the Gulf and East Coasts of the USA in the 2017 hurricane season. We observed a sudden and significant increase and decrease in GPS-IWV during and after the storms, respectively. We also observed enhanced GPS-IWV, particularly for landfall locations. The observed variation in the time series of IWV content can also be correlated with variations in the precipitation value constructed from GPM/IMERG satellite observations coinciding with the passage of the hurricane storm front. According to Sapucci et al. (2016), a sharp increase in GPS-IWV occurs before rainfall, and most of the maximum rainfall occurs near the GPS-IWV peak. As is evident from the comparison of GPS-IWV and precipitation values, stations around Houston and Corpus Christi (Hurricane Harvey’s first landfall) showed a sharp peak in GPS-IWV. In Houston, from August 25–31, 2017, GPS-IWV surged to 74 mm, far exceeding the typical value of 40 mm. An increase in GPS-IWV by an average of 28 mm (with a range of approximately 46–74 mm) is commonly observed in the general vicinity of the area crossed by a storm for four to five days. Essentially, all of the significant increases in GPS-IWV observed occurred due to the storm. Some variability in water vapor was found due to the local weather processes prestorm and poststorm, but nothing near the magnitude of the hurricane.

As the main advantage of GPS-IWV measures obtained from a dense continuous GPS network, it is possible to produce IWV distribution maps in close proximity to real-time or near-real-time. We observed a maximum IWV in distribution maps of IWV and were able to predict hurricane tracks. Our findings verify that temporal variation in GPS-IWV is closely related to the hurricane route, possibly tracking and monitoring hurricane operations typically at least several hours prior to the storm’s arrival. The capacity to base the maximum IWV field as an input for a spaghetti plot underlines the ability for a GPS to predict hurricane paths at least six hours before the arrival of a hurricane (see Fig. 4). Typically, the strongest convergence of water vapor occurs one hour before severe precipitation (Adams et al. 2013). This result is a main finding of the tropospheric GPS product, which offers complementary data for predicting a hurricane’s path. According to Trenberth et al. (2007), the atmospheric water budget of a hurricane is governed by the inflow of water vapor largely within the lowest 1 km area of the storm and within 100 km of the center of the storm. Moisture fluxes enhance the amount of heavy rainfall and hence provide latent heat that drives the hurricane. Similarly, GPS-derived IWV from a dense network of ground stations provides a high-resolution water vapor field. GPS-IWV represents the mean IWV in the atmosphere over an inverse-cone shape with a base radius of approximately 25 km at an altitude of 3 km with a 7-degree elevation cutoff angle.

When water vapor condenses into clouds or rain, it releases latent heat; thus, GPS measurements of water vapor can contribute to the forecasting of hurricane intensity. Kamineni et al. (2003) showed a reduction in tracking errors of 100 km and intensity errors of 20–25% for hurricanes when moisture profiles are assimilated into an operational NWP model. Iwabuchi et al. (2009) also showed significant improvements in forecasting severe storm intensity when GPS-derived precipitable water is assimilated into the NWP model, even for a single GPS station. More importantly, near-real-time GPS meteorology (hourly updates with 15-min estimates) is widely used in forecasting NWP models (Iwabuchi et al. 2006; Mahfouf et al. 2015). Generally, the contributions of GPS-IWV assimilated with NWP constitute approximately 1–2% of all data, such as data from radiosoundings, but they have a significant impact on the predicted precipitation pattern (Mahfouf et al. 2015).

Conclusions

This study examines the impact of the GPS-IWV during the 2017 hurricane season. Two different intense hurricane events in 2017, namely Hurricanes Harvey and Irma, were considered. We employed more than 360 (for Hurricane Harvey) and 562 (for Hurricane Irma) active GPS stations forming a dense ground network with interstation distances of approximately 20 km for a 1-month period. The GPS zenith total delay (ZTD) obtained from the NGL was converted to IWV using the observed surface pressure and mean atmospheric water vapor column temperature obtained from VMF, which in turn uses ECMWF Reanalysis-Interim. GPS-IWV values typically increased from 50 to 70 mm during the storm period. Virtually all of the significant increases in IWV are attributable to the storm. We investigated the spatial distribution of water vapor content during the storm’s passage and observed that the increase in IWV content occurred several hours before the storm arrived.

We analyzed the maximum possible GPS-IWV distribution map to examine whether PS can play an important role in tracking spatiotemporal variability in water vapor fields during hurricane storm events. In our IWV distribution maps, we observed large proportions of IWV for both hurricanes and were able to determine where the storm would travel next. The results confirm that temporal changes in GPS-IWV are strongly linked to a hurricane’s path and may facilitate the tracking and monitoring of hurricane activities at least several hours before the arrival of a storm. These findings are in agreement with Tang et al. (2016) in showing that TC development projections could be improved through the monitoring of water vapor in the lower and middle atmosphere. During Hurricanes Harvey and Irma, IWV values ranged mostly between 60 and 80 mm. The capacity to employ the maxima of the IWV field as inputs for spaghetti plots underscores the potential for GPS to predict hurricane paths at least six hours before a storm’s arrival (Fig. 4). The provision additional information for predicting the paths of hurricanes is a key outcome of this GPS tropospheric product.

Overall, our findings demonstrate that the use of GPS-IWV could advance the preliminary monitoring of hurricane activity in regions such as the Gulf and East Coasts of the USA and the Gulf of Mexico, where dense arrays of GPS stations are available. The monitoring results of this study were evaluated based on the locations of the areas affected by the studied hurricanes and the timing of such effects. The very high temporal resolution of GPS satellite data, typically of 5–15 min, enables a sequence of IWV distribution maps to support a better estimate of a hurricane’s path.