The radar network of the German Weather Service [Deutscher Wetterdienst (DWD)] provides 3D Doppler data in high spatial and temporal resolution, supporting the identification and tracking of dynamic small-scale weather phenomena. The software framework Polarimetric Radar Algorithms (POLARA) has been developed at DWD to better exploit the capabilities of the existing remote sensing data. The data processing and quality assurance implemented in POLARA include a dual-PRF dealiasing algorithm with error correction. Azimuthal shear information is derived and processed in the mesocyclone detection algorithm (MCD). Low- and midlevel azimuthal shear and track products are available as composite (multiradar) products. Azimuthal shear may be considered as a proxy for rotation. The MCD results and azimuthal shear products are part of the severe weather detection algorithms of DWD and are provided to the forecaster on the NinJo meteorological workstation system. The forecaster analyzes potentially severe cells by combining near-storm environment data with the MCD product as well as with the instantaneous azimuthal shear products (mid- and low level) and their tracks. These products and tracks are used to diagnose threat potential by means of azimuthal shear intensity and track longevity. Feedback from forecasters has shown the utility of the algorithms to analyze and diagnose severe convective cells in Germany and in adjacent Europe. In this paper, the abovementioned algorithms and products are presented in detail and case studies illustrating usability and performance are shown.
Severe weather associated with convective storms poses a serious threat to life and property. The presence of storm-scale rotation, or mesocyclones, increases the risk of heavy rain, severe hail, strong wind gusts, and tornadoes within these storms. Thus, the real-time detection of mesocyclonic shear gives valuable additional information for the meteorological warning management of a national weather service. The radar network of the German Weather Service [Deutscher Wetterdienst (DWD)] offers scanning of the lower atmosphere with high spatial and temporal resolution so that small-scale dynamic severe weather phenomena, such as mesocyclones, can be resolved.
In the context of radar meteorology, observations of mesocyclonic signatures in the form of storm-scale Rankine vortices have been analyzed (Donaldson 1970; Burgess 1976) and several mesocyclone detection algorithms have been developed in recent decades (Zrnić et al. 1985; Stumpf et al. 1998; Joe et al. 2004). The operational mesocyclone detection algorithm (MCD) at DWD uses a pattern vector approach as described in Zrnić et al. (1985) to identify the regions of high shear within the center of the mesocyclonic rotation. An important aspect of the DWD algorithm is the merging of 2D mesocyclone features from multiple network radars in overlapping regions to create a single 3D detection. A similar technique was demonstrated in the multiple-radar severe storm analysis program for three neighboring WSR-88D radar systems (Stumpf 2002).
The mesocyclone objects are ranked and assigned a severity index using mesocyclone object properties like shear and rotational velocity. This is provided to the forecaster, as guidance, through the NinJo meteorological workstation system (Koppert et al. 2004). Meteorologists then judge the significance of the mesocyclone detections using this severity metric, applying persistency, consistency checks (track of mesocyclone detections), and forecaster judgment (additional occurrences of severe weather features like hook echoes, interpretation of near-storm environment data, and large-scale forcing analysis). Further guidance is given by the azimuthal shear and corresponding track products.
In this paper the scan strategy and the data acquisition within DWD’s radar network are outlined in section 2. Doppler wind data quality assurance issues of great importance for downstream applications are discussed in section 3. The mesocyclone detection and the rotation compositing are introduced in sections 4 and 5, respectively, as network-based algorithms. The operational applicability of the MCD and the rotation products and the advantages of the network approach are demonstrated in section 6 using case studies. Comparative statistical analyses of the MCD (through comparison of the MCD with previous studies and through comparison of the MCD network and single-site modes) and studies showing correlations between the cell lifetime and the cell’s rotational characteristics are presented in section 7. The main conclusions are summarized in section 8.
DWD operates a network of 17 weather radar systems delivering 2D and 3D radar data in 5-min intervals (Helmert et al. 2014; see Fig. 1). Sixteen radar sites are equipped with dual-polarimetric systems. In Emden (northwestern Germany) a single-polarimetric radar system is in operation. The volume scan is composed of 10 elevation scans with elevation angles ranging from 0.5° to 25°. This scan is configured such that both an adequate maximum unambiguous Doppler velocity interval and a maximum unambiguous range are achieved (Seltmann et al. 2013). To obtain sufficient overlap between the radar sites located at a mean nearest-neighbor distance of approximately 140 km, a maximum range of 180 km was determined to be required for the six lowermost elevation scans. This is realized by implementing a 4:3 dual-PRF dealiasing scheme (Dazhang et al. 1984) with a so-called extended maximum unambiguous velocity of 32 m s−1. Since true Doppler velocities do not typically exceed 64 m s−1, extended aliasing effects will be at most a single fold. For the higher elevations, shorter ranges are sufficient for scanning the troposphere and less demanding dual-PRF ratios or single PRFs are used, always keeping the maximum unambiguous velocity at 32 m s−1 (see Table 1).
Strong quality control of the radial velocity data is needed for the mesocyclone detection and azimuthal shear algorithms. The first step in quality control is achieved by means of a zero Doppler velocity filter to remove stationary clutter and a set of data quality threshold filters, including clutter correction (CCOR) and signal quality index (SQI) (Seltmann 2000). It is notable that each pulse volume is composed of uniform PRF samples. The clutter filter process is not complicated by nonuniform samples, which is the case for the staggered PRT technique (Torres et al. 2004).
The volume scan data from all DWD radar stations are gathered at the central office of DWD in real time by means of an automated file distribution system (AFD; Deutscher Wetterdienst 2016a) and are available for downstream applications. To adequately exploit the possibilities of the polarimetric data, a software framework called Polarimetric Radar Algorithms (POLARA) was developed (Rathmann and Mott 2012). POLARA contains various algorithms, such as quality control, hydrometeor classification, and quantitative precipitation estimation. As part of quality control, the correction of dual-PRF dealiasing errors is performed to eliminate spurious velocity errors in high shear regions. This is an essential prerequisite for the mesocyclone detection and azimuthal shear algorithms. The quality-controlled radar elevation scan data are henceforth called base data. These base data are stored at the central office of the DWD in a long-term archive system (Deutscher Wetterdienst 2016b). The radar base data and the derived products are visualized in the radar and Storm Cell Identification and Tracking (SCIT) layers (Joe et al. 2005) of the NinJo workstation (Koppert et al. 2004).
3. Quality assurance of Doppler wind data at DWD
In this section one of the most critical problems in the quality assurance of Doppler wind data at DWD is addressed: the correction of dual-PRF dealiasing errors. The reader is referred to appendix A for details concerning the dual-PRF dealiasing technique, which is assumed to be known.
The two velocities and of azimuthally adjacent range bins, which enter into the dual-PRF dealiasing, should ideally be equal [see Eq. (A1)], that is, showing apart from possible aliasing no further deviations. However, small discrepancies are tolerable and the dual-PRF dealiasing works correctly as long as the (dealiased) absolute gate-to-gate velocity difference Δυ stays within the following limit (Joe and May 2003):
Here, and are the PRFs related to and , respectively; and and are the associated maximum unambiguous velocities. In case of the 800-/600-Hz scan mode, velocity differences up to 2.65 m s−1 will not cause dealiasing errors. However, this tolerance will be exceeded because of ubiquitous noise fluctuations. Corresponding errors typically appear as randomly scattered erroneous pixels. Additionally, in regions of high shear, such as in atmospheric vortices, the tolerance is naturally exceeded and erroneous data appear in small pixel groups. It shall be noted that dual-PRF dealiasing errors may occur, even if the data in the individual adjacent beams are not aliased (Joe and May 2003). The distribution of dual-PRF dealiasing errors is discrete, as the errors are multiples of 2 or 2 (May 2001). The dual-PRF dealiasing errors are corrected by means of a detection–correction filtering algorithm. Within the POLARA quality assurance module (Werner 2014), a Laplacian filter based on Joe and May (2003) is used. The Laplacian filter basically uses a pixel window, with k being an odd integer, to detect and then to calculate a correction value for the window’s central pixel’s radial velocity value . The correction uses the related maximum unambiguous velocity , which is either or . It is performed by comparing a phase-averaged mean value , calculated within the extended unambiguous velocity interval from the central pixel’s neighbors, to all possible corrections with n being the integer fold number,
The value , which differs least from , is used as the corrected value , which is therefore defined as
where is determined by minimizing with respect to n.
By using a phase-averaging procedure [see, e.g., section 2.3.8 on periodic variables in Bishop (2006)], acting within the extended velocity interval, aliasing within this interval does not impose problems, as and are treated as identical values (analogous to −π and π denoting the same actual phase angle). It should be noted that Eqs. (2) and (3) are like Eqs. (A3) and (A4), respectively, except that in Eq. (3) is determined more robustly.
The dual-PRF dealiasing error correction in POLARA is performed in two stages, which are described as follows:
In a first stage, a large window (11 × 11 pixels) Laplacian filter is applied to the data without correcting the data. Range bins with erroneous velocities are identified and marked. An identification is claimed when a corrected value different from the original value can be calculated according to Eqs. (2) and (3).
In a second stage, a small window (5 × 5 pixels) Laplacian filter is applied to the data, and the data are corrected. However, only those range bins not marked in the first stage (i.e., presumably valid, error-free pixels) are used as neighbors within the Laplacian filter’s averaging procedure.
It shall be noted that the dual-PRF dealiasing and subsequent error correction, even when working perfectly, provide velocity data that may show aliasing within the extended maximum unambiguous velocity interval. Figure 2 shows the performance of the dual-PRF dealiasing error correction by means of simulation. A vortex core situated at a distance of 62 km from the radar site (red-encircled area) embedded in a uniform wind field has been modeled. It is important that the erroneous data pixels are corrected but shear structures like Rankine-combined vortices are not eliminated.
Figure 3 shows observational radar data before and after applying the dual-PRF dealiasing error correction. At the radar site Munich (MUC), a test scan setup with a 250 m × 1° resolution was used. The 4-times-higher range resolution together with the 800-/600-Hz dual-PRF setup causes a lower number of samples per pulse volume (implying a degradation of data quality) in comparison to the current operational scan setup, thus resulting in a “worst-case scenario.” It can be seen that a vortex (white-encircled area) is recovered in the quality-controlled data.
Figure 4 shows the performance of the dual-PRF dealiasing error correction for different positions of a small (4-km diameter) and a large (10-km diameter) vortex with respect to the radar site. For the small-diameter simulation, an observed vortex with rotational velocity = 8.5 m s−1 was used. The large-diameter simulation is based on a vortex with = 15 m s−1. It can be seen that for distances below approximately 100 km, the quality assurance works almost perfectly. At distances between 100 and 180 km (the maximum unambiguous range), the number of pixels inside the vortex core showing erroneous values increases. However, this number remains below 28%.
The dual-PRF dealiasing error correction shows the obvious restriction of having enough pixels with valid radial velocity values within the Laplacian filter window.
Especially for contiguous radar echo areas consisting of only a few pixels, a correction is difficult. However, there have not yet been observed erroneous mesocyclone detections or significantly large artifacts in the azimuthal shear products (see the following sections) as a result of remaining dealiasing errors. This can probably be attributed to the use of geometrical thresholds, which inhibit such small areas from producing significant effects.
4. The mesocyclone detection algorithm of DWD
a. Outline of the algorithm
The mesocyclone detection algorithm can be configured for both the detection of cyclonic and anticyclonic rotation. In its current operational configuration, the MCD performs a search only for cyclonic rotation (positive azimuthal shear), which is the dominant mode of mesocyclonic rotation on the earth’s Northern Hemisphere. Mesoanticyclones are found less frequently, typically in conjunction with cyclonically rotating systems, such as in the case of bow echoes with “bookend” vortices and cyclone–anticyclone couplets, in the case of left movers with splitting storms and in the case of supercells with a cyclonically rotating updraft, where anticyclonic shear can be observed in the midlevel (Bluestein 2013).
The mesocyclone detection algorithm basically follows the pattern vector approach introduced by Zrnić et al. (1985). There are three parts to the approach
the detection of pattern vectors (1D),
the grouping of pattern vectors to features (2D), and
the combination of features to mesocyclone objects (3D).
1) The detection of pattern vectors
After correcting the radial velocity data for spurious velocities, the next step is the calculation of the azimuthal shear for each elevation scan, whereby the possibility of aliasing (folding within the extended maximum unambiguous velocity interval) must be considered. First, we define
where and are the quality assured (corrected) radial velocity values from azimuthally adjacent range bins. Let be that gate-to-gate velocity difference, which shows the smallest absolute value,
This defines the corrected gate-to-gate velocity difference, and the corresponding shear value is obtained as with b being the center-to-center distance between the two adjacent radial velocity bins.
Insignificant shear values are treated as missing values; that is, they are ignored by the subsequent pattern vector search. The original velocity values of all pixels are kept in memory to allow for exact calculations later on (of, e.g., shear, specific angular momentum, rotational velocity). The noise level is determined by the error of being 2 m s−1 with the single-pixel velocity error σ 1.5 m s−1, where σ is estimated from corrected observational radial velocity data by performing a 3 × 3 Laplacian filter analysis and investigating the deviations of the central pixel’s value from its mean value. Thus, shear values with related gate-to-gate shear < 2 m s−1 can be regarded as insignificant (further reasoning, discussed at the end of this subsubsection, leads to a weaker threshold used operationally). In cyclonic detection mode, negative shear is treated as noise, which reduces the influence of possible azimuthal shear artifacts (e.g., remaining clutter contamination) in regions of positive cyclonic shear.
A pattern vector is a sequence of positive azimuthal shear at constant range. Figure 5 shows an example of an observed Doppler rotation signature (typical velocity couplet) and the pattern vectors that are identified within the area of rotation. The pattern vector search accepts interruptions of one range bin within a pattern vector to account for data dropouts (data values below the reflectivity or azimuthal shear noise thresholds, or residual azimuthal shear caused by incorrect radial velocities). Each pattern vector is filtered with respect to Doppler (specific) angular momentum and azimuthal shear (Zrnić et al. 1985) so that only shear regions above the noise threshold are accepted for further analysis. The following set of thresholds is used:
In applying the latter set of thresholds (the OR condition), the idea is to identify high shear and small-scale rotations (more related to the late, mature stage of mesocyclone development) as well as high angular momentum and large-scale rotations (more related to the early stage of mesocyclone development). Thresholds as low as approximately one-third of those given in Zrnić et al. (1985) have proven suitable for the radar network of DWD. It shall be noted that threshold adjusting for pattern vectors is necessary so that environmental differences in mesocyclones and warning criteria can be considered (Crozier et al. 1991, p. 494, last paragraph).
Picking up noisy azimuthal shear during the pattern vector search is suppressed by the selection criteria for significant pattern vectors [Eqs. (7) and (8)] and by the conditions for clustering pattern vectors to features (and features to mesocyclone objects) as discussed in the sections 4a(2) and 4a(3) below. Using simulations we determined a threshold of significant gate-to-gate shear lower than that introduced above (approximately 0.8 m s−1).
For the specific angular momentum and azimuthal shear of a pattern vector, one can find that and with n being the number of azimuthal shear pixels that make up the pattern vector. Thus, the azimuthal shear of a pattern vector is in fact the mean azimuthal shear of the contributing pixels and the shear threshold of Eq. (7) could already be applied to the single pixels. Experience has shown that pattern vectors can possess an embedded higher shear sequence (with lower shear at the start and end points), likely since vortices observed in radar data are smoothed due to radar resolution effects (Wood and Brown 1997). As a result the mean shear drops below the threshold of Eq. (7) and an actually significant pattern vector can be lost. Applying the low shear threshold given in Eq. (7) to single pixels of the azimuthal shear data (prior to the pattern vector search) preselects higher shear and prevents such loss.
2) The grouping of pattern vectors to features
The mesocyclone detection proceeds by grouping the pattern vectors to features. Ideally, a feature corresponds to a complete rotation signature within an elevation scan and consists of several pattern vectors adjacent in range. However, in reality, there is always some noise present [missing data as a result of signal filters; corrupted data, such as remnants of clutter or wireless local area network (WLAN) interference spokes not fully removed by the quality assurance], and the pattern vectors may be noncontiguous. Thus, the following criteria are used to extract only meaningful features and, at the same time, allow for imperfect data:
First, the number of detected pattern vectors must exceed a certain threshold that depends on the range resolution and the minimum assumed vortex spatial extent. In Zrnić et al. (1985), for elevation scan data with a 150-m range resolution, a minimum of six pattern vectors is required for feature detection so that vortex diameters > 1 km can principally be resolved. Since at DWD the volume scan data have a 1-km range resolution, a smaller number of three pattern vectors is currently used (requesting just two adjacent pattern vectors is not enough for suppressing chance coincidences caused by noise, especially when considering the comparably low shear and angular momentum thresholds necessary for detecting significant rotational signatures). The adjacency criterion allows for a separation of one range bin between the pattern vectors associated with the same feature.
- Second, symmetry criteria must be passed: the spatial extension of features should be roughly equal in azimuth and range direction. The symmetry criteria, based on Zrnić et al. (1985), are as follows: is the feature diameter in the radial direction (distance between the pattern vector closest to the radar and the pattern vector farthest from the radar). The terms , , and are estimates for the feature diameter in the azimuthal direction: is the momentum weighted diameter, is a root-mean-square value, and is a simple average [for more details see Eqs. (A7)–(A9) in the appendix of Zrnić et al. (1985)]. However, only the pattern vectors with the largest specific angular momentum (typically found in the core region of the rotation) are chosen to enter into the symmetry calculations: in the case of a vortex field, these pattern vectors will cover the core region, whereas in the case of a linear wind event, such as a gust front, the pattern vectors are scattered along the gust front and will violate the symmetry conditions. The reason for introducing this core region restriction is that one frequently observes a vortex distorted by a convergence zone caused by the thunderstorm inflow. The pattern vectors extend into the convergence region and, without the large momentum restriction, would then violate the symmetry condition so that a possibly strong rotation would remain undetected. See Fig. 5 for an example of a feature detected during the occurrence of a tornadic supercell.
3) Combining features into 3D mesocyclone objects
Finally, the 2D georeferenced features from all radars in the network are combined into three-dimensional vertically aligned mesocyclone objects (see Fig. 6) by means of a connected components labeling algorithm (see, e.g., Shapiro and Stockman 2001, section 3.4). Vertical alignment (overlap) between two 2D features is given if the offset of the feature centers is less than the sum of the feature radii. In this way, tilted mesocyclones (physical and perhaps created as a result of the time offsets between the scans) are taken into account. Features of all elevations and of all radar sites are considered (network approach). These 3D objects are referred to as mesocyclone objects and are candidates for (detected) mesocyclones. Compare also to Stumpf (2002), in which a similar technique for the feature combination was demonstrated.
b. Mesocyclone object properties
Properties of the 3D mesocyclone object, such as maximum azimuthal shear, maximum rotational velocity, height of the lowest detected rotation signature above ground, total depth of the rotational column and mean and maximum reflectivity, vertically integrated liquid water (VIL), and VIL density, are calculated (see Table 2). These properties (also referred to as attributes) can be categorized into geographical coordinate, shape, Doppler velocity, and reflectivity attributes. For a more robust calculation of coordinates and reflectivity-related quantities, an expanded circular area of twice the diameter of the rotation features (i.e., mesocyclone object region) is considered. This is done in order to take into account the fact that the mesocyclone location and spatial extent do not usually coincide with the parent cell’s position and extent (as delineated by reflectivity). Calculated within this expanded area using the elevation scan base data, the reflectivity-related attributes reflect the parent cell’s properties. The expansion proved useful also for the calculation of geographical mesocyclone coordinates, giving more robust and stable mesocyclone positions especially in the case of gaps in the radial velocity data. In appendix B a detailed overview of all calculated attributes is given.
There are limitations to the accuracy of these attributes imposed by the geometry of the volume scan. Because of the spherical grid geometry of the elevation scans, a detected feature shows a range-dependent resolution in the directions of azimuth and elevation, the range resolution being constantly 1 km. For a range location at 60, 120, and 180 km, the horizontal resolutions are 1.05, 2.09, and 3.14 km, respectively. The related vertical resolutions show the same values because the radar main lobe has a symmetrical 1 width. The heights of the lowermost elevation scan at the given range locations are 0.21, 0.85, and 1.90 km, respectively. This implies a degradation of the accuracy of geometrical (diameter, height), reflectivity-based, and rotational attributes with increasing range. Furthermore, rotational signatures close to the ground may be below the lowermost elevation scan when situated at large distances from the radar.
The network approach is able to partly compensate for these degradations by merging all detected features to obtain a mesocyclone object with the best possible sampling. The estimation of the detected mesocyclone’s severity (see section 4c) is mainly based on “maximum properties” (maximum shear, maximum rotational velocity). This complies with the requirements of the warning service, where an overestimation is more acceptable than an underestimation. Furthermore, range degradation effects typically lead to underestimated parameters of rotational strength (too low shear, too low rotational velocity), so the “maximum picking” naturally selects the parameter with the best resolution available.
A further solution is to homogenize the feature’s parameters if possible. Concerning the azimuthal shear, this can be done by applying a range correction (Newman et al. 2013). This correction is also relevant for the rotation composites as introduced in section 5 and is one of the major steps to be addressed in a future development (see section 8).
c. The mesocyclone severity metric
To reject rotational signatures that were not associated with a thunderstorm, the thresholds given in Table 3 are applied to all mesocyclone objects. In particular, at least two features are needed for successful mesocyclone detection.
The severity of the mesocyclone is computed using the severity ranking shown in Table 4, which is based on geometry and the strength of rotation. Severity levels 1–5 imply mesocyclonic rotation with increasing strength and confidence. According to the Bureau of Meteorology (2016), the necessary conditions for deep mesocyclonic rotation are that the vortex should be of storm scale (2–10 km), have a reasonable vertical depth (3 km), have a sufficient duration, and meet a minimum strength criterion (a rotational velocity of 15 m s−1 or a vertical vorticity of 10−2 s−1). Severity level 3 is established to meet these criteria by applying the respective thresholds on , , , and . Sufficient duration of a mesocyclone detected by the MCD can be tested by considering the detection history and by comparison with rotation track composites (see section 5). A statistical analysis of the mesocyclone detection algorithm is provided in section 7.
5. The DWD rotation and rotation track products
The azimuthal shear composite products offer independent and additional visual insight. The preprocessing of azimuthal shear is inspired by the linear least squares derivative (LLSD) algorithm of the National Severe Storms Laboratory (Smith and Elmore 2004; Miller et al. 2013; Newman et al. 2013) in the sense that a noise-suppressing extraction of azimuthal shear from the elevation scan base data is attempted. The analysis of shear derived from radial velocity data is a specific challenge, since derivatives are particularly noisy as a result of the combination of variances. A single range bin with spurious velocity can cause a “hot spot” within the azimuthal shear composite. Thus, a series of smoothing and averaging techniques is used to obtain a more robust final product.
First, a smoothing 3 × 3 pixel window filter is applied to the radial velocity data. This filter essentially works identically to the one used for dual-PRF dealiasing error correction with the exception that the central pixel’s radial velocity value is replaced by the phase-averaged mean of the neighboring eight pixels. In this way, single noisy pixels are smoothed out, while regions where the radial velocities change linearly with increasing range and azimuth are less affected and thus emphasized. Such regions imply constant shear (as is typical for the core of mesocyclones). Therefore, the Laplacian filter smoothing mitigates single noisy pixels while affecting vortices of constant shear less.
Vertical averaging, according to Eqs. (12)–(14), is also performed so that deep rotation is emphasized in comparison to shallow rotation and noise. The algorithm works on a composite grid in polar-stereographic projection with 1100 × 1200 pixels of approximately 1 km × 1 km resolution. A given grid cell of the composite is assigned a sample containing (, , ) triples (values of shear , reflectivity , and height with , where n is the number of samples over the depth of the atmosphere). Each triple is extracted from a certain radar site’s elevation scan range bin (the projected center of the range bin must fall into the composite grid cell). All elevation scans are considered. In regions of radar overlap, the sample shows contributions from multiple radars. The situation is analog to that depicted in Fig. 6 with the orange patches now being interpreted as radar range bin data. The triples are sorted in increasing order of height and a rotation value ζ is calculated as a height–difference weighted mean,
The azimuthal shear product (further referred to as rotation composite) is created for two altitudes: The low-level rotation composite (ROT-LL) shows the mean shear between 0 and 3 km (AGL), while the midlevel rotation composite (ROT-ML) depicts the mean shear between 3 and 6 km (AGL). Both rotation composites show positive cyclonic and negative anticyclonic shear.
The composite rotation track products are created by accumulating the maximum shear so that in the operational setup, each grid cell shows the maximum (cyclonic) shear of the last 3 h (the accumulation time is configurable and can be set on demand for offline analyses). This product highlights persistent, significant, and moving cyclonic vortex areas, which is typical for organized deep convection, such as supercells or bow echoes. A separate product for anticyclonic shear is conceivable, but experience shows that it has limited value. False signatures caused by significant vertical shear near the surface, or outflow boundaries (Miller et al. 2013), can usually be distinguished by eye from real rotation tracks using these rotation products (false signatures typically do not produce line or stripe-like patterns, as is observed for “real” rotation tracks).
Any algorithm utilizing radar data suffers from the uneven polar grid spacing of the base radar data. Deterioration (loss of horizontal and vertical resolution) of the measurements of a single radar with increasing range is unavoidable. The network approach is better able to compensate for the loss in vertical sampling by combining all available radar measurements. However, problems remain in composite border regions (the edge of the radar network domain).
Rotation tracks are better defined in the central region of the composite and not too close to any of the radar sites, where vertical wind shear may cause high azimuthal shear signals (“hot spots”) that are not related to rotation. The forecaster must be aware of these limitations to carefully draw decisions in operational warning situations.
6. Case studies
In the following subsections, two case studies are presented to demonstrate the value of the network approach. Furthermore, the potential of the MCD and rotation composites to support forecasters in their warning decisions is shown. These products are visualized together with VIL track composites, which highlight cell tracks and give a general impression of the severe weather potential.
VIL values in excess of 10 kg m−2 indicate convection. For the most severe events, VIL can exceed 60 and approach 100 kg m−2. VIL and VIL track composites follow the network algorithm approach, as described for the rotation and rotation track composites. Analogous to the rotation track composites, the VIL track composites show a maximum value 3-h history in the operational configuration.
a. Tornadic supercells in southern Germany on 13 May 2015
Within a humid, unstable air mass, thunderstorms developed in the afternoon and reached high severity over southeastern Bavaria, where high CAPE and shear values were present and dynamical forcing supported the lifting processes. In the evening marginally severe thunderstorms followed from the west. They were associated with a surface low that propagated eastward over the most southern part of Germany.
The VIL track composite in Fig. 7 indicates severe weather in southwestern Germany. High VIL values in excess of 60 kg m−2 were measured. It shall be noted that the times given in Fig. 7 are reference times (start of volume scan). The volume scan starts at approximately 40 s after the full minute and takes approximately 4 min, 30 s to finish, so a time interval of approximately 4 h is covered (1700:00–2055:00 UTC reference time corresponds to 1700:40–2059:30 UTC actual observation time). An accumulation time of 4 h rather than 3 h was chosen for this case study to show most complete tracks.
Besides smaller signals three major supercells were detected by the MCD. These are also very visible in the rotation track composites (Fig. 8). One supercell moves directly over the radar at Feldberg, another on a parallel track displaced to the northeast, and the third supercell track can be seen close to the radar at Türkheim. All hail and tornado reports were located on or close to the VIL and rotation tracks. Most notably, two F3 tornadoes (see Fig. 8) can be associated with supercells detected by the MCD, which could be tracked constantly for the previous 2 h before the tornadoes were observed.
Figure 9 shows the mesocyclone detections for the same time range as in Fig. 8. Each detection is represented by a circle with a diameter corresponding to the equivalent diameter. In Figs. 9a and 9b mesocyclone detections obtained when running the MCD in the default network mode are depicted. From Fig. 9a one can see that up to six radar sites contributed to a detection. The two detections with the largest number of features show 22 and 25 detected rotation signatures, respectively. Figures 9b and 9c show the mesocyclone severity for the network and the single-site approach, respectively. In particular, differences can be seen for the northern mesocyclone track. Besides eliminating multiple simultaneous detections, it seems that the MCD in network mode generates better severity estimates. The combination of features detected by the radar at Türkheim and the neighboring radar sites results in detections with severity levels 4 and 5. Single-site detections alone fail to generate the highest severity levels as a result of insufficient rotational strength or vertical extent. The radar at Türkheim alone, though being close to the mesocyclone, is not able to capture the rotation in its full shape and characteristic. This can be understood from restrictions in the scan geometry (cone of silence) and the fact that the strongest rotation within a mesocyclone may occur at higher altitude.
Figure 10 illustrates the benefit of the network approach for the VIL and rotation track composites. The effects of the blanking sector (area labeled “1”) and the cone of silence visible within the area marked by ellipse 2 in Fig. 10a are mitigated in Fig. 10b. At far ranges the single-site VIL shows biased values [overestimation in the case of strong vertical reflectivity gradients (Brown and Wood 2000; Bellon and Zawadzki 2003) and underestimation in the case of overshooting precipitation]. From Figs. 10c and 10d, it can be seen that the (low level) rotation track composite shows less noise and shallow rotation signals when using the network mode. A single main rotation track with a bend toward the east becomes visible in the encircled area of Fig. 10d, showing a characteristic “fishbone pattern” (an effect of the periodic scanning) for the whole track.
b. A long-lived hail-producing supercell in northern Germany on 28 August 2016
On 28 August 2016, a midlevel trough was approaching western Europe accompanied by a strong southwesterly flow. Germany was dominated by an unstable, humid, and hot air mass. The first thunderstorms developed west of Stuttgart in southwestern Germany and in Emsland in the northwest. In the course of the afternoon, severe thunderstorms occurred in almost all parts of northwestern Germany. A single supercell showed a track of approximately 200-km length in an easterly direction starting south of Hamburg. Hail with diameters up to 8 cm was observed.
In Fig. 11 the VIL track of the single right-moving supercell is visible as a broad line of red, violet, and white extending from west to east (as opposed to the majority of smaller thunderstorms moving northeast). VIL values reached up to 100 kg m−2. The aforementioned hail observation (8-cm diameter) correlates with the easterly end of the VIL track. Severe weather events, reported in the European Severe Weather Database (ESWD), are overlaid on the VIL track in Fig. 11b. The hail-producing potential of the single supercell is evident.
The low- and midlevel rotation track composites (Fig. 12) indicate regions of stronger and weaker rotation, the midlevel rotation being generally stronger (especially for the easterly part of the track). The low-level rotation composite shows some isolated pixels with high azimuthal shear values. However, these are most likely artifacts, related to the fact that elevation scan base data at low elevations suffer more from “impurities,” such as clutter contamination, than elevation scan data at higher elevations. The mesocyclone was detected over the period from 1305 to 1650 UTC. Most of these detections are associated with a mesocyclone base > 3 km.
Figure 13 shows the effect of the network approach on the MCD. In Fig. 13a the MCD detections are depicted with a color coding that indicates the number of radar sites that contributed to a given detection. For almost all detections there are contributions from at least two sites and even for some detections from four sites. For an effective suppression of false detections, at least two features are needed for a successful mesocyclone object identification (see Table 3). Several MCD detections in Fig. 13a are composed of single features from single sites. If the MCD were used as a single-site algorithm (treating the base data of all radar sites independently), then these detections would be lost. For larger parts of its track, the mesocyclone is not adequately sampled by a single radar alone, since it has a higher base and is at some distance from the neighboring radars (approximately 80–200 km for the easterly part of the track).
Figures 13b and 13c show the detections obtained in the network and the single-site mode, respectively. Here, the color coding indicates the severity level. Ellipse 1 in Fig. 13c marks a region without any MCD single-site detections. The early detections visible in Fig. 13b are lost in Fig. 13c. Furthermore, in Fig. 13c there is a double detection of the same mesocyclone object within the area marked by ellipse 2. The MCD detection at 1410 UTC in Fig. 13b showing single-feature contributions from the radar sites Boostedt, Rostock, Hannover, and Ummendorf is lost in Fig. 13c. Ellipse 3 in Fig. 13c highlights an area where again many detections are lost when the MCD is operating in single-site mode.
7. Statistical analysis
a. Comparison of MCD statistics with previous studies
In this subsection mesocyclone detection statistics are presented and compared with previous studies. The characteristics of mesocyclones detected by the MCD are investigated in Wapler et al. (2016) based on a 3-yr analysis (April–September 2013–15). No restriction on the number of features per mesocyclone object was applied at that time. The same is true for the analysis presented in Bellon and Zawadzki (2003), who studied mesocyclones in southern Quebec, Canada, using a 9-yr dataset from one radar. In Wapler et al. (2016), characteristics of the MCD detections are compared to the results given in Bellon and Zawadzki (2003).
The annual and diurnal cycles of mesocyclone detections as obtained by the MCD show typical (broad) maxima in June–August and for 1600–1700 UTC, respectively (Wapler et al. 2016). It was found that for all six strong tornadoes (F2) in the years 2012–14 within Germany, a mesocyclone was detected by the MCD well before the tornado occurrence. Strong tornadoes were selected because, in the case of weaker tornadoes ranked as F0 or F1 on the Fujita scale, the possibility of nonmesocyclonic origin is high (Brooks and Dotzek 2008; Dotzek et al. 2005). In Wapler et al. (2016), the occurrence of mesocyclone detections for hailstorm events was investigated. This analysis revealed that about half of all hailstorms in Germany are associated with a mesocyclone detected by the MCD during the time of reported hail. However, newer studies that analyzed the entire life cycle of hailstorms showed that approximately three-quarters of all hail-producing cells in central Europe show rotational characteristics (Wapler 2017).
In addition to annual and diurnal cycles, the frequency plots for the mesocyclone attributes , , , , and the maximum and mean shear and moment (, , , ), as well as VIL, echo-top height , and are presented in Wapler et al. (2016).
The depth of 85% of the detected mesocyclones is below 3.5 km. Because of the thresholds applied in the detection algorithm, the minimum mesocyclone diameter is 3 km (at least three adjacent pattern vectors are required to form a feature). Half of the mesocyclones have a diameter between 3 and 6 km. Two-thirds of all mesocyclones have a maximum shear between 4 and 8 m s−1 km−1.
These shape and rotational characteristics of the MCD mesocyclone detections are in general agreement with those presented in Bellon and Zawadzki (2003).
b. Comparison of the MCD network with the MCD single-site approach
For a further investigation comparing the single-site and network approaches of the MCD, archived mesocyclone detections of the 4-yr period 2013–16 (April–September) were analyzed. Mesocyclone objects were restricted to be composed of at least two features.
The severity criteria are the same as described in section 4, except that the rotational velocity criterion is not used—this attribute is not available for all archived mesocyclone detections.
Figure 14a shows the frequency of mesocyclone detections as function of the severity level for the default network and the single-site approach.
Approximately 10%–15% more mesocyclone detections are gained when using the network mode. This gain is roughly evenly distributed across severity levels 1–4. For severity level 5, no detections were gained using the network approach, since mesocyclones detected at this level have a depth of at least 8 km, so multiple features are very likely detected even in single-site mode.
Figure 14b shows the spatial distribution of all mesocyclone detections for the years 2013–16 that are detected only within the network approach. Concerning the single-site mode, limitations of the scan geometry in the area close to a radar site will affect the accuracy of the detected mesocyclones properties (e.g., the vertical extension cannot be measured accurately), but a detection may still be accomplished. On the other hand, at a larger distance from a given radar site, a mesocyclone may be not optimally scanned by this radar alone (as a result of the decreasing resolution and the increasing beam height), so only a single feature or no feature is detected and the whole detection is lost in single-site mode. Thus, it is plausible that the detections gained in network mode are fairly uniformly distributed over Germany. Figure 14c shows the frequency of the gained mesocyclone detections as function of the distance to the nearest radar site. The frequency histogram has a linearly rising slope as is expected from a spatially uniform distribution of the additional detections. The decrease at distances larger than 65–70 km is due to diminishing radar overlap at the radar network border region. The peak at around 70 km is in line with the mean nearest-neighbor distance within the DWD radar network, which is 140 km. Applying the MCD in network mode increases the number of detections particularly for weak meocyclones and therefore is likely to extend the lead time in severe weather warning operation.
For cell lifetime studies, tracks of storm cells as detected by the reflectivity-based cell detection algorithm Konvektive Entwicklung in Radarprodukten (KONRAD; convective development in radar products; Lang 2001; Wapler et al. 2012) were investigated. Only those storms that show at least one associated mesocyclone detection are considered. In Figs. 15a and 15b, absolute and relative numbers of KONRAD tracks with associated mesocyclone detections are shown as a function of the KONRAD track duration (i.e., the estimated cell lifetime). A relative gain of 10%–15% in the number of “mesocyclonic” KONRAD tracks can be deduced for the network approach independent of cell lifetime.
Figure 15c shows the mesocyclone fraction, that is, the percentage of the KONRAD lifetime with attributed mesocyclone detections. The number of KONRAD tracks with very low mesocyclone fraction decreases, while the number of KONRAD tracks with significant mesocyclone fraction increases for the network mode. Figure 15d depicts the mean time between successive mesocyclone detections along a given KONRAD track. For the MCD in the standard network mode, the mean time decreases; that is, mesocyclone detections attributed to a KONRAD track are more closely spaced. Both findings suggest an improvement in the continuity of mesocyclone detections for the network approach, which is also demonstrated by the case studies (section 6).
Figure 16 shows the lifetime frequency distributions of mesocyclonic KONRAD cells. As the severity level of the associated mesocyclone detections increases, the frequency distributions become flatter; from all detections attributed to a given KONRAD track, the maximum severity was chosen. Thus, KONRAD cells with associated mesocyclone detection show a longer lifetime on average: the higher the mesocyclone severity, the longer the cell’s lifetime. The gain of additional mesocyclone detections for the network approach is distributed rather equally among the lifetimes.
Recent life cycle studies with a meteorological rather than an algorithmic aspect can be found in Wapler (2017), where the life cycles of hailstorms, including possible mesocyclone occurrence, are investigated.
8. Summary and discussion
In this paper, the operational mesocyclone detection algorithm (MCD) and the rotation composite algorithm of DWD were presented. These algorithms follow a network approach and process quality-controlled Doppler wind data.
The quality assurance successfully removes dual-PRF dealiasing artifacts while not eliminating mesocyclonic vortices. The MCD with its severity metric shows a good correlation with the rotation tracks.
Because of the dense radar network of DWD, the MCD shows significant additional detections in network mode. Therefore, the tracks of mesocyclone detections become more complete and a future tracking and nowcasting extension of the MCD will benefit from the network approach. In network mode, there are more samples available for a given mesocyclone detection to calculate the mesocyclone properties. This leads to an improvement in severity estimation. Also, the rotation tracks are improved using the network approach as a result of better suppression of shear artifacts and shallow rotation.
The shape and rotational characteristics of mesocyclones over central Europe detected by the presented MCD confirm the results reported in Bellon and Zawadzki (2003) for Quebec. An investigation of a correlation between a reflectivity-based cell detection (KONRAD) and the MCD revealed that KONRAD cells with an additional detected mesocyclone show a longer lifetime on average, whereby the lifetime increases with the severity of the mesocyclone detection. This finding is in line with mesocyclonic rotation being known to occur in highly organized, long-lived storms and confirms the practical significance of the MCD. Further studies suggest a correlation between the occurrence of hail-producing cells and the presence of mesocyclonic rotation (as detected by the MCD). Future cell lifetime studies will benefit from the detection of rotational cell characteristics.
A rigorous statistical analysis of the performance of the MCD requires a sufficiently large radar database of convective storms that are classified into mesocyclonic and nonmesocyclonic storm types (as well as further subdivisions; e.g., tornadic and nontornadic supercells). Using a recently introduced procedure, forecasters at DWD now regularly classify current weather phenomena so that a “case database with categorization” can gradually build up. However, the categorization of cases, done by human experts, is subjective to a certain degree and so are the derived performance measures. Therefore, a meaningful statistical analysis should have a comparative character (comparing performances of algorithms being applied to the same database and using the same criteria for hits and misses).
A proper setting for an algorithm intercomparison study can be found within WMO’s World Weather Research Programme as the example of the Sydney 2000 Forecast Demonstration Project (Keenan et al. 2001) shows (radar-based nowcasting played a major role within the project). In Sydney 2000, participants ran their own algorithms, properly adjusted to the available data (an important aspect to eliminate bias) on a common database. A mesocyclone detection algorithm intercomparison may well be feasible in the future within such a forecast demonstration project.
Quality assurance issues to be tackled by future development include a range correction of azimuthal shear so that rotation tracks and mesocyclone detections will be improved at the composite border regions. Also, azimuthal shear data in overlap areas (of different radar sites) naturally showing mixed resolution will be equalized by the range correction so that these data become better comparable (Newman et al. 2013) and that the MCD and rotation composites are improved.
In the near future, POLARA will include a 3D reflectivity-based cell detection algorithm following the network approach (Werner 2017). A generic tracking module will be available for both cell and mesocyclone detection algorithms so that tracking and nowcasting of 3D cell objects with detected reflectivity and rotation characteristics will be possible. The tracking of MCD objects opens new possibilities concerning the detection of weak mesocyclones. The thresholds as used in the severity calculation may be further lowered for the purpose of early detection. The resulting increase in false alarms can be reduced by requesting that a mesocyclone detected with low thresholds must be part of a track with a certain minimum duration.
In conclusion, the mesocyclone detection and rotation products provide a new perspective at looking at severe convective storms in Germany and adjacent European countries, and provide valuable information to improve the weather and warning service of DWD. The mesocyclone detections are integrated in the NowCastMIX system (James et al. 2013), which provides basic information for AutoWARN (Reichert 2013), the operational warning decision support system of DWD.
The support offered by other members of the radar meteorology R&D team at DWD is greatly appreciated. An evaluation of all DWD radar–based algorithms and products mentioned in this paper at the testbed organized by the European Severe Storms Laboratory (ESSL) in Wiener-Neustadt proved very useful. We are very grateful to Pieter Groenemeijer for the many fruitful discussions and recommendations, in particular concerning the estimation of the mesocyclone severity. Furthermore, this study has benefited from the authors’ discussions with members of the DWD forecasting department concerning the interpretation of the rotation and rotation track composites in particular. Last but not least, we thank the anonymous reviewers for the very valuable comments and recommendations.
Velocity Extension by Dual-PRF Dealiasing Technique
In the dual-PRF scan mode, the radar data are collected with azimuthally alternating PRFs on 1° azimuthally spaced rays (Dazhang et al. 1984). The indices 1 and 2 of the physical quantities given below relate to these two PRFs. For the lowermost six elevations within DWD’s volume scan, the PRFs are = 800 Hz and = 600 Hz, and the wavelength for the DWD radar systems is approximately 5.3 cm. The corresponding maximum unambiguous velocities are = 10.6 m s−1 and = 7.95 m s−1. The dual-PRF dealiasing technique relies on the basic assumption that the velocities (of azimuthally adjacent pixels) to be measured are identical. If this assumption is fulfilled, then an estimate with related extended maximum unambiguous velocity can be derived from the measured velocities υ1 and υ2,
Inserting numbers into Eq. (A1) yields 32 m s−1.
The dual-PRF dealiasing is performed at the radar site in the radar system’s signal processor in real time. Introducing i as azimuth and j as the range index, Eq. (A1) becomes
The estimate as given by Eq. (A2) is not used directly, since both measured velocities have errors that contribute to the error in . Instead, one analyses possible estimates , which represent possible aliasing corrections of the originally measured value and are given by
Here, n is an integer confined by the ratio of the extended maximum unambiguous velocity to the maximum unambiguous velocity of the given ray. Minimizing with respect to n yields and thus the sought for correction is
This procedure will result in data loss if is not available, since two velocity estimates at the same range are required. Thus, a dealiasing algorithm was developed running within POLARA at DWD’s Central Office, where the complete datasets are gathered. Then, is replaced by a phase-averaged mean calculated within from the six neighbor values with and . These six neighbor values with the given indices are related to the same maximum unambiguous velocity .
Complete List of the Properties of Mesocyclone Detections
An overview is given below of all calculated attributes for the mesocyclone objects detected by the MCD (see also Table 2).
a. Geographical coordinates and shape attributes
A digital elevation model with 1-km resolution is used for height above ground layer (AGL) calculations.
Latitude and longitude center coordinates and (°), respectively, are calculated as the shear weighted location of the lowest 3 km (or less if the whole mesocyclone does not extend over 3 km) of the 3D mesocyclone object. Here, all pixels and their shear values within the 3D mesocyclone object region are evaluated.
The mesocyclone base and top (km) are calculated as the height (AGL) of the lowest and highest 2D features, respectively.
The mesocyclone depth (km) is defined as the difference between the mesocyclone top and base heights.
For the mesocyclone diameter (km), the diameter of the largest 2D feature is used, defined by the maximum of the estimates of range and azimuth extensions as given in Zrnić et al. (1985).
The mesocyclone equivalent diameter (km) is calculated as the diameter of the circle that covers the same area as the group of pattern vectors composing the largest 2D feature. The term is used to make elongated features more comparable.
b. Doppler velocity attributes
Doppler velocity attributes are evaluated by considering all pattern vectors that belong to the 3D mesocyclone object. Each 1D pattern vector has an azimuthal shear and a specific angular momentum value. The maximum and mean values are the maximum and the arithmetic mean of all the pattern vectors, respectively. Shear values are given in units of 1 m s−1 km−1 = 10−3 s−1, and specific angular momentum values are measured in units of 1 m s−1 km = 103 m2 s−1.
Maximum and mean azimuthal shear and (10−3 s−1), respectively.
Maximum and mean specific angular momentum and (103 m2 s−1), respectively.
Maximum rotational velocity (m s−1).
c. Reflectivity attributes
Reflectivity attributes are evaluated using all reflectivity information within the 3D mesocyclone object region. For the echo-top height, a 18-dBZ threshold is used.
Maximum and mean reflectivity and (dBZ), respectively.
Cell-based VIL (kg m−2).
18-dBZ echo-top height (km).
VIL density (g cm−3).