Extratropical cyclones are a typical feature of the flow in the midlatitudes with a significant impact on the local weather. Their impact depends in particular on their intensity and their size. It is remarkable that despite their importance there is no accepted, universal definition of a cyclone (Neu et al., 2013). Furthermore, there is no unique definition of vortex properties such as intensity and size. Especially in sheared flow, the determination of these vortex properties is challenging and some common methods yield inconsistent results. In this work, we will apply a kinematic method in order to approach this problem. Kinematic methods have been used successfully in fluid mechanics, for example, in the identification of coherent structures in turbulent flows (e.g. Dubief and Delcayre, 2000). The advantage of such methods is that they distinguish between deformation and rotation in the flow field and integrate the additional information into the vortex identification procedure. Furthermore, this decomposition represents locally a complete description of the flow field. Hence, kinematic methods enable the determination of cyclone sizes and intensities (circulations) in a consistent way especially in regions of strong shear.
We will define the vortex intensity with the help of circulation which is a measure of the influence and importance of a vortex on the general atmospheric circulation (Sinclair, 1997). In addition, the circulation is an integral parameter taking into account the vortex area and therefore better represents the vortex as a whole. Hence, an accurate knowledge of the vortex size is necessary. Although extratropical cyclones have been analysed in numerous studies with emphasis on cyclone activity (for recent reviews, see e.g. Ulbrich et al., 2009; Neu et al., 2013), only a small fraction of these studies deal with the additional analysis of the geometric properties of cyclones. The latter concentrate on the analysis of the size evolution (Grotjahn et al., 1999; Simmonds, 2000; Rudeva and Gulev, 2007; Schneidereit et al., 2010), the interrelationships of extratropical cyclone properties (Nielsen and Dole, 1992) and their vertical structure (Gray and Dacre, 2006; Lim and Simmonds, 2007).
In these studies, cyclone sizes have been determined by traditional methods on the basis of well-known fields such as the pressure/geopotential height fields and the (geostrophic) vorticity fields. Schneidereit et al. (2010) give a detailed review over traditional methods mainly based on pressure/geopotential height fields. They divide the methods into three groups depending on the approach used: (1) based on the derivative of the pressure (e.g. searching for the nearest saddle-point (col) in the pressure field as in Nielsen and Dole, 1992; Rudeva and Gulev, 2007; Rudeva, 2008), (2) based on the determination of the enclosed area (e.g. definition of the outermost-closed isobar in Wernli and Schwierz, 2006) and (3) based on the application of functional fits (e.g. Gaussian fit in Schneidereit et al., 2010). Unfortunately, methods based on pressure fields give inconsistent results in strong ambient flow or in case of two cyclones that are close to each other (e.g. Grotjahn et al., 1999; Wernli and Schwierz, 2006). For these reasons, some studies concentrated on (geostrophic) vorticity fields as a basis of cyclone size determination. Even though the vorticity depends on the spatial resolution of the data (e.g. Ulbrich et al., 2009), the advantage of vorticity fields compared to pressure fields is that the locations of vortex centres are not affected by a strong background flow (Sinclair, 1994). Vortex sizes are determined by searching for the distance where either the vorticity falls to zero or the (radial) vorticity gradient changes its sign (Sinclair, 1997; Simmonds, 2000; Simmonds and Keay, 2000; Lim and Simmonds, 2007). Flaounas et al. (2014) used a fixed vorticity threshold (3·10^{−5}s^{−1}) for the identification of vortex sizes in low-level (850 hPa) vorticity fields. However, these traditional methods fail to capture vortex sizes properly in some flows (for a detailed discussion on the inadequacy of traditional methods, see e.g. Jeong and Hussain, 1995). While methods based on pressure fields are inadequate in the presence of strong background flows (e.g. Sinclair, 1994), vorticity is connected to the rotation in the flow field. The problem of vorticity-based methods is that vorticity alone cannot distinguish between sheared and curved flow. Jeong and Hussain (1995) state that the choice of fixed vorticity thresholds is therefore subjective and that in the presence of strong shear, a high threshold might misrepresent the vortex core.
In contrast to traditional methods, kinematic methods distinguish between deformation and rotation in the flow. Mathematical basis is the analysis of the velocity gradient tensor ∇v and its invariants. The velocity gradient tensor can be decomposed into a symmetric component that describes the deformations (strain-rate tensor S) and an antisymmetric component that describes the rotation of the flow (vorticity tensor Ω). The advantage of these variables is their invariance under rotations or translations of the coordinate system which is necessary in the definition of consistent vortex sizes. While vorticity is also an invariant of the flow field, the consideration of vorticity alone only describes the antisymmetric part of the flow while the symmetric part is not taken into account such as it is done in kinematic methods. Truesdell (1953) introduced the kinematic vorticity number ${W}_{k}=(W(/(\mathbf{S}$ as the ratio of the local rate-of-rotation ∣∣Ω∣∣ and the local rate of strain ∣∣S∣∣ considering both parts of the flow. Vortex areas are identified as regions where the rotation prevails over the deformation (W_{k}>1). In a similar manner with the difference that W_{k} is a dimensionless number, the Okubo-Weiss parameter (or Q-method) describes the absolute value of the excess/deficit of rotation over deformation for incompressible flows (∇·v=0) by $Q=(\mathrm{\Omega}(-(\mathbf{S}($. As one of the rare applications of kinematic methods to large-scale atmospheric flows, the Okubo-Weiss parameter has been applied successfully in tropical cyclone studies (e.g. Dunkerton et al., 2009; Tory et al., 2013). It should be noted that multiple more kinematic methods exist (for a comparison of different methods, see e.g. Jeong and Hussain, 1995; Chakraborty et al., 2005).
The main topic of this work is the introduction of a kinematic vortex size determination method based on the kinematic vorticity number W_{k} to atmospheric flows (Section 2) and its comparison to traditional methods regarding their ability to identify mid-latitude cyclone sizes (and volumes) in various idealised (Section 3) and real (Section 4) flow situations. It will be shown that traditional methods fail in identifying cyclone sizes in certain flow situations in a consistent manner, while the kinematic method based on W_{k} is capable of extracting vortex structures in a consistent way. This work also deals with the question what constitutes a vortex and what part of the vortex is extracted by the different methods. We will summarise our main findings and conclusions in Section 5.
Most of the vortex identification methods of turbulent flows are based on the local evaluation of the flow field with the help of a velocity gradient tensor ∇v and its invariants. The velocity gradient tensor calculated at a point r_{0} gives information about the structure of the infinitesimal flow field surrounding that point. This can be seen by a first-order Taylor series expansion of the velocity (e.g. Fortak, 1967; Batchelor, 2000):
The kinematic vorticity number W_{k} introduced by Truesdell (1953) is defined as the ratio of the tensor norms of Ω and S:
Considering the sign of vorticity will slightly modify the equation with the advantage to study vortices of positive and negative vorticity (or cyclones and anticyclones) likewise. This extended kinematic vorticity number W_{k} is given by:
There are multiple more kinematic vortex identification methods besides the W_{k}-method used in this work. However, in planar flows the methods coincide (Jeong and Hussain, 1995). Still, we decided in favour of the W_{k}-method over the other methods (e.g. Okubo-Weiss parameter) since it – as a dimensionless number – allows a comparison of vortex structures relative to the background deformation (or shear). We think this is the main advantage of the W_{k}-method compared to other kinematic methods. In addition, the ${\overline{W}}_{k}$number averaged in the region of a vortex represents a kinematic circulation and its value gives the resemblance of a vortex with a rotating solid object: The larger this value is the larger the equivalence to a solid body. This improves the conceptional understanding of atmospheric flows.
In this section, we will compare different size determination methods in idealised 2-D set-ups. After giving an overview over traditional cyclone size determination methods in Section 3.1, we specify the methods used for comparison reasons in the idealised cases in Section 3.2. The experimental set-ups are described in Section 3.3, followed by a presentation of the results (Section 3.4) and a discussion (Section 3.5).
Nielsen and Dole (1992) were probably the first who derived statistics on the sizes of synoptic cyclones in surface pressure data. In their work, Nielsen and Dole (1992) discussed different possible measures of cyclone sizes, namely the distances between the nearest (1) high, (2) low, or (3) col (saddle point) of sea level pressure (SLP) and (4) the horizontal area enclosed by the outermost-closed isobar around a low-pressure centre. As definitions (1) and (2) fail in the case of cyclone families and lee cyclogenesis, respectively, they concentrated on definition (3) in the analysis of surface weather maps. A col is given when the pressure gradient falls to zero. Rudeva and Gulev (2007) and Rudeva (2008) defined a col slightly different as the point where the radial pressure gradient falls to zero. They searched along radial lines starting from the cyclone centre for that point and defined the outermost-closed isobar as the pressure value at the nearest zero radial pressure gradient found. Their method also allows to study the asymmetry of cyclones. The outermost-closed isobar can be defined to include single centres (Wernli and Schwierz, 2006) or multiple centres (Hanley and Caballero, 2012). Functional fits to the surrounding field of a cyclone have been applied by Grotjahn et al. (1999), Grotjahn and Castello (2000) and Schneidereit et al. (2010). Grotjahn et al. (1999) applied Mexican hat fits to the longitudinal and latitudinal directions around cyclone centres. Grotjahn and Castello (2000) used a circular average of the geostrophic kinetic energy in order to determine cyclone sizes. Under the assumption of azimuthal symmetry, Schneidereit et al. (2010) applied a 1-D Gaussian fit to the radial geopotential height distribution surrounding the cyclone.
The (geostrophic) vorticity field as basis of cyclone size determination was favoured by, for example, Sinclair (1997); Simmonds (2000); Simmonds and Keay (2000); Lim and Simmonds (2007). Even though the vorticity depends on the spatial resolution of the data (e.g. Ulbrich et al., 2009), the advantage of vorticity fields compared to pressure fields is that it is possible to detect vortex centres even in a strong background flow (Sinclair, 1994). Vortex sizes are determined by searching for the distance where either the vorticity falls to zero or the radial vorticity gradient changes its sign similar to the nearest col definition in the pressure field. This search is done along radial lines (e.g. Sinclair, 1997; Lim and Simmonds, 2007) or along the directions of maximum (negative) gradient similar to the definition of a water catchment boundary (Simmonds and Keay, 2000). However, it is possible that a cyclone is embedded in an elongated vorticity streamer of a jet streak so that the vorticity as well as the radial (or tangential) vorticity gradient will neither fall to zero nor change its sign. Therefore, Sinclair (1997) restricted the change in distances between neighbouring search lines.
In a first step, local maxima (minima) in the 2-D vorticity (pressure) field are identified. A local maximum (minimum) is found when the eight surrounding points have lower (higher) values than the central point. In a second step, the following four size estimation methods are applied:
For methods (1), (3) and (4), contour lines are calculated by a standard contouring function. The area A is calculated by the sum over all triangle areas formed by two neighbouring contour points and the vortex centre. By assuming that the area A is circularly distributed around the centre, the system's radius R is calculated as $R=\sqrt{A/\pi}$. In method (2), the 2-D pressure field surrounding the low-pressure centre is mapped to a 1-D radial distribution: In a first step, the surrounding pressure distribution is interpolated on 36 radial lines (every 10°) starting from the vortex centre up to 1000 km (every 50 km). In a second step, the mean of the 36 pressure values for each distance r is determined. Finally, a Gaussian fit is applied to the resulting pressure distribution with a gnuplot fitting procedure, that fits the following function to the 1-D distribution: ${p}_{Gauss}(r)=a\xb7exp(-{(r-{r}_{0})}^{2}/2{b}^{2})$, where a gives the pressure drop and b represents the radius which is equal to the standard deviation of the Gaussian distribution. The vortex centre is located at r_{0}.
The four size estimation methods are tested and compared in different idealised set-ups. The aim of these tests is to find out, how well the different methods perform in re-extracting the predefined vortex sizes from various flow fields.
In the idealised test cases, the pressure field p will be predefined. Geostrophic wind v_{g} and geostrophic vorticity fields ζ_{g} are calculated from the pressure field by
A low-pressure disturbance defined by a 2-D Gaussian distribution with intensity Δp=5 hPa and radius R=250 km,
The superposition of two low-pressure systems on a flat pressure field of 1000 hPa is investigated for varied distances between the two centres. The pressure disturbances of the lows (${p}_{1}^{\mathrm{\star}}$, ${p}_{2}^{\mathrm{\star}}$) are given by two 2-D Gaussian distributions of different intensities (Δp_{1}=10 hPa, Δp_{2}=2.5 hPa) and different sizes (R_{1}=250 km, R_{2}=160 km) calculated by eq. (14). The total pressure field is given by $p=1000-{p}_{1}^{\mathrm{\star}}-{p}_{2}^{\mathrm{\star}}$. The first low indicated by index 1 is fixed at the location (x_{0}, y_{0})=(0,0). Low 2 changes its position stepwise in southwesterly direction starting at the location of low 1 (or a distance of 0 km) up to a distance of about 1400 km in 70.7 km steps (50 km to the south/50 km to the west; see Fig. 2 for set-up and two examples). The resolution of the calculated fields is 10 km.
In this test case, a low-pressure disturbance p* with Δp=5 hPa and R=250 km (see eq. 14) is superposed by a jet streak pressure-gradient p_{jet} on a flat pressure field of 1000 hPa. For changing distances between the low centre and the jet axis, the size of the low-pressure disturbance is determined by the different methods and compared to the original R. The jet streak's pressure gradient is calculated by a Gaussian error function (abbreviated by erf). The jet axis is oriented in the west–east direction. The (south to north) pressure profile is given by ${p}_{jet}=\mathrm{\Delta}{p}_{jet}\xb7\text{erf}((y-{y}_{Gauss})/(\sqrt{2}\sigma ))$ where ${y}_{Gauss}$ is the position of the jet axis and Δp_{jet}=7.5 hPa gives the pressure difference between the edges of the jet and the jet axis. As the geostrophic wind is proportional to the pressure gradient, the wind field associated with the jet is Gaussian distributed with a standard deviation of σ=350 km. The total pressure field in hPa is given by p=1000–p*–p_{jet}. The low is fixed at the location (x_{0}, y_{0})=(0,0). The position of the jet axis moves stepwise (50 km steps) from 1400 km south to 1400 km north of the low centre (see Fig. 3 for set-up and examples).
Note, subscripts of the radii R correspond to the method, for example, G for the Gaussfit-method, W for the W_{k}-method and so on. Multiple letters are used when radii coincide.
The vortex core radius identified by the W_{k}-method and by the Gaussian fit (not shown) coincides with the wind maximum at a radius of R_{W,G}=250 km (blue (wind), black(W_{k}) curves in Fig. 4) which is equal to the predefined radius. W_{k} equals 1 when rotation and deformation are of the same size: vorticity and deformation distributions cross at R_{WG}=250 km and converge far away from the vortex centre (red (ζ), yellow (deformation) curves in Fig. 4). The deformation peaks outside of the vortex area identified by the Wk- and Gaussfit-method with a lower, broader peak than the vorticity. Vorticity-dominated and deformation-dominated areas are adjacent regions in the vortex. The p-method identifies the largest radius (about 450 km, green line in Fig. 4) when the outermost-closed isobar is determined by increments of 1 hPa. A finer increment of 0.1 hPa determines a larger radius of about 700 km.
The splitting of the two systems, that is, the identification of two single instead of one system, occurs at a smaller distance in the vorticity field (for zeta-/W_{k}-method (red/blue curve in Fig. 5) at around 420 km which coincides approximately with the sum of the undisturbed radii R_{1}=250 km, R_{2}=160 km) compared to the pressure field (for p-/Gaussfit-method around 770 km; green/yellow curve in Fig. 5). With growing distances the radii of both systems increase until the values stabilise (around 900 km for the p-, Gaussfit-method and 700 km for the ζ-, W_{k}-method). The stepwise increase of R_{p} can be attributed to the coarse increment of 1 hPa for the contour lines since this behaviour is not observed for finer increments (not shown). While the p-method and ζ-method show strong variations in the vicinity of the splitting point, the W_{k}-method and the Gaussfit-method show only slight variations and otherwise coincide with the predefined radii.
When the jet axis is in the vicinity of the low centre, the methods based on pressure (p-method/Gaussian fit; green/yellow curve in Fig. 6) show strong variations and a lack of identification for distances between 0 and 500 km; while the W_{k}- and ζ-method are not strongly affected (blue/red curve in Fig. 6). For the application of the Gaussfit-method a local pressure minimum is needed^{4} . When the Gaussfit-method is modified such that the fit is applied to the pressure field surrounding the local vorticity maximum, radii can be identified over all distances (dashed yellow line in Fig. 6). However, the variations are strong when the jet axis is close to the vortex. The W_{k}-method reproduces the predefined radius of the low with slight variations. The radius identified by the ζ-method is proportional to that identified by the W_{k}-method, although the outermost-closed vorticity contour is not zero when the jet axis is south and near the vortex centre (not shown).
We have seen in the previous section that some methods are not capable in identifying the cyclone radii in particular flow situations. For example, near the splitting point of two lows, the ζ- and the p-method showed strong variations. In some asymmetric fields caused by the superposition of a low and a jet, the Gaussfit-method is strongly affected and the p-method partially fails to identify the cyclone. We will now discuss the reasons for the failure and what part of the vortex is seen by the different methods.
The part of the vortex that is identified by the outermost-closed isobar strongly depends on the flow situation and on the contour value/increment. This is in accordance with Wernli and Schwierz (2006) who observed an increase of 40% (decrease of 30%) of detected cyclones by a reduction (increase) of the contour increment from 2 to 1 hPa (4 hPa). Likewise to streamlines, the p-method only represents a snapshot of the flow at a certain timestep. This can be very different in various flow situations and from one timestep to another. As a result, the area of a cyclone is only poorly represented by pressure/geopotential height contours. This is especially important when investigating mobile and developing systems, respectively (Sinclair, 1994).
In the undisturbed case, the Gaussfit-method coincides with the wind maximum and the maximum of the radial pressure gradient, respectively (see also Schneidereit et al., 2010). Even in the case of the superposition of two lows, the Gaussfit-method nicely reproduces the predefined radii. This result was expected because the predefined low-pressure disturbances were already Gaussian distributed and the asymmetry of the pressure field surrounding a pressure minimum is only minor. On the other hand, the superposition of a jet and a low involves much more asymmetry (ideal test case 2 in Section 3.4.3). In this case, the Gaussfit-method fails in re-extracting the radius when all surrounding points are considered even though the predefined low-pressure disturbance was originally Gaussian distributed.
The vorticity can be split up into shear and curvature vorticity. In the undisturbed case, cyclonic curvature exists in the whole domain. At the wind maximum the shear vorticity changes its sign resembling the flow situation at a jet axis (Fig. 7), but curvature vorticity is still positive. Vorticity becomes zero when both parts are balanced. In the disturbed ideal test cases 1 and 2 the observed outermost-closed vorticity contour is partly different from zero. Hence, a fixed threshold would fail: either it would only identify strong vortices whose intensity might not be comparably strong since the background shear is misleading or it would only identify undisturbed systems, neglecting vortices embedded in shear. If no restriction to a fixed vorticity threshold is made, R_{ζ} changes approximately proportional to the W_{k}-method and it seems to be an alternative to that method. On the other hand, it is not easy to interpret which part of the vortex is then extracted. Here, an interpretation in terms of shear and curvature vorticity is difficult. It is even more complicated when the (contour) threshold changes along certain directions as is done in Sinclair (1997) and Lim and Simmonds (2007) who determined the boundary of a vortex when either the vorticity is zero or the radial gradient of vorticity changes its sign along a set of radial lines. That definition can lead to a zero contour in one direction and a different non-zero value in another direction for the same system.
The kinematic vorticity number is larger than one (W_{k}>1) when the rotation prevails over the deformation at a point and it is exactly one in the case of a pure shearing motion. In case of an idealised cyclone, it can be seen that for a point located at the radius of maximum wind its neighbouring wind field resembles a pure shearing motion and therefore the W_{k}=1 contour coincides with the radius of maximum wind (Fig. 7). In order to display the meaning of W_{k}>1 and W_{k}=1 in a more non-trivial case, we plotted the local flow field around a point (blue streamlines in Fig. 8) at the boundary (defined by W_{k}=1, thick black contours in Fig. 8) of a vortex and inside the vortex (defined as W_{k}>1) in case of the superposition of two lows (see Fig. 8). The local point at the boundary (point 1, Fig. 8a) is embedded in a shearing environment. Particles that at first are near to that point separate rapidly following the streamlines. In contrast, the local point inside the contour (point 2, Fig. 8b) is surrounded by particles that stay in the vicinity of that point moving in spirals or closed circles around that point.
Summarised, particles inside the W_{k}=1 contour stay close to each other, that is, mass is accumulated inside the vortex. Therefore, the part of the vortex identified by the W_{k}-method can be interpreted as a vortex core. This statement is also supported by a calculation of the positive vorticity concentrated inside the W_{k}=1 contour relative to the total positive vorticity. In the undisturbed reference case about 84% of the positive vorticity is concentrated inside the vortex core (inside the W_{k}=1contour). However, it should be noted that the area influenced by the vortex can be much larger than the area of its core.
After presenting details on the reanalysed data used, we will apply the W_{k}-method in a real storm case example and compare the results to traditional methods. Furthermore, we will present some statistics of midlatitudal cyclones of the northern hemisphere.
The data used for the analysis are the geopotential height and the horizontal wind fields of the National Centers for Environmental Prediction (NCEP)/National Center for Atmospheric Research (NCAR) Reanalysis provided by the Research Data Archive of the National Centers for Environmental Prediction National Weather Service NOAA US Department of Commerce (1994). These data are available four times per day on a regular 2.5°×2.5° grid (Kalnay et al., 1996). We analysed the geopotential height data on 12 pressure levels from 1000 to 100 hPa (100, 150, 200, 250, 300, 400, 500, 600, 700, 850, 925, 1000 hPa) for the northern hemisphere winter months (December, January, February) for the years 1999/2000 (abbreviated by DJF 1999/2000).
For every 6-hourly timestep in the period DJF 1999/2000, geostrophic ${W}_{k}^{\mid}$-fields were computed from the derivations of the geostrophic wind fields on each pressure level with the help of (11). The derivations were calculated as central differences omitting the poles. We further restricted the analysis to the northern hemisphere and a latitudinal band between 30°N and 80°N (including these latitudes). No terrain filtering was used in the lower levels. The geostrophic wind fields v_{g} were derived from the geopotential height fields Φ by ${\mathbf{v}}_{g}={f}^{-1}\mathbf{k}\times \nabla \mathrm{\Phi}$ with Coriolis parameter f=2Ω sinφ where $\mathrm{\Omega}=2\pi /\text{day}=7.2921\xb7{10}^{-5}$s^{−1} is the rotation rate of Earth and φ is the latitude. Every grid point that yield $\mid {W}_{k}^{\mid}\mid >1$ was set to 1, every point with $\mid {W}_{k}^{\mid}\mid \le 1$ was set to zero. In this way, we derive a vortex patch field which cuts out the vortex structures.
After calculating the ${W}_{k}^{\mid}$ fields as described above, simply connected regions^{5} of ${W}_{k}^{\mid}>1$ (positive circulations/lows) and of ${W}_{k}^{\mid}<1$ (negative circulations/highs) were separately identified in each field. We will call a single simply connected region of W_{k}>1 a W_{k} feature. Note that a W_{k} feature can include multiple vorticity centres and therefore rather represents a large-scale circulation area (or cyclone family) in such cases. Therefore, we additionally analysed single vortex centres including single vorticity extrema. Such single centres were determined by the outermost-closed vorticity contour enclosing only one vorticity centre. Note, that the area/circulations of the single centres are smaller or equal in total than that calculated for the W_{k} feature. In order to account for broad extrema we further restricted the minimum distance between isolated vorticity extrema inside the same W_{k}>1 region; if two systems are closer than 600 km (≈5° latitude = twice the resolution), they were considered as a single broad centre. Then the outermost-closed vorticity contour around both centres was calculated by a standard contouring method.
W_{k} features/single centres are composed of a set of grid points. Each grid point is associated with a grid box area of 2.5°×2.5°. Note that this area depends on the latitude. The sizes of the W_{k} features/single centres were determined by the sum of all grid box areas associated with the feature/centre, and corresponding circulations Γ_{total} were calculated by (12). The radius of a system was determined as $R=\sqrt{A/\pi}$ under the assumption that the area belongs to a circular system (likewise to the effective radius definition in Rudeva and Gulev, 2007). Furthermore, we calculated the centre of circulation of each feature/centre by $C=\sum _{i}{\mathrm{\Gamma}}_{i}{\mathbf{x}}_{i}/{\mathrm{\Gamma}}_{total}$ where the summation is done over each grid box included in the simply connected region/outermost-closed contour (cf. Müller et al., 2015). Γ_{i} is the circulation associated with the ith grid point and x_{i} is the coordinate vector of this grid point.
Only in the real storm case example (Section 4.3), we additionally did a temporal and vertical tracking for the explicit storm Anatol. Anatol was traced over its lifetime by manually connecting the appropriate storm centres of successive timesteps on the 1000 hPa level. The vertical tracking of the storm was done following the work of Lim and Simmonds (2007) by a numerical method that searches for the nearest vortex centre in superposed vertical levels starting from the 1000 hPa level. A vertical connection between two centres in neighbouring pressure levels was confirmed when the distance between those centres was less than about 340 km. This distance accounts for a diagonal vertical tilt (north–west/–east, south–west/–east) from about 50° latitude polewards (the diagonal distance between grid points in 50° latitude is about 330 km which further decreases polewards).
For the statistics, all identified W_{k} features/single centres of positive circulation per timestep are taken into account irrespective of their temporal evolution. We will analyse the frequency distributions concerning the radius with a box width of 50 km for systems with radii larger than 200 km. This allows a comparison to the existing literature such as Golitsyn et al. (2007) and Schneidereit et al. (2010). However, the identified absolute circulations cover several orders in magnitude complicating the definition of linear box widths. Therefore, we will compute complementary cumulative distributions for the analysis of the circulations of lows and highs. These distributions are statistically more stable and were already successfully used in the analysis of cyclone/anticyclone kinetic energies in Golitsyn et al. (2007).
The capability of the W_{k}-method to identify cyclones even in the upper troposphere and in high-shearing situations is tested exemplarily in a real winter storm case. The investigated example storm – known as storm Anatol in Germany – occurred from 2 to 4 December 1999 (see Fig. 9 for Anatol's track). The lowest observed pressure was 953 hPa recorded at 3 December 1999 18 UTC near the north-east coast of Denmark (Jutland; see Ulbrich et al., 2001). Anatol hit Denmark and northern Germany at the afternoon and evening of 3 December 1999 with gusts up to 50 m/s. It was one of three extreme storm events affecting Europe in December 1999 (Ulbrich et al., 2001), and it was among the costliest European winter storms between 1980 and 2013 (NatCatSERVICE of Munich Re, 2014). It caused a record storm surge at the Danish and German North Sea coast (Ulbrich et al., 2001). Furthermore, the storm was associated with a strong jet in the middle and upper troposphere and is therefore a challenging situation for size estimation methods.
The temporal development of Anatol reveals its rapid intensification from a wave-like structure over the North Atlantic to a mature cyclone in less than 24 hours (Fig. 10). Note, that the contours of the geopotential in 1000 and 500 hPa show a rather wave-like pattern and the low-level (850 hPa) vorticity is rather weak in the beginning (2 December 1999 12 UTC, see Fig. 10a). In addition, the vorticity centre is embedded in a large-scale cyclonic (positive) vorticity field. The black arrow and the white cross in Fig. 10 indicate the positions of the vortex centres near the surface and in the upper troposphere, respectively. During the intensification period, the system is strongly baroclinic (see Fig. 10b and c). Lower- and upper-level centres become aligned during maturity (see Fig. 10d at 4 December 1999 00 UTC). In comparison, the extended W_{k}-field considering the sign of vorticity isolates nicely the cyclonic and anticyclonic vortices from the continuous (vorticity) fields (Fig. 11). The storm centre can easily be detected (see black arrows). The intensification of Anatol is also mirrored in the maximum W_{k}-value inside the storm, which is 2 in the beginning (2 December 1999 12 UTC, Fig. 11a), 12 hours later W_{k}≈10 and 24 hours later W_{k}≈15 implying that the rotation is 15 times larger than the deformation inside the storm. Note, however, that the maximum of the W_{k}-value does not always coincide with the vorticity maximum or minimum.
The vertical structure and development of storm Anatol as seen from the W_{k}-method's perspective (Fig. 12a–d, top row) support that the system is only shallow in the beginning (black arrow, Fig. 12a) but rapidly intensifies due to the interaction with an upper-level vortex (white cross, Fig. 12) that leads to strong stretching of the vortex (Fig. 12b and c). The strongest baroclinic tilt is observed at 3 December 1999 00 and 12 UTC (Fig. 12b and c). In addition, horizontal interactions between vortex centres are observed, for example, see the interaction of the Icelandic low (centre located at about −20°W, 65° N) with storm Anatol at the 1000 hPa level. The area of the Icelandic low is deformed over time so that it appears to rotate around storm Anatol (Fig. 12b and c) and later it follows Anatol (Fig. 12d).
The evolution of Anatol's circulation over its lifetime shows a rapid intensification over the first six timesteps (36 hours) by one order in magnitude from about 10^{7} to 10^{8} (Fig. 13a). Compared to that rapid intensification at the beginning, the circulation dissipated much slower by a nearly constant gradient of about 264 m^{2}/s^{2} after reaching its maximum at 3 December 1999 18 UTC/4 December 1999 00 UTC (Fig. 13a). Simultaneously, the area on the 1000 hPa level broadens over Anatol's lifetime nearly constantly (see Fig. 13b). After about 13 timesteps at 5 December 1999 18 UTC (78 hours after initiation), the vertical connection between lower- and upper-level vortex becomes less organised as can be seen by the drop of the vertical means relative to the vortex characteristics identified in 1000 hPa. At the end of Anatol's lifetime the connection between the vertical levels is less pronounced.
In order to compare the W_{k}-method's view on Anatol with the ζ-method's perspective, vorticity isosurfaces of 1, 3, 5·10^{−5}s^{−1} of the geostrophic vorticity field are plotted in Fig. 12e–h (bottom row). The main difference between the fields is obvious at upper levels; due to the stronger shear in the upper levels, the vorticity centres are rather embedded in regions of positive vorticity than clearly separated. This complicates a study of upper- and lower-level vortices by means of fixed threshold values of vorticity alone. By fixing the threshold to a value (e.g. 3·10^{−5}s^{−1}), lower-level features – especially during formation – would not be detected since the vorticity is too small as in Fig. 12a and e. Flaounas et al. (2014) use a fixed threshold of ζ=3·10^{−5}s^{−1} applied to the 850-hPa level vorticity fields of the ERA-Interim data set which has a higher horizontal resolution of 1.5°×1.5° than NCEP. Flaounas et al. (2014) reasoned that this value is adequate even in the initial stage of the cyclone development at this specific level. However, in the coarsely resolved NCEP data set used in our analysis, the formation of Anatol would have been missed especially near the ground. A fixed vorticity threshold needs to be carefully chosen for each height level because vorticity magnitudes generally increase with height due to an increase in shear with height. Still, it is not clear if this mixture of thresholds gives a consistent measure of the extent of a cyclone and if different thresholds in different levels lead to comparable sizes. However, an adjustment of thresholds is not necessary when the W_{k}-method is used since it relates the rotation to the background deformation (and shear) and therefore separates the relevant parts of the vortices (i.e. the vortex cores) from the rest of the flow field.
Campa and Wernli (2012) used a different approach in order to study the vertical distribution (and interaction) of potential vorticity (maxima) in cyclones. They used a fixed radius of 200 km around a surface cyclone centre (SLP field) and analysed the vertical column of air limited by this area. They restricted their analysis to surface cyclone centres at the moment of maximum intensity (lowest core pressure during lifetime) and up to 24 hours before that maximum was reached. Campa and Wernli (2012) found the radius of 200 km to fit best their needs as a compromise between a radius that is too small and possibly misses upper-level features (100 km) and a radius that leads to too much averaging (300 km). Especially during the development of a cyclone, the system is usually (strongly) tilted and a fixed radius might miss the upper-level features (more on the relationship between tilt, forcing and cyclogenesis, as well as a classification of cyclones concerning those parameters can be found in Gray and Dacre, 2006). An advantage of the W_{k}-method in such an analysis is that it allows to account for the vertical tilt when the area would be limited to the vortex tube surrounding the cyclone axis. Compared to the method of Campa and Wernli (2012) who use a more Eulerian perspective connected to the Lagrangian tracing of surface centres, the W_{k}-method would describe a rather Lagrangian perspective following the whole vortex tube. A vortex tube is defined as a surface composed of vortex lines that has a constant circulation for every cross-section at an instant of time. However, the circulation can change over time due to, for example, baroclinic production. That the W_{k}-method detects approximately a vortex tube can be seen in Fig. 13a where the circulation at the 1000 hPa level is compared with the vertical mean circulation over the identified vortex height. The identified circulations are approximately similar like in a vortex tube.
Furthermore, we have seen that the W_{k}-method can visualise horizontal vortex interactions (e.g. the interaction of Icelandic low and storm Anatol). Although we did not focus on them in this work, we have studied successfully the horizontal interaction of low- and high-pressure systems at the 500 hPa level in Müller et al. (2015) where we introduced a pattern recognition technique based on the W_{k}-method in order to determine the circulations and locations of vortices in omega-blocking situations using point vortex equilibria.
We have seen that the W_{k}-method is able to extract vortex structures even in upper levels of the atmosphere in a real winter storm case in the previous section. In order to gain even more confidence in the results obtained by the W_{k}-method, we will compare the results of the identified cyclones with existing statistics.
On average about 41 ±6 W_{k} features occur at the 1000 hPa level and a smaller number of about 30 ±5 W_{k} features at the 600 hPa level. The number of single centres is only a bit larger (≈46 ±7 in 1000 hPa vs. ≈39 ±7 in 600 hPa). The numbers in lower and upper levels seem to be correlated over long time periods (Fig. 14). That more systems are detected at the 1000 hPa level might be related to the fact that near the surface more disturbances of the geopotential height field are initiated due to topography and friction. Interestingly, the occurrence of Anatol at the beginning of December is connected with a rather low number of cyclones compared to the rest of the plotted period (Fig. 14). Likewise, a rather low number of W_{k} features is observed at the end of December where two intensive storms hit Europe (storms Lothar and Martin; see Ulbrich et al., 2001, for more details on the storms). With this low number of events it cannot be clarified, if this connection between intense storms and low total numbers of cyclones is only random. Future work is necessary.
For the relative frequency distributions concerning the radii, two general observations can be given (see Fig. 15): (1) The majority of the systems is subsynoptic with radii smaller than 1000 km in both levels and (2) systems in the upper level tend to be larger than the systems at the lowest level. At the 1000 hPa level a broad peak occurs at radii around 300–500 km for W_{k} features as well as for single centres (solid lines in Fig. 15). This peak is shifted and sharpened to larger radii at the upper levels (sharper peak around 400–700 km; dashed lines in Fig. 15). Especially, the W_{k} features can be very large at the upper level reaching synoptic scale, while only a small number of single centres reach radii larger than 1000 km at the 600 hPa level. The observation that the majority of the radii are subsynoptic is in accordance with the literature, that is, Schneidereit et al. (2010) observed the highest frequency between 300 and 500 km at the 1000 hPa level (Gaussfit-method). While methods based on pressure usually show larger radii and less systems per timestep (e.g. Rudeva and Gulev, 2007, observe around 14–20 cyclones with an effective radius of about 600 km). However, it should be noted that the W_{k} method identifies vortex cores. Therefore, the total area influenced by the vortex can be considerably larger as was seen in the idealised cases in Section 3.
W_{k} features that are, in general, larger than single centres also have higher circulation magnitudes than single centres (Fig. 16). Furthermore, W_{k} features in the upper level are more intense reaching higher values of circulations. However, the circulation distributions of single centres in both levels are nearly equal (yellow lines in Fig. 16). Note that the curves decrease nearly exponentially indicating the existence of a characteristic scale of circulation of the vortices. The majority of the systems has circulations of the order of 10^{7} m ^{2}/s which is in accordance with the results in Sinclair (1997). Only about 1% of the single centres at the 1000 hPa level reach circulations of more than 1·10^{8} m^{2}/s. At its maximum intensity, Anatol reached a circulation of about 0.9·10^{8} m^{2}/s. An estimate from synoptic-scale characteristic values of velocity (U=10 m/s) and radius (R=1000 km) leads to a circulation of approximately Γ≈2πRU=6·10^{7} m^{2}/s which is in accordance with our observations, too.
In this work we used a kinematic method (W_{k}-method) to extract cyclone sizes from different (idealised and real) flow situations and compared the results to traditional methods. The W_{k}-method relates the rotation to the (background) deformation of the flow field and therefore has several advantages compared to the traditional methods in the real and idealised data. Precisely, the kinematic vorticity number W_{k} is defined as the ratio of the rates of rotation and deformation. In contrast to absolute measures, W_{k} – as a dimensionless number – describes the excess or deficit of rotation relative to the deformation (including divergence). A W_{k} value of 1 is related to the balance of rotation and strain rate (i.e. a pure shearing motion). We identified a vortex area as a connected region where the rotation prevails over the deformation, that is, where W_{k}>1. The main findings of our work are as follows:
In summary, the W_{k}-method seems to be a promising tool for the determination of vortex properties and the study of vortex interactions. So far, we applied successfully the 2-D definition of the W_{k} number in different vertical atmospheric layers which is sufficient for (quasigeostrophic) synoptic-scale systems. How well the method applies to smaller scale vortices where it may be necessary to use the 3-D version of W_{k} based on the 3-D velocity gradient tensor will be seen in future work. Furthermore, it is possible to reduce or increase the W_{k} =1 threshold in order to study either early circulations in the genesis state or very intense ones compared to the background deformation. Further investigations in this direction are planned.
^{2}^{2}The isotropic expansion in a 2-D flow equals a change of the area in the coordinate directions with a rate of one half of the total divergence. It is zero in case of an incompressible fluid.
^{3}^{3}which is given as $(A(=\sqrt{Tr(A{A}^{T})}$ for any tensor A (e.g. Kunnen et al., 2010)
^{4}^{4}That is slightly different from the method of Schneidereit et al. (2010) which only needs the average gradient of geopotential height to exceed a certain threshold.
^{5}^{5}Two points are simply connected if they are direct neighbours in either north, south, east or west direction.
Lisa Schielicke was funded by the Helmholtz graduate research school GeoSim. We would like to thank two anonymous reviewers for their helpful comments.
BatchelorG. K. An Introduction to Fluid Dynamics. 2000; Cambridge Mathematical Library. Cambridge: Cambridge University Press. 14th ed.
CampaJ., WernliH. A PV perspective on the vertical structure of mature midlatitude cyclones in the northern hemisphere. J. Atmos. Sci. 2012; 69(2): 725–740.
ChakrabortyP., BalachandarS., AdrianR. J. On the relationships between local vortex identification schemes. J. Fluid Mech. 2005; 535: 189–214.
DubiefY., DelcayreF. On coherent-vortex identification in turbulence. J. Turbul. 2000; 1(1): 11.
DunkertonT., MontgomeryM., WangZ. Tropical cyclogenesis in a tropical wave critical layer: easterly waves. Atmos. Chem. Phys. 2009; 9: 5587–5646.
FlaounasE., KotroniV., LagouvardosK., FlaounasI. CycloTRACK (v1.0) tracking winter extratropical cyclones based on relative vorticity: sensitivity to data filtering and other relevant parameters. Geosci. Model Dev. 2014; 7(4): 1841–1853.
FortakH. Vorlesungen über Theoretische Meteorologie – Kinematik der Atmosphäre. 1967; Berlin, Germany: Institut für Theoretische Meteorologie der Freien Universität Berlin.
GolitsynG. S., MokhovI. I., AkperovM. G., BardinM. Y. Distribution functions of probabilities of cyclones and anticyclones from 1952 to 2000: an instrument for the determination of global climate variations. Doklady Earth Sci. 2007; 413(1): 324–326.
GrayS. L., DacreH. F. Classifying dynamical forcing mechanisms using a climatology of extratropical cyclones. Q. J. R. Meteorol. Soc. 2006; 132(617): 1119–1137.
GrotjahnR., CastelloC. A study of frontal cyclone surface and 300-hPa geostrophic kinetic energy distribution and scale change. Mon. Wea. Rev. 2000; 128: 2865–2874.
GrotjahnR., HodyssD., CastelloC. Do frontal cyclones change size? Observed widths of north pacific lows. Mon. Wea. Rev. 1999; 127: 1089–1095.
HanleyJ., CaballeroR. Objective identification and tracking of multicentre cyclones in the ERA-Interim reanalysis dataset. Q. J. R. Meteorol. Soc. 2012; 138(664): 612–625.
JeongJ., HussainF. On the identification of a vortex. J. Fluid Mech. 1995; 285: 69–94.
KalnayE., KanamitsuM., KistlerR., CollinsW., DeavenD., co-authors. The NCEP/NCAR 40-year reanalysis project. Bull. Am. Meteorol. Soc. 1996; 77(3): 437–471.
KunnenR., ClercxH., GeurtsB. Vortex statistics in turbulent rotating convection. Phys. Rev. E. 2010; 82(3): 036306.
LimE.-P., SimmondsI. Southern hemisphere winter extratropical cyclone characteristics and vertical organization observed with the ERA-40 data in 1979–2001. J. Clim. 2007; 20(11): 2675–2690.
MüllerA., NévirP., SchielickeL., HirtM., PueltzJ., co-authors. Applications of point vortex equilibria: blocking events and the stability of the polar vortex. Tellus A. 2015; 67: 29184. DOI: http://dx.doi.org/10.3402/tellusa.v67.29184.
MunichRe. Significant Natural Disasters Since 1980 – Costliest Winter Storms in Europe (Overall Losses). 2014. NatCatSERVICE. Online at: http://www.munichre.com/en/reinsurance/business/non-life/natcatservice/significant-natural-catastrophes.
National Centers for Environmental Prediction/National Weather Service/NOAA/U.S. Department of Commerce. 1994. NCEP/NCAR Global Reanalysis Products, 1948-Continuing. Research Data Archive at the National Center for Atmospheric Research, Computational and Information Systems Laboratory. Online at: http://rda.ucar.edu/datasets/ds090.0/.
NeuU., AkperovM. G., BellenbaumN., BenestadR., BlenderR., co-authors. IMILAST: a community effort to intercompare extratropical cyclone detection and tracking algorithms. Bull. Am. Meteorol. Soc. 2013; 94(4): 529–547.
NielsenJ. W., DoleR. M. A survey of extratropical cyclone characteristics during GALE. Mont. Wea. Rev. 1992; 120: 1156–1168.
RudevaI. A. On the relation of the number of extratropical cyclones to their sizes. Izv. Atmos. Ocean. Phys. 2008; 44: 273–278.
RudevaI. A., GulevS. K. Climatology of cyclone size characteristics and their changes during the cyclone life cycle. Mon. Wea. Rev. 2007; 135: 2568–2587.
SchneidereitA., BlenderR., FraedrichK. A radius-depth model for midlatitude cyclones in reanalysis data and simulations. Q. J. R. Meteorol. Soc. 2010; 136: 50–60.
SimmondsI. Size changes over the life of sea level cyclones in the NCEP reanalysis. Mon. Wea. Rev. 2000; 128(12): 4118–4125.
SimmondsI., KeayK. Mean southern hemisphere extratropical cyclone behavior in the 40-year NCEP-NCAR reanalysis. J. Clim. 2000; 13(5): 873–885.
SinclairM. R. An objective cyclone climatology for the southern hemisphere. Mon. Wea. Rev. 1994; 122(10): 2239–2256.
SinclairM. R. Objective identification of cyclones and their circulation intensity, and climatology. Wea. Forecas. 1997; 12(3): 595–612.
ToryK. J., DareR., DavidsonN., McBrideJ., ChandS. The importance of low-deformation vorticity in tropical cyclone formation. Atmos. Chem. Phys. 2013; 13(4): 2115–2132.
TruesdellC. Two measures of vorticity. Indiana Univ. Math. J. 1953; 2: 173–217.
UlbrichU., FinkA., KlawaM., PintoJ. Three extreme storms over Europe in December 1999. Weather. 2001; 56(3): 70–80.
UlbrichU., LeckebuschG., PintoJ. Extra-tropical cyclones in the present and future climate: a review. Theor. Appl. Climatol. 2009; 96(1–2): 117–131.
WernliH., SchwierzC. Surface cyclones in the ERA-40 dataset (1958–2001). Part I: novel identification method and global climatology. J. Atmos. Sci. 2006; 63: 2486–2507.