# Data Assimilation of Tropopause Height Using Dry Intrusion Observations

Yann Michel

► To cite this version:

Yann Michel. Data Assimilation of Tropopause Height Using Dry Intrusion Observations. Monthly Weather Review, 2010, 138, pp.PP. 102-122. 10.1175/2009MWR2912.1. meteo-00460309

HAL Id: meteo-00460309

https://meteofrance.hal.science/meteo-00460309

Submitted on 2 Nov 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Data Assimilation of Tropopause Height Using Dry Intrusion Observations

YANN MICHEL

CNRM-GAME and CNRS, Météo-France, Toulouse, France

(Manuscript received 15 December 2008, in final form 22 May 2009)

## ABSTRACT

This article investigates the problem of initializing upper-level potential vorticity by using the detection of dry intrusions that can be seen in water vapor images. First, a satellite image processing technique has been developed for the identification and tracking of dry intrusions on geostationary satellite images. This technique can also be applied to images derived from model fields through a radiative transfer model. A linking algorithm automatically compares the trajectories of the dry intrusions in the model and in the satellite images. Differences of brightness temperatures are then converted to differences of tropopause height through a simple linear model, which is based on the correlation found in the background. As the scheme is likely to provide observations of the tropopause height, it also suggests that a space-alignment representation of the errors be used. A simple one-dimensional study provides a depiction of the background error covariance in alignment space, which is compared to the traditional approach of background error covariance in amplitude space. An approximate form of the Ertel potential vorticity operator is then used to incorporate pseudo-observations inside a global four-dimensional variational assimilation scheme. A case study of cyclogenesis highlights the difference between the amplitude-space assimilation of potential vorticity values and the alignment space assimilation of the height of the tropopause.

## 1. Introduction

Data assimilation combines imperfect model short-range forecasts with uncertain data, acknowledging their respective uncertainties. Modern data assimilation algorithms are able to cope with a high flux of observations and correct the background with time-consistent and flow-dependent increments. However, the current representation of errors is based on an additive representation (Lawson and Hansen 2005) and is fundamentally restricted by the truncation of the probability distributions of errors to their first two moments (this is known as the Gaussian assumption). Such an approach can perform poorly when locations of coherent structures in the atmospheric flow are horizontally or vertically displaced from their observations (Ravela et al. 2007).

Global models sometimes fail to predict midlatitude cyclogenesis, even if suspected upper- or lower-level precursors are visible in the images from geostationary satellites. Different operational procedures, often based on potential vorticity inversion, have been developed to

exploit the link that can be seen between water vapor (WV) images and the initial state. Cyclogenesis has been shown to be sensitive to potential vorticity structures from different spatial scales (Beare et al. 2003). At the same time, potential vorticity (PV) can be inverted to mass and wind fields, provided a balance relation and a boundary condition (Davis 1992; Arbogast et al. 2008). The correction of the initial state has often been performed through the following:

* a comparison between satellite WV images and PV background fields,

* a definition of a correction in PV space, and

* an initialization through PV inversion.

This methodology was first described and applied by Demirtas and Thorpe (1999) and further studied by Swarbrick (2001).

The comparison between satellite WV images and PV background fields is based on the PV-WV relationship. The common statement (Røsting et al. 1996; Mansfield 1996; Browning 1997; Demirtas and Thorpe 1999; Swarbrick 2001, among others) is based on a thermal relationship interpretation. A stratospheric air mass is characterized by low humidity, warm potential temperature, and high potential vorticity content. As a result,

Corresponding author address: Yann Michel, CNRM/GMAP, Météo-France 42 av. G. Coriolis, 31057 Toulouse CEDEX 1, France. E-mail: yann.michel@meteo.fr

DOI: 10.1175/2009MWR2912.1
these first two characteristics cause PVs to appear dark in WV images, and the last characteristic results in a local maximum in the PV fields. Therefore, the imprecise nature of the PV-WV relationship is mainly attributed to the occurrence of tropopause folds, air mixing, and sometimes upper-level clouds such as cirrus (Swarbrick 2001). Demirtas and Thorpe (1999) suggested that one has to avoid situations where there is a cutoff low, and polar or tropical regions, where the WV radiance is more likely to be determined by temperature. At this point, an objective determination of the situations where the PV-WV relationship holds (or not) is still lacking, which is one of the main drawbacks of this methodology (Swarbrick 2001).

Wirth et al. (1997) have developed a conceptual model for this relationship, based on the quasigeostrophic dynamics of Juckes (1994). They show that the thermal aspect of the PV-WV relationship may be roughly parameterized as a link between the tropopause height perturbation $\delta_z$ and the brightness temperature perturbation $\Delta T_{BB}$ [appendix, Eq. (A3)]:

$$ \Delta T_{BB} \approx \gamma_1 \delta_z, \quad \gamma_1 = -4 \text{ K km}^{-1}. \quad (1) $$

However, Wirth et al. (1997) explain that considering only the thermal aspect of the PV-WV relationship is insufficient. The WV brightness temperature can be closely related to an isosteric surface (Ramond et al. 1981), which is seriously impacted by vertical motion. Therefore, dry intrusions are likely to be the consequence of vertical motions that arise from the past. This would imply the following two consequences on the PV-WV methodology:

* it is crucial to consider a whole sequence of images to try to relate darkening in the image to PV increase through vertical motion;

* local displacements between PV fields and WV images should be interpreted with great care as they can be the consequence of differences between the motions of the isosteric surface and tropopause height.

As shown by Arbogast (1998), an isolated PV anomaly is associated with vertical motion whose sign depends on the larger-scale wind environment. Therefore, the WV image darkening would depend on the larger-scale wind structure. Demirtas and Thorpe (1999) and Santurette and Georgiev (2005) mention that the PV-WV method does not apply to cutoff systems, because of air mixing, but we suggest that it could also be because of the dependence on vertical motion over the deformation of the larger-scale wind field. A diagnostic approach that statistically compares tropopause height and brightness temperatures over the trajectory of the dry intrusion will

be used in this study. Another caveat of the PV-WV methodology relates to the tropopause fold. The bidimensional WV image cannot represent the vertical tilt of PV anomalies, which are sometimes overlaid by clouds. The approach described here does not handle these important cases. In particular, in the view of Wirth et al. (1997), the darkening signal in WV images is linked with vertical motion. Therefore, cases involving both upper-level and lower-level anomalies, or tropopause folds, may have to be discarded in a first step despite their significance.

Referring to the second point, Hello and Arbogast (2004) pointed out that WV and PV fields may be sufficiently different such that the first comparison should be done in image space, using a radiative transfer model to generate synthetic images. Swarbrick (2001) found the use of synthetic images of little practical use, but pointed to the inadequacy of the radiative transfer code as a possible cause for his statement.

The problem of initializing the PV structures with the help of images can be split into these three updated points:

* to automatically compare the model and satellite images and to observe differences in image space;

* to convert this measure to observations of the characteristics of the PV structures (e.g., the tropopause localization on the horizontal and on the vertical);

* and to incorporate the modifications to a new analysis, taking into account the error made in both the observation and the background.

Sometimes, however, forecasters are able to provide an estimate of the position of the tropopause (Santurette and Georgiev 2005; Guérin et al. 2006). Even in this case, there is still a lack of knowledge on the structure of the increments to provide a new initial state. Authors often rely on the imagery to specify the horizontal shape of the PV increments, sometimes with the help of data assimilation structure functions (Guérin et al. 2006). Verkley et al. (2005) use instead a so-called influence function, which is isotropic and nearly Gaussian-shaped, to incorporate PV modifications.

In the vertical, Verkley et al. (2005) used a nearly uniform, or barotropic, structure function. Manders et al. (2007) used the same formulation in the horizontal, but further constrained the vertical increment to a set of levels centered around the tropopause. Swarbrick (2001) extends the modification to a set of levels around the isobaric level where the PV distribution matches the WV images (typically between 300 and 500 hPa). Røsting and Kristjánsson (2006) suggest the use of singular vectors to constrain the modifications in the vertical. There is clearly an emerging point of view that PV increments—at least
in the vertical—should be defined using more adequate structure functions.

Guérin et al. (2006) suggested that the PV inversion procedure could be replaced by the assimilation of pseudo-observations. He therefore relied on the structure functions of the data assimilation scheme by specifying pseudo-observations only at the dynamical tropopause (e.g., the 1.5-PVU<fnref n="1" /> surface). The covariances of observation errors, known as the R matrix in data assimilation, would provide the appropriate framework to specify the uncertainty associated with PV modifications. Many methods indeed provide estimates of observation error variances (Desroziers and Ivanov 2001; Chapnik et al. 2004, 2006; among others). When using three-dimensional PV inversion, one may rely on ensemble methods to incorporate this uncertainty (Plu and Arbogast 2005; Manders et al. 2007).

The goal of this work is to provide additional improvement to the PV-WV methodology by trying to make it more objective. Wirth et al. (1997) suggested that the WV image is a tracer of the midtroposphere. Therefore, it is likely that the PV-WV relationship requires a Lagrangian viewpoint, as also suggested in the case study by Georgiev (1999). The first section introduces the tracking algorithm of dry intrusions on water vapor images. Only the recent improvements over the version fully described in Michel and Bouttier (2006) will be reported. There is a lack of knowledge on the structure of PV increments, and especially the vertical structure. We use the tracking algorithm to study the vertical structure of potential vorticity errors in dry intrusion regions. As WV images are likely to provide measurements of the height of the tropopause, or of its variation, rather than measurements of values of PV on levels, we also develop a representation of errors in position space, following the ideas of Lawson and Hansen (2005). The last section is a case study of cyclogenesis applying the preceding concepts, and mainly compares the vertical structure of PV increments using an additive representation of errors that assimilates PV values to the alignment model of errors that assimilates the height of the tropopause.

## 2. The tracking algorithm of dry intrusions

The tracking algorithm is very similar to the description given by Michel and Bouttier (2006). It is based on an adaptive thresholding technique from Morel and Sénési (2002), and additional information can be found in these references.

### a. A short description

We use WV images from the spinning enhanced visible and infrared imager (SEVIRI) onboard the Meteosat Second Generation (MSG) satellite. The algorithm iteratively applies thresholds for the brightness temperature of WV images from $T_{min}$ to $T_{max}$ at every $\Delta T = 0.5$ K (e.g., the temperature resolution of SEVIRI in WV channels). Values for $T_{min}$ and $T_{max}$ are given for several images by Michel and Bouttier (2006). Each connected set of pixels, referred to as a cell, is selected only if it has a greater surface S than a threshold (3000 km²) and if it has a relative maximum of brightness temperature:

$$T_{Cell} - T_{i} \geq \Delta T_{d}, \quad (2)$$

where $T_i \in [T_{min}, T_{max}]$ is the iterative threshold temperature, and $T_{Cell}$ is the maximum temperature of the pixels covering the detected cell. The parameter $\Delta T_d = 2.5$ K is close to the one of Michel and Bouttier (2006) and well suited for MSG images. Cells are required to be in a $\Delta T_d$ deep temperature tower. They are then grouped into trajectories using instantaneous estimates of their speeds. An overlapping criterion between cells at successive times links cells with primary or secondary links.

### b. Improvements in the selection procedure

Michel and Bouttier (2006) showed that the tracking algorithm was likely to provide good probability of detection of the relevant dry intrusions, but that additional filters were necessary to screen spurious features. They relied on the following characteristics of the tracking algorithm:

* the lifetime,

* the temperature evolution in trajectories, and

* the position relative to the jet stream.

The latter criterion employed the wind at upper levels from background fields. This wind was then transformed using a mathematical morphology skeleton operator, which extracts the line of the jet streaks. The filter uses the relative position of a cell with respect to this jet streak. External dynamical information indeed proved to be necessary to reduce the false alarm rate. If background fields happen to be inaccurate enough to misrepresent the jet streak, the dry intrusions may be inappropriately filtered out. Unfortunately, this procedure is necessary, and we may reconsider a purely observation-based representation of the jet streak when high quality wind data become available near the tropopause. This could happen with the development of active Doppler Wind lidars, such as Atmospheric Dynamics Mission (ADM) Aeolus mission (Stoffelen et al. 2005). Two new selection criteria involve a morphological test and a test based on the distribution of temperature around the cell.

<footnote n="1">¹ Following Hoskins et al. (1985), a convenient unit for PV is 1 PVU = $10^{-6}$ m² s⁻¹ K kg⁻¹.</footnote>
<figure><img src="image_2.png" /><figcaption>FIG. 1. The distribution of the minimum temperature, $T_{min}$, within disk B (shaded). The closed solid contour delineates the cell as detected by the tracking algorithm. The closed dashed contour delineates the convective system overhead. Also shown are the warm intrusion area (generally larger than the cell) and the position of the jet streak (arrows).</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>FIG. 2. FAR and POD (%) for different tunings of the filter based on the distribution of temperature around the cell. By increasing FAR, the curves correspond to reference cloudy temperatures every 2 K (starting with dots and $T_{cld}$ = -48°C and ending with crosses and $T_{cld}$ = -58°C). Each curve is generated by having different values for the threshold quantile (between 0 and 100% every 10%, all starting at the same point for $Q_{th}$ = 0, and in decreasing POD for increasing values of $Q_{th}$).</figcaption></figure>

## 1) A MORPHOLOGICAL TEST

This test insures that cells are close enough to ellipticity. This is implemented as a test on the ratio between the surface S of the cell and the surface of the approaching ellipse (from axes a and b):

$$ \frac{S}{\pi ab} \geq r_{th} = 0.05. \qquad (3) $$

This filter proved necessary to eliminate some cells that have very particularly elongated shapes, and which are not suited to the data assimilation procedure described in section 3. This filter has a minor impact. It ensures that detected cells have aspect ratios that are reasonable and that will permit linking between model and satellite cells (see below). The filter typically screens out only 3% of the cells.

## 2) THE DISTRIBUTION OF TEMPERATURE AROUND THE CELL

This kind of filter was mentioned by Michel and Bouttier (2006) as necessary, but not straightforward to implement as it typically requires information outside of the cell. The environment of the cell is checked to detect cells that are false alarms because they are “holes” in cloudy areas. For each point of the contour, we define a disk “B” surrounding the contour point, and compute minimum temperature $T_{min}$ within this disk as depicted in Fig. 1. Defining $T_{cld}$ as a reference cloudy temperature, the test is implemented on the quantile of distribution of minimum temperature $Q_{T_{min}}$ to detect cells embedded in cloudy environments:

$$ Q_{T_{min}} < Q_{th}. \qquad (4) $$

Of course, there may be different choices for the threshold $Q_{th}$ and of the cloud temperature $T_{cld}$ that may yield similar results. One therefore needs to compute the probability of detection (POD) and the false

alarm rate (FAR) for each sample value of ($T_{cld}$, $Q_{th}$). This is done in Fig. 2. The curves do not generally cross over, which means that for a given POD, the lowest FAR is obtained with the warmer value for the reference cloudy temperature. Using $T_{cld}$ = -50°C and $Q_{th}$ = 0.1 is considered a good trade-off between FAR (reaching 0.77) and POD (near 0.13). This value of POD means that typically a dry intrusion will only be selected after a few hours of life. Strong cyclogenesis however is associated with a more intense signal in the WV images, such that a POD computed only over strong cyclogenesis events would be higher as shown by Michel and Bouttier (2006). The FAR is lowered as well by the additional constraint that a dry intrusion is required to be detected in both satellite and model images (see below).

Figure 3 shows the visualization tools developed for the tracking of dry intrusions on satellite images, and the effect of the filtering procedure for a particular case.

### c. Comparing model and satellite image processing results

The methodology intends to make a direct comparison of the characteristics of the tracking algorithm applied to satellite and model images.<fnref n="2" /> One shall have, in the data assimilation terminology, innovations of positions, or of brightness temperatures, of the detected

<footnote n="2">² This should be understood as images produced by Radiative Transfer for Tiros Vertical Sounder (RTTOV 8; Saunders et al. 1999) from background fields.</footnote>
<figure><img src="image_4.png" /><figcaption>FIG. 3. Tracking of dry intrusions on WV MSG images (6.2-µm SEVIRI channel). The contours of the cells are green solid lines; the trajectories of the gravity center of the cells are purple lines. The black arrow shows the estimated displacement speed. The black line indicates the projection on the jet streak, as it is derived from the wind background, see Michel and Bouttier (2006). The red arrows show the wind at this projection. (top) Raw results and (bottom) with relevant filtering procedures. Case at 1200 UTC 1 May 2007 (the tracking algorithm has cycled for more than a week).</figcaption></figure>
cells. An automatic linking between model ($C_m$) and satellite cells ($C_s$) is needed. The one used here relies on a distance between two cells $d(C_m, C_s)$. One first builds a link for all the cells in a model image with their closest cells in the satellite image. Then one builds a second link for all the cells in a satellite image with their closest cells in the model image. The final link is created if and only if these last two links prove to be the same, which means that the satellite and model cells are the closest cells one to another. An additional time constraint is imposed, requiring the link to last for at least three hours. Two kinds of distances have been tested; the first one is a classic geometrical distance on the sphere. The second one, analogous to the Mahalanobis distance, uses the relative dimensions of the cells in the direction of displacement as a weight of the geometrical distance. This latter formulation proved to be more robust and provides reasonable linking between model and satellite image processing results. Objective evaluation of the performance of the linking has not however been performed at this stage of the study.

### 3. The vertical structure of potential vorticity errors in amplitude and position spaces

#### a. Alternative error models

##### 1) AMPLITUDE AND ALIGNMENT ERRORS

In the classical framework of amplitude errors, the background $\mathbf{x}_b$ is taken as being an additive perturbation of the truth $\mathbf{x}_i$; for example,

$$ \mathbf{x}_b = \mathbf{x}_t + \varepsilon_b, \qquad (5) $$

where $\varepsilon_b$ is the background error. The analyzed state $\mathbf{x}_a$ is naturally defined as an additive combination of the background and the innovation $\mathbf{d}$:

$$ \mathbf{x}_a = \mathbf{x}_b + \mathbf{Kd}, \qquad (6) $$

where the innovation vector $\mathbf{d}$ is the difference between the observations and the background projected into observations space. However, one could imagine searching for an analyzed state that would be a spatial distortion of the background; for example,

$$ \mathbf{x}_a = \mathbf{x}_b(\mathbf{s} + \mathbf{Kd}). \qquad (7) $$

In this latter equation, $\mathbf{s}$ is the spatial grid, and $\mathbf{d}$ represents innovations in alignment (or position) space. There is a need for a spatial interpolation in this representation. This is known as an alignment error model, following the terminology of Lawson and Hansen (2005). It is less general than the additive error model, as a bound on the

background would imply a similarly bounded analysis. Classical data assimilation schemes rely on the general additive error model of Eq. (6) and the hypothesis of Gaussian distribution of errors (or more generally they truncate the distribution of errors to their first two moments). The background error $\varepsilon_b$ in Eq. (5) generally has a Gaussian distribution when the error is the sum of numerous independent random errors, referring to the Central Limit Theorem. The alignment error model may be written

$$ \mathbf{x}_b = \mathbf{x}_t(\mathbf{s} + \varepsilon_D), \qquad (8) $$

where $\varepsilon_D$ is the alignment error. It is likely that the additive error $\mathbf{x}_b - \mathbf{x}_t$ has a non-Gaussian distribution (Ravela et al. 2007). Data assimilation algorithms that rely on the Gaussian hypothesis may exhibit unacceptably distorted analysis states when the background exhibits large scale spatial errors, as shown by Lawson and Hansen (2005), Ravela et al. (2007), Chen and Snyder (2007), and Beezley and Mandel (2008).

##### 2) LINEARITY CONDITION

Considering the mixed error model:

$$ \mathbf{x}_b = \mathbf{x}_t(\mathbf{s} + \varepsilon_D) + \varepsilon_A, \qquad (9) $$

where $\varepsilon_D$ and $\varepsilon_A$ are respectively the alignment and amplitude errors, it is possible to derive a linearity condition (Lawson and Hansen 2005). The linearization of (9) yields

$$ \mathbf{x}_b = \mathbf{x}_t + \frac{d\mathbf{x}_t}{ds}\varepsilon_D + \varepsilon_A + O(\|\varepsilon_D\|), \qquad (10) $$

such that if the alignment error $\varepsilon_D$ is Gaussian, the additive error $\varepsilon_b$ is Gaussian as well, and its covariance matrix is

$$ \mathbf{B} = \mathbf{E}(\varepsilon_b \varepsilon_b^T) \qquad (11) $$

$$ \begin{aligned} \mathbf{E}(\varepsilon_A \varepsilon_A^T) &= \varepsilon_A \varepsilon_D^T \frac{d\mathbf{x}_t^T}{ds} + \frac{d\mathbf{x}_t}{ds} \mathbf{E}(\varepsilon_D \varepsilon_A^T) \\ &\quad + \frac{d\mathbf{x}_t}{ds} \mathbf{E}(\varepsilon_D \varepsilon_D^T) \frac{d\mathbf{x}_t^T}{ds}. \end{aligned} \qquad (12) $$

The linearization is valid if the second term in the Taylor expansion in Eq. (10) can be neglected. Using the scaling $s = L\tilde{s}$ and $x_t = X\tilde{x}$ and scalar-valued fields, one can write the following:

$$ \frac{\varepsilon_D}{L} \ll 2 \left( \frac{d\tilde{x}/d\tilde{s}}{d^2\tilde{x}/d\tilde{s}^2} \right), \qquad (13) $$

which shows that the alignment error relative to the length scale of the field has to be small with respect to relative spatial variations of the field.
### 3) ASSIMILATION IN POSITION SPACE

When the distribution for $ε_D$ is thought to be Gaussian, it is possible to use the framework of the Kalman filter in alignment space, for example,

$$ \mathbf{x}_a = \mathbf{x}_b (\mathbf{s} + \delta \mathbf{s}_a), \qquad (14) $$

$$ \delta \mathbf{s}_a = \mathbf{DH}_d^T (\mathbf{H}_d \mathbf{DH}_d^T + \mathbf{R}_d)^{-1} (\mathbf{y}_d - \mathbf{H}_d \mathbf{x}_b), \quad (15) $$

where $δs_a$, called the analyzed vector displacement, is the solution of a least squares problem where the weights are defined by:

* **D**, the matrix of covariances of errors in position space;

* **R_d**, the matrix of covariances of observation errors with errors defined in position space;

* **H_d**, the (observation) operator that observes the position in the model state; and

* **y_d**, the vector of observations (of the positions of a structure).

## b. A simple study of the vertical structure of PV errors

### 1) METHODOLOGY

This section is designed to improve our knowledge of PV errors in the regions of dry intrusions. Outside of an idealized framework, samples of background error realizations are lacking. Two widespread techniques have been used to estimate background error covariance matrix, the so-called National Meteorological Center [NMC, now known as the National Centers for Environmental Prediction (NCEP)] method (Parrish and Derber 1992), and the use of an ensemble of data assimilation and prediction cycles (Belo-Pereira and Berre 2006). The first method uses differences between forecasts valid at the same time, but for different ranges, and is known to be deficient, for example, in data rich areas. The second method provides adequate samples of background errors, when the perturbations introduced in the ensemble have the proper statistics (Belo-Pereira and Berre 2006).

We shall in any case rely on an ergodic assumption, as in the NMC method, mainly because it is convenient. This approximation is believed to be relevant given the preliminary nature of this work. As the focus is on the vertical structure of PV errors within dry intrusions, we extract PV profiles of the cells detected by the tracking algorithm on a sample of model images. We then make the hypothesis that the real background errors in PV have a similar error distribution to 6-hourly differences of PV profiles along the trajectories of the tracking algorithm. Figure 4 shows a sample of $ε_b$ and $ε_D$ background errors. On the left, the two PV profiles are extracted along the same trajectory of the tracking algorithm but are 6 h apart. During this time, the tropopause has lowered, but the lower-level anomaly is

<figure><img src="image_5.png" /><figcaption>FIG. 4. (a) Two PV profiles whose differences are used to compute the background error covariances. (b) The difference between these two profiles in amplitude space. (c) The difference between these two profiles in position space.</figcaption></figure>

a little bit higher. This is interpreted in amplitude space as a strong increase of PV at 200 hPa, and a dipolar structure aimed to move the lower-level anomaly around 600 hPa. The error can be fully described in position space using an interpolating operator. Thus, the tropopause displacement is approximately 40 hPa in the upper levels. The alignment error is negative at 600 hPa, in agreement with the structure of the PV profiles. Note that this methodology uses the same number of degrees of freedom for the background errors in amplitude and position space.

### 2) THE STRUCTURE OF ERROR COVARIANCES

The background error covariance matrices are then simply sampled through

$$ \mathbf{B} = (\overline{\varepsilon_b - \bar{\varepsilon}_b})(\overline{\varepsilon_b - \bar{\varepsilon}_b})^T, \qquad (16) $$

$$ \mathbf{D} = (\overline{\varepsilon_d - \bar{\varepsilon}_d})(\overline{\varepsilon_d - \bar{\varepsilon}_d})^T, \qquad (17) $$
<figure><img src="image_6.png" /><figcaption>FIG. 5. Background error standard deviation in (a) amplitude and (b) alignment as a function of pressure (hPa).</figcaption></figure>

where the overbar refers to the average over the 1825 samples of errors (which is about 10 times larger than the dimension of the matrices). Figure 5 shows the standard deviations of background errors in amplitude and alignment spaces. A striking feature of the two profiles of the standard deviations is that they exhibit opposite behavior: in amplitude, $\sigma_b$ increases with height (except a local minimum at 200 hPa, which may be due to the occurrence of tropopause folds within the samples, with weaker gradients of PV around this altitude), whereas $\sigma_d$ generally decreases with height. Large vertical displacements errors ($\sigma_d \sim 200$ hPa) seem to occur in the troposphere. Looking at individual profiles associated some cases of large displacement errors, it is apparent that the alignment error model sometimes tries to convert small amplitude errors to excessively large alignment errors. This could be understood as the reverse version of statements by Lawson and Hansen (2005) and Ravela et al. (2007) for the specific problem of the PV vertical structure: amplitude errors (alignment errors) following a Gaussian distribution to induce a non-Gaussian distribution in

alignment (amplitude) space. This will be discussed further after the description of correlation structures.

Figure 6 shows the structure of background error correlations in amplitude and alignment spaces. The structure functions for PV, as deduced from the sample, exhibit a rather classic small anticorrelation between the upper levels (100 to 500 hPa) and the lower levels (600 to 800 hPa; Fig. 6a). The correlation length scale looks smaller in the troposphere than in the stratosphere (but this may be a consequence of the choice of pressure as a vertical coordinate). In alignment space, one finds, again on the contrary, broader correlations in the lower troposphere. There is a localized maximum of correlation length scale toward the average position of the tropopause (350 hPa).

### 3) AN EXAMPLE OF ASSIMILATION

The lines of the **B** matrix provide direct meaningful interpretation, as the increment when assimilating an observation is proportional to it. It is more complicated for the **D** matrix, as its lines describe the covariances in

<figure><img src="image_7.png" /><figcaption>FIG. 6. Correlation matrices for the background errors (a) amplitude and (b) alignment space.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>FIG. 7. Assimilation of a perfect observation (circle). Background, and analyzed state amplitude and alignments shown as solid, and dotted and dashed lines, respectively: (a) 6 PVU at 300 hPa and (b) 3 PVU at 550 hPa.</figcaption></figure>

alignment space. Therefore, the increment when assimilating an observation (of position) has its shape determined by a regularized displacement of the background [e.g., Eqs. (14)–(15)]. The consequence of the assimilation of an observation is shown in Fig. 7. The observation is supposed to be perfect; for example, $\mathbf{R} = 0$ (in amplitude space) and $\mathbf{R}_d = 0$ (in alignment space), and consist of a PV value for the amplitude space (and the position of the observation is taken into account through the observation operator), and of a position value for the alignment space (where the PV value of the observation is taken into account through the observation operator).

When the displacement is small (Fig. 7a), the two methods provide comparable results. When assimilating in amplitude space, the analysis is the sum of the background and of an increment, which can be considered an anomaly yielding to an inversion of PV toward 250 hPa. In alignment space, one can see that the monotonic behavior of PV with height is conserved in the analysis, providing maybe a more physically appealing solution, as the general structure of the background is conserved. Going to larger displacements (Fig. 7b), it becomes clear that the two algorithms diverge. The amplitude error model builds a strong lower-level PV anomaly, whereas the alignment error model lowers the PV profile. This is thought to be an indication that the alignment error model should be used only with care in the troposphere. It is very likely that the PV total error is the combination of the tropopause alignment error and of an amplitude error, and that neither of the two error models is fully able to represent the statistics of the total error well. When the alignment and amplitude error are separable; for example, $\mathbf{E}(\varepsilon_D \varepsilon_A^T) = 0$, Eq. (12) reads

$$ \mathbf{E}(\varepsilon_b \varepsilon_b^T) = \mathbf{E}(\varepsilon_A \varepsilon_A^T) + \frac{d\mathbf{x}_l}{ds} \mathbf{E}(\varepsilon_D \varepsilon_D^T) \frac{d\mathbf{x}_l^T}{ds}, \quad (18) $$

which states that the alignment error is likely to inflate the variances and extend the correlation, and this is re-

ciprocal. As a result, one may interpret the above description as the amplitude error model (the alignment error model) being likely to be more consistent with the real distribution of errors in the troposphere (the stratosphere).

### 4) GAUSSIANITY OF THE ERRORS

We can give some more weight to the previous statement by examining the total errors at different levels, shown in Fig. 8, and overlaid with their Gaussian approximation in amplitude and alignment spaces. An obvious feature of these distributions is their relative departure from the Gaussian case. In particular, the distribution of $\varepsilon_D - \bar{\varepsilon}_D$ at 550 hPa shown in Fig. 8d has a tripolar structure, with an overoccurrence of strong departures of positions, as already mentioned. This can be quantified, of course, through standard measurements of non-Gaussianity, such as kurtosis and skewness. This does not indicate however which error model is most suited for our purpose. The figure suggests that the alignment error model may fail to depict fairly the total error at lower levels, and that the amplitude error may be inappropriate to fairly depict the total error in the upper levels, such that a mixed—also called “two steps”—error model might be useful. These conclusions are preliminary and could be improved by the use of background error samples derived from ensemble methods.

### c. Conclusion

There is an apparent mismatch between the vertical structure of PV modifications made on one hand by Guérin et al. (2006) and on the other hand by Verkley et al. (2005), Røsting et al. (2003), Hello and Arbogast (2004), and Røsting and Kristjánsson (2006). Guérin et al. (2006) use 4D-Var assimilation of the 1.5-PVU value at the estimated tropopause pressure levels. The latter authors use PV tridimensional inversion with a much deeper vertical structure. The main drawback of
<figure><img src="image_9.png" /><figcaption>FIG. 8. Distributions of background errors in (left) amplitude and (right) alignment spaces at (top) 300- and (bottom) 550-hPa levels. The solid curve shows the approximate Gaussian distribution with obtained from the sample and used in the B or D error covariance matrices.</figcaption></figure>

PV inversion is that the observation error is not taken into account, and that ensemble methods have to be defined to evaluate this uncertainty (Plu and Arbogast 2005; Manders et al. 2007). We suggest that these differences may be explained through the particular nature of the tropopause, which is nearly a materially conserved surface (Wirth et al. 1997). An alignment model may therefore be more suited for the assimilation of the tropopause height. A preliminary study using PV profiles extracted from a dry intrusion shows that the vertical structure of both models is indeed quite different, and generally deeper for the alignment model than for the additive model.

## 4. Tropopause pseudo-observations and their error standard deviation

### a. Automatic linking between PV structures and WV cells

As shown by Wirth et al. (1997), the WV image is below the tropopause (e.g., their Figs. 3 and 6), yet it evolves under a different advection wind than the one at the tropopause level. This induces a possible shift between WV and PV fields. This horizontal shift is sometimes combined with a vertical shift, when PV structures are tilted in the vertical, as it often happens in so-called tropopause folds. Therefore, an explicit operator linking

PV structures and WV cells was needed. This operator is taken as a simple horizontal match between a WV cell and the closest structure of PV, defined as being a local maximum by an amount of at least 10% with respect to the environmental value of PV (defined as the mean value of PV along the contour of the displaced WV cell). The operator is computed for each vertical level, without checking for time or vertical consistency. This is a very basic first approach, using standard methods from image processing tools.

### b. The potential vorticity operator within 4D-Var

The PV observation operator, based on a simplified form of Ertel PV, has been implemented with its tangent-linear and adjoint versions by Guérin et al. (2006) into the Météo-France global spectral model, Action de Recherche Petite Echelle Grande Echelle (ARPEGE; see Courtier et al. 1991).

The expression for Ertel PV is

$$ \varepsilon_{\text{PV}} = \frac{1}{\rho} \zeta_a \cdot \nabla \theta, \quad (19) $$

where $\rho$ is density, $\zeta_a$ the vertical component of absolute vorticity vector, and $\theta$ potential temperature. The operator is based on the low Rossby number approximation with the hydrostatic assumption:
$$ \mathcal{H}_{\text{PV}} = -g\xi_a \frac{\partial\theta}{\partial p} - g\frac{fp}{R} \left(\frac{p_0}{p}\right)^{R/C_p} \left[ \left(\frac{\partial U}{\partial p}\right)^2 + \left(\frac{\partial V}{\partial p}\right)^2 \right], \quad (20) $$

where g is gravity, f the Coriolis parameter, R the gas constant, C<sub>p</sub> the specific heat of air at constant pressure, p the pressure, and p<sub>0</sub> a reference pressure. This later expression has the advantage of only involving vertical derivatives of θ, U, and V, and is thus easier to implement in the present ARPEGE formulation of observation operators as stated by Guérin et al. (2006). We used the ARPEGE configuration that was operational at the moment of the case study; the resolution reaches T358 C2.4 for the forecasts, and T107 C1.0, T149 C1.0 for the resolution of the increments in the two loops of the minimization. The “C” denotes the stretching factor of the ARPEGE grid (Yessad and Bénard 1995), which allows the grid resolution to reach T850 over France (but only T150 over New Zealand). The assimilation is performed in the unstretched space. The model has 46 vertical levels.

### c. The PV–WV relationship and PV observations

There have been several attempts to derive a statistical relationship between PV and WV. For example, Georgiev (1999) performed regressions between PV on isobaric levels and WV radiances for two trajectories of dry intrusions. The correlation is maximized near 500 hPa, but only reaches 0.6. Juckes (1994) introduced an idealized framework for the dynamics of the tropopause. The tropopause is taken as a material surface between two volumes of homogeneous quasigeostrophic PV, with one value in the troposphere and a different value in the stratosphere. The tropopause displacement δ<sub>z</sub> from its unperturbed position is then linearly related by Juckes (1994) to the potential temperature anomaly at the tropopause θ'<sub>TP</sub> through

$$ \delta_z(x) = \frac{g}{N_s N_t \theta_{00}} \theta'_{\text{TP}}(x), \quad (21) $$

where N<sub>t</sub> and N<sub>s</sub> are the Brunt–Väisälä frequencies of the troposphere and stratosphere, respectively; g is the gravitational acceleration; and θ<sub>00</sub> a reference value of the potential temperature. This is derived in the case when the boundary influence can be neglected.

Wirth et al. (1997) further extended this conceptual model to WV images. With the broad assumption that anomalies of brightness temperature represent tropopause temperature anomalies θ'<sub>TP</sub>, one could linearly relate δ<sub>z</sub> to the variation of brightness temperature ΔT<sub>BB</sub> defined in Eq. (1). It may be more convenient to write it in its inverse form and in pressure coordinates:

$$ \delta_p \approx \Gamma_1 \Delta T_{\text{BB}}, \quad (22) $$

where Γ<sub>1</sub> = 13 hPa K<sup>-1</sup>, as shown in the appendix, using the same parameters as in Juckes (1994) and Wirth et al. (1997). However, Wirth et al. (1997) actually showed that the PV–WV relationship is much more complex, even under the quasigeostrophic assumption. It may be scale-dependent if there is a lower-level PV anomaly (Juckes 1994).

By running the tracking algorithm on both satellite and model images, PV profiles can be extracted at the location of the cells that have been detected on model images. It is therefore possible to study linear regressions between background samples of brightness temperatures T<sub>BB</sub> = [T<sup>(1)</sup><sub>cell</sub>, ..., T<sup>(N)</sup><sub>cell</sub>] and tropopause pressure levels P = [P<sup>(1)</sup>, ..., P<sup>(N)</sup>]:

$$ \mathbf{P} = \Gamma \mathbf{T}_{\text{BB}} + \mathbf{\Lambda}, \quad (23) $$

where Γ and Λ are the parameters of the linear regression. Climatological values of Juckes (1994) yield Γ ≈ Γ<sub>1</sub> = 13 hPa K<sup>-1</sup>. As reported by Juckes (1994), the occurrence of large-scale surface anomalies may increase γ<sub>1</sub> by a factor of 1.5, thus showing that there is quite an uncertainty on Γ ∝ (1/γ<sub>1</sub>) and that the relationship is case-dependent.

An extended dataset of PV data and model images has been used to study the statistical relationship between P and T<sub>BB</sub>, again with very poor results (Michel 2008). However, the correlation was sometimes found to be very high for some specific cases, as will be shown in the next section. When the correlation is large, the regression coefficient Γ may be used together with the brightness temperatures T<sub>sat</sub> and T<sub>cell</sub> to provide estimates of the pressure level of the tropopause. It is important to note, however, that the frequent lack of good correlation in the general case between PV and WV may prevent the described procedure from being applied in any operational setting.

### d. Observation error standard deviations

#### 1) ALIGNMENT MODEL ERROR

The residual of the regression from Eq. (23) can be written as

$$ \sigma_p = \text{stdev}(\mathbf{P} - \Gamma \mathbf{T}_{\text{BB}}). \quad (24) $$

We then make the approximate assumption that the error standard deviation of the tropopause pressure level observation is

$$ \sigma_0^D \approx \sigma_p. \quad (25) $$
<figure><img src="image_10.png" /><figcaption>FIG. 9. Operational synoptic analysis of the Landes storm, following the graphical summary of Santurette and Joly (2002). The storm is noted lond on the map.</figcaption></figure>

Clearly, this is only an approximation of the true error standard deviation in position space. Indeed, errors on the determination of Γ are not taken into account; for example, one makes a hypothesis that the background is able to produce the good PV-WV relationship. Taking into account this uncertainty would require an ensemble of forecasts, which would be associated with different values for Γ and the corresponding uncertainty. In fact, this would be very close to the framework of the Ensemble Kalman Filter. One may define the analyzed position for the tropopause with the Best Linear Unbiased Estimate equation, where the increment in position (pressure) space is taken as

$$d\mathbf{P}_{\text{tropo}}^{a} = \mathbf{D}\mathbf{H}_{d}^{\Gamma} \frac{1}{1 + (\sigma_{0}^{D}/\sigma_{d})^{2}} \Delta_{P}^{0}, \quad (26)$$

where D is the previously described background error correlation in position space, σd is the background error in position space (cf. Fig. 5), and Hd the observation operator that measures the background tropopause pressure level. However, we do not yet have the ability to produce the analysis through Eq. (14), which requires an interpolating operator and PV inversion. Moreover, horizontal structure functions needed to spatially spread the increment are lacking in this 1D study.

We decided to rely on the 4D-Var assimilation of PV observations. They sample the background at each level where the PV anomaly is detected. In the position space formulation, the observations are the pressure levels rather than the PV values. To mimic this within traditional 4D-Var assimilation, we add the vector dP<sup>a</sup><sub>tropo</sub> to the pressure levels of the sample observations. As the uncertainty on the determination is already taken into

account by Eq. (26), the 4D-Var assimilation is used as an inversion tool, with a high confidence to fit the observations

$$\sigma_b = 0.1 \hat{\sigma}_b, \quad (27)$$

## 2) ADDITIVE MODEL ERROR

Observations y of the 1.5-PVU value at the pressure level P<sub>y</sub> are available as follows:

$$y = 1.5, \quad (28)$$

$$\mathbf{P}_{\text{y}} = \mathbf{P}_{\text{tropo}}^{b} + \Delta_{P}^{0}. \quad (29)$$

The 4D-Var assimilation directly incorporates them, probably with a quite different B matrix than the one described in Figs. 5 and 6. Here, the problem mainly consists of estimating the standard deviation of the errors of these observations. One has a measure of uncertainty in pressure space, which can be converted to PV space by making the rough assumption that the background has a typical PV profile

$$\sigma_{0}^{\text{PV}} \approx \sigma_{P} \frac{\partial \mathcal{H}_{\text{PV}}(x_{b})}{\partial p} \Big|_{p=P_{\text{tropo}}^{b}}. \quad (30)$$

The assimilation of PV observations is then directly performed with ARPEGE 4D-Var, as in Guérin et al. (2006). The expected differences lie in the way that the uncertainty is taken into account and in the shape of the vertical increment, which may be wider when assimilating the observations in position space.
<figure><img src="image_11.png" /><figcaption>(a) 2006100112</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>(c) 2006100212</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>(b) 2006100200</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>(d) 2006100300</figcaption></figure>

FIG. 10. Some coherent structures visible in the ARPEGE analysis fields up to 48 h before the storm reaches France, valid every 12 h. The tropopause (1.5-PVU surface) height is in solid black lines (contour every 1 km, below 10 km). Relative vorticity at 850 hPa is shaded for values above $5 \times 10^{-5} \text{ s}^{-1}$, with shading intervals of $5 \times 10^{-5} \text{ s}^{-1}$. The 850-hPa winds in excess of $15 \text{ m s}^{-1}$ are plotted with barbs (half barbs) representing $10 (5) \text{ m s}^{-1}$. Also noted on the graph are the tropopause anomalies $A_1, A_2, A_3, A_4$, the lower-level cyclonic system $B_1$ that will become the storm, and the larger scale cutoff system $C$ behind the storm.

**5. A case study**

*a. Synoptic description*

The Landes<fnref n="3" /> storm reached France at 0600 UTC 3 October 2006, yielding measured winds above $40 \text{ m s}^{-1}$.

The lower-level vortex was advected by a powerful jet stream, and reached Germany around 1200 UTC. The operational analysis of the situation by the forecasters at Météo-France is depicted in Fig. 9, following the graphical summary introduced by Santurette and Joly (2002). The cyclone has rather a short length scale of around 1000 km.

We describe now the temporal behavior of the most important coherent structures that may play a role in the deepening of the low. The lower-level vortex structure

<footnote n="3">³ The Landes is a region of France, located near the southwest Atlantic coast.</footnote>
<figure><img src="image_15.png" /><figcaption>FIG. 11. ARPEGE forecast and ECMWF verifying analysis of the storm of Landes valid at 0600 UTC 3 Oct 2006. Tropopause (2-PVU surface) height (solid black line, every 1000 m), relative vorticity at 850 hPa (gray shades, every 5 × 10⁻⁵ s⁻¹, for values above 5 × 10⁻⁵ s⁻¹), and wind at 850 hPa (for values above 15 m s⁻¹) are shown.</figcaption></figure>

that will become the storm is visible even a few days before (noted $B_1$ in Fig. 11) and has a very good spatial and temporal coherence. The relative vorticity in the analysis is about $2 \times 10^{-4} \text{ s}^{-1}$ until 1200 UTC 2 October when it suddenly deepens to $3.5 \times 10^{-4} \text{ s}^{-1}$. This vortex is associated with strong convective activity associated with warm tropical air (not shown). As shown by Moore and Montgomery (2004, 2005), this structure may evolve by itself, relying on the diabatic heating to move eastward. The interaction with an anomaly of temperature at the surface could even be the cause of deepening, following the concept of a diabatic Rossby vortex. In our case, however, the presence of strong dry intrusion in the images suggests that dynamic features at the upper levels may play a key role as well.

The upper-level PV vortex that seems to be directly involved in the storm is marked as $A_1$ in Fig. 10. It is spatially very close to the lower-level vortex $B_1$. Based on the operational analysis, the pressure level of the tropopause at the $A_1$ location is around 400 hPa at 1200 UTC 1 October. Between 0000 and 1800 UTC 2 October, $A_1$ deepens strongly to reach about 550 hPa. During this cyclogenesis, $A_1$ is stretched to a much more compact shape by the surrounding upper-level wind. Other upper-level PV vortices are depicted in Fig. 10c. Among them, the anomaly $A_2$ deepens as well, starting from 0000 UTC 2 October.

Finally, one can see a larger-scale cutoff system behind the storm of interest, marked $C$ in Fig. 10. It splits

from the ambient dynamics toward 1200 UTC 2 October. Several upper-level PV vortices, including $A_3$ and $A_4$, are trapped around $C$ and circulate very quickly. A lower-level vortex can also be noticed to appear within $C$, and is associated with strong winds at 850 hPa.

### b. Operational forecast

The forecast of the cyclone by the operational version of ARPEGE significantly varied from cycle to cycle. For forecast ranges beyond 24 h, ARPEGE performed poorly in forecasting both the position and the intensity of the lower level winds. As shown in Fig. 11, the 36-h forecast predicts a stretched, lower intensity vortex (Fig. 11a) than the more compact one that occurred (Fig. 11b). Winds at 850 hPa, which are taken as representative of 10-min wind gusts at the surface, reach $30 \text{ m s}^{-1}$ for the European Centre for Medium-Range Weather Forecasts (ECMWF) analysis, in quite good agreement with the observations (up to $40 \text{ m s}^{-1}$), but only $24.5 \text{ m s}^{-1}$ for the 36-h forecasts. Also, the upper-level PV anomaly has a different shape and location. Differences are very visible at 0900 and 1200 UTC, as well, when the low is crossing over France (not shown).

### c. Generation of pseudo-observations

We first run the tracking algorithm on both satellite and model images. Dry intrusions that seem to be directly linked with the storm are detected and labeled $I_0$ and $I_1$ in Fig. 12. Surprisingly, there is a split of cell $I_0$
<figure><img src="image_16.png" /><figcaption>(a) Tracking on WV MSG satellite images, valid on 01st October 2006, 18 TU</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>(b) Tracking on WV MSG satellite images, valid on 03rd October 2006, 00 TU</figcaption></figure>

FIG. 12. Tracking of dry intrusions in satellite images. Legend identical to that of Fig. 3, except that all lines are in black. Relevant dry intrusions are labeled *I*₀, *I*₁, *I*₂, and *I*₃.
<figure><img src="image_18.png" /></figure>

(a) Time evolution of $I_0$

<figure><img src="image_19.png" /></figure>

(c) Time evolution of $I_1$

<figure><img src="image_20.png" /></figure>

(e) Time evolution of $I_2$

<figure><img src="image_21.png" /></figure>

(g) Time evolution of $I_3$

<figure><img src="image_22.png" /></figure>

(b) Linear regression for $I_0$

<figure><img src="image_23.png" /></figure>

(d) Linear regression for $I_1$

<figure><img src="image_24.png" /></figure>

(f) Linear regression for $I_2$

<figure><img src="image_25.png" /></figure>

(h) Linear regression for $I_3$
around 0200 UTC 2 October (not shown), which yields two different trajectories: the cell $I_0$ is caught by the cutoff *C* and quickly disappears, whereas the other cell (denoted $I_1$) seems linked with the storm (Fig. 12b).

We use the same, fully automatic extraction of PV profiles at the position of the cells from the model images for the study of the PV–WV relationship for the dry intrusions. Results are shown in Fig. 13. The dry intrusion $I_0$ seems to show a precise relationship, thus allowing us to introduce observations by comparing model and satellite cells. This is interesting, as $I_0$ is directly linked with PV anomaly $A_1$ associated with the cyclogenesis (Fig. 11). The parameters of the linear regressions are summarized in Table 1, following the methodology described in section 4.

Figure 14 shows the comparison of the brightness temperatures of the linked satellite and model trajectories for the two dry intrusions $I_0$ and $I_1$. For the latter dry intrusion, there is a very good agreement in the time evolution despite an offset of about 7 K, which is partly due to the global radiance bias (about 1.7 K for the 6.2-μm channel on MSG-1). For the dry intrusion $I_0$, there is a pronounced discrepancy between the model and satellite trajectory, roughly corresponding to a difference of 4 K in warming. Using the linear regression of Eq. (23), together with parameters of Table 1, yields a difference of about 80 hPa in the vertical position of the dynamical tropopause.

### d. PV initial increment

A cross section of the potential vorticity of the background at 1800 UTC 1 October 2006 is shown Fig. 15b. A rough estimate gives $\partial\mathcal{H}_{\text{PV}}(x_b)/\partial p|_{p=500hPa} \simeq 1$ PVU $\times (100 \text{ hPa})^{-1}$. We therefore use the value $\sigma_0^{\text{PV}} = 0.155$ PVU for the assimilation in amplitude space. On the contrary, $\sigma_0^D = 15.5$ hPa is used in alignment space.

Both the amplitude and alignment assimilation procedures use the same operational configuration of ARPEGE and start with the same background, a 6-h forecast valid at 1800 UTC 1 October 2006. Observations are generated and assimilated every hour of the 6-h time window of ARPEGE 4D-Var.

Figure 16 illustrates the difference in the PV analysis increment between the assimilation of the 1.5-PVU value with an additive model (Fig. 16a) and the assimilation of the estimated tropopause pressure level with an approximate alignment model (Fig. 16b). As pointed out in

section 3, the vertical structure is broader for the alignment model, which also vertically displaces the upper-level PV fields (2-PVU and 3-PVU surfaces). Moreover, the amplitude of the increment is larger. This effect is likely to be more important when the background exhibits stronger vertical gradients of PV. The amplitude of the analysis increment in the additive model is likely to be bounded by the 1.5-PVU value (as the observation value is equal to 1.5 PVU and PV is mainly positive). In contrast, by vertically displacing PV surfaces, the alignment model may yield larger increments.

Figure 15 shows the difference between the background at 1200 UTC 1 October 2006 and at 1800 UTC, using a horizontal displacement to roughly remove the effect of advection. There is a clear deepening of the tropopause (cf. Figs. 15a and 15b), such that the temporal difference exhibits a strong positive increment in PV space (Fig. 15b). The vertical extension of this increment is in better agreement with the alignment model (Fig. 16b) than with the additive one (Fig. 16a). To the degree that the analysis error is well approximated by a difference of the deepening of the tropopause under the vertical velocity, the alignment model may be more consistent. These broader increments are also in better agreement with the large vertical extent of singular vectors (Røsting and Kristjánsson 2006).

### e. Impact on forecast

The effect of assimilation in these “nearly single” observation experiments is rather limited, but positive, as depicted in Fig. 17. The deepening of the tropopause by the assimilation procedure leads to a partial correction for the position error of the surface cyclone especially at 1200 UTC (cf. Figs. 17d,f and 17b). The impact is more pronounced in the alignment experiment, where the PV increments were stronger (Fig. 16) than with the amplitude experiment. The wind analysis is also improved: the maximum wind over the Landes regions at a 36-h range increases from the value 23.6 m s⁻¹ (background forecast) to 24.9 m s⁻¹ (alignment model) and 23.8 m s⁻¹ (additive model), respectively.

## 6. Conclusions

Dry intrusions appear as small-scale, darkening features in WV images that are often associated with cyclogenesis. Even if modern data assimilation techniques

←

FIG. 13. Tropopause (1.5-PVU surface) pressure levels and brightness temperatures along the trajectories of dry intrusions $I_0$, $I_1$, $I_2$, and $I_3$. Shown are the time evolution of brightness temperature (dashed line) and tropopause pressure level (solid line) for (a), (c), (e), (g) background detected cells and (b), (f), (h) the linear regressions between the two data sets.
<table><caption>TABLE 1. Parameters of the linear regressions.</caption><thead><tr><th>Dry intrusion</th><th>I<sub>0</sub></th><th>I<sub>1</sub></th><th>I<sub>2</sub></th><th>I<sub>3</sub></th></tr></thead><tbody><tr><td>Γ (hPa K<sup>-1</sup>)</td><td>20.91</td><td>33.52</td><td>11.41</td><td>1.195</td></tr><tr><td>R<sub>corr</sub></td><td>0.885</td><td>0.803</td><td>0.3672</td><td>0.609</td></tr><tr><td>σ<sub>P</sub> (hPa)</td><td>15.50</td><td>63.84</td><td>15.48</td><td>2.153</td></tr></tbody></table>

are able to extract part of this information through direct radiance assimilation (Köpken et al. 2004), PV initialization methods have been developed in many operational centers since the pioneer work of Demirtas and Thorpe (1999). Several case studies have highlighted the potential of modifying the upper-level potential vorticity according to the comparison with water vapor images. Other studies are more reserved on the usefulness of the methodology (Swarbrick 2001). The PV–WV relationship is very complex, such that it is difficult to make the PV modifications objectively. The main issues associated with the PV–WV problem are to make the modifications more objective and to properly take into account the covariance structure of observation and background errors.

Wirth et al. (1997) have introduced a qualitative model for the joint evolution of brightness temperatures in dry intrusions and tropopause height. First, they confirm that great care is necessary when trying to detect analysis errors, as discrepancies between the image and the tropopause may be likely. Second, the tropopause height is not a natural variable, and assimilating it requires either a conversion to PV values or a position–space formulation of the problem. To understand the potential differences between these two approaches, a 1D-Var scheme in alignment space is formulated. It incorporates observations of the tropopause pressure

levels, whereas traditional amplitude-based data assimilation uses observations of PV values at this pressure level. We show that the alignment model may be more in agreement with some “empirical” vertical structure functions developed for PV inversion and used by some authors, in the sense that the vertical spread seems to be larger than for direct additive assimilation. Alignment and additive data assimilations differ in the way the uncertainty is taken into account, and on the vertical extent of the increment. Moreover, as PV is mainly positive in the Northern Hemisphere, the increment of the additive amplitude scheme is likely to be bounded by the value chosen for the tropopause (typically 1.5 PVU), whereas the increment of the alignment scheme is not.

We compare both methods on a case study of the Landes storm. The dry intrusion that is linked with the cyclogenesis is detected on both satellite and model imageries, and show a quantitative difference in their warming. The background shows a good correlation between the pressure level of the tropopause and the brightness temperature of the dry intrusion. This regression is used to provide observations that are assimilated with the preceding amplitude and alignment space assimilation procedures. The work described in this paper may bring two improvements to the PV–WV methodology.

1) We show that the automated use of the tracking algorithm together with tropopause height estimates can yield high correlations on some cases. This permits to put an objective constraint on the tropopause amplitude corrections performed by forecasters.

2) We provide a new way to study and define the structure functions associated with this kind of observations, and highlight that the two data assimilation models differ in specifying error standard deviation.

<figure><img src="image_26.png" /><figcaption>FIG. 14. Comparison of the temperatures of the satellite cell (solid line) and of the model cell (dashed line).</figcaption></figure>
<figure><img src="image_27.png" /><figcaption>(a) Advected background at 20061001, 12TU</figcaption></figure>

<figure><img src="image_28.png" /><figcaption>(b) Advected background at 20061001, 18TU</figcaption></figure>

<figure><img src="image_29.png" /><figcaption>(c) Six hours difference</figcaption></figure>

FIG. 15. Potential vorticity cross sections for (a) background at 1200 UTC 1 Oct 2006 horizontally advected such that the anomaly roughly matches (b) the background at 1800 UTC 1 Oct 2006, with PV plotted every 1 PVU, starting with 1 PVU; (c) differences between these two fields, with a contour interval of 0.5 PVU and solid (dashed) lines representing positive (negative) values.

However, the lack of reliability of the PV–WV relationship prevents a large number of assimilation cases from being successful (Swarbrick 2001). Two approaches might be followed. The first is to try to further understand the PV–WV relationship from a deterministic point of view, in order to figure out when the PV–WV relationship

<figure><img src="image_30.png" /><figcaption>(a) Additive model</figcaption></figure>

<figure><img src="image_31.png" /><figcaption>(b) Alignment model</figcaption></figure>

FIG. 16. Potential vorticity cross sections for the analysis at 1800 UTC 1 Oct 2006. The contour interval is 0.2 PVU, with solid (dashed) for positive (negative) values.

is working. In particular, the work of Wirth et al. (1997) appears to give some clues by putting the vertical velocity in the center of the relationship.

Another alternative approach is to use the PV–WV relationship in a completely statistical manner. The Ensemble Kalman Filter directly samples the covariance matrices over an ensemble of forecasts, so that any kind of observation operator can be included (e.g., Chen and Snyder (2007)). This could be applied to assimilate the brightness temperatures of the cells, or other characteristics. The size of the required ensemble is a significant unknown. The sampling error of the covariance matrices could be reduced by adopting a hybrid approach where the increments are projected into PV space, and then inverted back, thus ensuring that the modifications lie in PV space. Future work may consider the use of this powerful data assimilation algorithm to further exploit the potential of WV images related to PV structures.

*Acknowledgments*. We thank F. Bouttier, F. Rabier, and two anonymous reviewers for their careful reading
<figure><img src="image_32.png" /><figcaption>FIG. 17. Comparison of the ARPEGE forecasts of relative vorticity at 850 hPa (shaded for values above $5 \times 10^{-5}$ s$^{-1}$, with a shading interval of $5 \times 10^{-5}$ s$^{-1}$) with observed 10-m winds, full (half) barbs representing 10 (5) m s$^{-1}$: (a),(b) background forecasts for 0600 and 1200 UTC, respectively; (c),(d) amplitude; and (e),(f) alignment forecasts for the same times.</figcaption></figure>
of a previous version of this paper. Financial support has been provided by the French National Agency ANR with the project ADDISA.

APPENDIX

Inverse of the Linear Coefficient in Pressure
Coordinates

Taking the same parameters as in Juckes (1994) gives
a Coriolis parameter $f = 1 \times 10^{-4} \text{ s}^{-1}$, the gravity con-
stant $g = 10 \text{ m s}^{-2}$, a reference potential temperature
$\theta_{00} = 300 \text{ K}$, and Brunt–Väisälä frequencies in the tro-
posphere and the stratosphere $N_t = 10^{-2} \text{ s}^{-2}$, $N_s = 2N_t$.
Using pressure as the vertical coordinate requires
specification of the relative variations of $z$ and $P$. As in
Juckes (1994), we use the combination of the hydrostatic
relationship and an exponential profile for the density $\rho$

$$
\frac{\partial P}{\partial z} = -\rho g = -g\rho_{00}e^{-z/H_{\rho}}, \quad (A1)
$$

$$
P(z) - P(H_{\rho}) = g\rho_{00}H_{\rho}(e^{-z/H_{\rho}} - e), \quad (A2)
$$

$$
\frac{z}{H_{\rho}} = 1 - \ln \left[ 1 + \frac{P(z) - P(H_{\rho})}{P_{r}} \right], \quad (A3)
$$

where the reference pressure is taken as $P_r = g\rho_{00}H_\rho/e \approx 433$ hPa. The pressure-altitude relationship for pressure and altitude perturbations $\delta P$ and $\delta z$ follows:

$$
\delta P = - \frac{P_r}{H_\rho} \left[ 1 + \frac{P(z) - P(H_\rho)}{P_r} \right] \delta z. \quad (A4)
$$

The simple first-order model derived in the appendix
A of Wirth et al. (1997) assumes that the isosteric surface
is located near the tropopause. Here, the first order
writes, assuming $z \sim H_\rho$

$$
\delta P \approx \Gamma_1 \Delta T_{BB} = - \frac{P_r}{\gamma_1 H_\rho} \Delta T_{BB}, \quad (A5)
$$

$$
\Gamma_1 = 13 \text{ hPa K}^{-1}. \tag{A6}
$$

REFERENCES

Arbogast, P., 1998: L’inversion du tourbillon potentiel: méthod-
ologie, application l’étude des interactions non-linéaires dans
la formation des dépressions météorologiques. Ph.D. thesis,
Université Paul Sabatier, Toulouse, 227 pp.

—, K. Maynard, and F. Crepin, 2008: Ertel potential vorticity inversion using a digital filter initialization method. *Quart. J. Roy. Meteor. Soc.*, **134**, 1287–1296.

Beare, R., A. Thorpe, and A. White, 2003: The predictability of ex-
tratropical cyclones: Nonlinear sensitivity to localized potential-
vorticity perturbations. *Quart. J. Roy. Meteor. Soc.*, **129**, 219–237.

Beezley, J. D., and J. Mandel, 2008: Morphing ensemble Kalman filters. *Tellus*, **60**, 131–140.

Belo-Pereira, M., and L. Berre, 2006: The use of an ensemble ap-
proach to study the background error covariances in a global
NWP model. *Mon. Wea. Rev.*, **134**, 2466–2489.

Browning, K., 1997: The dry intrusion perspective of extra-tropical cyclone development. *Meteor. Appl.*, **4**, 317–324.

Chapnik, B., G. Desroziers, F. Rabier, and O. Talagrand, 2004: Properties and first application of an error-statistics tuning method in variational assimilation. *Quart. J. Roy. Meteor. Soc.*, **130**, 2253–2275.

—, —, —, and —, 2006: Diagnosis and tuning of obser-
vational error in a quasi-operational data assimilation setting.
*Quart. J. Roy. Meteor. Soc.*, **132**, 543–565.

Chen, Y., and C. Snyder, 2007: Assimilating vortex position with an
ensemble Kalman filter. *Mon. Wea. Rev.*, **135**, 1828–1845.

Courtier, P., C. Freydier, J. Geleyn, F. Rabier, and M. Rochas, 1991: The Arpege project at Météo-France. *Proc. ECMWF Seminar on Numerical Methods in Atmospheric Models*, Reading, United Kingdom, ECMWF, 193–232.

Davis, C., 1992: Piecewise potential vorticity inversion. *J. Atmos. Sci.*, **49**, 1397–1411.

Demirtas, M., and A. Thorpe, 1999: Sensitivity of short-range weather forecasts to local potential vorticity modifications. *Mon. Wea. Rev.*, **127**, 922–939.

Desroziers, G., and S. Ivanov, 2001: Diagnosis and adaptive tuning of observation-error parameters in a variational assimilation. *Quart. J. Roy. Meteor. Soc.*, **127**, 1433–1452.

Georgiev, C., 1999: Quantitative relationship between METEOSAT water vapour data and positive vorticity anomalies: A case study over the Mediterranean. *Meteor. Appl.*, **6**, 97–109.

Guérin, R., G. Desroziers, and P. Arbogast, 2006: 4D-Var analysis of pseudo potential-vorticity observations. *Quart. J. Roy. Meteor. Soc.*, **132**, 1283–1298.

Hello, G., and P. Arbogast, 2004: Two different methods to correct the initial conditions applied to the storm of 27 December 1999 over southern France. *Meteor. Appl.*, **11**, 41–57.

Hoskins, B., M. McIntyre, and A. Roberson, 1985: On the use and
significance of isentropic potential vorticity maps. *Quart. J.
Roy. Meteor. Soc.*, **111**, 877–946.

Juckes, M., 1994: Quasigeostrophic dynamics of the tropopause. *J. Atmos. Sci.*, **51**, 2756–2768.

Köpken, C., G. Kelly, and J.-N. Thépaut, 2004: Assimilation of METEOSAT radiance data within the 4D-Var system at ECMWF: Assimilation experiments and forecast impact. *Quart. J. Roy. Meteor. Soc.*, **130**, 2277–2292.

Lawson, W., and J. Hansen, 2005: Alignment error models and ensemble-based data assimilation. *Mon. Wea. Rev.*, **133**, 1687–1709.

Manders, A., W. Verkley, J. Dieveveen, and A. Moene, 2007: Application of a potential vorticity modification method to a case of rapid cyclogenesis over the Atlantic Ocean. *Quart. J. Roy. Meteor. Soc.*, **133**, 1755–1770.

Mansfield, D., 1996: The use of potential vorticity as an operational forecast tool. *Meteor. Appl.*, **3**, 195–210.

Michel, Y., 2008: Assimilation de données d'images télédétectées en météorologie. Ph.D. thesis, Université Paris Est, Paris, 195 pp.
—, and F. Bouttier, 2006: Automated tracking of dry intrusions on satellite water vapour imagery and model output. *Quart. J. Roy. Meteor. Soc.*, **132**, 2257–2276.
Moore, R., and M. Montgomery, 2004: Reexamining the dynamics of short-scale, diabatic Rossby waves and their role in mid-latitude moist cyclogenesis. *J. Atmos. Sci.*, **61**, 754–768.

—, and —, 2005: Analysis of an idealized, three-dimensional diabatic Rossby vortex: A coherent structure of the moist baroclinic atmosphere. *J. Atmos. Sci.*, **62**, 2703–2725.

Morel, C., and S. Sénési, 2002: A climatology of mesoscale convective systems over Europe using satellite infrared imagery. Part I: Methodology. *Quart. J. Roy. Meteor. Soc.*, **128**, 1953–1971.

Parrish, D., and J. Derber, 1992: The National Meteorological Center's spectral statistical interpolation analysis system. *Mon. Wea. Rev.*, **120**, 1747–1763.

Plu, M., and P. Arbogast, 2005: A cyclogenesis evolving into two distinct scenarios and its implications for short-term ensemble forecasting. *Mon. Wea. Rev.*, **133**, 2016–2029.

Ramond, D., H. Corbin, M. Desbois, G. Szejwach, and P. Waldteufel, 1981: The dynamics of polar jet streams as depicted by the Meteosat WV channel radiance. *Mon. Wea. Rev.*, **109**, 2164–2176.

Ravela, S., K. Emanuel, and D. McLaughlin, 2007: Data assimilation by field alignment. *Physica D*, **230**, 127–145.

Røsting, B., and J. Kristjánsson, 2006: Improving simulations of severe winter storms by initial modification of potential vorticity in sensitive regions. *Quart. J. Roy. Meteor. Soc.*, **132**, 2625–2652.

—, J. Sunde, and K. Midtbo, 1996: Monitoring of NWP models by use of satellite data. *Meteor. Appl.*, **3**, 331–340.

—, J. Kristjánsson, and J. Sunde, 2003: The sensitivity of numerical simulations to initial modifications of potential vorticity—A case-study. *Quart. J. Roy. Meteor. Soc.*, **129**, 2697–2718.

Santurette, P., and A. Joly, 2002: ANASYG/PRESYG, Météo-France's new graphical summary of the synoptic situation. *Meteor. Appl.*, **9**, 129–154.

—, and C. Georgiev, 2005: Weather Analysis and Water Vapor Imagery and Potential Vorticity Analysis. Academic Press, 179 pp.

Saunders, R., M. Matricardi, and P. Brunel, 1999: An improved fast radiative transfer model for assimilation of satellite radiance observations. *Quart. J. Roy. Meteor. Soc.*, **125**, 1407–1425.

Stoffelen, A., and Coauthors, 2005: The atmospheric dynamics mission for global wind field measurement. *Bull. Amer. Meteor. Soc.*, **86**, 73–87.

Swarbrick, S., 2001: Applying the relationship between potential vorticity fields and water vapour imagery to adjust initial conditions in numerical weather prediction. *Meteor. Appl.*, **8**, 221–228.

Verkley, W., P. Vosbeek, and A. Moene, 2005: Manually adjusting a numerical weather analysis in terms of potential vorticity using three-dimensional variational data assimilation. *Quart. J. Roy. Meteor. Soc.*, **131**, 1713–1736.

Wirth, V., C. Appenzeller, and M. Juckes, 1997: Signatures of induced vertical air motion accompanying quasi-horizontal roll-up of stratospheric intrusions. *Mon. Wea. Rev.*, **125**, 2504–2519.

Yessad, K., and P. Bénard, 1995: Introduction of a local mapping factor in the spectral part of the Météo-France global variable mesh numerical forecast model. *Quart. J. Roy. Meteor. Soc.*, **122**, 1701–1719.