During the Sydney 2000 Forecast Demonstration Project (FDP) a four-dimensional variational assimilation (4DVAR) scheme was run to analyze the low-level wind field with high spatial and temporal resolution. The 4DVAR scheme finds an optimal fit to the data and a background field under the constraints of a dry boundary layer model. During the FDP, the system assimilated data from two Doppler radars, a surface mesonet, and a boundary layer profiler, and provided low-level analyses every 10 min. After the FDP, a number of experiments have been performed to test the ability of the system to provide short-term forecasts (0–60 min) of the low-level wind and convergence. Herein, the performance of the system during the FDP and the forecast experiment's performed after the FDP are described. Two strong gust front cases and one sea-breeze case that occurred during the FDP are also examined. It is found that for the strong gust front cases, the numerical forecasts improve over persistence in the 1-h time frame, whereas for the slower-moving sea-breeze case, it is difficult to improve over a persistence forecast.
During the Sydney 2000 (S2000) Forecast Demonstration Project (FDP) a local-scale analysis system was run for the assimilation of radar, surface, and profiler data. The system is called the Variational Doppler Radar Analysis System (VDRAS) and has been under development for approximately five years. VDRAS was run operationally for two years (1999 and 2000) at the Weather Forecast Office in Sterling, Virginia, using data from a single Weather Surveillance Radar-1998 Doppler (WSR-88D) and provided analyses of low-level wind and convergence in real time. The single-Doppler assimilation system is described in detail in Sun and Crook (2001). Recently we enhanced the system to assimilate data from multiple Doppler radars. For the Sydney 2000 FDP, the variational analysis system assimilated data from two Doppler radars spaced approximately 50 km apart as well as observations from an upper-level profiler and a surface mesonet. The system was run locally at the Bureau of Meteorology in Sydney, Australia, and provided analyses for the forecasters every 10 min.
During the FDP, forecasts of the low-level wind were not performed, as the reliability of such forecasts had not yet been determined. However, after the FDP we have run a series of experiments to test the ability of the system to provide short-term (60 min) forecasts. As a practical matter, these forecasts are straightforward to perform since the initial conditions for the numerical forecast are retrieved by the analysis system. However, the reliability of such forecasts depends on a number of factors, which we seek to address in this paper.
One of these factors is the use of a dry model in the operational version of VDRAS. Although a moist version of VDRAS exists (Sun and Crook 1997), it is computationally too expensive to run in real time. In Sun and Crook (2001) it was shown that the low-level wind could be retrieved reasonably well using a dry model even in cases where moist processes were important. However, it remains an open question whether a dry model can demonstrate skill at forecasting the low-level wind in cases where moist processes are involved. We seek to address this question herein by examining two gust front cases that are driven by the moist process of evaporative cooling. As discussed in Crook (1994), the skill in forecasting gust front motion depends on an accurate estimation of the low-level buoyancy field. With a dry model it is possible to retrieve the buoyancy field given accurate estimates of the convergence field and its temporal evolution. Furthermore, the near-surface buoyancy can be estimated from surface temperature observations. Thus, we argue that given accurate radar and surface temperature measurements, the initial buoyancy field can be estimated using a dry model. The question that remains is whether the motion of the gust front is substantially affected by additional cooling that occurs during the forecast period, a forcing which is not included in the dry model.
During the FDP, the analyzed wind and convergence fields were used as input to a rule-based thunderstorm forecasting system, the Auto-nowcaster (see Wilson et al. 2004, in this issue). The horizontal convergence from VDRAS is one of the fields used to indicate areas of convection initiation as well as the sustenance of preexisting convection. The Auto-nowcaster was designed to provide short-term forecasts of thunderstorms; hence, it is important to provide not just wind analyses from VDRAS, but also short-term wind forecasts. Furthermore, since moist convection is often strongly tied to boundary layer convergence (Wilson and Schreiber 1986), the skill of explicit numerical forecasts of convection is largely dependent on the skill of numerical forecasts of boundary layer convergence. The goal of the present paper is to examine the reliability of such boundary layer forecasts.
Herein, we describe a system for variational assimilation of multiple-Doppler observations (section 2). In section 3 we describe the data sources that were available for the Sydney 2000 FDP and present three case studies from the FDP in section 4. In section 5 we examine short-term forecasts of these cases that were performed after the FDP and then verify the results. In the next section a number of experiments are performed with single-Doppler data. Although many large metropolitan areas in the United States are covered by both WSR-88D and Terminal Doppler Weather Radar (TDWR), there are extensive regions that have less coverage. Hence, it is important to determine the forecast skill given just single-Doppler information. Finally, the results are summarized and conclusions given in section 7.
2. Analysis method
The real-time system for assimilation of single-Doppler information is described in detail in Sun and Crook (2001) while the enhancements that were made to assimilate multiple-Doppler data are described in Crook and Sun (2002). For this reason we will only briefly describe the assimilation system that was used for the Sydney 2000 FDP.
The assimilation scheme finds a solution to the equations of motion that fits the data and a background field as closely as possible, in a least squares sense. The numerical model used in the real-time system consists of the dry Boussinesq equations of motion (Sun et al. 1991). The analyses fit the numerical model exactly; that is, we assume a perfect model. The equations are discretized on a Cartesian grid (i.e., with flat terrain)1 that has constant grid spacing in the vertical and horizontal. The optimal fit is achieved by reduction of the following cost function:
where Jo is the discrepancy of the analysis from the radar data (radial velocity), Ja the discrepancy from a background field, and Jp are penalty terms. Although reflectivity can be assimilated into the dry version of VDRAS by assuming that it acts as a tracer, we have found that the subsequent analyses can be quite noisy. Hence, for the FDP, we only assimilated radial velocity data.
The term Jo is defined by
where yo is the data, x are the model variables, 𝗢 is the observational error matrix (assumed diagonal in this work), and 𝗛 is the forward operator which maps the model variables on the model grid to the data variables on the data grid. (The structure of the data grid is described in more detail in section 3.) In our radar data assimilation system, 𝗛 has two components. The first component maps model variables onto the data grid by taking a weighted average of the model variables in a radar beam; that is,
where G is the power gain in the radar beam, υr,e is the observed radial velocity, υr is the model radial velocity, and Δz is the model's vertical grid spacing. The second component of 𝗛 maps the model Cartesian velocity components (u, υ, w) onto the model radial velocity υr and is given by
where r denotes the distance between a grid point (x, y, z) and the radar location (xo, yo, zo).
The second term in the cost function, Ja, is the discrepancy from a background field. The background field is important for providing the lateral boundary conditions, filling regions of missing radar data, and providing a first guess for the minimization problem. Typically, the background field is some form of large-scale field (analysis or forecast) in which the errors between grid points are correlated. For S2000, a number of options for the background field presented themselves. These included 1) using a forecast field from an operational model such as the Australian Limited Area Prediction System (LAPS) model, 2) using a forecast initialized from the previous VDRAS analysis, 3) using the previous VDRAS analysis, or 4) performing our own analysis of mesoscale datasets that were available during S2000. Option 1 would have required a significant amount of setup time that our short deadline did not allow and options 2 and 3 were not possible because of problems with data quality (discussed in the next section). We thus chose to perform our own analysis of some of the special datasets that were available during S2000. The analysis was obtained by combining mesonet wind observations at the surface with profiler data above. The surface wind observations were first interpolated to the model grid using a Barnes analysis scheme (Barnes 1964) and then merged with the profiler data above by using a local linear least squares fit. The background error covariance model that is used is pseudo Gaussian in the horizontal (Sun and Crook 2001). As will be shown, by treating the surface and profiler data separately from the higher-resolution radar data we are able to maintain the details contained in the radar data while fitting to the background field in regions without radar data. Because of the irregular nature of radar data, particular attention needs to be paid to the boundaries between regions with and without radar data. The assimilation system reduces the discontinuities in these regions in two ways: by fitting a numerical model to the data and background fields and through the use of a background error covariance model.
The last term Jp in the cost function is the temporal and spatial smoothness term and is described in Sun and Crook (2001). The weighting coefficients in the penalty terms are determined by trial and error and depend to a certain degree on the amount of smoothing that is required.
The cost function in (1) is minimized using a limited-memory quasi-Newton iterative procedure (Liu and Nocedal 1989). For each iteration, the numerical model is integrated forward and the cost function computed. The adjoint of the forward model (which gives the gradient of the cost function with respect to the initial conditions) is then integrated backward. We have found that when using real data the cost function reduction typically levels out after 30 iterations. We thus decided to terminate the minimization procedure after 40 iterations.
As noted in Crook and Sun (2002), the resultant wind fields do not fit the radial wind data exactly but are the best fit to the data and a background field under the constraints of a numerical model. This contrasts with a dual-Doppler analysis where the radial wind data are fit exactly. This means that the VDRAS wind fields are less subject to noise in the data and are dynamically consistent. Also estimates of the wind field can be obtained in regions not accessible to dual-Doppler analysis, such as outside of the dual-Doppler lobes and in regions that have only single-Doppler coverage.
3. Data sources
Figure 1 shows the observational network that was available for the Sydney 2000 FDP. The following data sources were used in the analysis scheme:
C-band Doppler radar at Kurnell near the Sydney airport—beamwidth, 1.0°; 11 elevation angles (0.7°, 1.5°, 2.5°, 3.5°, 4.5°, 5.5°, 6.9°, 9.2°, 12.0°, 15.6°, 20.0°); and volume scan rate, 5 min.
C-band polarimetric Doppler radar (C-Pol), located approximately 48 km west northwest of Kurnell—beamwidth, 1.0°; 11 elevation angles (0.6°, 1.5°, 2.4°, 3.3°, 4.2°, 5.1°, 6.3°, 7.8°, 9.2°, 10.6°, 13.4°); and volume scan rate, 10 min.
Twenty-nine mesonet stations within a domain of 150 km × 150 km, locations shown in Fig. 1—wind measurement height, 10 m; most stations had an update rate between 1 and 5 min; however, four stations updated every 30 min.
A 54.1-MHz profiler located at Sydney Airport for specification of the upper-level flow—update rate, 15 min.
The analysis procedure is conducted in the following steps:
Data collection—All data sources that are valid within the assimilation time window of 10 min are collected.
Data preparation—Some of the quality control (QC) that is performed within VDRAS is described in Sun and Crook (2001). During the Sydney 2000 FDP, two additional problems arose with data quality: (i) Unfolding of C-Pol radial velocity—C-Pol was run with a Nyquist velocity of 12.5 m s−1. Folded velocities occurred often during the FDP. If the region of aliased velocities was small (less than ∼10 km in width), then the analysis system was relatively unaffected by the contamination. However larger areas of folding did cause significant problems. An unfolding algorithm was run on both radar datasets; however, there were times when it failed to unfold aliased velocities and other times when it unfolded correct velocities. (We estimate that over the course of the FDP 5% of the radar volumes were affected by incorrectly unfolded data. At those times, the analyses in the affected regions showed winds that were in the opposite direction to the winds in the unaffected regions.) In order to mitigate the effects of folded velocity data on the analysis, we decided to turn off the cycling procedure and use the mesoscale analysis as a background field (rather than the analysis from the previous cycle). (ii) Sea clutter—The Kurnell radar, which is located close to the coast, was often contaminated by sea clutter. One of the signatures of sea clutter is a rapid decrease with height of the reflectivity field. To remove sea clutter contamination, data over the ocean that showed a decrease of 10 dBZ between the first and second tilts were removed.
After QC, the data from both radars are interpolated to the model Cartesian grid in the horizontal but are left on their original constant elevation surfaces in the vertical. The justification for this processing procedure is that the poorest azimuthal resolution of the radar data (approximately 2 km at 120 km, which is the farthest range distance in the analysis domain) is still better than the horizontal resolution of the model grid (2.5 km). In this way we reduce the size of the data volume by an order of magnitude without introducing significant interpolation error. In the vertical, however, the data resolution varies from being much finer than the model grid (375 m) near the radar to much coarser far from the radar near the domain boundaries.
The mesoscale analysis—The next step in the analysis procedure is to calculate the background wind field from the surface wind observations and upper-level profiler data. The surface wind observations are first analyzed to the model grid using a Barnes interpolation scheme. Surface temperature observations are also analyzed to the model grid and used in the buoyancy retrieval (discussed in the next step). The radius of influence in the scheme was set at the average station spacing (∼30 km). Above the surface, the background wind field is obtained by applying a local-linear least squares fit (with a 1500-m radius of influence) to the surface value and the profiler data above. In Sydney, the lowest level of usable profiler data was typically at 600 m or above.
Assimilation of the radar data—The final step in the analysis is the assimilation of the Doppler radar data. As already noted, this is achieved by minimizing the cost function (2) by running the numerical model forward and the adjoint of the model backward. An estimate of the background error variance was obtained by comparing typical rms errors in our analyses (the background field) with typical observational errors in radial velocity. A comparison of our VDRAS wind analyses in Sterling, Virginia, against independent aircraft wind data (which admittedly have their own sources of error) indicated an average rms difference of 3–4 m s−1. Observational error in radial velocity for typical values of spectrum width and signal-to-noise ratio in the boundary layer are of the order of 1 m s−1 (see Fig. 6.6 in Doviak and Zrnić 1984). Hence, we estimate that the ratio of background to observational variance (rms difference squared) is of the order of 10. As noted earlier, during the FDP the cycling procedure had to be turned off, so that the mesoscale analysis was used as the background field. At these times, the ratio of background to observational variance was set to a factor of 50. Clearly, more research is needed in the area of background and observational error estimation at these scales and we are actively working in this area.
By fitting the equations of motion to a time series of radar data, VDRAS has the capability of retrieving the buoyancy field (see, e.g., Sun and Crook 1996). For S2000, the buoyancy retrieval was enhanced by using the surface temperature analysis as a background field. In this way, the retrieved buoyancy field fits closely to the surface temperature analysis at the surface as well as being dynamically consistent above the surface.
For the field test, the system was run on a Digital Equipment Corporation (DEC) Alpha workstation at the Sydney Bureau of Meteorology. In order to run in real time, it is essential that the algorithm complete its calculations faster than the time window used for assimilation. For the Sydney system, we used a time window of 10 min, which typically captured two volumes from C-Pol and three volumes from Kurnell. Consequently, we tuned the size and resolution of the analysis domain so that 40 iterations could be completed in just less than 10 min of CPU time. This resulted in a domain size of 120 km × 120 km × 3 km with a horizontal resolution of 2.5 km and vertical resolution of 375 m.
With a few exceptions, the analysis system ran continuously from mid-August 2000 until February 2001. The exceptions were primarily caused by interruptions in the data flow and with problems with folded velocity data from C-Pol. In this section we will describe three cases that occurred during the FDP and show the VDRAS analyses that were calculated in real time. In the following section we will present short-term forecasts of these events that were performed after the FDP. As discussed previously, in real time it was necessary to turn off the cycling procedure in order to prevent folded radial velocity data from propagating from one analysis to the next. In the rerun cases, the radial velocity data have been correctly unfolded so it was possible to include cycling in the assimilation procedure. Because of the different handling of folded velocities, we will not make a quantitative comparison of the real-time and rerun analyses, but rather concentrate on the forecasts that have been performed from the rerun analyses.
a. The 3 November hailstorm and tornado case
On 3 November 2000, a tornadic hailstorm developed south of Sydney and moved north-northeastward. The propagation direction of the hailstorm was more northerly than most of the storms on that day as it followed the boundary layer forcing produced at the intersection of the sea breeze with the storm's outflow. The events of this day are discussed in more detail in a companion paper in this special issue (Sills et al. 2004, in this issue). The position of the sea breeze and outflow are shown in Fig. 2 by thick solid lines at two times: (a) 1428 and (b) 1540 LST. (Although boundaries were inserted in real time by forecasters at the Sydney Bureau of Meteorology, the fine lines shown throughout this paper were drawn after the FDP based on low-level reflectivity and radial velocity scans.) The hailstorm can be seen just to the south of this intersection point in the 0.6° reflectivity scan from C-Pol (grayscale) in both Figs. 2a and 2b. The analyzed wind vectors and horizontal convergence at the lowest VDRAS level (z = 180 m) are also overlaid.2 As can be seen, strong low-level convergence has been analyzed at the intersection point of the gust front with the sea breeze. This convergence “bull's-eye,” which reached a maximum of 1.1 × 10−3 s−1, successfully tracked the motion of the gust front–sea-breeze intersection point as it moved northward.
b. The 30 November gust front case
On 30 November, two strong gust fronts developed and moved into the analysis domain, the first in the southwest quadrant (between 1300 and 1400 LST) and the second in the northwest quadrant (around 1430 LST). These two gust fronts then converged with the sea breeze, which was active on this day, to produce a large area of convective rain around 1500 LST. Figure 3 shows the analyzed low-level winds and convergence field at 1330 LST. The low-level (0.6°) reflectivity from C-Pol (in dBZ) is shown by the gray shade. Although a convergence signature along the sea-breeze front is not analyzed, the sea-breeze flow is shown by the onshore winds in the eastern half of the domain. Note, however, that offshore winds have been analyzed along the coast south of Kurnell. This is a result of the fact that the only data in this region come from a mesonet station that is just to the west of the sea breeze and shows offshore flow at this time. The strong gust front in the southwest is evident in the analysis at 1330 LST (the gust front in the northwest has not yet developed at this time). The maximum convergence along the southwest gust front is 0.35 × 10−3 s−1. After 1400 LST, an increasing proportion of the C-Pol radial velocity data became folded, which in turn corrupted the velocity and convergence fields in the real-time analyses. In the next section we will present analyses (and forecasts) of this case in which the radial velocity data has been unfolded.
c. The 18 September sea-breeze case
Figure 4 shows the VDRAS winds and convergence fields at the lowest level (z = 180 m) on 18 September 2000, at (a) 1634 and (b) 1734 LST. During this event, C-Pol data were available only intermittently. The low-level (0.6°) reflectivity from C-Pol is also plotted and shows a fine line associated with the sea breeze over Sydney. Typically the sea breeze would first become evident in C-Pol data along the coastline south of Sydney. In many cases a bulge would form in the sea breeze as it moved inland over the lower terrain of Sydney but was held back by the higher terrain north and south of Sydney. On this day, the leading edge of the sea breeze moved inland at approximately 2.5 m s−1.
The analysis system has captured the sea-breeze flow along with the horizontal convergence at the leading edge of the breeze. The convergence field is shown (with a contour interval of 2 × 10−4 s−1) and clearly delineates the position of the sea breeze at both analysis times. The convergence reaches a maximum of 1.5 × 10−3 s−1 along the leading edge of the sea breeze at 1734 LST.
The C-Pol reflectivity data in Fig. 4 also show an east–west-oriented convergence line in the westerly flow west of the sea breeze. This convergence line, which was occasionally observed during the Sydney 2000 FDP, appears to be a result of the interaction of a westerly flow with the higher terrain to the west of Sydney. On this day, the westerly flow ahead of the sea breeze was quite strong (6–7 m s−1 at z = 180 m). Note that the analysis system has picked up some of the convergence associated with the fine line ahead of the sea breeze although the analyzed convergence is approximately 10 km south of the radar-observed fine line.
d. Verification of real-time analyses
In Crook and Sun (2001), verification statistics for the real-time analyses were presented. (We have not attempted to calculate statistics for the 30 November case because of the corruption by folded velocities.) The rms difference between the observed and analyzed radial velocity was calculated for each analysis volume. This rms difference was then averaged over four consecutive days (29 October–1 November) when there was no contamination from folded velocities. The average rms difference between the analyzed radial wind and radial wind observed by Kurnell (C-Pol) was 1.54 (1.91) m s−1.
Although these statistics indicate that the fit of the observed radial velocity to the analyzed radial velocity is quite close, it is not an independent test of the accuracy of the analysis. As an independent test, wind observations from aircraft taking off and landing at Sydney airport have been used. (There was only two aircraft wind observations that were valid within 10 min of the analyses shown in Figs. 2–4. These wind observations are plotted in Figs. 4a and 4b as a thick vector, labeled with wind speed and direction.) Table 1 shows the mean vector difference (MVD) between the VDRAS and aircraft-observed wind vector for flights on the 4 days described in Crook and Sun (2001) (18 September, 3 and 8 October, and 3 November). The mean vector difference between the two horizontal vectors (u1, υ1) and (u2, υ2) is given by
where N is the total number of grid points summed over. Table 1 shows that the mean vector difference ranges from 2.1 to 3.5 m s−1 while the mean aircraft-observed wind speed ranges from 4.5 to 9.0 m s−1. Calculated over the 4 days, the MVD is 2.6 m s−1 while the mean wind speed is 6.3 m s−1. We consider this to be a reasonable agreement given that the aircraft-observed wind also contains observational error and is a line average along the aircraft's flight track, whereas the analysis wind is an average over a model grid volume.
5. Forecasts and verification
In this section we present short-term forecasts of the three cases described above. In these reruns, a cycling procedure was used in which the VDRAS analysis from the previous cycle is used as a background field and as a first guess for the minimization procedure. A schematic of the cycling procedure is shown in Fig. 5. Each assimilation cycle spans 10 min and incorporates three volumes of Kurnell data and two volumes of C-Pol data. A 60-min forecast is performed every 10 min. Since wallclock time was not such an issue in these rerun cases we used a larger domain of 150 km × 150 km × 6 km. The horizontal and vertical resolutions were kept the same (2500 and 375 m, respectively).
a. The 3 November case
Figures 6 and 7 show forecasts of the 3 November case initialized at 1425 and 1515 LST, respectively. The first panel shows the initial conditions for horizontal velocity and convergence. The second and third panels depict the 30- and 60-min forecasts, respectively. In these panels, the forecast convergence is shown by the line contours, while the analyzed convergence at the verifying time is plotted in grayscale.
By comparing Fig. 2a with Fig. 6a and Fig. 2b with Fig. 7a the difference between the rerun analyses and the real-time analyses can be seen.3 The convergence analyzed along the gust front in these rerun cases is more compact and has a stronger amplitude than in the real-time analyses. (Note that the contour interval in Figs. 6 and 7 is 2.0 × 10−4 s−1.) The maximum convergence analyzed along the gust front is 1.6 × 10−3 s−1 at 1425 LST and 1.4 × 10−3 s−1 at 1515 LST. These values are approximately 1.5 times larger than the maximum values in the real-time run. The primary reason for this increase in amplitude is that the background field in these rerun cases already contains the convergence signature along the gust front (from the previous analysis) whereas the background field in the real-time analyses was from a large-scale analysis without significant horizontal variability.
Figures 6b and 6c indicate that for the forecast initialized at 1425 LST, the simulated gust front (line contours) does not propagate as quickly as the analyzed gust front (grayscale contours). The analyzed gust front moves toward the north at approximately 6 m s−1 whereas the forecast gust front moves at only 2 m s−1. The maximum convergence along the gust front also decreases significantly over the 60-min forecast period to a value of 0.6 × 10−3 s−1, a reduction of a factor of 3 compared with the initial convergence.
The forecast initialized at 1515 LST (Figs. 7b and 7c) does a much better job at tracking the location of the analyzed gust front. In this latter simulation, the forecast gust front moves at ∼5 m s−1 compared to the analyzed speed of 6 m s−1.
To examine the reason for the forecast failure at 1425 LST we first plot the temperature perturbation field at the lowest VDRAS level (z = 180 m). As discussed in section 3, this is a combination of the retrieved buoyancy field and the surface temperature analysis. Figure 8a shows the temperature perturbation at the initial time, while Fig. 8c shows the 60-min forecast field (at 1525 LST). The first point to note is that at the initial time, the temperature perturbation behind the gust front is approximately 1.5°C. This is most likely an underestimate of the true temperature perturbation across the gust front. As the hailstorm moved over mesonet station G11 (marked in Fig. 1) the observed temperature dropped by 5°C in 5 min.4 There are a number of reasons why the temperature deficit in the initial conditions at 1425 LST is significantly less than that measured at this one station. First, the model temperature deficit of 1.5°C is an average over the lowest grid volume and not a point measurement. Second, the model is not able to retrieve the full buoyancy difference across the gust front due to the fact that the convergence at the leading edge of the gust front is underestimated by the radar observations. Third, the temperature deficit across the gust front in the surface analysis is reduced, since it takes into account the surrounding stations that were not impacted by the hailstorm. In Fig. 8b, we plot the surface analysis at the initial time of 1425 LST. The first point to note is that the temperature gradient in the surface analysis is broader than in the initial conditions (Fig. 8a) since it is based only on observations that are spaced approximately 30 km apart. Although the analyzed temperature decreases to the south it does not delineate the temperature gradient across the gust front as well as in the VDRAS analysis. The difference between the surface temperature analysis and the VDRAS analysis is between 0.5° and 1.0°C. As already discussed, this difference is due to the fact that the VDRAS temperature field is the average over the lowest grid volume and is a combination of both the surface analysis and the dynamically retrieved buoyancy field.
As discussed above, the initial conditions are not able to capture the full temperature deficit across the gust front. This is one of the reasons why the model underpredicts the motion of the gust front (especially at 1425 LST). However, this is not the only reason, since the forecast does significantly better at 1515 LST, and the analyzed temperature deficit across the gust front is a similar 1.5°C at this time. We speculate that the improved forecast at 1515 LST is due to the fact that there is a larger extent of southerlies analyzed behind the gust front at this time than at 1425 LST. Therefore, even though the full buoyancy is not retrieved in the initial conditions, there is sufficient integrated momentum, or inertia, in the gust front at the latter time to allow the front to propagate. To test this possibility, a series of forecasts was run with the initial buoyancy set to zero. In these experiments, the gust front forecasts were significantly degraded indicating the importance of the buoyancy field in gust front propagation. However, due to the gust front's initial intertia, the northerly flow did persist for some time during the forecast period. The persistence of the gust front flow was significantly greater in the forecast initialized at 1515 LST compared to the 1425 LST forecast.
A third possible reason for the underprediction of the gust front motion at 1425 LST is that additional cooling may have occurred during the forecast period that is not predicted by the dry forecast model. To examine this possibility, we have plotted in Fig. 8d the temperature perturbation analysis (from VDRAS) at the verifying time (1525 LST). By comparing Fig. 8c with Fig. 8d it can be seen that there is an additional cooling of ∼0.5°C that is not predicted by the model. This additional cooling may have been caused by the effects of evaporation that are not incorporated in the dry model. However, we would argue that this effect is not the primary reason for the underprediction of the gust front motion since the additional cooling is only 1/3 of the total temperature deficit across the gust front. Furthermore, this effect also occurs in the forecast at 1515 LST, which does a much better job at predicting the motion of the gust front.
Verification statistics for the forecasts on 3 November are plotted in Fig. 9a, b. The first panel (a) shows the MVD averaged over all of the forecasts between 1400 and 1600 LST. The second panel (b) gives the rms difference in the convergence fields averaged over all the forecasts. Also plotted are the statistics for a persistence forecast.5 As can be seen, the numerical forecast of the wind and divergence fields improves over persistence at both 30 and 60 min. The improvement over persistence is in the range of 20%–30% for both fields.
Also plotted in Fig. 9 are forecast statistics for analyses in which no radar data are assimilated (i.e., from the mesoscale analysis field). As can be seen, the MVD and convergence error are significantly less for these forecasts, which may cast some doubt on the utility of the VDRAS forecasts. However, the background analyses and forecasts completely miss the strong wind and convergence of the gust front. For example, the strongest convergence in the background analyses is 0.06 × 10−3 s−1 compared to 1.2 × 10−3 s−1 in the analyses using radar data. The wind and convergence field in the analysis at 1425 LST without radar data is shown in Fig. 10. As can be seen, the maximum convergence falls below the minimum contour interval of 1.0 × 10−4 s−1. Since the horizontal variability in the background field is low, forecasts from that field verify reasonably well. However, the utility of such analyses and forecasts, which do not capture the gust front structure at all, is somewhat questionable.
b. The 30 November case
Figure 11 shows a 1-h forecast of the gust front case on 30 November 2000. The forecast shown in Fig. 11 is initialized with the analysis at 1335 LST, which is the same time shown in the real-time analysis in Fig. 3. As in the previous case, the convergence along the gust front in this rerun analysis is significantly stronger than in the real-time analysis (1.0 × 10−3 s−1 compared to 0.3 × 10−3 s−1). Figure 11b shows that the forecast position of the convergence line at 30 min agrees quite well with the analyzed position. At 60 min, however, the agreement is not as good. The forecast convergence has decayed significantly (to a maximum value of 0.1 × 10−3 s−1) and lags behind the analyzed convergence by approximately 10 km.
The verification statistics for this case are plotted in Fig. 12. The statistics are averaged over 10 forecasts between 1300 and 1500 LST. As in the previous case, the numerical forecasts of wind and convergence improve over persistence at both 30 and 60 min. (Even though the convergence field is decaying in the forecast shown in Fig. 11b, it still improves over persistence.) Again the improvement is in the range of 20%–30%.
c. The 18 September case
Figure 13 shows a 1-h forecast of the 18 September sea-breeze case initialized at 1535 LST. A comparison of Figs. 4a and 13a shows that the convergence fields from the real-time and rerun analyses are similar. However, since C-Pol data were assimilated only sporadically into the real-time system it is difficult to draw too many conclusions from this similarity.
The 30- and 60-min forecasts of this event (Figs. 13b and 13c, respectively) are not very successful. The observed sea breeze moves very slowly inland at a speed of approximately 1.5 m s−1 while the east–west convergence line ahead of the sea breeze remains stationary. However, in the numerical forecast, the sea-breeze convergence decays during this 60-min period, while the east–west convergence is pushed northward. A likely reason for the forecast failure is that the sea-breeze evolution is strongly controlled by land–sea temperature contrasts that are not included in the current analysis and forecast system.
The verification statistics for this event are plotted in Fig. 14. Persistence does better than the numerical forecast for wind throughout the 60-min period and for convergence at 30 min. At 60 min, the error in the persistence and numerical forecast convergence fields are approximately the same. The results from this case indicate that for a system that evolves only slowly with time, it is difficult to improve over persistence with the current system.
6. Experiments with a single radar
For S2000, data from two closely spaced Doppler radars were available. Although many large metropolitan areas in the United States are covered by both WSR-88D and TDWR radars, there are large regions that have less coverage. In this section, we examine how well the analysis system would do given only single-Doppler information.
We first examine forecasts produced using just C-Pol data. Qualitatively, the analyses and forecasts (not shown) are similar to those produced by the two-radar system. The difference between the initial analyses is small (0.1 × 10−4 s−1 in rms divergence and 0.16 m s−1 in mean vector wind difference) and much smaller than in the Kurnell-only analysis (shown next). This is not surprising since the radar coverage from C-Pol was significantly greater than that from Kurnell. The forecasts also show a similar behavior. At the earlier time, 1425 LST, the forecast convergence line decays over the 1-h period and does not keep apace with the analyzed line. At the latter time, the forecast is significantly improved and the strength of the line stays roughly constant.
Statistics on the wind and convergence forecasts are plotted in Fig. 15. Statistics are calculated from the 10 forecasts between 1400 and 1600 LST. To calculate these statistics we have verified against the C-Pol-only analyses since we are interested in how well the forecasting system would do if only single-Doppler information was available. Figure 15 shows that the forecasts using just C-Pol data have only mixed success. The divergence forecasts improve over persistence (by 27% at 60 min, compared with 25% from the two-radar forecasts). However, the wind forecasts are slightly worse than persistence at both 30 and 60 min. To compare the two-radar statistics (Fig. 9) against the C-Pol-only statistics (Fig. 15) it should first be noted that these statistics use different verification analyses (two-radar analyses and C-Pol-only analyses, respectively). Comparing Fig. 15 with Fig. 9 shows that the wind statistics from the two-radar forecasts improve over the C-Pol-only forecasts. However, for the convergence field, the two radar forecast statistics are actually larger than the single-radar statistics. Again, this is for the same reason that the no-radar forecasts show lower verification scores; the single-radar analyses miss some of the convergence features in the two-radar analyses and thus the single-radar forecasts can show improved verification scores.
The analysis and forecast using just Kurnell data at 1515 LST are plotted are Fig. 16. The convergence field from the two-radar analysis system is plotted in grayscale. As can be seen, the convergence line is only barely detected at this time. This is not surprising since 1) the detection of the gust front from Kurnell was limited and 2) during the period between 1400 and 1600 LST, the gust front was roughly parallel to the Kurnell radar beam. The difference between the Kurnell-only analysis and the two-radar analysis is 0.23 × 10−4 s−1 in rms divergence and 0.7 m s−1 in mean vector wind difference, which is significantly larger than the difference between the C-Pol-only and two-radar analyses.
Since the gust front is only barely captured in the Kurnell-only analysis, the forecast produced from this analysis is poor as well (Figs. 16b and 16c). Statistics for the 10 forecasts between 1400 and 1600 LST are plotted in Fig. 17. As can be seen there is little difference between persistence (based on the Kurnell-only analyses) and the numerical forecasts, since neither the analysis nor the forecasts have been able to capture the gust front convergence.
Although the numerical model has shown some skill (compared to persistence) when data from two radars are assimilated, the single-Doppler experiments have been less encouraging. It should be recognized, however, that these are very preliminary attempts at forecasting gust front motion. Improvements in skill should be gained by improvements in model physics (especially the inclusion of moist physics) and the data assimilation procedure (especially improved background and observational error statistics). However, cases where the gust front is parallel to the radar beam will continue to pose a challenge to numerical forecast schemes initialized with radar data.
7. Summary and conclusions
In this paper we have described a real-time variational wind analysis scheme (VDRAS) and shown its performance during the Sydney 2000 Forecast Demonstration Project. We have also presented short-term (1 h) forecasts of three events that were performed after the FDP. The main findings from this study are listed below.
With a few exceptions, the analysis system ran reliably from mid-August 2000 until February 2001 providing wind fields every 10 min. The exceptions were primarily caused by problems with incorrectly unfolded data from the C-Pol radar and interruptions in the data flow.
The analysis system produced realistic wind fields under a variety of flow situations. In the present study and that of Crook and Sun (2002) we have described the performance for two strong gust front events, two sea-breeze cases, and one cold front event.
The analyzed wind fields have been verified using both dependent radar data and independent aircraft-observed data. Results showed that the fit to the radial wind data was in the range of 1.5–1.9 m s−1. Comparison with independent aircraft-observed data for four events indicated a mean vector difference of 2.6 m s−1.
A series of 1-h forecasts of three events (two gust fronts and one sea breeze) have been examined. For the two gust front events, which propagated rapidly with time, the numerical model showed consistent skill in forecasting wind and convergence compared with persistence. However, persistence was a better forecast for the slower-moving sea-breeze event.
For the two gust front cases, it was found that the numerical model did not propagate the gust front as quickly as observed and that the convergence along the leading edge of the gust front tended to decay with forecast time. We offered two explanations for this underprediction. 1) It is difficult to obtain the full temperature difference across the gust front using surface observations and the retrieved buoyancy. 2) Additional cooling can occur during the forecast period that is not predicted by the dry model. However, the 3 November case showed that even with a poor retrieval of the buoyancy field, if there is a sufficient fetch of flow analyzed behind the gust front, then the numerical model can propagate the gust front at a speed close to that observed.
Analyses and forecasts using single-Doppler information have also been examined. It was found that the analyses using just C-Pol data were qualitatively similar to the two-radar analysis. However, in terms of forecasting skill, the results were more mixed. Although the convergence forecasts did show skill compared to persistence, the wind forecasts did not. Finally, the analyses and forecasts using just Kurnell data were poor, because of the limited coverage from Kurnell and the fact that the radar beam roughly paralleled the gust front.
In a previous study by Crook and Tuttle (1994) 1-h forecasts of two gust front cases in northeastern Colorado were examined. These forecasts were initialized with a more traditional analysis scheme using a Cressman filter applied to radial velocity, Tracking Radar Echoes by Correlation (TREC) data, and retrieved buoyancy fields. The TREC method estimates the cross-beam component of velocity by tracking reflectivity features using a cross-correlation technique (Tuttle and Foote 1990). Data from a 50-station surface mesonet were also used. For the wind forecasts, the improvement over persistence was in the range of 20%–30%. However, the convergence forecasts were not as successful, with very little improvement at 30 min and only a 10% improvement at 60 min. Caution should be taken when comparing the two studies since the analysis techniques, data sources, and events are different. However, possible reasons for the improvement in the convergence forecasts with the current system include the more consistent fitting of the forecast model to the data, the filtering achieved by the numerical model and smoothness penalty terms, and finally the fact that the TREC winds often contained considerable noise.
The improvement in the convergence forecasts with the current system has a number of important implications. First, it suggests that convergence forecasts from the current system could add value to automated algorithms that forecast the growth, propagation, and decay of moist convection. One of these algorithms is the Auto-nowcaster, which is described in another paper in this special issue (Wilson et al. 2004). Currently, the VDRAS analyses (but not the forecasts) are being used as one of the input fields to the Auto-nowcaster for making 1-h forecasts of convection. In these forecasts it is assumed that the convergence field is constant (i.e., persistent) over the 1-h forecast period. The results from the current study suggest that in some cases (e.g., rapidly moving gust front cases) it would be better to use the short-term numerical forecasts rather than assuming persistence. However, this is not a general result since in some cases persistence is a better forecast. Hence, before the numerical forecasts can be considered for use in the Auto-nowcaster, it is necessary to perform a long-term study over many convective cases to determine if the forecasts on average improve over persistence.
A second implication of the improved convergence forecasts is for explicit numerical forecasts of convection. The growth, propagation, and decay of moist convection are strongly tied to boundary layer convergence (Wilson and Schreiber 1986). Hence, the skill of explicit numerical forecasts of convection is largely dependent on the skill of numerical forecasts of boundary layer convergence.
The success of the wind forecasts produced by the current system has been more mixed. For the two gust front cases, the wind forecasts using data from both radars did improve over persistence. However, the single-Doppler forecasts of the 3 November case did worse than persistence as did the two-radar forecasts of the 18 September sea-breeze case. Short-term forecasts of the low-level wind are useful in a number of circumstances, for example, outdoor sporting events, recreation/leisure activities, and wildfire fighting. In fact, during the period 18–19 September, there were a number of wildfire events in the region of Hornsby, 20 km northwest of central Sydney (R. Webb 2002, personal communication). Wind forecasts, particularly forecasts of wind shift lines, would have been very useful during this event. For this reason we are paying particular attention to improving the skill of wind forecasting using the current system.
This paper represents a status report on our development of a real-time variational wind analysis and forecasting scheme. We are currently exploring a number of enhancements to the present system. As discussed in Crook and Sun (2002), we are examining ways to improve the background analysis by incorporating additional data sources and improving the fit to the surface observations. Studies are being conducted to improve the estimation of the observational and background error statistics. We are also examining methods to reduce the anomalous generation of horizontal convergence in the conversion of the radar data from spherical to Cartesian coordinates.
In order to improve the short-term forecasts we have been focusing on improvements in the initialization of the buoyancy field. To reliably predict the motion of gust fronts and convergence lines it is necessary to retrieve the low-level buoyancy field that drives these phenomena. As already noted, in a dry model the retrieval of the buoyancy field depends on an accurate estimation of the convergence field and its temporal evolution. Hence, improvements in the estimation of the convergence field should lead to improvements in the skill at forecasting gust fronts and convergence lines. In a moist model, evaporative cooling can significantly affect the low-level buoyancy field. Thus, as a first step toward the implementation of a full cloud model, we are planning to include an evaporative cooling term in the current real-time system, which should lead to an improvement in the estimation of the buoyancy field.
This project would not have been possible without the assistance of the following people: Terri Betancourt, Jaimi Yee, Frank Hage, Rita Roberts, Jim Wilson, and Cindy Mueller. Reviews of an earlier version of this manuscript by Rita Roberts and Jim Wilson led to substantial improvements. Thanks are also given to Tom Keenan for supplying Fig. 1. This research is sponsored by the National Science Foundation through an interagency agreement in response to requirements and funding by the Federal Aviation Administration's Aviation Weather Development Program. The views expressed are those of the authors and do not necessarily represent the official policy or position of the U.S. government.
Corresponding author address: Dr. N. Andrew Crook, NCAR MMM/RAP, P.O. Box 3000, Boulder, CO 80307-3000. Email: firstname.lastname@example.org
The numerical model and adjoint do not include the effects of variable terrain. For this reason, we have only run the system in regions where terrain effects are small. Over most of the S2000 VDRAS domain, the variation in terrain is less than 200 m. There is a small region with terrain that rises 500 m above the rest of the domain; however, this region is limited to the western edge of the domain.
Throughout the paper, all analyzed fields are shown at the lowest VDRAS level of z = 180 m. The difference between the analyzed wind fields and the surface observations, which is evident in some of the plots, is due to both the difference in height as well as the smoothing inherent in the surface analysis.
The real-time analyses for this case were not significantly impacted by folded velocity data, which means that the primary difference between the real-time and rerun analyses is the inclusion of cycling in the latter analyses. Note there is a 3–7-min time difference between the valid times of the rerun and real-time analyses; however, this time difference does not have a significant effect on the analyses.
G11 was one of the mesonet stations with a 30-min update rate. Although the 1425 LST analyses did include the temperature decrease of 5°C observed at G11, it does indicate the importance of a fast update rate (<5 min) in forecasting these rapidly evolving systems.
In this study, we have only compared our numerical forecasts with a persistence forecast. A second possibility would be an extrapolation forecast; however, it is difficult to make such comparisons objective. For example, to produce an extrapolation forecast, a motion vector must be chosen. The choice of such a motion vector is often very subjective. (During Sydney 2000, motion vectors for convergence lines had to be hand drawn by a forecaster.) Once a motion vector has been chosen, it is then necessary to decide what portion of the velocity or convergence field to extrapolate. This means that it is possible to obtain a wide variety of forecasts and verification statistics, depending on the subjective choices that are made.