Digital_banner_NL

Forecasting clear-air turbulence

Peter Bechtold (ECMWF), Martina Bramberger, Andreas Dörnbrack (both DLR, Germany), Martin Leutbecher, Lars Isaksen (both ECMWF)

 

Forecasting severe turbulence in the free troposphere and stratosphere is challenging. The turbulence is generated by processes such as shear instabilities (Kelvin–Helmholtz instabilities), upper-level fronto-genesis, large-amplitude mountain waves and breaking of convectively generated gravity waves (Ellrod & Knapp, 1992). In particular, turbulent eddies and waves with length scales of a hundred metres to several kilometres pose an important hazard to civil aviation. There is thus a great demand in aviation forecasting for reliable turbulence estimates (Sharman et al., 2014). The eddy dissipation rate (EDR), which is the cube root of the dissipation rate of turbulent kinetic energy, and hence has units of m2/3 s–1, has become the International Civil Aviation Organization (ICAO) standard for aircraft reporting and therefore the standard measure for clear-air turbulence (CAT).

A recent survey among the ECMWF forecast user community has revealed that a CAT parameter has a high priority on the user wish list, just below convection. Furthermore, we expect that an observable measure of turbulence intensity computed in the forecast will be beneficial for research purposes, in particular the development of turbulence parametrizations in the free troposphere and stratosphere. Our objective was therefore to develop a calibrated CAT parameter for the ECMWF Integrated Forecasting System (IFS) that serves the purpose of aviation forecasting as well as research in turbulence. The adopted solution was developed in a strong and useful collaboration with scientists from the German Aerospace Center (DLR). It is planned to become available later this year in IFS Cycle 47r3.​

Method

We essentially follow the method presented in Sharman & Pearson (2017). This means that a certain number of turbulence predictors are projected onto the observed climatological distribution of the EDR, and then a suitable CAT parameter is defined. The method has the advantage that it can be applied to widely used predictors of CAT in aviation forecasting and to the output of the turbulence scheme of the IFS. As the IFS currently does not use a prognostic turbulence scheme for the turbulent kinetic energy, the EDR is not directly available from the IFS. The following three predictors have been included to obtain a measure of the EDR:

  • Ellrod1 The positive definite Ellrod1 index has been developed by Ellrod & Knapp (1992). It is the product of the vertical wind shear with the total deformation. The deformation terms, involving horizontal gradients, are still small at 9 km resolution compared to the vertical gradients. However, using their product with the vertical shear term has been shown to improve the correlations compared to using the vertical shear term only (Ellrod & Knapp, 1992). The Ellrod1 index provides highly valued guidance in aviation forecasting and was also shown to be the best performing index in Sharman & Pearson (2017).
  • GWD While the Ellrod1 accounts for resolved flow features only, we also include a subgrid contribution from the drag by breaking convectively generated gravity waves. A simple approach is to scale the dissipation rate from the non-orographic gravity wave scheme (assuming a globally uniform departure wave spectrum) with the normalised vertically integrated convective heating between 500 hPa and the cloud top. The cube root of the dissipation is then in units of the EDR.
  • DISS The total turbulent dissipation rate of the IFS in units of EDR is derived from the model’s physical tendencies for horizontal momentum. It includes tendency contributions from the vertical diffusion scheme and the convective momentum transport. The largest contribution to the total turbulent dissipation is from the vertical diffusion scheme that includes dissipation due to turbulent mixing, orographic wave drag and orographic blocking (Beljaars et al., 2004).

The numerical expressions for Ellrod1, GWD and DISS are presented in detail in Bechtold et al. (2021). The computation of the Ellrod1 index increases the numerical cost of an IFS forecast by 4–5% as it requires additional inverse spectral transforms of the vorticity and the divergence fields to compute the required horizontal wind gradients.

%3Cstrong%3EFIGURE%201%3C/strong%3E%20Probability%20density%20distributions%20of%20the%20natural%20logarithm%20of%20the%20Ellrod1%20index,%20gravity%20wave%20drag%20(GWD)%20and%20turbulent%20dissipation%20(DISS)%20for%20the%20atmospheric%20layer%20between%20600%20and%20150%20hPa%20(4.5%E2%80%9315%20km)%20as%20obtained%20from%20six%20months%20of%20IFS%20simulations%20for%202019%20at%20the%20TCo1279%20resolution%20(about%209%20km%20horizontal%20grid%20spacing).%20The%20solid%20lines%20show%20the%20fitted%20log-normal%20distributions.%20The%20climatological%20EDR%20distribution%20is%20denoted%20by%20the%20dashed%20grey%20curve.
FIGURE 1 Probability density distributions of the natural logarithm of the Ellrod1 index, gravity wave drag (GWD) and turbulent dissipation (DISS) for the atmospheric layer between 600 and 150 hPa (4.5–15 km) as obtained from six months of IFS simulations for 2019 at the TCo1279 resolution (about 9 km horizontal grid spacing). The solid lines show the fitted log-normal distributions. The climatological EDR distribution is denoted by the dashed grey curve.

Figure 1 shows the climatological distributions of the natural logarithm of Ellrod1, GWD and DISS as obtained from daily 24-hour hindcasts with the IFS during the first half of 2019. The distributions are valid for the 4.5–15 km atmospheric layer. Also included in Figure 1 (solid lines) are the log-normal fits to the individual distributions, as well as the climatological distribution of the EDR derived from aircraft data by Sharman & Pearson (2017). We note that all parameters closely follow a log-normal law and that the distributions of GWD and DISS, which naturally have units of the EDR, are reasonably close to the climatological distribution of the EDR.

The aim of the extensive pre-computations in Figure 1 is to fit log-normal distributions for each index in order to be able to project indices on the climatological EDR. This computation has only to be done once. Then, the final step and actual forecasting of the EDR consists in computing in each forecast Ellrod1, GWD and DISS and projecting the values on the climatological EDR. The simple projection procedure, together with a table of the fitting parameters (the mean and variance of the distributions as used in Figure 1), is given in Bechtold et al. (2021).

Observations

It can be difficult to obtain calibrated EDR measurements onboard civil aircraft as the data is generally the property of the airline. We have retrieved 'max EDR', that is the peak EDR of 12 individual EDR measurements over a one-minute interval, from the US National Oceanic and Atmospheric Administration (NOAA) Meteorological Assimilation Data Ingest System (MADIS) public archive for aircraft data. The dataset and algorithm onboard civil aircraft to compute the EDR is thoroughly described in Sharman et al. (2014). In summary, the aircraft response to turbulent eddies with wavelengths of 10 m to 1 km is felt as bumpiness. The aircraft is more sensitive to vertical gusts than to lateral gusts, i.e. the steering is more sensitive in the vertical as aircraft are built to be stable in the lateral axis. The EDR is proportional to the root-mean-square vertical acceleration experienced by an aircraft under specific flight conditions. The EDR algorithm uses either the measured vertical accelerations or the aircraft vertical winds that are determined using the aircraft's angles of attack and cruising speed. A Fourier transform is performed on the 8 Hz sampled vertical velocity time series (with a cruise speed of 250 m s–1, this retains eddies of 30 m that can significantly affect the aircraft wings), and a von Kármán spectrum is fitted to the retrieved vertical velocity spectrum in the inertial turbulent subrange. When the algorithm detects a turbulent event, a report is generated and downlinked at a 1 min (~15 km) interval. Follow-up reports are then also generated.

In order to compare the IFS with observations, we have performed daily 24-hour forecasts for January to March 2019 and used hourly model output on the full vertical model resolution of 137 levels, but on a reduced horizontal grid of 0.3º x 0.3º to make the data volume more manageable. The comparison focuses on the height range 5–12 km, i.e. the cruising altitude, where the IFS vertical resolution is roughly 300 m. The projection of the forecasts onto the observations is done by retaining all observations 15 minutes before and after the full hour and allowing a maximum height difference between observations and model data of 160 m.

We had to account for the fact that, out of the more than 4 million observations onto which the IFS data has been projected, a large majority has zero value, while the model is producing a quasi-continuous field. We therefore included an EDR threshold of 0.005 m2/3 s–1 for the observations. We were thus able to retain only 197,000 observations for the statistics. The statistics are robust, however, even for a single month.​

Definitions of CAT

We have defined CAT parameters/products based on the EDR projections of the three indices, Ellrod1, GWD and DISS, as follows:

CAT1 = 0.7 x Ellrod1* + GWD*

CAT2 = 0.66 x DISS* + GWD*

CAT12 = 0.5 x (CAT1 + CAT2)

where the asterisk denotes the value of the index after the EDR projection. In addition, scaling factors have been added to better fit the aircraft observations we have used during the 3-month evaluation period. While CAT1 is essentially based on the Ellrod1 index (GWD* is relatively small), CAT2 represents the total dissipation rate of the IFS, and CAT12 is the arithmetic mean of CAT1 and CAT2.

%3Cstrong%3EFIGURE%202%3C/strong%3E%20Probability%20density%20distributions%20of%20the%20EDR%20for%20the%20period%20January%20to%20March%202019%20and%20for%20heights%20between%205%20and%2012%20km,%20as%20obtained%20from%20the%20NOAA/MADIS%20observational%20dataset%20(black)%20and%20from%20daily%200%20to%2024-hour%20IFS%20forecasts%20projected%20onto%20the%20observation%20locations,%20where%20CAT1,%20CAT2%20and%20CAT12%20are%20as%20defined%20in%20the%20text.
FIGURE 2 Probability density distributions of the EDR for the period January to March 2019 and for heights between 5 and 12 km, as obtained from the NOAA/MADIS observational dataset (black) and from daily 0 to 24-hour IFS forecasts projected onto the observation locations, where CAT1, CAT2 and CAT12 are as defined in the text.

Figure 2 shows, for the whole 3-month period, the probability distribution functions of the observations and the corresponding projected IFS data for the different EDR estimates. Both CAT1 and CAT2 underestimate the relative occurrence of weak turbulence intensities with EDR < 0.05 m2/3 s–1 compared to the observations, which deviate from the log-normal law for small values. CAT2 overestimates the observed distribution for medium to strong turbulence, while CAT1 underestimates the occurrence of strong to severe turbulence with EDR > 0.3 m2/3 s–1. The best overall match with the observed distribution is obtained by a linear combination of CAT1 and CAT2.

Case study 2 March 2019

Figure 3 shows global maps of two EDR-based diagnostics for CAT on 2 March 2019 as obtained from a 0–24 h forecast with a horizontal resolution of about 9 km. The values represent averages over the whole day and the 10–12 km atmospheric layer. To illustrate the role of wind shear, isotachs of the 250 hPa wind speed have also been included in Figure 3. We notice that the global distributions of EDR are very similar between the independent products CAT1 and CAT2. Maximum daily and layer mean values of between 0.1 and 0.18 are present near the flanks of the subtropical jets and in the storm tracks, but also over orography such as the Rocky Mountains and the Himalayas. We also notice significant EDR in tropical regions, especially near the equator and in convective regions over land. The convective contribution to CAT1 stems from the convective gravity wave drag, while for CAT2 it also includes a contribution from convective momentum transport.

%3Cstrong%3EFIGURE%203%3C/strong%3E%20Global%20daily%20mean%20distributions%20of%20CAT1%20and%20CAT2%20for%202%20March%202019%20averaged%20over%20the%2010%E2%80%9312%20km%20layer%20from%20a%2024-hour%20IFS%20hindcast%20at%20the%20TCo1279%20resolution%20(about%209%20km%20horizontal%20grid%20spacing).%20The%20isotachs%20of%20250%20hPa%20wind%20speed%20(m/s)%20are%20overlaid%20by%20red%20contours.
FIGURE 3 Global daily mean distributions of CAT1 and CAT2 for 2 March 2019 averaged over the 10–12 km layer from a 24-hour IFS hindcast at the TCo1279 resolution (about 9 km horizontal grid spacing). The isotachs of 250 hPa wind speed (m/s) are overlaid by red contours.

A comparison with all non-zero EDR aircraft reports between 5 and 12 km altitude for 2 March 2019 is presented in Figure 4. The different CAT products predict the turbulence regions over the southern Rockies and northern Florida, and also the observed elevated turbulence over the central US, southern Greenland and the north-eastern Atlantic. However, while CAT1 tends to underpredict the high turbulence regions with observed EDR above 0.3 m2/3 s–1, CAT2 tends to overestimate the horizontal extension of high turbulence regions. The linear combination of CAT1 and CAT2 seems to perform best. This is also consistent with the probability density functions in Figure 2. The representativity of the observations has to be treated with some caution, however, as aircraft try to avoid and/or rapidly quit regions with strong turbulence.

%3Cstrong%3EFIGURE%204%3C/strong%3E%20EDR%20for%202%20March%202019%20including%20(a)%20all%20available%20observations%20between%205%20and%2012%20km,%20(b)%20CAT1,%20(c)%20CAT2%20and%20(d)%20a%20linear%20combination%20of%20CAT1%20and%20CAT2%20as%20defined%20in%20the%20text.%20The%20CAT%20values%20are%20taken%20from%20a%200%20to%2024-hour%20hindcast%20at%20the%20TCo1279%20resolution%20(about%209%20km%20horizontal%20grid%20spacing).%20Colours%20(see%20legend)%20and%20symbols%20are%20used%20to%20denote%20the%20strength%20of%20the%20EDR,%20with%20different%20symbols%20for%20each%200.1%20m%3Csup%3E2/3%3C/sup%3E%20s%3Csup%3E%E2%80%931%3C/sup%3E%20EDR%20interval.
FIGURE 4 EDR for 2 March 2019 including (a) all available observations between 5 and 12 km, (b) CAT1, (c) CAT2 and (d) a linear combination of CAT1 and CAT2 as defined in the text. The CAT values are taken from a 0 to 24-hour hindcast at the TCo1279 resolution (about 9 km horizontal grid spacing). Colours (see legend) and symbols are used to denote the strength of the EDR, with different symbols for each 0.1 m2/3 s–1 EDR interval.

Need for ensembles

We have verified the daily high-resolution forecasts of CAT during the January to March 2019 test period. However, turbulence is inherently intermittent and reliable estimates of turbulence require an ensemble approach. Therefore, we have also evaluated ensemble forecasts (ENS) of CAT, comprising a 15-member ensemble at a resolution of about 18 km that was run daily during the first two weeks of January 2019.

 
EDR parameter Corr HRES Jan–Mar MAE HRES Jan–Mar Corr HRES 1–14 Jan Corr ENS 1–14 Jan MAE HRES 1–14 Jan CPRS ENS 1–14 Jan
CAT1 0.33 0.050 0.33 0.38 0.049 0.030
CAT2 0.30 0.057 0.32 0.37 0.054 0.034
CAT12 0.35 0.045 0.36 0.40 0.049 0.029

TABLE 1 Verification of different EDR parameters against observations for the high-resolution forecasts (HRES - grid spacing of about 9 km) for January–March 2019 and for HRES and the ensemble forecasts (ENS - grid spacing of about 18 km) for the period 1–14 January 2019. Verification statistics are correlation (Corr), mean absolute error (MAE) and continuous ranked probability score (CRPS).

 

Table 1 shows the point correlations and mean absolute error (MAE) against the observations of the high- resolution forecasts for January to March 2019. For the period 1–14 January 2019, which comprises 19,600 observations, we compare the high-resolution forecasts to the ensemble forecasts. The latter are evaluated using the ensemble mean correlation and the continuous ranked probability score (CRPS). The CRPS of the ensemble against observations directly compares to the MAE of the high-resolution forecasts. All data has been reduced to a 0.3º x 0.3º output grid on model levels.

The 3-month average point correlations from the high-resolution forecasts with the observations is ​0.33 for CAT1, 0.30 for CAT2 and 0.35 for CAT12, while the MAE is around 0.05 m2/3 s–1. The errors are slightly higher for CAT2 compared to CAT1, which is expected, as CAT2 is based on model tendencies and therefore more variable than CAT1, which is based on state variables. For the shorter January test period, the correlations attain values between 0.32 for CAT2 and 0.36 for CAT12, while the MAE remains around 0.05 m2/3 s–1. Comparing these values to the ensemble, we see that the ensemble performs significantly better, with ensemble mean correlations attaining 0.37 for CAT2 and up to 0.40 for the combined product CAT12, while the CRPS is around 0.03 m2/3 s–1. We therefore think that the ensemble forecasts are of reasonable reliability, given that turbulence intensity is typically classified in EDR intervals of 0.1 m2/3 s–1. The obtained point correlations might still appear relatively low, but they have to be put into perspective with point correlations of around 0.53 obtained for 10-metre wind speed forecasts over land, and point correlations of 0.2–0.4 obtained when verifying forecasts of daily rainfall accumulations over tropical land against synoptic observations.

Outlook

We think there is now sufficient evidence that a clear-air turbulence (CAT) diagnostic based on the eddy dissipation rate (EDR) will be a useful addition to the standard IFS output for both forecasters and research in turbulence, particularly when used in the context of ensembles. Such a turbulence diagnostic will also be interesting for the next climate reanalysis, ERA6.

We have put CAT based on the total dissipation rate (referred to as CAT2 in the article) into IFS Cycle 47r3, which is expected to become operational in the autumn of 2021. The diagnostic will be available on model levels for the high-resolution forecast and the ensemble, therefore roughly satisfying the 0.1º horizontal and 300 feet vertical resolution requirements of the turbulence guidance product that is in development for the International Civil Aviation Organization by 2030 (Kim et al., 2018). In the future we might also access real-time EDR data as provided by the International Air Transport Association (IATA) for forecast verification and possibly data assimilation.

The total dissipation rate is a more useful CAT predictor in the planetary boundary layer than the Ellrod1 index. It also avoids the increased computational burden of 4% that is associated with the computation of the horizontal gradients in the latter. However, both CAT1 and CAT12 can be computed online by internal users of the IFS from Cycle 47r3 onward under optional computational flags and may be turned on operationally as part of a later upgrade. Finally, the detailed information provided in the Technical Memorandum by Bechtold et al. (2021) should enable users to also compute offline an additional CAT product based on Ellrod1 and enable advanced postprocessing using non-linear regression and/or machine learning.


Further reading

Bechtold, P., M. Bramberger, A. Dörnbrack, L. Isaksen & M. Leutbecher, 2021: Experimenting with a Clear Air Turbulence (CAT) index from the IFS. ECMWF Technical Memorandum No. 874.

Beljaars, A.C.M., A.R. Brown & N. Wood, 2004: A new parametrization of turbulent orographic form drag. Q. J. R. Meteorol. Soc., 130, 1327–1347.

Ellrod, G.P. & D.I. Knapp, 1992: An objective clear-air turbulence forecasting technique: verification and operational use. Weather and Forecasting, 7, 150–165.

Kim, J.-H., R. Sharman, M. Strahan, J. ScheckC. Bartholomew, J.C.H. Cheung et al., 2018: Improvements in nonconvective aviation turbulence prediction for the world area forecast system. Bull. Am. Meteor. Soc., 99, 2295–2311.

Sharman, R., L.B. Cornman, G. Meymaris, J. Pearson & T. Farrar, 2014: Description and derived climatologies of automated in situ eddy dissipation rate reports of atmospheric turbulence. J. Appl. Meteor. Climatol., 53, 1416–1432.

Sharman, R.D. & J. Pearson, 2017: Prediction of energy dissipation rates for aviation turbulence. Part I: Forecasting nonconvective turbulence. J. Appl. Meteor. Climatol., 56, 317–337.