Abstract

In operational applications lateral-boundary coupling data are provided to one-way nested limited-area models with time intervals of more than an order of magnitude larger than the time step of the coupled model. In practice, these fixed coupling-update frequencies are established by common-sense guesswork and by technical restrictions rather than by rigorous methods. As a result, situation-dependent failures are never completely excluded when coupling data enter the domain more rapidly than can be sampled by the a priori fixed frequency. To avoid misinterpreting such failures, the data transfer between the coupling and the coupled model should be monitored.

The present paper approaches this as a problem of undersampling. It investigates how the coupling-update frequency can be monitored by using a digital recursive filter in the coupling model. A response function for such a filter is derived. Its implementation in a NWP model is discussed and some extensive tests are presented. A possible application is discussed in which this monitoring is used for assessing the data transfer to the coupled model and additionally for adapting the coupling updates to the actual meteorological content of the coupling-model output.

1. Introduction

Warner et al. (1997) gave an extensive overview of the limitations to be critically taken into account when using limited-area models (LAMs). One of them is the use of lateral-boundary data of a much poorer temporal resolution than the time step of the coupled model. Surprisingly almost no efforts have been made to invent methods to evaluate these coupling-update frequencies in a well-founded way. Termonia (2003) showed in a specific setup that it is possible to monitor the efficiency of the coupling-update frequency, allowing the data transfer between coupling and coupled model to be improved. This paper proposes a more general approach.

The problem of using large coupling-update time intervals can be stated as follows: In the process of preparing and providing coupling data from the coupling model to an operational one-way nested coupled model, the time series of the coupling data of the coupling model is resampled with a typical time resolution of, say, 3, 6, or 12 h. This is at least an order of magnitude larger than the current model time steps. As a result, the original time series may be severely undersampled. Any information potentially present in those frequencies in the temporal spectral decomposition larger than the Nyquist sampling frequency of the new time resolution is irrevocably lost. Moreover, the frequencies lower than the Nyquist frequency may be corrupted because of aliasing. These undersampled data are then temporally interpolated to get data for each time step of the coupled model (which is not necessarily the same as for the coupling model).

In operational applications, the frequency of these coupling updates are fixed by common-sense guesswork and pragmatical considerations to deal with technical limitations in data transfer. Once this coupling-update frequency is established, the user of the LAM has no choice but to trust it; no guarantee is given that situation-dependent cases in which crucial information is lost in the higher frequencies will never occur.

In the present paper, the problem is approached from the viewpoint of signal sampling. It is investigated to what extent the information loss due to the undersampling of the coupling data in the coupling model can be estimated before its use in the coupled model. This is done by filtering the time series of some pertinent model variable in the coupling model to estimate the information present in the high frequencies, that is, higher than the Nyquist frequency of the resampling.

A detailed discussion of an implementation of a recursive digital filter in the source code of the ALADIN limited-area model (ALADIN International Team 1997) is discussed. It is argued that this implementation is sufficiently independent from the specifications of the ALADIN model to be adapted to other models. Tests are presented for December 1999. During this month some of the severest storms passed over France. It is explained how the coupling-data information loss can be monitored for such storms.

The present paper shows that the data transfer can be monitored in a better way than proposed in Termonia (2003), providing a useful tool for the person using some particular LAM output. When crucial information turns out to be missing in the data transfer, the user of the coupled model has to consider relying on the model output of the coupling model, instead of misinterpreting the corrupted forecast of the coupled model while mistakenly believing that it improves the coupling-model forecast.

In a second stage, one may try to optimize the data transfer between the coupling and the coupled model. Ideally, coupling data should be provided with the time resolution of the coupled-model time step. In operational applications this can never be accomplished. However, it is explicitly shown in a specific example that the occasional damage caused by the resampling can be kept within predetermined limits by adapting the coupling updates according to the estimated information loss.

A detailed discussion of an operational failure due to this problem can be found in Termonia (2003).

2. The information loss in interpolated data

This section investigates the resampling and the linear interpolation of the time evolution Φn of a model variable Φ, where n labels the time step; that is, Φn = Φ(nδt) is the value of Φ at t = nδt, with δt the model time step.

To consider the effect of the resampling, the Fourier transform of Φn has to be considered:

 
formula

where γ(ω) are the Fourier components and where the integral is restricted to the interval [−ωN, ωN], with ωN = π/δt the Nyquist frequency of the original time resolution.

After resampling with a new time interval Δt > δt, any information in the modes with frequencies larger than the new Nyquist frequency ΩN = πt is lost; that is, the information contained in the coefficients γ(ω) with |ω| > πt is not passed from the coupling model to the coupled model. Given the time series of a variable Φn this information loss can be computed in the time domain by means of a high-pass filter with a pass band for all the frequencies |ω| > πt, having a frequency response function of 1.

There is some arbitrariness in the choice of such a response function in the lower part of the frequency domain |ω| < πt. However, a response function can be chosen in order to mimic the output of the linear interpolation in this frequency interval. Denoting the linear interpolation between time step n0 and time step n1 by

 
formula

with Δt = (n1n0)δt the time interval in which the interpolation is taken, the linear interpolation of Φ can be written as

 
formula

The interpolation can thus be completely characterized by the linear interpolation of each mode eiωt,

 
eiωt = (ℒeiωτ) × eiωt0,
(4)

with t0 = n0δt and τ = tt0 taken in the interval [0, Δt]. In fact, the interpolation of

 
c(τ) = exp(iωτ),
(5)

at τ0.5 = Δt/2 in the middle of the interval Δt is given by

 
formula

This is geometrically illustrated in Fig. 1. The mode c(τ) moves in time along the unit circle from 1 to ct). Figure 1a shows a mode with frequency satisfying |ω| < πt. In this case the interpolation dampens the mode with a factor cos(ωΔt/2) at τ = Δt/2, losing 1 − cos(ωΔ/2) of the amplitude. Figure 1b shows a mode with πt < ω < 2πt, where the interpolation generates the opposite phase (dashed line), completely corrupting the mode, as expected since it exceeds the Nyquist frequency.

Fig. 1.

Illustration of the linear interpolation of an oscillating mode c(τ) = exp(iωτ) with given frequency ω in a time interval Δt. The mode c(τ) moves in time along the unit circle from 1 to ct). A distinction can be made between (a) a mode with frequency satisfying |ω| < πt, in which case the interpolation dampens the mode with a factor cos(ωΔt/2) at τ = Δt/2, loosing 1 − cos(ωΔ/2) of the amplitude, and (b) a mode with πt < ω < 2πt, where the interpolation generates the opposite phase (dashed line), completely corrupting the mode

Fig. 1.

Illustration of the linear interpolation of an oscillating mode c(τ) = exp(iωτ) with given frequency ω in a time interval Δt. The mode c(τ) moves in time along the unit circle from 1 to ct). A distinction can be made between (a) a mode with frequency satisfying |ω| < πt, in which case the interpolation dampens the mode with a factor cos(ωΔt/2) at τ = Δt/2, loosing 1 − cos(ωΔ/2) of the amplitude, and (b) a mode with πt < ω < 2πt, where the interpolation generates the opposite phase (dashed line), completely corrupting the mode

So the loss of the information in a mode with frequency ω at time τ0.5 in the middle of the time interval [t0, t0 + Δt] can be estimated by filtering with a frequency response function,

 
formula

In this way, the arbitrariness in the choice of the response function in the low frequencies |ω| < πt is used to have a response that mimicks the effect of the linear interpolation in this part of the spectrum.

In the conclusions of Termonia (2003), it was proposed to try to use the root-mean-square difference between the original forecast and the interpolated value. It is now interesting to see what this yields on modes of the form in Eq. (5). In particular, the difference between the interpolated and the original field becomes

 
formula

The components with frequency ω = 4πt do not contribute to this integral, whereas they are 4 times bigger than the Nyquist frequency πt and are lost after the resampling. Physically, this corresponds to time scales of the order of Δt/4. For a coupling-update interval of Δt = 3 h this corresponds to completely missing the information loss in phenomena with typical time scales of about 45 min. For a coupling-update interval of Δt = 12 h this corresponds to 3 h. So using measures based on some norms of the differences ΔΦ(t0 + Δt/2) in (8) may deceptively underestimate the loss of phenomena of such time scales.

In fact, it could be argued that the usefulness of a measure based on differences of the kind in (8) depends on whether frequencies of about ω ≈ 4πt occur in operational forecasts or not. To exhaustively answer this question, the frequency analysis should be performed explicitly. Moreover these 4πt modes depend on the used time interval Δt. So such an analysis should be repeated when changing this coupling-update interval in other setups. Spectral decompositions in the time domain are not part of NWP models. To tackle this question, some techniques to perform such a spectral analysis should be investigated or developed first. In contrast, by implementing a digital filter in the coupling model having a response function resembling Λ in (7) close enough, the signal can be filtered in the time domain to obtain an estimate of the information loss. Such an implementation in the ALADIN model will be explained in the next section.

3. Implementation in the ALADIN model

To estimate the information loss in the linearly interpolated undersampled coupling data, it has been investigated whether model variables can be filtered during the integration of a coupling model.

Although the filter should in principle be applied to all model variables, it is applied here with the aim to monitor the coupling of storms only. The best variable to detect a passing storm is the pressure. In the hydrostatic ALADIN limited-area model the logarithm of the surface pressure Π = lnPs is the prognostic variable. Hence the filter was implemented for Π.

The tests presented in this paper were performed in an implementation of ALADIN on a domain having France in its center, henceforth referred to as ALADIN-France. ALADIN-France receives its coupling data from the global model Action de Recherche Petite Echelle Grande Echelle (ARPEGE) of Météo France, but provides the coupling data for the ALADIN implementation over Belgium used at the Royal Meteorological Institute of Belgium, that is, ALADIN-Belgium. Hence, the tests in the ALADIN-France model pertain to the monitoring of the data transfer from ALADIN-France to ALADIN-Belgium. It is explained later that the results obtained here do not depend on this specific setup of the ALADIN model but can also be applied to the data transfer between, say, a global model and a LAM.

To apply a filter operationally during the forecast run of the coupling model a filter is needed that only uses a marginal amount of the computing resources. Indeed, the monitoring cannot delay the operational run of the coupling model in a substantial way. For this reason, a second-order recursive filter was chosen, since recursive filters are known to be computationally cheap, see, for example, Shanks (1967). Moreover it was argued by Lynch and Huang (1994), in the context of digital filtering initialization (DFI), that they may be applied during a forecast run (e.g., for noise control in climate modeling). To approximate the function Λ in (7), Butterworth high-pass filters (Holtz and Leondes 1966) were tested, having a response function H(ω) whose modulus square has the form

 
formula

where | · | denotes the modulus of the complex number, δt the resolution of the digital signal (in this case the coupling-model time step), N the order of the filter, and ωc the critical frequency defining the high-pass band. A filter with such a response can be realized by the recursive formula

 
formula

where xn is the incoming signal and yk the filtered output. This expression is particularly convenient for use during the forecast of the coupling model. Indeed, k corresponds to model time step t = kδt, and km to tmδt. To monitor the coupling-update frequency for coupling storms, x was taken to be Π and y yields the filtered Πfilt at time t. For given ωc, N, and δt, the coefficients am and bn can be computed as explained in the appendix.

The recursive expression (10) only uses values of the fields computed at previous time steps. These can be kept in computer memory during the coupling-model run. For the two-dimensional field Π this constitutes very little memory. The implementation of Eq. (10) in a model code is easy. The computation of Eq. (10) itself for one two-dimensional field is very cheap compared to the rest of the model integration.

Figure 2 shows the response |H(ω)| function of a Butterworth filter with ωc = 0.9πt, with δt = 415.385 s and Δt = 3 h as a function of ωΔt, in comparison with the loss Λ in (7). This Butterworth filter is maximally flat in the low-frequency range, providing an excellent fit of Λ in the range [0, π/2] on Fig. 2. This avoids erroneously detecting a slowly changing background. For this reason, this filter is used throughout the rest of this paper.

Fig. 2.

The loss Λ(ωΔt) (solid line) and the modulus of the response function of the second-order Butterworth filter with δt = 415.385 s, Δt = 3 h, and ωc = 0.9πt as a function of ωΔt (dashed line)

Fig. 2.

The loss Λ(ωΔt) (solid line) and the modulus of the response function of the second-order Butterworth filter with δt = 415.385 s, Δt = 3 h, and ωc = 0.9πt as a function of ωΔt (dashed line)

It is important to note that recursive filters of this kind generate a phase shift causing a time delay of the filtered signal. This is also explained in the appendix. Figure 3 shows the group delay α(ωΔt) as a fraction of the coupling-update frequency Δt, for the recursive filters in (9) with ωc = 0.9πt and δt = 415.385 s, for Δt = 12 h and Δt = 30 min. For Δt = 6 h, 3 h, 1 h, the curves lie between these two curves (not shown on the figure). So for the second-order Butterworth filters used in this paper, this causes a group delay of about half the coupling-update frequency of the modes with frequencies close to ωπt. For higher frequencies it is less than 0.5. The filter slightly distorts the form of the feature representing the information loss. Nevertheless, the variability is the same as would be obtained by a filter with an equal real response function, but with a constant time delay.

Fig. 3.

Group delay of the second-order Butterworth filter with ωc = 0.9πt and δt = 415.385 s as a fraction of the coupling-update frequency Δt, for Δt = 12 h (solid line) and Δt = 30 min (dashed line). For Δt = 6 h, 3 h, 1 h, the curves lie between these two curves (not shown)

Fig. 3.

Group delay of the second-order Butterworth filter with ωc = 0.9πt and δt = 415.385 s as a fraction of the coupling-update frequency Δt, for Δt = 12 h (solid line) and Δt = 30 min (dashed line). For Δt = 6 h, 3 h, 1 h, the curves lie between these two curves (not shown)

In fact, such phase shifts can be avoided by applying the recursive filter twice: once in the forward time direction and once in the backward time direction, see, for example, Lynch and Huang (1994). The aim of the present paper is to keep the implementation as minimalistic as possible, in order to keep the used resources in the coupling model as small as possible. So results are only shown for the forwardly applied recursive filter, which is the most straightforward to implement in a model code. The price to be paid is the phase shift, but it will be shown in the examples later that this does not constitute a serious problem for the monitoring.

The monitoring of the coupling-update frequency is entirely done in the coupling model ALADIN-France at the end of each time step. This model was run with a time step of δt = 415.385 s. In fact, ALADIN is a spectral model (Haugen and Machenhauer 1993) so that at the end of each time step, Π is available by its spectral coefficients Π̃ defined by

 
formula

where ΠIJ are the values of Π in the grid points with spatial indices I = 0, … , M − 1, and J = 0, … , N − 1; n is the time index; and Π̃kl the complex spectral coefficients with k = 0, … , M − 1 and l = 0, … , N − 1, where the complex conjugate of Πkl satisfies the relation Π̃*kl = Π̃Mk,Nl to have reality of ΠIJ. It can be seen from this expression that the recursive filter can be applied to the gridpoint values ΠIJ as well as on the spectral coefficients Π̃kl. Indeed, filtering the spectral coefficients first and transforming them to gridpoint space gives the same result as applying them on gridpoint space directly, and vice versa. In this spectra model the spectral coefficients Π̃kl are filtered by Eq. (10) to obtain Π̃filtkl.

The model is run on a domain with M = N = 300, including an artificial extension zone to make the fields periodic (see Haugen and Machenhauer 1993) and a spatial resolution of Δx = Δy = 9.51 km.

For the coupling of storms only the large-scale part of Π is of interest. In spectral space this is obtained by using only those components having wavelengths larger than a certain scale. As will be shown later it is useful to retain lengths larger than 158 km, corresponding to a truncation in wavenumber space to the modes with wavenumber k and l satisfying

 
(k2 + l2)1/2 ≤ 10.
(12)

As will be shown later, it is useful in the spectral model to compute the sum of the moduli of these spectral coefficients Π̃filtkl for diagnostic purposes. This yields an upper bound for the values of the gridpoint field ΠfiltIJ.

 
formula

which is obtained from the Schwarz inequality applied to (11). When μ is small, ΠfiltIJ is small in all grid points. In that case there will be no information loss larger than μ after the resampling anywhere on the ALADIN-France domain.

The ALADIN code was adapted to use five filters of the type (10) on the components of Π̃, reserving some extra memory to hold Π̃filt,kn and the previous time steps of Π̃km, and to compute Π̃filt,k at the end of each time step. Additionally, at forecast times in the middle of the predefined 3-h coupling time intervals 0130, 0430, 0730 UTC, etc., the norm μ is computed. If it exceeds a predetermined value μc, the spectral field Π̃filt is transformed to gridpoint space and Πfilt is written to the output files.

4. Tests

This section discusses some tests with the implemented filters as described in the previous section.

a. Details about a single model run

The case of an ALADIN-France forecast of one of the French Christmas storms on 26 December 1999 is discussed here. The results for the rest of December 1999 will be discussed in the next subsection.

Figure 4 shows the ALADIN-France forecast of the mean sea level pressure of one of the French Christmas storms at 0300, 0600, 0900, and 1200 UTC 26 December 1999, based on the 0000 UTC analysis. The monitoring of the coupling-update frequency is illustrated on the indicated rectangular domain. The storm enters this domain between 0300 and 0600 UTC and leaves the domain at 1200 UTC.

Fig. 4.

ALADIN-France forecast of the mean sea level pressure at (a) 0300, (b) 0600, (c) 0900, and (d) 1200 UTC 26 Dec 1999. A rectangular boundary of an embedded domain is indicated. The contour interval is 2 hPa

Fig. 4.

ALADIN-France forecast of the mean sea level pressure at (a) 0300, (b) 0600, (c) 0900, and (d) 1200 UTC 26 Dec 1999. A rectangular boundary of an embedded domain is indicated. The contour interval is 2 hPa

Please note that a limited-area model on such a small domain is not suited to improve the forecasts of such a storm since it crosses the entire domain in less than 12 h. The forecast is in that case completely dominated by the lateral-boundary conditions. However, this domain is used here solely for the sake of the presentation; the actual coupling itself is not studied. The problems of the limitations in temporal resolution of the coupling data treated in the present paper are also present on the boundaries of larger domains. On such domains the impact on the forecast quality in a central region of interest is estimated based on knowledge of the error propagation from the lateral boundaries to the interior; see, for example, Errico et al. (1993), Warner et al. (1997), and Vannitsem (2003) for a recent paper. In contrast to this standard practice, the monitoring of the information loss studied in the present paper anticipates the errors at the lateral boundaries.

Figure 5 shows the filtered field Πfilt obtained after spectrally truncating Π by Eq. (12), and filtering it according to Eq. (10), with Δt = 3 h, for the forecast in Fig. 4, at 0430, 0730, 1030, and 1330 UTC, in the middle of the 3-h coupling-update time intervals [0300, 0600], [0600, 0900], [0900, 1200], and [1200, 1500], respectively. The estimated information loss of the storm is clearly distinguishable. The time delay can be noticed by comparing this figure to Fig. 4. Indeed, the information loss at, for instance, 1030 UTC on Fig. 5c almost coincides with position of the storm at 0900 UTC on Fig. 4c.

Fig. 5.

Filtered Π = lnPs by the second-order Butterworth filter, after truncating the spectral decomposition to scales larger than 158 km, at (a) 0430, (b) 0730, (c) 1030, and (d) 1330 UTC 26 Dec 1999. The contour interval is 0.001 hPa. The time delay of about 1 h 30 min can be noticed by comparing this figure to Fig. 4 

Fig. 5.

Filtered Π = lnPs by the second-order Butterworth filter, after truncating the spectral decomposition to scales larger than 158 km, at (a) 0430, (b) 0730, (c) 1030, and (d) 1330 UTC 26 Dec 1999. The contour interval is 0.001 hPa. The time delay of about 1 h 30 min can be noticed by comparing this figure to Fig. 4 

In fact, when using the gridpoint fields only (e.g., in a gridpoint model), it turns out not to be necessary to truncate the spectrum as in (12). However, in that case, some tests showed that the filter became unstable during the first few hours of the forecast because of some small residual noise in the small scales left over after the DFI1 [this has been verified by computing ΣIJ |∂Ps/∂t|IJ as in Lynch and Huang (1992), but is not shown here]. This instability could be avoided by applying the filter from 0300 UTC onward.

Figure 6 shows the untruncated filtered surface pressure field Πfilt at 1030 UTC. It can be verified that Fig. 5c shows the same signal as Fig. 6 but without the small-scale details. If the aim is to monitor the coupling update frequency for coupling storms between the two models, the untruncated field Πfilt can be used as well as the truncated Πfilt by (12).

Fig. 6.

Filtered Π = lnPs by the second-order Butterworth filter, without truncating the spectral decomposition, at 1030 UTC 26 Dec 1999, to be compared to Fig. 5c. The contour interval is 0.001 hPa

Fig. 6.

Filtered Π = lnPs by the second-order Butterworth filter, without truncating the spectral decomposition, at 1030 UTC 26 Dec 1999, to be compared to Fig. 5c. The contour interval is 0.001 hPa

Figure 7a shows the norm (13) of the truncated Πfilt by Eq. (12). After some initial period during the first 3 h of the forecast, the storm can be recognized. The maximum takes place between 0900 and 1200 UTC 26 December 1999. In the afternoon of 27 December 1999 the norm starts to increase again, indicating the creation of a second storm (as is explained in the next subsection).

Fig. 7.

Time series of the norm of the filtered field Π = lnPs as a function of the time of (a) the truncated field for scales larger than 158 km only and (b) summing over all length scales

Fig. 7.

Time series of the norm of the filtered field Π = lnPs as a function of the time of (a) the truncated field for scales larger than 158 km only and (b) summing over all length scales

Figure 7b similarly shows this norm computed without truncation. Please note that the filter starts at 0300 UTC, as explained above. Considering the interval between 1200 UTC 26 December and 0000 UTC 27 December, it turns out that the information loss in the large scales of the storm is hidden in the small-scale noise. So the norm (13) has to be computed after truncating the spectral coefficients.

The information in Fig. 5 (or Fig. 6) is useful for correctly interpreting the ALADIN forecast on the exemplifying domain indicated on the maps. Indeed, it can be checked whether some values of Πfilt are exceeded in the domain of the coupled model. By testing the values of Πfilt restricted to the coupled-model domain, the person using the coupled model is capable of assessing the information loss in the coupling data before interpreting the forecast of the coupled model. Establishing a critical information loss of, for instance,

 
Πfiltc = 0.003,
(14)

it is decided that such loss propagates through the domain in the model output from 0600 until 1200 UTC. Based on this the user, keeping in mind the time delay due to the phase shift, is warned about a potential lack in the data transfer between the coupling and the coupled model, and more precisely between the coupling time 0300 and 0600 UTC. In that case it is advisable to take the model output of the coupling model more seriously than the higher-resolution coupled model.

b. Systematic tests and operational implications

This subsection shows how the user can be warned automatically in the case of some predetermined information loss and how the coupling updates can be adapted to the meteorological situation to improve the data transfer. This is illustrated by means of the ALADIN-France runs during the month December 1999.

Figure 8 shows the norm (13) computed for the ALADIN-France forecasts during the month December 1999 at the forecast times 0430, 0730, 1030, … , 2530 UTC of the filtered surface pressure with the second-order Butterworth filter, with ωcΔt = 0.9π for Δt = 3 h, Δt = 1 h, and Δt = 30 min.

Fig. 8.

Spectral norms of the filtered field Π = lnPs by the second-order high-pass filters with ωc = 0.9πt and δt = 415.385 s for Δt = 3 h (diamonds), Δt = 1 h (triangles), and Δt = 30 min (squares) at forecast time ranges 0430, 0730, 1030, … , 2530 UTC during Dec 1999

Fig. 8.

Spectral norms of the filtered field Π = lnPs by the second-order high-pass filters with ωc = 0.9πt and δt = 415.385 s for Δt = 3 h (diamonds), Δt = 1 h (triangles), and Δt = 30 min (squares) at forecast time ranges 0430, 0730, 1030, … , 2530 UTC during Dec 1999

Taking a critical value of

 
μc = 0.004
(15)

for the coupling-update interval of Δt = 3 h, exceedance is found 19 times in the five runs listed in Table 1. Of a total of 232 times considered in Fig. 8, this constitutes 8% of the total 3-h model output of the ALADIN-France model. The spectral representation Πfilt is transformed to the gridpoint Πfilt field at these times only.

Table 1.

The five ALADIN-France forecast dates in Dec 1999 during which the critical value μc = 0.004 for the norm of the filtered and truncated surface pressure field Π was exceeded, as given in Fig. 8. The 19 times when this value was exceeded are given in the second column, and the maximum values of μ during these dates are given in the third column

The five ALADIN-France forecast dates in Dec 1999 during which the critical value μc = 0.004 for the norm of the filtered and truncated surface pressure field Π was exceeded, as given in Fig. 8. The 19 times when this value was exceeded are given in the second column, and the maximum values of μ during these dates are given in the third column
The five ALADIN-France forecast dates in Dec 1999 during which the critical value μc = 0.004 for the norm of the filtered and truncated surface pressure field Π was exceeded, as given in Fig. 8. The 19 times when this value was exceeded are given in the second column, and the maximum values of μ during these dates are given in the third column

These five dates, of course, include the forecast based on the 0000 UTC analysis on 26 December. Figure 9 shows the mean sea level pressure of the other runs. It can be seen that all these forecasts contained a storm propagating through the ALADIN-France domain. Note that this also includes a storm during the night between 27 and 28 December 1999, for which increased values of μ can be seen in Figs. 7a and 8. Figure 10 shows the filtered truncated fields Πfilt at the corresponding times. The user of the coupled model on the exemplifying domain, considering the different storms in Figs. 4 and 9, needs to be only concerned about a failing data transfer in the 0000 UTC forecast on 26 December. Indeed, in the other four cases the information loss takes place outside the coupled model domain, as can be seen in Fig. 10. So using the information in Figs. 5, 8, and 10, the user has to consider taking the output of ALADIN-France instead of the coupled model on the indicated coupled domain in the run on 26 December 1999.

Fig. 9.

ALADIN-France forecast of the mean sea level pressure at (a) 1500 UTC 3 Dec 1999, (b) 0900 UTC 17 Dec 1999, (c) 1500 UTC 25 Dec 1999, and (d) 2100 UTC 27 Dec 1999. The contour interval is 2 hPa

Fig. 9.

ALADIN-France forecast of the mean sea level pressure at (a) 1500 UTC 3 Dec 1999, (b) 0900 UTC 17 Dec 1999, (c) 1500 UTC 25 Dec 1999, and (d) 2100 UTC 27 Dec 1999. The contour interval is 2 hPa

Fig. 10.

Filtered Π = lnPs by the second-order Butterworth filter, after truncating the spectral decomposition to scales larger than 158 km, at (a) 1630 UTC 3 Dec 1999, (b) 1030 UTC 17 Dec 1999, (c) 1630 UTC 25 Dec 1999, and (d) 2230 UTC 27 Dec 1999. The contour interval is 0.001 hPa

Fig. 10.

Filtered Π = lnPs by the second-order Butterworth filter, after truncating the spectral decomposition to scales larger than 158 km, at (a) 1630 UTC 3 Dec 1999, (b) 1030 UTC 17 Dec 1999, (c) 1630 UTC 25 Dec 1999, and (d) 2230 UTC 27 Dec 1999. The contour interval is 0.001 hPa

Additionally, as was shown in Termonia (2003), the forecast of the French Christmas storm on 26 December can be improved by using a 1-h coupling-update interval instead of a 3-h interval to give qualitatively comparable results as the forecast of the coupling model ALADIN-France. In fact, the reduction of the upper bound μ when using Δt = 1 h or even Δt = 30 min, can be also seen in Fig. 8. Even if the coupled model cannot improve the ALADIN-France forecast of the storm on 26 December 1999, the information in this figure can be used to keep the information loss below the predefined limit μc by providing the coupling data with Δt = 1 h in the time intervals when this value is exceeded. This avoids some serious underestimation by an unaware interpreter of the coupled-model output.

The field Πfilt in Figs. 5 and 10 can be treated as part of the postprocessing of the ALADIN-France output. It can then be decided operationally at script level to provide coupling data by ALADIN-France to the coupled model with Δt = 1 h for the forecast ranges with exceedance only. For the other forecast ranges, a default 3-h coupling update can still be used. For the illustrated domain on Figs. 5 and 10, this only increases the data transfer for December marginally, that is, for the two or three coupling intervals when the storm crossed the boundary of the exemplifying coupled domain. This represents about 1% of the total data transfer of that month and is technically feasible.

In fact, other technical solutions to deal with the information loss could be investigated. For instance, in the case of the storm entering the exemplifying domain on 26 December 1999, Fig. 5 shows that it was inside the domain at 0600 UTC. The standard practice is to take LAM domains as large as possible hoping the cyclogenesis takes place inside the domain. Based on the operational monitoring it could be decided to provide the 0600 UTC coupling-model state as an “analysis” for the coupled model instead of the analysis at 0000 UTC, to have the storm present in the initial conditions of the coupled forecast.

In this respect it might also be useful to investigate spectral nudging techniques in spectral models of the kind used by von Storch et al. (2000) for downscaling analyses in a regional climate model. Specifically, it could be investigated whether the large scales of the coupled model could be spectrally nudged to the large scales of the coupling model (having the storm inside the coupled area) after the information loss has been detected.2

5. Discussion and conclusions

A strategy was proposed to monitor the data transfer between coupling and coupled models by filtering the forecast of the coupling model during the forecast run. Some critical values can be established to decide whether the data transfer is insufficient. Specifically, this paper studied the values for Πfiltc and μc taken as in (14) and (15). Lower values imply more severely assessing information loss and more overall data transfer, as can be seen from Fig. 8. Note that taking a value of, for instance, μc = 0.001, Δt = 1 h could be taken as a default coupling interval, and Δt = 30 min has to be used five times. Whether this is feasible is a technical question and not a scientific one. The strategy proposed in the present paper allows us to do this in a quantitative way.

The values studied here were fixed with the aim to improve storms of the type in Figs. 4 and 9. To improve or to motivate the choice for Πfiltc and μc more profoundly, the sensitivity of the forecast to lateral-boundary coupling errors corresponding to the data-transfer information loss should be studied as a function of μc and Πfiltc. This lies outside the scope of this paper.

A first attempt to monitor the coupling frequency was done in the paper of Termonia (2003). Compared to this first trial the monitoring strategy proposed here has a number of serious advantages:

  • Flexibility. In Termonia (2003), Δt = 3 h was taken and it was difficult to change such a setup to monitor other coupling-update intervals. In contrast, the coefficients am and bn can be specified as arguments of the model, and the information loss can be estimated at any candidate time step.

  • Robustness. This strategy does not depend on assumptions about the cases to be encountered operationally, as would be needed to apply norms of the differences between interpolation and coupling-model output in Eq. (8) and also for the method proposed in Termonia (2003).

  • Maintenance. The filter is easy to implement—see Eq. (10)—in way that is easy to understand. This is also an advantage for the maintenance of the code.

  • Usage. It is possible to understand the output of the filter in terms of signal handling, whereas norms of differences as in (8) are potentially deceptive depending on the chosen coupling-update interval and the meteorological situation in the forecasts.

  • Independence of model structure. It was argued in the paper that the recursive filter can be implemented on spectral coefficients as well as on gridpoint values.

Applying the filter to gridpoint fields in a gridpoint model has the advantage that the spatial distribution of the information loss is obtained immediately, in contrast to the need of a spectral transform in a spectral model. Removing the small scales is less straightforward in this case. However, looking at Fig. 6 it can be seen that the minimum value of Πfilt for the storm is smaller than −0.008. Elsewhere it does not become smaller than −0.001, except west of Brittany, where a small low present on Fig. 4c is detected. This suggests that for a solution in a gridpoint model, this spatial truncation is not really necessary.

There are no constraints to implement the filter in a global model. For spectral global models the filtering can be done on the spectral coefficients of the decomposition in spherical harmonics. However, in that case computing the norm μ as in the ALADIN-France model may become useless, since μ could depend on storms at the other side of the globe that are of no interest for the regional forecasts under interest. So, the filtered fields should be written to the output files at regular times as part of the model output and be used in the postprocessing to get spatial details about the information loss to be encountered in the created coupling data.

The obtained results can also help establishing the coupling-update frequencies for LAMs in other contexts than operational NWP, in particular, in regional downscaling studies. It could be investigated whether filtering coupling-model fields are useful for establishing the coupling-update intervals for phenomena other than severe storms. This paper did not profoundly treat the choice of the digital filter. The recursive filter was explicitly chosen here to have a minimum of computing costs in the coupling model and to have a good fit of the information loss derived in section 2. Extra research of the benefits and disadvantages of alternative time-symmetric nonrecursive digital filters or recursive filters with nearly constant time delay could be beneficial, even if this might be of little practical importance when it comes to operational implementation details.

Acknowledgments

This research benefited from many useful discussions with F. De Meyer. This work was supported by the Belgian Federal Office for Scientific, Technical and Cultural Affairs (OSTC/DWTC/SSTC) and the ALATNET training network. ALATNET is supported by the TMR/IHP Programme of the European Community, but the information provided here is the sole responsibility of the ALATNET team and does not reflect the Community's opinion.

REFERENCES

REFERENCES
ALADIN International Team
,
1997
:
The ALADIN project: Mesoscale modelling seen as a basic tool for weather forecasting and atmospheric research.
WMO Bull
,
46
,
317
324
.
Errico
,
R. M.
,
T.
Vukićević
, and
K.
Raeder
,
1993
:
Comparison of initial and lateral boundary condition sensitivity for a limited-area model.
Tellus
,
45A
,
539
557
.
Haugen
,
J. E.
, and
B.
Machenhauer
,
1993
:
A spectral limited-area formulation with time-dependent boundary conditions applied to the shallow-water equations.
Mon. Wea. Rev
,
121
,
2618
2630
.
Holtz
,
H.
, and
C. T.
Leondes
,
1966
:
The synthesis of recursive digital filters.
J. Assoc. Comput. Mach
,
13
,
262
280
.
Kulhánek
,
O.
,
1976
:
Introduction to Digital Filtering in Geophysics.
Elsevier, 168 pp
.
Lynch
,
P.
,
1997
:
The Dolph–Chebyshev window: A simple optimal filter.
Mon. Wea. Rev
,
125
,
655
660
.
Lynch
,
P.
, and
X-Y.
Huang
,
1992
:
Initialization of the HIRLAM model using a digital filter.
Mon. Wea. Rev
,
120
,
1019
1034
.
Lynch
,
P.
, and
X-Y.
Huang
,
1994
:
Diabatic initialization using recursive filters.
Tellus
,
46A
,
583
597
.
Lynch
,
P.
,
D.
Giard
, and
V.
Ivanovici
,
1997
:
Improving the efficiency of a digital filtering scheme for diabatic initialization.
Mon. Wea. Rev
,
125
,
1976
1982
.
Shanks
,
J. L.
,
1967
:
Recursion filters for digital processing.
Geophysics
,
32
,
33
51
.
Termonia
,
P.
,
2003
:
Monitoring and improving the temporal interpolation of lateral-boundary coupling data for limited-area models.
Mon. Wea. Rev
,
131
,
2450
2463
.
Vannitsem
,
S.
,
2003
:
Intrinsic error growth in a large-domain Eta regional model.
Mon. Wea. Rev
,
131
,
2697
2704
.
von Storch
,
H.
,
H.
Langenberg
, and
F.
Feser
,
2000
:
A spectral nudging technique for dynamical downscaling purposes.
Mon. Wea. Rev
,
128
,
3664
3673
.
Warner
,
T. T.
,
R. A.
Peterson
, and
R. E.
Treadon
,
1997
:
A tutorial on lateral boundary conditions as a basic and potentially serious limitation to regional numerical weather prediction.
Bull. Amer. Meteor. Soc
,
78
,
2599
2617
.

APPENDIX

The Coefficient of the Recursive Filter

This appendix shows how the coefficients am and bn in Eq. (10) were obtained. For details about computing such coefficients see, for instance, Kulhánek (1976).

The coefficients am and bn of the filter in Eq. (10) are obtained by first writing the response function H(ω) in Eq. (9) as a quotient,

 
formula

where z = eiωδt. To obtain the coefficients am and bn it is useful to consider the poles of the function |H(z)|2,

 
formula

and with

 
formula

Choosing the N poles zn with n = 1, … , N lying outside the unit circle in the complex plane to get a stable filter, the coefficients bn are determined by the relation

 
formula

The coefficients am with m = 0, … , N are obtained by

 
formula

with Q a normalization factor such that |H(z = −1)| = |H(ω = π/δt)| = 1.

It is important to know that recursive filters generate a phase shift. By construction H(z) has a frequency-dependent phase ϕ(ω)

 
formula

Taking the convolution with the signal, this means that eiωt is multiplied by this complex factor. Such a phase shift causes a time delay ϕ(ω)/ω. It is more useful to consider the group delay

 
formula

expressed as a fraction of the time interval Δt, showing the backward time shift of the signal as a function of the frequency ω.

Footnotes

Corresponding author address: Piet Termonia, Royal Meteorological Institute of Belgium, Ringlaan 3, B-1180 Brussels, Belgium. Email: Piet.Termonia@oma.be

1

See Lynch and Huang (1992), Lynch (1997), and Lynch et al. (1997) for details in the ALADIN context.

2

In the operational ALADIN setup, the coupling files contain the spectral coefficients. So the information is present for the complete domain of the coupled model.