# Correcting Fast-Mode Pressure Errors in Storm-Scale Ensemble Kalman Filter Analyses

Authors: Corey K. Potvin; Louis J. Wicker
Journal: Advances in Meteorology
 (2013)

---

## Abstract

A typical storm-scale ensemble Kalman filter (EnKF) analysis/forecast system is shown to introduce imbalances into the ensemble posteriors that generate acoustic waves in subsequent integrations. When the EnKF is used to research storm-scale dynamics, the resulting spurious pressure oscillations are large enough to impact investigation of processes driven by nonhydrostatic pressure gradient forces. Fortunately, thermodynamic retrieval techniques traditionally applied to dual-Doppler wind analyses can be adapted to diagnose the balanced portion of an EnKF pressure analysis, thereby eliminating the fast-mode pressure oscillations. The efficacy of this approach is demonstrated using a high-resolution supercell thunderstorm simulation as well as EnKF analyses of a simulated and a real supercell.


---

## Body

## 1. Introduction

The EnKF [1] has become a popular and valuable tool for storm-scale research [2–12]. Particularly when dual-Doppler radar data are available, EnKF data assimilation can provide reliable analyses of wind and, to a lesser degree, temperature and microphysical variables in convective storms. EnKF analyses of pressure, on the other hand, are subject to severe errors, at least with some compressible model configurations (the first tests of the EnKF with a compressible model were performed by Tong and Xue [5]). This problem is illustrated in Figure 1 using output from the National Severe Storms Laboratory Collaborative Model for Multiscale Atmospheric Simulation (NCOMMAS; [13, 14]) ensemble square root filter [15]. Similar behavior occurs using the Data Assimilation Research Testbed [16] EnKF with the Advanced Research Weather Research and Forecasting (WRF-ARW; [17]) model (James Marquis and Thomas Jones, personal communication 2013). The pressure analysis errors severely impede investigation of critical storm processes that are, in part, driven by dynamic pressure gradient forces, including supercell occlusion downdrafts [18], lifting of negatively buoyant air [19], supercell propagation [20], the descending rear inflow and ascending front-to-rear flow in mesoscale convective systems [21], and possibly descending reflectivity cores [22].Figure 1
True (left column) and EnKF mean posterior (right column)p′ (shading; hPa) and radar reflectivity factor (contoured at 10, 30, and 50 dBZ) at z=0.9 km. Fields are valid after (top row) zero, (middle row) one, and (bottom row) fifteen 2 min forecast cycles. The true p′ were filtered and averaged as in Potvin et al. [23] to permit more direct comparisons with the (coarser) EnKF p′.The pressure oscillations shown in Figure1 are associated with acoustic waves generated within the data assimilation region. The waves are presumably excited as each ensemble member adjusts to an updated initial condition that is dynamically inconsistent with the model (i.e., unbalanced). This hypothesis is supported by two observations. First, the acoustic waves occur whether or not pressure is updated during the data assimilation and therefore cannot be attributed to erroneous ensemble covariances between the pressure and other variables (though the latter could conceivably exacerbate the problem in cases where pressure is updated). Second, spurious waves are not evident in perfect-model EnKF observing system simulation experiments (OSSEs) with the NCOMMAS (not shown), indicating that the waves arise only when analysis increments are substantially unbalanced. The generation of high-amplitude fast modes due to unbalanced initial conditions is a long-recognized problem in numerical weather prediction, and many approaches have been used to improve dynamical balance during the data assimilation procedure [24, 25]. Given the small influence of the pressure field on the remaining state variables in certain compressible cloud models on the relevant spatiotemporal scales, errors due to the acoustic waves are largely confined to the pressure field in storm-scale EnKF analyses. This makes it possible to retrieve the portion of the pressure field that is in balance with the remaining model fields (hereafter, the “balanced” pressure). This obviates the need to mitigate the acoustic waves during the data assimilation procedure, at least when the analysis is not being used to initialize a numerical forecast. Instead, the balanced pressure can be retrieved after the data assimilation is complete.Several methods are available for diagnosing perturbation (from the hydrostatic base state) pressure, typically in its nondimensional Exner function formπ (the prime symbol is omitted herein to simplify the notation), from the equations of motion. Early methods satisfied the horizontal equations of motion on individual horizontal planes [26, 27]. A major limitation of that approach is that the analyzed π is offset from the true π by a vertically varying constant, precluding unique solution of ∂π/∂z unless independent pressure measurements are available at each analysis level. This problem is avoided when all three equations of motion are satisfied [28–30], in which case the retrieved π is offset by a volume-wide (rather than vertically varying) constant, thereby permitting the impact of vertical pressure gradients to be quantitatively considered. The tradeoff is that errors in the analyzed local derivative and buoyancy terms in the vertical momentum equation typically result in ∂π/∂x and ∂π/∂y analyses that are inferior to those obtained using two-dimensional retrieval. As will be shown, however, there is a very simple procedure for obtaining the advantages of both the 2D and 3D approaches. This combined method permits useful pressure retrievals to be obtained from acoustic wave-contaminated EnKF analyses.


## 2. Pressure Retrieval Method

We adopt a variational framework for our pressure retrieval scheme. The momentum equation constraints use the Klemp and Wilhelmson [31] formulation of the equations of motion except with the Coriolis term (which is negligible) omitted. A horizontal smoothness constraint is imposed to filter noise. The cost function we seek to minimize can therefore be expressed as
(1)J=12∑i∑j∑kC1J12+C2J22+C3J32+C4J42,
where
(2)J1=(∂π∂x-F),J2=(∂π∂y-G),(3)J3=(∂π∂z-gθ′θ02-H),J4=(∂2π∂x2+∂2π∂y2);i, j, and k are the model grid indices; F, G, and H are the sums of the local time derivative, advection, and turbulent mixing terms for the u-, v-, and w-equations, respectively; g is the gravitational acceleration; and θ0 and θ′are the horizontally homogeneous base state and perturbation potential temperature, respectively. The subscripted C’s represent the weighting coefficients for each respective constraint, computed similarly to the coefficients in [30]:
(4)C1=[∑i∑j∑k(F2)ΔxΔyΔz∑i∑j∑kΔxΔyΔz]-1,C2=[∑i∑j∑k(G2)ΔxΔyΔz∑i∑j∑kΔxΔyΔz]-1,C3=[∑i∑j∑k(H2)ΔxΔyΔz∑i∑j∑kΔxΔyΔz]-1,C4=k0-4,
where Δx, Δy, and Δz are the analysis grid spacings and k0 is the cutoff wavenumber, for which the theoretical filter response is 0.5. The present study uses a 3Δx cutoff wavelength.In the first step of the pressure retrieval scheme, a 2D retrieval is performed (i.e.,C3 is set to zero) over each analysis level. The resulting analysis, π2D, is stored, and then a 3D retrieval is performed, yielding π3D. The vertically varying constant A(k) by which π2Dis offset from the balanced π, πbal, is then estimated by computing the mean difference between π3D and π2D for each k. To see why A(k) can be estimated in this way, consider that the 2D and 3D analyses can be written as
(5)π2D=πbal+A(k)+π2Derror,π3D=πbal+π3Derror,
where π2Derrorand π3Derrorare the retrieval errors, apart from the vertically varying constant. Assuming that the 2D and 3D retrievals are approximately unbiased relative to each other (apart from the vertically varying constant), averaging the differences between the retrievals over each horizontal analysis plane (symbolized by brackets) yields
(6)〈π2D-π3D〉=〈π2Derror-π3Derror〉+A(k)≈A(k).
The final analysis πfinal is then obtained by subtracting the A(k) estimate from π2D:
(7)πfinal=π2D-〈π2D-π3D〉.
The πfinal is still subject to a volume-wide constant, which can be estimated as the volume-mean difference between πfinal and the π from the simulation or EnKF analysis to which the retrieval method is applied. This estimate is reasonable to the extent that the errors in the input π are unbiased. Errors in the estimate of this constant will not impact most applications since it is typically only the spatial derivatives of π that are sought. The estimated volume-wide constant has been subtracted from the πfinal in the experiments below. To further facilitate interpretation, π has been converted to dimensional perturbation pressure p.


## 3. Experiments

### 3.1. Verification of Retrieval Procedure

An NCOMMAS supercell simulation was used to test the robustness of the retrieval procedure and the integrity of the computer code. The simulation was performed on a stationary 102.6 km × 102.6 km × 20.4 km domain with 600 m grid spacing in all three dimensions and a large (small) time step of 4 (2/3) s. A fully dual-moment version of the Ziegler [32] microphysical parameterization (Ziegler Variable Density or ZVD) scheme [33] was used. The remaining model settings were identical to those of Potvin et al. [23]. The supercell in the present simulation appears reasonably realistic, and its evolution qualitatively resembles that of Potvin et al. [23].Pressure retrievals were performed over the entire simulation domain. The local derivatives ofu, v, and w were computed using centered finite differences with default Δt=30 s (e.g.,  ∂u/∂t=(u(t+30)-u(t-30))/60). In general, spatial derivatives were computed using centered differencing on the Arakawa A (unstaggered) grid, necessitating averaging of the model u, v, and w from the Arakawa C (staggered) grid. At the model boundaries, the π derivatives in J1-J3 were computed using one-sided differences. In preliminary retrievals, setting the turbulence terms in the equations of motion to zero generally had negligible or mildly positive impact, presumably due to discretization errors being of similar magnitude to the turbulence terms themselves, which are typically much smaller than the remaining momentum equation terms. The turbulence terms were consequently omitted in the experiments shown.Visual comparisons of the model and retrievedpfinal, ∂pfinal/∂x, and ∂pfinal/∂z reveal high fidelity in the retrieval technique (a representative example is shown in Figure 2), as do vertical profiles of root mean square error (RMSE) and relative (i.e., as a percentage of RMS ptrue) RMSE (RRMSE) computed within the storm (Figure 3). (Both the model and retrieved pressure gradients were computed using centered differences on the Arakawa A grid.) Calculations of the horizontal domain-wide p bias at each level indicate that the 3D retrieval largely eliminates the vertically varying constant present in the 2D retrieval (Figure 4). To assess the sensitivity of the technique to hydrometeor mixing ratio errors (which are often large in EnKF analyses due to gross imperfections in current microphysical parameterization schemes), the retrieval was repeated with the hydrometeor drag term (in the vertical equation of motion) omitted. While ignoring water loading substantially degraded p3D within the storm (Figure 5), p3D was barely impacted in precipitation-free regions (not shown), resulting in relatively small domain-wide p3D biases at each level and, thus, only minor error increases in the A(k) estimates (Figure 4). The spurious vertical trend introduced to pfinal (generally < 0.05 hPa km−1) by the degraded A(k) was practically negligible. The relative insensitivity of the proposed retrieval method to hydrometeor errors is a major advantage over existing 3D approaches.Figure 2
True (left) and retrieved (right)p (top; hPa), ∂p/∂x (middle; hPa km−1), and ∂p/∂z (bottom; hPa km−1) at z=0.9 km. Radar reflectivity factor is contoured at 10, 30, and 50 dBZ.Figure 3
Vertical profiles of RMSE (left column) and RRMSE (right column)pfinal (top row), ∂pfinal/∂x (middle row), and ∂pfinal/∂z (bottom row) for default retrieval (Δt=30 s; solid) and retrieval with Δt=120 s (dashed). All quantities valid where the true dBZ > 0.Figure 4
Vertical profiles of horizontal domain-wide bias inp2D (dotted), pfinal (solid), and pfinalretrieved with the hydrometeor drag term omitted in the cost function (dashed).As in Figure3 (top row) but for p3D (solid) and p3Dretrieved with hydrometeor drag term omitted (dashed).
(a)
(b)The retrievals were also relatively insensitive to temporal discretization errors in the local wind derivatives. Computing the latter usingΔt=120 s (rather than 30 s) to match typical storm-scale EnKF analysis cycle periods substantially increased the RMSE ∂pfinal/∂x and ∂pfinal/∂z (Figure 3) but did not seriously qualitatively degrade the retrieval (Figure 6). While temporal discretization errors will be much larger in cases of rapidly moving storms (the present storm travels eastward at ~10 m s−1), such errors may be substantially reduced using advection correction methods [34–36].Figure 6
As in Figure2 except for Δt=120s in the pressure retrieval.


### 3.2. Retrievals from a Simulated EnKF Mean Analysis

Having verified the robustness of the pressure retrieval scheme, we then applied the technique to mean posteriors from an EnKF OSSE conducted by Potvin and Wicker [11] (pictured in Figure 1). In that experiment (labeled “2-LFO” in [11]), mobile dual-Doppler velocity pseudoobservations were generated from the “truth” simulation of Potvin et al. [23], perturbed with random errors, and assimilated at two-minute intervals. Model error was introduced by using coarser ensemble grid spacing (600 m versus 200 m) and different microphysical parameterization (the Gilmore et al. [37] version of the Lin et al. [38] scheme versus the ZVD scheme) from those in the “truth” simulation. The use of an imperfect model in the EnKF resulted in dynamically unbalanced analysis increments that instigated acoustic waves during the data assimilation. The ensemble was integrated for 30 min prior to the first analysis update to develop physically realistic covariance structures.Local wind derivatives for the pressure retrieval scheme were computed using mean posteriors separated by 4 min (i.e.,Δt=2 min). This interval is similar to the volume scan periods for Weather Surveillance Radar—1988 Doppler (WSR-88D) convection sampling patterns. All fields required for the momentum equation constraints were obtained from the EnKF posteriors (i.e., treated as known quantities). The pressure retrievals were performed offline from the EnKF analysis; that is, the retrieved posterior pressure field at each time was not used to initialize the subsequent forecast cycle in the data assimilation experiment. Due to the very weak impact of the pressure field on the remaining state variables in compressible models (the property exploited by the proposed method to remove the fast-mode pressure oscillations), rebalancing the pressure field at each analysis cycle should have generally negligible impact on the remaining model fields throughout the data assimilation. We have verified this with separate, real data experiments (not shown), in which the differences due to the pressure rebalancing produced only tiny differences in the final ensemble analysis and in subsequent 30 min forecasts.Given the insensitivity of the pressure retrievals in Section3.1 to discretization and hydrometeor mixing ratio errors, it was expected that pfinal obtained from the EnKF analyses would well capture the balanced portion of the EnKF p (pbal) and thereby improve upon the EnKF p  and generally resemble ptrue (keeping in mind that errors in the other EnKF fields will create some discrepancies between the EnKF pbal and ptrue). This was indeed the case throughout the assimilation period, except at higher altitudes. After 20 min of data assimilation (t=50 min), the EnKF wind analyses are sufficiently accurate that pfinal reasonably resembles ptrue at lower and middle levels, much more so than does the EnKF p (Figures 7, 8, and 9). The improvement of pfinal over the EnKF p, however, decreases with height (Figure 9), presumably due in large part to the shallower fast-mode pressure oscillations aloft (not shown). Within the top half of the storm, pfinal is generally mildly inferior to the EnKF p, indicating that retrieval errors dominate improvements from removing the (small) fast-mode pressure errors. This suggests that application of the proposed pressure retrieval method be restricted to altitudes where spurious pressure oscillations are large.Figure 7
True (left column), retrieved (middle column), and EnKF mean posterior (right column)p (top row), ∂p/∂x (middle row; hPa km−1), and ∂p/∂z (bottom row; hPa km−1) at z=0.9 km, t=50 min. The true p were filtered and averaged as in Potvin et al. [23] to permit more direct comparisons with the (coarser) EnKF p. Radar reflectivity factor is contoured at 10, 30, and 50 dBZ.Figure 8
As in Figure7 except for z=3.3 km.Figure 9
As in Figure3 but for pfinal retrieved from mean EnKF posterior (solid) and mean EnKF posterior p (dashed).While the retrieval results from this OSSE are likely somewhat optimistic despite our efforts to mitigate the “identical twin” problem, larger EnKF analysis errors in real cases should not substantially impede the technique’s ability to recover the EnKFpbal and, thus, to remove the fast-mode pressure errors (which is the objective of the proposed method). Furthermore, the previous results, combined with the resilience of dual-mobile-Doppler EnKF wind analyses to microphysical and background wind errors [11, 39], suggest that pressure retrievals obtained from high-quality radar datasets such as those collected during the second Verification of the Origins of Rotation in Tornadoes Experiment (VORTEX2; [40]) may be useful for investigating pressure-gradient-driven storm processes. Rigorous testing of this hypothesis, unfortunately, is currently precluded by the lack of dense 3D pressure observations within storms.


### 3.3. Retrievals from a Real Data EnKF Mean Analysis

While our ability to verify the real-world performance of the proposed technique is quite limited, qualitative evaluation of pressure retrievals from EnKF analyses of the May 29-30, 2004, Geary, OK, USA, tornadic supercell supports cautious optimism. The EnKF analyses were obtained by assimilating quality-controlled Doppler velocity data from two Shared Mobile Atmospheric Research and Teaching (SMART; [41]) radars using the NCOMMAS EnKF with the ZVD microphysics scheme (see Potvin et al. [39] for further details). Pressure retrievals were performed in the same manner as with the simulated EnKF analyses in Section 3.2.A representative comparison of thepfinal and EnKF p, valid after 44 min of data assimilation (0037 UTC 30 May), is shown in Figure 10. As in the EnKF OSSE (Section 3.2), high-amplitude pressure oscillations generated during the data assimilation are not evident in the retrieved pressure fields. Furthermore, the pfinal comports with present understanding of the pressure distribution within supercells, which suggests retrieval errors are not unduly large. For example, both the pfinal and EnKF p exhibit an inflow low and a deep pressure minimum associated with the strong mesocyclonic rotation (Figures 10(a) and 10(b)), as well as upward-directed perturbation pressure gradient force along the gust front (Figures 10(c) and 10(d)). The pfinal, but not the EnKF p, persistently indicates a region of substantial downward-directed perturbation pressure gradient force near the analyzed occlusion downdraft (Figures 10(c)–10(f)). This is consistent with current understanding of occlusion downdraft formation [42], suggesting that the EnKF pressure analysis near the occlusion downdraft is contaminated with fast-mode errors that the retrieval successfully corrects.Figure 10
EnKF mean posterior (left column) and retrieved (right column)p at z=250 m (hPa; top row) and ∂p/∂z at z=1250 m (hPa km−1; middle and bottom rows). The black squares in the middle row outline the domain (centered near the occlusion updraft) in the bottom row. EnKF mean posterior vertical velocity (red; dashed = negative) is contoured at −5, −2, 2, and 5 m s−1 and radar reflectivity factor (black) at 10, 30, and 50 dBZ. Arrows represent EnKF mean posterior horizontal winds.


## 3.1. Verification of Retrieval Procedure

An NCOMMAS supercell simulation was used to test the robustness of the retrieval procedure and the integrity of the computer code. The simulation was performed on a stationary 102.6 km × 102.6 km × 20.4 km domain with 600 m grid spacing in all three dimensions and a large (small) time step of 4 (2/3) s. A fully dual-moment version of the Ziegler [32] microphysical parameterization (Ziegler Variable Density or ZVD) scheme [33] was used. The remaining model settings were identical to those of Potvin et al. [23]. The supercell in the present simulation appears reasonably realistic, and its evolution qualitatively resembles that of Potvin et al. [23].Pressure retrievals were performed over the entire simulation domain. The local derivatives ofu, v, and w were computed using centered finite differences with default Δt=30 s (e.g.,  ∂u/∂t=(u(t+30)-u(t-30))/60). In general, spatial derivatives were computed using centered differencing on the Arakawa A (unstaggered) grid, necessitating averaging of the model u, v, and w from the Arakawa C (staggered) grid. At the model boundaries, the π derivatives in J1-J3 were computed using one-sided differences. In preliminary retrievals, setting the turbulence terms in the equations of motion to zero generally had negligible or mildly positive impact, presumably due to discretization errors being of similar magnitude to the turbulence terms themselves, which are typically much smaller than the remaining momentum equation terms. The turbulence terms were consequently omitted in the experiments shown.Visual comparisons of the model and retrievedpfinal, ∂pfinal/∂x, and ∂pfinal/∂z reveal high fidelity in the retrieval technique (a representative example is shown in Figure 2), as do vertical profiles of root mean square error (RMSE) and relative (i.e., as a percentage of RMS ptrue) RMSE (RRMSE) computed within the storm (Figure 3). (Both the model and retrieved pressure gradients were computed using centered differences on the Arakawa A grid.) Calculations of the horizontal domain-wide p bias at each level indicate that the 3D retrieval largely eliminates the vertically varying constant present in the 2D retrieval (Figure 4). To assess the sensitivity of the technique to hydrometeor mixing ratio errors (which are often large in EnKF analyses due to gross imperfections in current microphysical parameterization schemes), the retrieval was repeated with the hydrometeor drag term (in the vertical equation of motion) omitted. While ignoring water loading substantially degraded p3D within the storm (Figure 5), p3D was barely impacted in precipitation-free regions (not shown), resulting in relatively small domain-wide p3D biases at each level and, thus, only minor error increases in the A(k) estimates (Figure 4). The spurious vertical trend introduced to pfinal (generally < 0.05 hPa km−1) by the degraded A(k) was practically negligible. The relative insensitivity of the proposed retrieval method to hydrometeor errors is a major advantage over existing 3D approaches.Figure 2
True (left) and retrieved (right)p (top; hPa), ∂p/∂x (middle; hPa km−1), and ∂p/∂z (bottom; hPa km−1) at z=0.9 km. Radar reflectivity factor is contoured at 10, 30, and 50 dBZ.Figure 3
Vertical profiles of RMSE (left column) and RRMSE (right column)pfinal (top row), ∂pfinal/∂x (middle row), and ∂pfinal/∂z (bottom row) for default retrieval (Δt=30 s; solid) and retrieval with Δt=120 s (dashed). All quantities valid where the true dBZ > 0.Figure 4
Vertical profiles of horizontal domain-wide bias inp2D (dotted), pfinal (solid), and pfinalretrieved with the hydrometeor drag term omitted in the cost function (dashed).As in Figure3 (top row) but for p3D (solid) and p3Dretrieved with hydrometeor drag term omitted (dashed).
(a)
(b)The retrievals were also relatively insensitive to temporal discretization errors in the local wind derivatives. Computing the latter usingΔt=120 s (rather than 30 s) to match typical storm-scale EnKF analysis cycle periods substantially increased the RMSE ∂pfinal/∂x and ∂pfinal/∂z (Figure 3) but did not seriously qualitatively degrade the retrieval (Figure 6). While temporal discretization errors will be much larger in cases of rapidly moving storms (the present storm travels eastward at ~10 m s−1), such errors may be substantially reduced using advection correction methods [34–36].Figure 6
As in Figure2 except for Δt=120s in the pressure retrieval.


## 3.2. Retrievals from a Simulated EnKF Mean Analysis

Having verified the robustness of the pressure retrieval scheme, we then applied the technique to mean posteriors from an EnKF OSSE conducted by Potvin and Wicker [11] (pictured in Figure 1). In that experiment (labeled “2-LFO” in [11]), mobile dual-Doppler velocity pseudoobservations were generated from the “truth” simulation of Potvin et al. [23], perturbed with random errors, and assimilated at two-minute intervals. Model error was introduced by using coarser ensemble grid spacing (600 m versus 200 m) and different microphysical parameterization (the Gilmore et al. [37] version of the Lin et al. [38] scheme versus the ZVD scheme) from those in the “truth” simulation. The use of an imperfect model in the EnKF resulted in dynamically unbalanced analysis increments that instigated acoustic waves during the data assimilation. The ensemble was integrated for 30 min prior to the first analysis update to develop physically realistic covariance structures.Local wind derivatives for the pressure retrieval scheme were computed using mean posteriors separated by 4 min (i.e.,Δt=2 min). This interval is similar to the volume scan periods for Weather Surveillance Radar—1988 Doppler (WSR-88D) convection sampling patterns. All fields required for the momentum equation constraints were obtained from the EnKF posteriors (i.e., treated as known quantities). The pressure retrievals were performed offline from the EnKF analysis; that is, the retrieved posterior pressure field at each time was not used to initialize the subsequent forecast cycle in the data assimilation experiment. Due to the very weak impact of the pressure field on the remaining state variables in compressible models (the property exploited by the proposed method to remove the fast-mode pressure oscillations), rebalancing the pressure field at each analysis cycle should have generally negligible impact on the remaining model fields throughout the data assimilation. We have verified this with separate, real data experiments (not shown), in which the differences due to the pressure rebalancing produced only tiny differences in the final ensemble analysis and in subsequent 30 min forecasts.Given the insensitivity of the pressure retrievals in Section3.1 to discretization and hydrometeor mixing ratio errors, it was expected that pfinal obtained from the EnKF analyses would well capture the balanced portion of the EnKF p (pbal) and thereby improve upon the EnKF p  and generally resemble ptrue (keeping in mind that errors in the other EnKF fields will create some discrepancies between the EnKF pbal and ptrue). This was indeed the case throughout the assimilation period, except at higher altitudes. After 20 min of data assimilation (t=50 min), the EnKF wind analyses are sufficiently accurate that pfinal reasonably resembles ptrue at lower and middle levels, much more so than does the EnKF p (Figures 7, 8, and 9). The improvement of pfinal over the EnKF p, however, decreases with height (Figure 9), presumably due in large part to the shallower fast-mode pressure oscillations aloft (not shown). Within the top half of the storm, pfinal is generally mildly inferior to the EnKF p, indicating that retrieval errors dominate improvements from removing the (small) fast-mode pressure errors. This suggests that application of the proposed pressure retrieval method be restricted to altitudes where spurious pressure oscillations are large.Figure 7
True (left column), retrieved (middle column), and EnKF mean posterior (right column)p (top row), ∂p/∂x (middle row; hPa km−1), and ∂p/∂z (bottom row; hPa km−1) at z=0.9 km, t=50 min. The true p were filtered and averaged as in Potvin et al. [23] to permit more direct comparisons with the (coarser) EnKF p. Radar reflectivity factor is contoured at 10, 30, and 50 dBZ.Figure 8
As in Figure7 except for z=3.3 km.Figure 9
As in Figure3 but for pfinal retrieved from mean EnKF posterior (solid) and mean EnKF posterior p (dashed).While the retrieval results from this OSSE are likely somewhat optimistic despite our efforts to mitigate the “identical twin” problem, larger EnKF analysis errors in real cases should not substantially impede the technique’s ability to recover the EnKFpbal and, thus, to remove the fast-mode pressure errors (which is the objective of the proposed method). Furthermore, the previous results, combined with the resilience of dual-mobile-Doppler EnKF wind analyses to microphysical and background wind errors [11, 39], suggest that pressure retrievals obtained from high-quality radar datasets such as those collected during the second Verification of the Origins of Rotation in Tornadoes Experiment (VORTEX2; [40]) may be useful for investigating pressure-gradient-driven storm processes. Rigorous testing of this hypothesis, unfortunately, is currently precluded by the lack of dense 3D pressure observations within storms.


## 3.3. Retrievals from a Real Data EnKF Mean Analysis

While our ability to verify the real-world performance of the proposed technique is quite limited, qualitative evaluation of pressure retrievals from EnKF analyses of the May 29-30, 2004, Geary, OK, USA, tornadic supercell supports cautious optimism. The EnKF analyses were obtained by assimilating quality-controlled Doppler velocity data from two Shared Mobile Atmospheric Research and Teaching (SMART; [41]) radars using the NCOMMAS EnKF with the ZVD microphysics scheme (see Potvin et al. [39] for further details). Pressure retrievals were performed in the same manner as with the simulated EnKF analyses in Section 3.2.A representative comparison of thepfinal and EnKF p, valid after 44 min of data assimilation (0037 UTC 30 May), is shown in Figure 10. As in the EnKF OSSE (Section 3.2), high-amplitude pressure oscillations generated during the data assimilation are not evident in the retrieved pressure fields. Furthermore, the pfinal comports with present understanding of the pressure distribution within supercells, which suggests retrieval errors are not unduly large. For example, both the pfinal and EnKF p exhibit an inflow low and a deep pressure minimum associated with the strong mesocyclonic rotation (Figures 10(a) and 10(b)), as well as upward-directed perturbation pressure gradient force along the gust front (Figures 10(c) and 10(d)). The pfinal, but not the EnKF p, persistently indicates a region of substantial downward-directed perturbation pressure gradient force near the analyzed occlusion downdraft (Figures 10(c)–10(f)). This is consistent with current understanding of occlusion downdraft formation [42], suggesting that the EnKF pressure analysis near the occlusion downdraft is contaminated with fast-mode errors that the retrieval successfully corrects.Figure 10
EnKF mean posterior (left column) and retrieved (right column)p at z=250 m (hPa; top row) and ∂p/∂z at z=1250 m (hPa km−1; middle and bottom rows). The black squares in the middle row outline the domain (centered near the occlusion updraft) in the bottom row. EnKF mean posterior vertical velocity (red; dashed = negative) is contoured at −5, −2, 2, and 5 m s−1 and radar reflectivity factor (black) at 10, 30, and 50 dBZ. Arrows represent EnKF mean posterior horizontal winds.


## 4. Conclusion

As in other data assimilation frameworks, storm-scale EnKF radar data assimilation can introduce dynamical imbalances that induce severe pressure errors in analyses, inhibiting investigation of important storm processes. Retrieval procedures that have traditionally been used to retrieve the pressure and buoyancy fields from dual-Doppler wind analyses can be used to rebalance the pressure field in such cases. Using a combination of the 2D and 3D pressure retrieval approaches that mitigates the deficiencies of both, we have demonstrated that useful analyses of the perturbation pressure field and its derivatives can be obtained from EnKF analyses despite model and discretization errors. The described procedure could be modified to separately retrieve, for example, the linear and nonlinear components of the dynamic perturbation pressure [43]. The method could also be extended to other nonhydrostatic equation sets including that used by the WRF-ARW model.It should be noted that some of the existing methods for suppressing acoustic and gravity waves in model initial conditions could be adapted to removing pressure oscillations in EnKF analyses. For example, during data assimilation, digital filter initialization (DFI; [44]) could be applied to the analysis increments at each cycle prior to integrating the ensemble forward to the next time. All such approaches, however, would increase computational cost, and it is not clear whether storm-scale analyses would improve given the model approximations used in these methods (e.g., the adiabatic backward model integration in DFI). Alternatively, a diagnostic pressure equation (e.g., [45]) could be solved (either exactly, or weakly as in the proposed method), thereby avoiding temporal discretization errors, but at the cost of increased spatial discretization errors. While we recommend evaluation of the relative efficacy of alternative methods for recovering the balanced pressure component, we also emphasize that the proposed method provides an existing, effective, simple, and low-cost way to remove pressure oscillations from storm-scale EnKF analyses and does not require data assimilation experiments to be rerun.



---

*Source: 624931-2013-10-07.xml*