Multiple Doppler Radar Analysis for Retrieving the Three-Dimensional Wind Field Within Thunderstorms

Multiple Doppler radar analysis has been widely used to retrieve three-dimensional wind fields within thunderstorms and meso-scale convective systems (MCS) since the late 1960s. A number of countries have constructed dense operational radar networks, such as the Operational Programme for the Exchange of Weather Radar Information (OPERA; Kock et al., 2000), to monitor and forecast severe weather in metropolitan regions. Multiple Doppler radar analysis using such operational radar networks improves 1) understanding of the physical mechanisms behind heavy rainfall and severe wind, 2) detection and forecasting of hazardous weather phenomena, and 3) planning for mitigation of human and socioeconomic losses in metropolitan regions.


Introduction
Multiple Doppler radar analysis has been widely used to retrieve three-dimensional wind fields within thunderstorms and meso-scale convective systems (MCS) since the late 1960s. A number of countries have constructed dense operational radar networks, such as the Operational Programme for the Exchange of Weather Radar Information (OPERA; Köck et al., 2000), to monitor and forecast severe weather in metropolitan regions. Multiple Doppler radar analysis using such operational radar networks improves 1) understanding of the physical mechanisms behind heavy rainfall and severe wind, 2) detection and forecasting of hazardous weather phenomena, and 3) planning for mitigation of human and socioeconomic losses in metropolitan regions.
Early single Doppler radar measurements provided a basic understanding of storm morphologies and their three-dimensional structures, including concepts for single-cell, multicell, and supercell storms (Browning, 1964(Browning, , 1965. Single Doppler radar observations can only provide information on the radial component of wind (i.e., velocity which is directed toward or away from the radar), rather than the full three-dimensional structure. Armijo (1969) formulated a method that allowed the deduction of the three-dimensional wind structure by combining the data from several Doppler radars. Improvements in this multiple Doppler radar analysis method were reported during the 1970s and 1980s, including the design of optimal radar networks (Ray et al., 1979(Ray et al., , 1983, the development of alternative analysis schemes for solving the mass continuity equation (Ray et al., 1980), and the introduction of floating boundary conditions (Chong & Testud, 1983). These improvements were primarily motivated by the need to overcome errors in the estimation of vertical velocity using upward integration of the mass continuity equation (Doviak et al., 1976). Errors in estimates of vertical velocity tend to amplify during such upward integration because of the stratification of density in the atmosphere (Doviak et al., 1976;Ray et al., 1980). Theoretical demonstrations indicate that downward integration of the mass continuity equation could yield more accurate estimates of vertical velocity than those that can be obtained from upward integration (Ray et al., 1980). Many subsequent studies have therefore applied downward integration schemes to determine the three-dimensional structure of winds within severe storms (Kessinger et al., 1987;Biggerstaff & Houze, 1991;Dowell & Bluestein, 1997). Downward integration requires observations at the storm top; however, typical radar scan geometries are configured for operational monitoring of low-level precipitation and severe low-level wind phenomena (e.g., downburst and tornado). Such configurations do not often provide detailed observations of the storm top, where the vertical velocities may be significantly different from zero during storm development (Mewes & Shapiro, 2002). Several studies have used low pass filters, such as the Leise filter (Leise, 1981), to apply high wavenumber adjustments to the lower boundary conditions prior to upward integration of the mass continuity equation. Such adjustments reduce noise in estimates of upper-level winds (Parsons & Kropfli, 1990, Wakimoto et al., 2003. One alternative is to apply the anelastic mass conservation equation as a weak constraint (Gao et al., 1999, hereafter G99;Gao et al., 2004, hereafter G04). This method is based on a three-dimensional variational approach, and removes the need to explicitly integrate the anelastic continuity equation. This prevents the accumulation of severe errors in the vertical velocity and ensures that uncertainties in the upper and lower boundary conditions do not propagate vertically. Furthermore, multiple Doppler radar analysis is usually performed in a Cartesian coordinate system; Doppler velocity data are often interpolated into this Cartesian coordinate system using a Cressman filter (Cressman, 1959). The scheme introduced by G99 bypasses this step by allowing reverse linear interpolation (from the regularly spaced Cartesian grid to the irregularly spaced radar observation points) during calculation of the cost function. This reverse interpolation procedure preserves the radial nature of radar observations; however, as noted above, operational radar networks are often incapable of providing dense observations, especially at upper levels. In such cases, the G99 scheme requires accurate background information, such as sounding data, to fill in the datavoid regions between successive elevation angles. It is frequently difficult to obtain accurate background information in these cases, due in part to the coarse temporal resolution of sounding data. If spatially continuous Doppler velocity data could be obtained in Cartesian coordinates through the careful use of Cressman filters, accurate vertical velocity could be obtained from operational radar scans without the need for additional information. Otherwise, additional information regarding upper-level winds is necessary to reduce errors in estimates of vertical velocity near the storm top.
This chapter presents a simplified version of the G99 scheme that applies a threedimensional variational approach on a regular Cartesian grid. The accuracy of calculated winds and the dependence of this accuracy on the density of upper-level radar observations are investigated using a set of idealized data sampled from a simulated supercell storm. A detailed description of the structure of this simulated supercell has been provided by Shimizu et al. (2008). The objective of this chapter is to propose an optimal method for analyzing severe thunderstorms using typical configurations of current operational radar data (less than 20 Plan Position Indicators, or PPIs, within 5-6 minutes).

Analysis method and variational scheme
This section briefly reviews the variational scheme for multiple Doppler radar analysis; a detailed description has been provided by Gao et al. (1999). The variational technique minimizes a cost function (J), which is defined as the sum of squared errors due to discrepancies between observations and analyses and additional constraint terms: Vr is the analysed radial velocity on a specified Cartesian grid, where i, j, and k indicate spatial location in the x, y, and z directions, respectively, and m indicates the m th radar in the network. N is the total number of observations, which is equal to the product of the number of grid cells and the number of radars in the network. C o is the reciprocal of the mean squared error in the observations. Vr_obs is the radial velocity interpolated to the regular Cartesian grid. The cost function is evaluated at each grid point in the Cartesian coordinates, rather than in spherical coordinates.
Each constraint is weighted by a factor that accounts for its respective proportion of the reciprocal of the mean squared error. As noted by G99, it is usually difficult to obtain appropriate values for the weighting coefficients. In particular, the value of the weighting coefficient for the anelastic mass conservation constraint plays an important role in determining the vertical wind component. This study uses the coefficient value used by G99, although G04 introduced a more objective method for estimating this coefficient.
The variational method uses the derivative of J with respect to the analysis variables to obtain an optimal solution. The gradient of the cost function is derived with respect to the control variables, namely the two horizontal wind components (u, v) and the vertical wind component (w). The form of the gradient used here differs slightly from that used by G99 because the form of the observational constraint differs. The gradient of the observational constraint with respect to u is given by where r is the distance between the radar and the grid point and x is the component of r in the x direction.
After the cost function is evaluated and its gradients are obtained, a quasi-Newton-type optimization scheme is used to update the control variables. This analysis uses a limitedmemory Broyden-Fletcher-Goldfarb-Shanno (L-BFGS) method (Liu & Nocedal, 1989). For most meteorological applications, the L-BFGS method is more efficient than the conjugate gradient method (Navon & Legler, 1987). L-BFGS uses an approximation of the secondorder derivative, so that an iteration of the L-BFGS method typically requires less computation than an iteration of the CG method. L-BFGS is therefore a better choice for optimizing a computationally expensive cost function.

Observational system simulation experiment and model description
The performance of the variational technique is evaluated in the context of an observational system simulation experiment (OSSE). This OSSE is conducted using numerical simulations of a supercell thunderstorm observed near Tokyo on 24 May 2000 (Shimizu et al., 2008). The numerical simulations are generated using the Cloud-Resolving Storm Simulator (CReSS; Tsuboki & Sakakibara, 2002). CReSS is a three-dimensional nonhydrostatic model. The microphysical and other parameterization schemes used in CReSS have been described in detail by Tsuboki & Sakakibara (2002). The model grid comprises 300 x 300 x 70 grid points, with grid intervals of 1 km in the horizontal directions. The vertical grid interval increases with height from 0.2 km near the surface to 0.37 km at the model top. The OSSE focuses on the three-dimensional distribution of wind within a 50 × 50 km domain around the simulated supercell at 1206 local standard time (LST), assuming that four Doppler radars are observing the storm (Fig. 1). Updrafts associated with the simulated storm reached approximately 12 km above sea level (ASL). The locations of the four radars are chosen so that the distance between each radar and the storm is approximately 30-40 km (the supercell is located in the center of the shaded domain shown in Fig. 1). High elevation angle (~20°) radar scans are required to observe the storm top. Figure 2 shows the heights of the simulated radar beams. Three different volume scan strategies are considered in this chapter. The first strategy assumes that one volume scan consists of 17 PPIs. This strategy corresponds to X-band radar surveillance of a thunderstorm with echo top below 10 km ASL during an interval of 5-6 minutes. The elevation angles used in this strategy are 0.7°, 1.2°, 1.7°, 2.2°, 2.8°, 3.3°, 3.9°, 4.7°, 5.6°, 6.5°, 7.4°, 8.3°, 9.3°, 10.3°, 11.8°, 13.5°, and 15.6°. The second strategy adds three high elevation angles (16.7°, 17.8°, and 18.9°) to the previous volume scan (blue lines in Fig. 2). This strategy corresponds to a sector-or adaptive-scanning mode (Junyent et al., 2010) for a tall thunderstorm located near the radar. The third strategy adds nine additional high elevation angles (20.0°-32.0° spaced at 1.5° intervals) to the preceding volume scan (red lines in Fig. 2).

www.intechopen.com
This corresponds to an ideal observing mode. The third scanning strategy appears to be impossible to implement at a rate of one volume scan every 5-6 minutes using current technology, but it may be possible using the high temporal resolution capabilities of nextgeneration phased array radars (Heinselman & Torres, 2011 Figure 3 shows the mixing ratio of rain and the distribution of horizontal and vertical winds at an altitude of 4 km ASL. The Weak Echo Region (WER) within the strong updraft region was well simulated (the maximum updraft was 25 m s -1 ). These strong updrafts were fed by southeasterly inflow below 1.5 km ASL (data not shown). Northwesterly wind was dominant at 4 km ASL, and advected the area of heavy precipitation toward the southeast (Fig. 3). Three downdraft cores were simulated at 4 km ASL. The first of these was located in the heavy rain region to the east of updraft, and was associated with precipitation loading. The second downdraft core was located in the light rain region to the southeast of the updraft, and was related to the melting and sublimation cooling of ice-phase precipitation (Shimizu et al., 2008). The third downdraft core was located in the non-precipitating region to the south of the updraft, and was associated with compensation for the nearby strong updraft. The strong updraft, first downdraft, and second downdraft cores were also simulated at 2 km ASL (Fig. 4). The maximum updraft speed exceeded 18 m s -1 at 2 km ASL. Anticlockwise wind rotation with a vertical vorticity of 0.08 s -1 was simulated along with this strong vertical velocity (Fig. 4). At 2 km ASL, the downdraft in the heavy precipitation region (maximum velocity -3 m s -1 ) covered a broad area to the northeast of the updraft at 2 km ASL.     Fig. 4. The strong updraft extended upward to 12 km ASL, with a maximum updraft speed of 24 m s -1 at 4 km ASL. Three downdraft cores flanked the strong updraft core. A weak compensating downdraft was located to the west of the strong updraft at approximately 6 km ASL. A second strong downdraft (6 m s -1 ) formed at 8 km ASL, and was associated with a high graupel mixing ratio (data not shown). A third strong downdraft was simulated at heights above 11 km ASL in the non-precipitating region. This downdraft was likely related to gravity wave dynamics.

Results
Four pseudo-radars are assumed to observe this simulated three-dimensional wind field (the locations of these pseudo-radars are shown in Fig. 1). The wind components are bilinearly interpolated from the model grid to sampling locations along the radar beams. Radial velocity is calculated from the three wind components interpolated to each radar grid point. The maximum range of detection is set to 70 km for all four pseudo-radars, as shown in Fig. 1. Each radar sweep observes a total of 90 azimuthal angles, with a gate spacing of 100 m and an azimuthal resolution of 1°. The radial velocities along each radar beam in the volume scan (see Fig. 2) are interpolated back onto the Cartesian coordinate system using a Cressman scheme with an influence radius R of 1.0 km. Some upper-level velocities will be lost because the highest elevation angle was less than 30° (Fig. 2). The robustness of the variational analysis method to noise is shown by adding random errors (mean 0 m s -1 , variance 1 m 2 s -2 ) to the radial velocities after interpolation back to the Cartesian coordinate system.
The accuracy of the multiple Doppler radar analysis is evaluated using root-mean-square errors in retrieved radial velocity (RMSE_VR) and vertical velocity (RMSE_W), defined as follows: The L-BFGS optimization scheme is able to successfully minimize the cost function. Figure 6 shows the relationship between the value of the cost function and the number of iterations performed. The value of the cost function is effectively constant after 20 iterations. The results of the variational multiple Doppler radar analysis using volume scans with 30 PPIs (experiment name: EL30) is shown in Fig. 7 and Fig. 8. Figure 7 shows that the EL30 experimental setup successfully retrieves the strong updraft core shown in the original model output (Fig. 4). The EL30 results indicate a maximum updraft speed at 2 km ASL of 12 m s -1 , and successfully reproduce the anticlockwise rotation at this level. The downdraft region located in the area of heavy precipitation to the east of the strong updraft is retrieved by EL30, although the size of this downdraft region is too small (Fig. 4). The horizontal wind field retrieved by EL30 is similar to the original model output (Fig. 4), but there are two major discrepancies in vertical velocity. First, the strength of the updraft speed at 2 km ASL is underestimated by 3 m s -1 . Second, a spurious downdraft is identified to the southwest of the strong updraft. The estimated northwesterly winds are much stronger in this spurious downdraft region than in the original model output (Fig. 4). This stronger northwesterly wind causes stronger divergence, which in turn induces the spurious strong downdraft. The strength of the enhanced northwesterly winds at 2 km ASL is similar to the strong northwesterly winds at 3-4 km ASL (cf. Fig. 2). Vertical smoothing of the radial velocity field by the Cressman interpolation procedure likely plays a major role in the erroneous vertical velocity field. Errors in radial velocity from vertical interpolation tend to occur near the boundaries of the storm, as shown in Fig. 7. This concentration of errors near the storm boundaries occurs because the number of radial velocity samples measured at neighboring grid cells is limited, so that the relative influence of radial velocities measured at distant grid cells grows. Fig. 7. Horizontal (vectors) and vertical wind velocity (shading) deduced by variational multiple Doppler radar analysis using 30 PPIs (EL30) at 2 km ASL. Solid contours indicate updrafts at 3 m s -1 contour intervals and downdrafts at 1 m s -1 contour intervals. The blue contour line outlines the region where the mixing ratio of rain exceeds 0.1 g kg -1 . A vertical cross-section along the thick horizontal line is shown in Fig. 8. Figure 8 shows a vertical cross-section of the wind fields retrieved by EL30 along the y = 26 km transect in Fig. 7. Compared with the original model output (Fig. 5), the strong updraft and downdraft cores in the heavy precipitation region (at 8 km ASL) and the downdraft core associated with gravity wave dynamics (at 11 km ASL) are well retrieved. The maximum retrieved updraft speed is 21 m s -1 , and occurred at approximately 4 km ASL. The maximum retrieved downdraft speed is 7 m s -1 . The maximum updraft speed was underestimated by 3 m s -1 , while the maximum downdraft speed was overestimated by 1 m s -1 . The speed of the downdraft associated with gravity wave dynamics was underestimated by 1 m s -1 relative to the original model output (Fig. 5). Several spurious updrafts and downdrafts can be identified in the EL30 retrieval, especially at 7-12 km ASL. Figure 7 and Fig. 8 indicate that the EL30 pseudo-radar configuration provides a good estimation of vertical velocity. Errors in the retrieved vertical velocity are uniformly less than 3 m s -1 .
www.intechopen.com Contours are as in Fig. 5. Fig. 9. As in Fig. 8 but with radial velocity observations contaminated by random noise.
As noted above, the use of the variational scheme removes the need to explicitly integrate the mass continuity equation, thus eradicating the vertical propagation of errors in the estimated wind field. Figure 9 shows the same vertical cross-section as in Fig. 8, but with random errors added to the input radial velocities. The distribution of retrieved vertical velocity is nearly identical in the cases with and without random errors. By contrast, the use of an upward integration scheme without random noise generates an erroneous vertical www.intechopen.com velocity field (too strong), especially at upper levels (Fig. 10). The addition of random noise to the upward integration method exacerbated these errors (data not shown).
These experiments reveal 1) that the variational approach provides more realistic estimates of vertical velocity than the upward integration method, and 2) that observational noise does not propagate upward when the variational method is used. These results have been obtained using a volume scan with 30 PPIs (experiment EL30). This choice of volume scan is currently unrealistic because it would take more than 10 minutes using typical current antenna rotation speeds (a few rotations per minute). For operational use, volume scans with 15-20 PPIs are realistic given the need for rapid updates (less than 400 seconds). Two further experiments are performed to mimic volume scan data with 20 PPIs (experiment name: EL20) and 17 PPIs (experiment name: EL17). The analysis method is the same for the EL20 and EL17 experiments as for the EL30 experiment. In EL20, the highest elevation angle is 18.9° (Fig. 2). A beam at this elevation angle reaches 14 km ASL 40 km away from radar. This is sufficient to ensure a valid upper boundary condition (w = 0) in this case, because the storm top is located at 12 km ASL. In EL17, however, the highest elevation angle is 15.7° (Fig. 2). This beam does not reach 12 km ASL within the detection range of the radar. Fig. 10. As in Fig. 8 but using an upward integration scheme. Figure 11 is the same as Fig. 8 but for EL20 rather than EL30. Vertical velocities retrieved by EL20 are similar to those retrieved by EL30 (Fig. 8) below 8 km ASL; however, the downdraft associated with gravity wave dynamics is not retrieved by EL20 because the pseudo-radar configuration does not observe that location. The maximum updraft was 20 m s -1 , a 1 m s -1 underestimate of the EL30 retrieval. As with EL30 (Fig. 8), a spurious updraft and downdraft are retrieved near 8-12 km ASL to the west of the strong updraft (Fig. 11). The area of this spurious downdraft is larger in EL20 than in EL30. This implies that EL20 is not capable of fully observing the storm top, so that erroneous upper boundary conditions induced the spurious downdraft.  Figure 12 is the same as Fig. 11 but for EL17 rather than EL20. The downdraft region associated with heavy precipitation to the east of the strong updraft is extended further downward in EL17. The maximum retrieved downdraft speed exceeds 9 m s -1 , a 3 m s -1 overestimate of the original model output (Fig. 5). These results indicate that even with the use of the variational approach, incomplete upper boundary conditions may lead to incorrect estimates of vertical velocity, especially near the storm top. www.intechopen.com Figure 13 shows root-mean-square errors in radial velocity (RMSE_VR) for the three experiments (EL30, EL20, and EL17). The EL30 RMSE_VR (black line in Fig. 13) was the smallest among the three experiments. The EL20 (red line) and EL17 (blue line) RMSE_VR both increased with height, while the RMSE_VR for all three experiments was small below 5 km ASL. The RMSE_VR and the root-mean-square error of the retrieved vertical velocity (RMSE_W) are closely related (Fig. 14). Differences in RMSE_W among the three experiments are relatively small below 5 km ASL and relatively large above 7 km ASL. Fig. 13. Vertical profile of RMSE_VR for the EL30 (black), EL20 (red), and EL17 (blue) scanning strategies.
For analyses of real thunderstorms the true value of vertical velocity is unknown, rendering the calculation of RMSE_W impossible; however, RMSE_VR can be calculated. Therefore, users could apply the variational multiple Doppler radar analysis approach to operational scans if the vertical profile of RMSE_VR is used to validate upper-level vertical velocities. The variational approach is a useful tool even for operational volume scans (less than 20 PPIs) because accurate three-dimensional winds can be retrieved at lower levels without contamination from sampling error in upper level or uncertainty in the upper boundary condition. These errors do not propagate vertically under the variational approach. Future work should examine the dependence of the retrieved wind field on scan strategy using observational datasets generated from real radar networks.

Conclusions
This chapter has introduced a variational multiple Doppler radar analysis for retrieving three-dimensional wind fields in a severe thunderstorm. A simplified version of the method presented by Gao et al. (1999) has been used to investigate the dependence of retr ie ve d ver ti ca l v e lo ci ty o n s can str at egy. Thr ee vo lu me s can s trat egi es h av e b e en considered in this chapter: 1) a typical operational volume scan (17 PPIs), 2) a dense operational volume scan (20 PPIs), and 3) an extremely dense volume scan (30 PPIs). The variational approach has notable advantages over the upward integration method, particularly the avoidance of error accumulation during the upward integration; however, incomplete observations of the upper boundary condition can cause errors in estimates of vertical velocity near the storm top even when the variational approach is used. Users should limit their use of upper-level wind retrievals according to the root-mean-square error of radial wind, as described in this chapter. The variational method provides accurate estimates of the three-dimensional wind field at lower altitudes regardless of the upper boundary conditions. The density of operational radar networks in metropolitan regions has been increased in recent years to better monitor and forecast severe weather. Together with this increase in operational radar network density, the variational analysis method presented in this chapter will provide new information on the three-dimensional structure of wind within thunderstorms, and advance understanding of the physical mechanisms underlying heavy rainfall and severe winds. Doppler radar systems have been instrumental to improve our understanding and monitoring capabilities of phenomena taking place in the low, middle, and upper atmosphere. Weather radars, wind profilers, and incoherent and coherent scatter radars implementing Doppler techniques are now used routinely both in research and operational applications by scientists and practitioners. This book brings together a collection of eighteen essays by international leading authors devoted to different applications of ground based Doppler radars. Topics covered include, among others, severe weather surveillance, precipitation estimation and nowcasting, wind and turbulence retrievals, ionospheric radar and volcanological applications of Doppler radar. The book is ideally suited for graduate students looking for an introduction to the field or professionals intending to refresh or update their knowledge on Doppler radar applications.