**3. Data assimilation of SPD at storm scale**

#### **3.1. Data assimilation system**

Data assimilation (DA) techniques optimize differences between first guess fields provided by numerical models and observations. One advanced DA theory is a variational method (Sasaki [33, 34]), and the NHM-4DVAR is a variational DA system for predicting thunderstorms (Kawabata et al. [35, 36]). The NHM-4DVAR is based on the Japan Meteorological Agency Nonhydrostatic Model (JMANHM; Saito et al. [37]), which includes three-ice bulk cloud microphysics, as the forward model. In the first version of NHM-4DVAR, the adjoint model considered only dry dynamics and advection of water vapor. Kawabata et al. [36] implemented an additional warm rain process for assimilating radar reflectivity. The horizontal resolution of NHM-4DVAR is 2 km. This 2-km spacing is known as "storm scale" because in this range it is possible to explicitly represent cumulonimbus and thunderstorms by numerical models. The control variables are the three wind components, potential temperature, surface pressure, nonhydrostatic pressure, total water (water vapor and cloud water), the relative mixing ratio of rainwater, and the pseudo-relative humidity (only for lateral boundary conditions).

This method optimizes the difference between observations and numerical models using their error statistics. To measure the difference, we need observation operators to convert the numerical model fields to observation spaces. For SPD DA, Eq. (2) is used as the operator and provides delays of GNSS radio waves from GNSS antennas to the model top boundary as

$$
\Delta \mathbf{L} = \int\_{\text{Model\\_top}}^{\text{Model\\_top}} (n - 1) d\mathbf{s}.\tag{7}
$$

Here, ∆L denotes the atmospheric delay (m), *n* denotes the refractive index along the path in a grid cell, and *ds* is the path length (m) in each model grid cell. Because the top boundary of numerical models is limited to a certain height level (usually 20–40 km), we need an assumption to add a delay from the top boundary to a GNSS satellite.

The authors assumed that the delay decreased to 1/e every 10 km from the top boundary to 200 km height (Bean and Thayer [38]) and that above 200 km the amount of delay was zero. In addition to this assumption, a straight-line assumption was adopted; the bending effect was thus eliminated. To account for these assumptions, a relatively large observational error was implemented. To calculate SPD in the model, first, the linear path is determined from the receiver to the GNSS satellite. Then, the middle point of the path within the model grid cells are set down. Delays are interpolated and averaged from eight points surrounding the middle point, with each weight set according to the distance. Finally, the slant path delays are calculated by integrating each delay in model grid cells from the receiver to the top boundary of the model.

The observations over elevation angles of 5° or more were assimilated in their study. The operator works on the World Geodetic System 1984 (WGS84; National Imagery and Mapping Agency 1997) [39].

#### **3.2. Assimilation experiment**

The new method demonstrated the capability to capture a strong PWV gradient associated with the parent storm of the F3 tornado that struck Tsukuba City in Ibaraki Prefecture, Japan, on May 6, 2012, with a numerical model simulation at a super-high resolution of 50 m (NHM50m; **Figure 6**). An area of large PWV contrast centered on strong precipitation (**Figure 6a**) implies a strong upward wind in front of and a strong downdraft behind the parent storm. The PWV contrast toward the tornado is also seen in **Figure 6b**, whereas, no such PWV gradient is seen at all in **Figure 6c**.

Data assimilation (DA) techniques optimize differences between first guess fields provided by numerical models and observations. One advanced DA theory is a variational method (Sasaki [33, 34]), and the NHM-4DVAR is a variational DA system for predicting thunderstorms (Kawabata et al. [35, 36]). The NHM-4DVAR is based on the Japan Meteorological Agency Nonhydrostatic Model (JMANHM; Saito et al. [37]), which includes three-ice bulk cloud microphysics, as the forward model. In the first version of NHM-4DVAR, the adjoint model considered only dry dynamics and advection of water vapor. Kawabata et al. [36] implemented an additional warm rain process for assimilating radar reflectivity. The horizontal resolution of NHM-4DVAR is 2 km. This 2-km spacing is known as "storm scale" because in this range it is possible to explicitly represent cumulonimbus and thunderstorms by numerical models. The control variables are the three wind components, potential temperature, surface pressure, nonhydrostatic pressure, total water (water vapor and cloud water), the relative mixing ratio

of rainwater, and the pseudo-relative humidity (only for lateral boundary conditions).

∫

tion to add a delay from the top boundary to a GNSS satellite.

*Model*\_*surface*

Here, ∆L denotes the atmospheric delay (m), *n* denotes the refractive index along the path in a grid cell, and *ds* is the path length (m) in each model grid cell. Because the top boundary of numerical models is limited to a certain height level (usually 20–40 km), we need an assump-

The authors assumed that the delay decreased to 1/e every 10 km from the top boundary to 200 km height (Bean and Thayer [38]) and that above 200 km the amount of delay was zero. In addition to this assumption, a straight-line assumption was adopted; the bending effect was thus eliminated. To account for these assumptions, a relatively large observational error was implemented. To calculate SPD in the model, first, the linear path is determined from the receiver to the GNSS satellite. Then, the middle point of the path within the model grid cells are set down. Delays are interpolated and averaged from eight points surrounding the middle point,

*Model*\_*top* (*n* − 1)*ds*. (7)

This method optimizes the difference between observations and numerical models using their error statistics. To measure the difference, we need observation operators to convert the numerical model fields to observation spaces. For SPD DA, Eq. (2) is used as the operator and provides delays of GNSS radio waves from GNSS antennas to the model top boundary as

**3. Data assimilation of SPD at storm scale**

**3.1. Data assimilation system**

152 Multifunctional Operation and Application of GPS

∆L =

#### *3.2.1. Single set of SPD observations at a single site*

First, three experiments were performed using NHM-4DVAR with a 2-km horizontal grid spacing in which three SPD observations (SO\_SPD), one ZTD observation (SO\_ZTD), or one PWV observation (SO\_PWV) from a single observation site were assimilated. Note that the ZTD and PWV observations were derived from the SPD observations originally. These experiments were conducted to confirm the effects of SPD assimilation on a single analysis and to examine the differences between SPD, ZTD, and PWV assimilations. The assimilation window was set at 10 min, and the observations were assimilated every 5 min (at 0, 5, and 10 min in the assimilation window). The observational data set was chosen by considering the horizontal distributions, elevation angles, and the first guess field from an experimental data set.

The propagation paths of GNSS signals from three satellites to a receiver in the model atmosphere in both the horizontal (**Figure 7a**) and the vertical plane (projected from the south; **Figure 7b**) are illustrated. Path I with the smallest elevation angle is also the longest, while path II at a near 90° elevation angle is the shortest. The large amount of delays is illustrated mostly in low altitudes area by warm colors.

The distributions of the increments (analysis minus first guess) of PWV in SO\_SPD (**Figure 8a**) are different from that in SO\_ZTD (**Figure 8b**) at the end of the assimilation window. In SO\_ZTD, the increment distribution (**Figure 8b**) is axisymmetric and elliptical mostly (i.e., isotropic). Although 4D-Var captures analysis increments shaped flow-dependent (anisotropic) in general, it is not

**Figure 7.** Propagation paths of GNSS signals from satellites (I, II, and III) to a receiver in the simulated atmosphere, (a) viewed in the horizontal plane and (b) the vertical plane. Colors illustrate the values of delay in each model grid cell (each cell is shown by one pixel). The open circle in (a) displays the GNSS receiver (the observational site). After Kawabata et al. [4].

0.0 1.0 2.0 3.0

0 10km

(c)

(a) (b)

Qv g/kg

top at the upper right corner. Modified Figure 11 of Kawabata et al. [4].

(a) (b)

domain, topography and GNSS stations. After Kawabata et al. [4].

0.0 1.0 2.0 3.0

Applications of GNSS Slant Path Delay Data on Meteorology at Storm Scales

0 10km

0 10km

http://dx.doi.org/10.5772/intechopen.75101

155

Qv g/kg

0.0 1.0 2.0 3.0

**Figure 9.** Vertical cross sections along path III. Analysis increments of the mixing ratio of water vapor (g kg−1). (a) SO\_SPD, (b) SO\_ZTD and (c) SO\_PWV. The observation site is at the lower left corner of each panel, and path III leaves the model

**Figure 10.** (a) Propagation paths of radio waves from GNSS satellites to receivers in the model atmosphere, viewed in the horizontal plane. Actual SPD observations obtained over Okinawa Island at 1130 JST. (b) the assimilation and forecast

Qv g/kg

**Figure 8.** Analysis increments of precipitable water vapor (mm) at the end of the assimilation window for (a) SO\_STD and (b) SO\_ZTD. After Kawabata et al. [4].

long enough for the 10-min assimilation window to capture such flow dependency. As a result, the shape of the analysis increment seen in SO\_ZTD is isotropic mostly. In contrast, the inhomogeneous distribution of the analysis increment in SO\_SPD (**Figure 8a**) is given by the various distribution of the slant paths. Moreover, the maximum value of the increments is much larger in SO\_SPD (10 mm; **Figure 8a**) than that in SO\_ZTD (6 mm, **Figure 8b**). The increment distributions in SO\_PWV and SO\_ZTD were close to each other (not shown).

Seeing vertical cross sections of the mixing ratio of water vapor (Qv) along path III (**Figure 9**), the distributions of analysis increment in SO\_ZTD and SO\_PWV are close, and the magnitudes of the increment are quite similar: The distributions of Qv reach at high altitudes from 1 to 5 km in vertical and extend to 5–8 km away in horizontal. The increment along path III in SO\_SPD is distributed along a distance of 15 km in horizontal and over 8 km in height; furthermore, the magnitudes of the increment are much larger than others, at low altitude especially (around 3 km), because all of the slant paths are within a narrow area of the lower troposphere above the observation site. The total weights in the cost function show good agreement among the assimilations of SO\_SPD, SO\_ZTD, and SO\_PWV, but their effects are seen at different places.

#### *3.2.2. Actual observational data assimilation*

Next, an actual **observation assimilation** and forecast experiments were conducted using NHM-4DVAR with the assimilation window of 30 min. Three-hour forecasts were performed after the 4D-Var analysis with multiple observations. The assimilation was started at 1100 JST on August 19, and the forecast was performed from 1100 to 1400 JST. The SPD, ZTD, or PWV data were assimilated every 10 min (hereafter, A\_SPD, A\_ZTD, and A\_PWV, respectively). Only the observations listed above, and no other kinds, were assimilated. No cumulus parameterization was applied in the experiments. The experimental domain and the distribution of GNSS receivers are illustrated in **Figure 10**. Unlike the PWV and ZTD data, it can be said that SPDs in the model atmosphere (**Figure 10a**) contain a great deal of horizontal information.

Applications of GNSS Slant Path Delay Data on Meteorology at Storm Scales http://dx.doi.org/10.5772/intechopen.75101 155

**Figure 9.** Vertical cross sections along path III. Analysis increments of the mixing ratio of water vapor (g kg−1). (a) SO\_SPD, (b) SO\_ZTD and (c) SO\_PWV. The observation site is at the lower left corner of each panel, and path III leaves the model top at the upper right corner. Modified Figure 11 of Kawabata et al. [4].

long enough for the 10-min assimilation window to capture such flow dependency. As a result, the shape of the analysis increment seen in SO\_ZTD is isotropic mostly. In contrast, the inhomogeneous distribution of the analysis increment in SO\_SPD (**Figure 8a**) is given by the various distribution of the slant paths. Moreover, the maximum value of the increments is much larger in SO\_SPD (10 mm; **Figure 8a**) than that in SO\_ZTD (6 mm, **Figure 8b**). The increment distributions

**Figure 8.** Analysis increments of precipitable water vapor (mm) at the end of the assimilation window for (a) SO\_STD

Seeing vertical cross sections of the mixing ratio of water vapor (Qv) along path III (**Figure 9**), the distributions of analysis increment in SO\_ZTD and SO\_PWV are close, and the magnitudes of the increment are quite similar: The distributions of Qv reach at high altitudes from 1 to 5 km in vertical and extend to 5–8 km away in horizontal. The increment along path III in SO\_SPD is distributed along a distance of 15 km in horizontal and over 8 km in height; furthermore, the magnitudes of the increment are much larger than others, at low altitude especially (around 3 km), because all of the slant paths are within a narrow area of the lower troposphere above the observation site. The total weights in the cost function show good agreement among the assimilations of SO\_SPD, SO\_ZTD, and SO\_PWV, but their effects are seen at different places.

Next, an actual **observation assimilation** and forecast experiments were conducted using NHM-4DVAR with the assimilation window of 30 min. Three-hour forecasts were performed after the 4D-Var analysis with multiple observations. The assimilation was started at 1100 JST on August 19, and the forecast was performed from 1100 to 1400 JST. The SPD, ZTD, or PWV data were assimilated every 10 min (hereafter, A\_SPD, A\_ZTD, and A\_PWV, respectively). Only the observations listed above, and no other kinds, were assimilated. No cumulus parameterization was applied in the experiments. The experimental domain and the distribution of GNSS receivers are illustrated in **Figure 10**. Unlike the PWV and ZTD data, it can be said that SPDs in the model atmosphere (**Figure 10a**) contain a great deal of horizontal information.

in SO\_PWV and SO\_ZTD were close to each other (not shown).

(a) (b)

*3.2.2. Actual observational data assimilation*

and (b) SO\_ZTD. After Kawabata et al. [4].

154 Multifunctional Operation and Application of GPS

**Figure 10.** (a) Propagation paths of radio waves from GNSS satellites to receivers in the model atmosphere, viewed in the horizontal plane. Actual SPD observations obtained over Okinawa Island at 1130 JST. (b) the assimilation and forecast domain, topography and GNSS stations. After Kawabata et al. [4].

The horizontal distributions of rainfall in A\_ZTD and A\_PWV (**Figure 11c** and **d**) were close to each other. In A\_SPD (**Figure 11b**), the initiation of intense rainfall enhanced over southern Okinawa at 1200 JST, and the maximum intensity reached 47 mm h−1 (not shown). By 1300 JST, a rainband had formed over the island, and the distribution and intensity resembled the radar observation (**Figure 11a**). Therefore, it can be concluded that the assimilation of SPD data showed improvement of the rainfall forecast with respect to both timing and intensity compared with the assimilation of PWV and ZTD data. This improvement was obtained through modifications in the atmospheric profile in the simulation after the SPD assimilation (not shown).

**4. Summary**

toring of cumulus convection.

over Okinawa Island.

taking advantages.

should speed the reality of new usage of GNSS-derived SPD.

GNSS analysis is performed with a number of models for solid earth tide, ocean tide, earth rotation, and atmospheric delay. However, conventional GNSS analysis inevitably has errors caused by local-scale variations in the atmosphere, even if estimates of the gradient parameters are applied. These errors affect positioning as well as atmospheric delays. Spatial correlation in screened post-fit phase residuals (i.e., the higher order inhomogeneity part) must be caused by local-scale variations in the atmosphere, which are not resolved by ordinary ZTD and gradient parameter analysis. Thus, GNSS postfit residuals can be used to detect the inhomogeneous distribution of water vapor on a local scale of less than 10 km caused by weather disturbances like convective thunderstorms. In this chapter, we introduced two new indices and one experimental method to retrieve a several-kilometer-scale PWV distribution. These approaches may be a preliminary way to serve operational moni-

Applications of GNSS Slant Path Delay Data on Meteorology at Storm Scales

http://dx.doi.org/10.5772/intechopen.75101

157

Currently, multiple GNSSs (e.g., the Russian GLONASS, EU's GALILEO, China's BeiDou, and Japan's Quasi-Zenith Satellite System (QZSS)) can serve to more precisely resolve localscale water vapor variation. In addition, the advancement of real-time GNSS analysis has been progressing rapidly. Real-time orbit and clock corrections have offered officially by the International GNSS Service (IGS) since April 2013. Moreover, a Multi-GNSS orbit and clock estimator called MADOCA (Multi-GNSS Advanced Demonstration tool for Orbit and Clock Analysis) has been developed by the Japan Aerospace Exploration Agency (JAXA) and they started to provide real-time ephemerides via the Internet (https://ssl.tksc.jaxa.jp/madoca/public/public\_index\_en.html) in September 2015. These rapid advancements in GNSS technology

Another application of SPD data at a local scale is data assimilation (DA). As shown in this chapter, a promising avenue for DA applications is to improve the initial conditions of a numerical weather model in high-resolution simulations. NHM-4DVAR with a 2-km horizontal grid spacing was used for assimilating GNSS slant path delay data. Compared with simulations after assimilating PWV and ZTD data, the assimilation of SPD showed inhomogeneous increments. Moreover, the analysis increments in the assimilation of SPD at low altitudes were larger than in the others. Conducting actual observation assimilations, the SPD assimilation clearly improved the forecast of both the timing and intensity of the rainband. This improved forecast was given through the decreased atmospheric stability

The SPD data include information on temperature, dry atmospheric pressure, and water vapor pressure. In addition, the data contain both horizontal and vertical information on those atmospheric parameters. Because atmospheric inhomogeneity is greatest in the lower troposphere, assimilating SPD data is a promising way to improve forecast of a rainband. A GNSS signal along a path with a 15° elevation angle propagates about 11 km in horizontal and travels 3 km in vertical. Therefore, it is concluded that data assimilation systems at storm scales like used in this study are necessary to assimilate SPD data for

**Figure 11.** Horizontal distribution of 1-h accumulated rainfall amount between 1300 and 1400 (FT = 02) JST. (a) Radar observation, (b) A\_SPD, (c) A\_ZTD and (d) A\_PWV. Modified Figures 2 and 15 of Kawabata et al. [4].
