Natural Hazards and Earth System Sciences A severe blizzard event in Romania – a case study

During winter cold strong winds associated with snowfalls are not unusual for South and Southeastern Romania. The episode of 2–4 January 2008 was less usual due to its intensity and persistence. It happened after a long period (autumn 2006–autumn 2007) of mainly southerly circulations inducing warm weather, when the absolute record of the maximum temperature was registered. The important snowfalls and snowdrifts, leading to a consistent snow layer (up to 100 cm), produced serious transport and electricity supply perturbations. Since this atypical local weather event was not correctly represented by the operational numerical forecasts, several cross-comparison numerical simulations were performed to analyze the relative role of the coupler/coupling models and to compare two ways of process-scale uncertainties mitigation: optimizing the forecast range and performing ensemble forecast through the perturbation of the lateral boundary conditions. The results underline, for this case, the importance of physical parametrization package on the first place and secondary, the importance of the model horizontal resolution. The resolution increase is beneficial only in the local process representation; on larger scale it may either improve or decrease the accuracy effect, depending on the specified nudging between large-scale and small-scale information. The event capture is likely to be favored by two elements: a more appropriate time-scale of the event’s physics and the quality of the transmitted large-scale information. Concerning the time scale, the statistics on skill as a function of forecast range are shown to be a useful tool in order to increase the accuracy of the numerical simulations. Ensembles forecasting versus resolution increase experiments indicate, for such atypical events, an interesting supply in the forecast accuracy through the ensemble method when applied to correct the minimum skill of the deterministic forecast. Correspondence to: S. Tascu (simona.tascu@meteoromania.ro)


Introduction
At polar and temperate latitudes the blizzard is a severe weather event of highest risk, often leading to important damages.Every year Romania experiences at least one blizzard episode.Consequently the conditions when it is produced and its characteristics were studied by Romanian meteorologists, many times in connection with the local wind crivatz dynamics (Draghici, 1983(Draghici, , 1984;;Cordoneanu et al., 1997;Popa and Soci, 2002).
The blizzard of 1-4 of January 2008 was chosen for the present study due to its severity and to the fact that it occurred after a long period of warm weather.The blizzard affected the southern part of Romania and the north of Bulgaria by considerable snowfalls and snowdrifts that were associated with low visibility; it seriously perturbed the socio-economic activity.For instance, in Bucharest just in twelve hours the snow layer reached seventy centimeters.
An essential condition for generating blizzards is the presence of intense gradients of pressure, and often, temperature.In Western, Northern and Central Europe they are produced at the extremity of very active cyclones, usually developed in the North Atlantic basin, afterward engaged on continental trajectories (Pinto et al., 2007).In these situations the severity of the weather is due to the energetic input of the surface heat fluxes (latent and sensible) contrast between land and ocean.
In the South-East Europe the winter storms appear in specific conditions in the coupling zone of an anticyclone (i.e.East-European, Azores or Scandinavian origin) and a Mediterranean cyclone.Especially in this season, the largescale lows and troughs evolving along the polar front jet seem to directly influence the cyclogenesis in the Mediterranean Sea (Trigo et al., 2002;Flocas and Karacostas, 1996).In the presence of a cut-off structure the mass penetration from the middle latitudes has a decisive influence on precipitation (Alpert et al., 1990;Alpert and Reisin, 1986).The climatological studies (Maheras et al., 2001) reveal the presence of three cyclogenetic regions in Western, Central and Eastern Mediterranean.The Eastern Mediterranean is in generally a weaker cyclogenetic region (Radinovic, 1965;Alpert et al., 1990) in comparison with the cyclogenetic regions of the Western Mediterranean.Generally, the surface fluxes due to the sea-land contrast contribute to the cyclogenesis process.On the other hand, in winter the orography plays a major role in the Genoa Gulf cyclogenesis (where the highest number of intense cyclones is registered) while in the Eastern Mediterranean the upper level dynamics seems to constitute a dominant feature (Maheras et al., 2002).
When the generated atmospheric flow in the coupling zone of an anticyclone and a Mediterranean cyclone is a blocking one, the duration of the blizzard may be of several days (Tayanc ¸et al., 1998).In Romania, besides the large-scale conditions, an important role in the evolution of the phenomenon is the orographic forcing induced by the presence and the shape of the Romanian Carpathians (see Fig. 1) and the diabatic one due to the proximity of the Black Sea (Draghici, 1983(Draghici, , 1984)).
The data and models used within the paper are shortly presented in Sect. 2 while the analysis of the specific characteristics of the blizzard of 1-4 January 2008 is the subject of Sect. 3. Since this atypical local weather event was not correctly represented by the operational numerical prediction models (usually having a good skill), several numerical sim- ulation have been carried out in order to find a solution for a better representation of such an extreme event.They allow us to analyze the sensitivity to coupler/coupled models, forecast range, initial and boundary conditions and to compare alternate solutions for uncertainty mitigation through resolution increase and ensemble prediction systems (EPS).The experiments organization corresponding to their aim and the obtained results are presented in Sect. 4. Finally, the conclusions of this study are summarized in Sect. 5.

Data and methods
For the large-scale conditions analysis of this blizzard episode the analysis of the ECMWF (European Center for Medium Range Weather Forecast) deterministic model and the satellite water vapor channel imagery (METEOSAT 9, WV6.2 µm channel) have been used.For a more detailed analysis we benefited from the numerical limited area model (LAM) ALADIN (Bubnova et al., 1995;Horanyi et al., 2006) operationally integrated in Romania at a 10 km horizontal resolution (the model domain and orography are presented in Fig. 1).The surface observations (with a quite high density over Romania) complement the information for the case analysis.
For the numerical simulations (explained in detail in Sect.4) the following models were used: -The LAM ALADIN coupled either with the IFS/ARPEGE or with the IFS/ECMWF deterministic global models for short range forecasts.The integration domain (1300×1300 km) covers Romania and its surroundings (Fig. 1).
-The regional model RegCM3 at 10 km resolution (the integration domain is presented in Fig. 2) was coupled with ECMWF deterministic model.For longer range forecasts up to 10 days it was used for downscaling (at 50 km resolution) a reduced number of the global ECMWF ensemble prediction system (EPS).
The development of IFS (ECMWF's Integrated Forecasting System) was begun in 1987, in collaboration with Météo -France.This is the basis of all numerical weather applications at global level in both organizations.Sharing the dynamical and data assimilation structure, the IFS/ECMWF and IFS/ARPEGE models use different physical parametrization packages.
The ALADIN model, developed as a limited area counter part of the IFS/ARPEGE model within an international project (initiated in 1990 by Météo -France which now encompasses 16 countries, including Romania), was especially designed for mesoscale phenomena simulation.It was built entirely on the notion of compatibility with his "mother" system ARPEGE (Courtier and Geleyn, 1988), but in a flexible manner, being possible to integrate it at different resolutions and geographical positions.Thus, the high compatibility between the models offers a special framework for the development and research activities.
The system based on RegCM3 (Giorgi, 1990;Giorgi et al., 1993) coupled with IFS/ECMWF deterministic model or with its EPS (Buizza, 2007) has been developed at the Romanian National Meteorological Administration during the last years for medium and long-term applications.It was tuned and tested in research studies over Romania.The system is able to represent multi time-scale physical processes interactions, especially those involving the sea surface temperature time-variations with strong impact in winter for Romania (R. Bojariu, personal communication) and in deep-soil hydrology.Those characteristics are absolutely necessary for the representation of the analyzed event in medium range simulations.
For the validation of the models precipitation forecasts for the entire integration domain there were used data from the Romanian observation network (168 synop and 150 rain gauges) and from the GTS -Global Telecommunication System (about 115 synop stations).For the graphical representation the precipitation data were interpolated in a geographical grid (0.2 • ×0.2 • ).

The case analysis
In the evolution of the analyzed blizzard episode four stages were revealed: -a precursory stage (31 December 2007, 18:00 UTC -1 January 2008, 06:00 UTC) when the upper air trough extended towards the Central Mediterranean enhancing the continental polar air advection; -a beginning stage (1 January 2008, 06:00 UTC -2 January 2008, 00:00 UTC) when the snowfalls and wind gusts started to occur more and more frequently; -a maturity stage (2 January 2008, 00:00 UTC -3 January 2008, 12:00 UTC) when the area coverage and the intensity of the phenomenon reached its maximum; -the final stage (3 January 2008, 12:00 UTC -4 January 2008, 06:00 UTC) when the phenomenon gradually extinct.

The precursory stage
On 31 December 2007, in the lower troposphere a high pressure system dominated the most part of Europe.Over the central Mediterranean Sea the pressure had relatively low values (minimum mean sea level pressure values were below 1012.5 hPa; see the green curve south of Italy on Fig. 3a).
In the upper troposphere, a deep trough with a closed low was extended from the Baltic Sea towards the Adriatic Sea, enabling the cold air and vorticity advection towards Western Romania (Figs. 4a and 5a).The wind vertical cross section (Fig. 6) performed along the 12.5 • E, between 36 • and 60 • N (the red line in Fig. 7) indicates the position of the polar front jet around the 300 hPa height.The jet streak (Fig. 7) was situated on the upstream side of the trough (see also Fig. 4a), conditioning its further deepening.

The beginning stage
In the next thirty hours the trough continued to extend towards the east of the Mediterranean Sea and the upper level geopotential low, positioned in the previous stage over the Baltic Countries descended over Serbia (Fig. 4b).In this stage high upper level baroclinicity was noticed in the southwest of Romania, where the relative vorticity increased to 12×10 −5 s −1 (Fig. 5b).The cold air advection on the upstream side of the trough enhanced the low level cyclonic vorticity as well.Over Romania and Bulgaria there was a significant humidity supply: in the lower layers from the Black Sea area (Fig. 3b) and in the upper ones, from the Mediterranean Sea on the downstream side of the trough (Fig. 4b).The snowfalls appeared in the southwestern part of Romania, where the orographic forcing was important.

The maturity stage
On 2 January 2008 the snowfalls and the wind gusts extended towards the eastern part of the country.The cyclogenesis process, initiated in the first part of the night over the Southwestern Black Sea basin, was responsible for the triggering of the blizzard on the Southeastern Romania, including Bucharest.
The analysis of the upper and lower troposphere structure together with the water vapor satellite image patterns led to the recognition of the elements of the cyclogenesis conceptual model within the baroclinic trough (Santurette and Georgiev, 2005).The continuous cold air and vorticity advection of the last two days in the upper troposphere determined the generation of a potential vorticity anomaly (values of 6-8 PVU at 300 hPa level -Figs.4c and 8) over Bulgaria, at the trough base.In the satellite image this anomaly corresponds well to the gray zone south of Romania (Fig. 9) while the black band, surrounding it, indicates the polar front jet position.In the same image the moist/cloud area (the nearly white zone southeast of Romania) has the specific "baroclinic leaf" shape and is located just downstream of the potential vorticity anomaly.In the lower troposphere the distinct baroclinic zone (high gradient of the wet bulb potential temperature over the Southwestern Black Sea -Fig.10) is overrun by the polar front jet, the jet streak being located upstream of this zone.
On the other hand, in the lower layers, the latent heat fluxes from warmer and more humid maritime air mass, in cyclonic twist towards the northwest generate a thermal positive anomaly ahead (Flocas and Karacostas, 1996), deepening the surface vortex (see in Fig. 10 the closed low of 1014 hPa over the Black Sea) and sustaining the altitude one.This process of mutual conditioning continued until the end of the day of 3 January, when the advection of polar air mass from the Russian Plain became dominant, breaking the surface-altitude coupling.
The cold air advected by the northeast air flow (Fig. 10) following the Romanian Carpathians intensified the pressure and temperature gradients on the southeastern part of the country, generating a low level jet.The Fig. 11 shows a vertical cross section of the ALADIN model analysis on 3 January 2008, 00:00 UTC along a line (the red one in Fig. 1) through the northeasterly low-level jet.Accordingly to this section the low level jet developed within a 1-2 km layer close to the ground, where the wind speed reached 16 m/s (the green area at the figure bottom).The intensity of the blizzard was related to the presence of the low-level jet, which is expressed climatologically as a regional wind called crivatz.

Final stage
This stage corresponds to the interval when the upper air geopotential low moved eastwards, this time decoupled from the lower troposphere circulation where the cold air was continuously advected over Romania from the Russian Plain (not shown).

Numerical experiments
The operational LAM ALADIN could represent this weather event only with a six hour anticipation, while the ECMWF global model could indicate it 24 h before, being contradictory to the results of both operational ALADIN and of its coupler (ARPEGE) initialized at the same moment.For this reason few numerical experiments have been conducted to assess: i) The effect of the coupling model for the same limited area model For this purpose the ALADIN model was coupled with ARPEGE and ECMWF models over the same domain (the operational one) at 10 km.Due to the three models' configuration (see Sect. 2), these experiments allow the analysis of the relative role of physics and resolution.
ii.1)The effect of the forecast range (from 1 to 10 days) The RegCM3 model was integrated for 10 days at the same 10 km resolution, being coupled with the ECMWF deterministic model.These experiments allow the large time-scale processes to develop and impact on event generation and thus, to indicate the best range that reproduces more appropriately the process. ii.
2) The effects of using a LAM ensemble prediction system The first 10 members of the ECMWF/EPS were dynamically downscaled by using RegCM3 model over the same domain used in ii.1) but at 50 km resolution (chosen to keep the same computation cost).The experiments were carried out in order to compare two ways of mitigating the sub-grid uncertainty: the direct resolution increase (10 km) and averaging a coarser resolution of a LAM EPS (Palmer, 1999).d), the ALADIN model's resolution was 10 km, like in the operational system, with the physical parametrization set-up close to that of ARPEGE; some differences concern the radiation scheme and the diagnostic/prognostic cloud condensed water.Both the global model ARPEGE (experiment a) and the LAM ALADIN model coupled with ARPEGE (experiment c), having as main difference the horizontal resolution, underestimate the observed precipitation amount (Fig. 13).The ALADIN model additionally generates two local maximums but they are positioned too southerly.The ECMWF model (experiment b) forecasts a more realistic precipitation amount.However, ALADIN coupled with the ECMWF model could not generate, despite of its higher resolution, the local main maximum over central Romanian Southern Plain.This shows that the resolution increase is not sufficient; it only creates/enhances local features when appropriate information transmission function between large-small scale is put to work (as in experiment d) and dump it in the opposite case.It also turns out that this function has to depend on the main physical parametrization hypothesis in order to let the high resolution model extract crucial large-scale information for the process development.

Uncertainty mitigation
The high-resolution solution accuracy was poor, while different coupling model physics showed an improvement of the global solution over the target limited area domain.However this improvement was not transmitted further to the coupled LAM solution.This sub-section is looking for a better solution for uncertainty mitigation both in physical process development representation and in the transmission of coupling information towards the limited area.
Because the non-systematic failures of physical process representation are often connected to non-linear interaction and model feedback' strengths rather than to a given parametrization the first set (ii.1) of experiments aimed to determine a better forecast range for capturing those timescales needed by the physical process to correctly develop.
On the other hand, the previous results (i) showed the role of the initial and lateral boundary conditions (LBC) in reproducing the event and in conditioning the efficiency of model resolution increase.In the second set -LAM/EPS -the unique LBC was replaced by the different solution of ECMWF/EPS.An additional input expected from this experiment would be to see if the approach of ensemble forecast (standing for another way of parametrization of both initial

Sensitivity to forecast range
A first series of experiments has been performed in order to analyze the error growth as a function of forecast range and its possible application to uncertainty estimation for operational and research activities.In all simulations the RegCM3 model was integrated at a 10 km horizontal resolution, over the same domain (see Fig. 2) but started on different days before the interval of interest (i.e.2-3 January 2008, the blizzard maturity stage).The initial state and boundary conditions were provided by the ECMWF deterministic model.The results were compared for different forecast ranges, from 30 to 240 h, all valid for 3 January, 06:00 UTC.
Figure 14 shows the last 24 h cumulated precipitation simulated by the model for 30, 54, 126 and 174 h forecast ranges (corresponding start day: 2 and 1 January 2008, 29 and 27 December 2007, 00:00 UTC), covering the 2 January 2008, 06:00 UTC-3 January 2008, 06:00 UTC interval.As one can see, the 54 h forecast range solution (Fig. 14b) is comparable to that of 30 h (Fig. 14a).With respect to the observations (Fig. 13) it has better dynamics of the precipitation system (moving south-eastwards), with better precipitation pattern, both solutions being able to represent the system phase speed and the good position of the isolated non-precipitating area observed around Bucharest before midnight.Still the 126 h forecast range solution (Fig. 14c) shows a good system position, with the precipitating system center over Western Black Sea (as it was observed, before midnight) while according to the 174 h forecast range solution (Fig. 14d) the system remains completely blocked in the western part of the domain.In the case of 126 h forecast the system developed backwards with respect to the Black Sea, over Mid-Southern Romanian Plain, due to the latent heat supply provided by the sea surface diabatic source.This source enhanced more the eddy energy conversion compared to the 54 h forecast range solution, which, with less history behind, produced a shallower system over the seacoast.
Since the longest RegCM3 forecast was limited to 3 January 2008, 00:00 UTC due to the availability of LBC (ECMWF deterministic model was operationally integrated up to 240 h) for the longer simulations the precipitation was cumulated on a slightly different interval (2 January 2008, 00:00 UTC-3 January 2008, 00:00 UTC) and was compared with the corresponding observation (Fig. 15d).For an easier comparison, Fig. 15 presents only the southern part of the integration domain, where the maximum precipitation amount was registered at 00:00 UTC.One can notice an unexpected improvement of the solution for the 240 h forecast range (Fig. 15c -initial state 24 December 2007, 00:00 UTC) compared to those for 216 and 168 h, respectively (Fig. 15b -initial state 25 December 2007, Fig. 15a -initial state 27 December 2007).
To have more insights into this behavior and to understand if this may contain some predictability information connected to the model-intrinsic response, we analyzed the precipitation prediction skill as a function of forecast range for the past winter: 2006-2007.For that cold season, daily 10-day forecasts were already available.The cumulated precipitation over 24 h intervals were compared, for the nearest grid points (10 km-spaced) to observations (A. Enculescu, personal communication).
The normalized (to the maximum value) bias, averaged for the observation points in the southern half of the domain (Fig. 15d), is displayed in Fig. 16.This error growth function shows a minimum skill for the 192 h forecast range (8 days) and a qualitative improvement for the 240 h (10 days) range, in agreement with our case study results (Fig. 15c in comparison with observations in Fig. 15d).It also indicates a steep slope of predictability loss for ranges between 144 and 192 h (six to eight days), which is reproduced by our case study as well.A similar agreement between the 2006-2007 winter bias and the model simulations for the blizzard case is noticed for shorter ranges.
Accordingly to this result the uncertainty of the numerical solution studied here seems be strongly related to the modeling-system internal variability, and if so, this might be used when estimating the likely probability of a forecast event development.
This result should be further studied in connection with the evolution of the energy spectrum contained in all motionscales represented by the model, as well as in connection with observed-spectral distribution (Peixoto and Oort, 1992), in order to understand the predictability "gaps" in that variability.The result should be analyzed for longer time-series in order to assess its statistical significance in both operational work and research.

Ensemble simulation
Since statistically the predictability of the deterministic simulations (10 km resolution) has a rapid decrease for timescales between six and eight days (in our case the seven days   forecast range -168 h), an ensemble simulation was built by using the REGCM3 model coupled with the first 10 members of the ECMWF/EPS initiated on 27 December 2007 (like in the experiment presented in Fig. 14d).The aim was to analyze the performance of the ensemble solution in filling this longer inertia predictability gap.The resolution of the ensemble simulation was however decreased (by a factor of 5 with respect to the deterministic simulation: x=50 km for RegCM3 and x=125 km for ECMWF) in order to be able to eliminate the models resolution when attributing the result to different contributions.The reason was that, if assuming ensemble technique as a sub-grid parametrization concept (Palmer, 1999), the effect of model and LBC resolution increase would be otherwise over-accounted for.
Comparing the deterministic and ensemble simulations for the same forecast range, thus having the same time-scales enabled for physical processes, the differences between them might primarily be attributed to the method itself.
Figure 17 shows the precipitation field averaged for the 10 members of the ensemble: one can notice, despite the coarser resolution, the precipitation system dynamics improvement (see for comparison in Fig. 14d the deterministic simulation and in Fig. 13 the observed precipitation for comparison).The ensemble solution allows a more realistic diabatic supply from the sea, leading to the system's further development.
This turned out to be a more complete technique to resolve spatial scale uncertainty, at a similar computational cost.This solution may be applied regardless the physical "gaps" in the energy spectra of the motions we looked for, at ii.1).

Conclusions
The blizzard from 1-4 January 2008 was the meteorological phenomenon with the highest degree of risk throughout the 2007-2008 winter in Romania.The damages and the disturbance of social and economic activities justify the analysis of this episode.
In the precursory stage and the beginning one, the motion towards the southeast of the deep trough generated in the middle and upper troposphere represented the initiating factor of the phenomenon.
The relative vorticity advection and the potential vorticity anomaly from the highest layers triggered, in the maximum development phase, the cyclogenesis process in the Black Sea.The ground vortex, generating a thermal anomaly, reactivated at its turn the cyclogenesis in the Black Sea and the anomaly of vorticity from the highest levels.The pressure and temperature gradients from the Southeastern Romania, generated at the contact between the East-European anticyclone and the cyclonic low above the Black Sea, subsequently favored a typical evolution of the blizzard episode, with the formation of low level jet and abundant snowfalls.
The numerical simulations performed in order to assess the main element leading to a better solution indicate that the model physics parametrization appropriateness is, in such mesoscale winter events, more important than the resolution gain.The model resolution might be an element of improvement only after the large-scale information is allowed to pass in an optimal nudged way.To sustain this idea, we intend to further perform a series of experiments on the large-scale information transmission for a few coupling models and to establish its relation with the physical packages and processes involved.The study also emphasized the usefulness of model-specific error growth function, in the prediction of severe weather events.The ensemble prediction system succeeds, for this case, to cover the lack of accuracy in the less predictable time-scales indicated by the error growth statistical function of the deterministic forecast.Again, for a same physical package, the main role of LBC was proved.This result finally emphasizes the need to better understand the physical basis of the underlying mechanism of scale-predictability and to statistically assess error growth behavior for severe mesoscale weather forecast.

Fig. 1 .
Fig. 1.The ALADIN model domain ( x=10 km, Lambert conform projection) and orography (the same color palette as in Fig. 2); the red line indicates the position of the vertical cross section (see Sect. 3.3).

F
. Georgescu et al.: A severe blizzard event in Romania -a case study

Figure 12
Figure 12 shows the 24 h cumulated precipitation forecast by the global models ARPEGE (a) ECMWF (b) and by the ALADIN limited area model coupled with ARPEGE (c) and with the ECMWF model (d), respectively.The global models have almost the same horizontal resolution (about 25 km) over Romania (ARPEGE having a variable mesh), but different physical parametrization packages.For both experiments (c) and (d), the ALADIN model's resolution was 10 km, like in the operational system, with the physical parametrization set-up close to that of ARPEGE; some differences concern the radiation scheme and the diagnostic/prognostic cloud condensed water.Both the global model ARPEGE (experiment a) and the LAM ALADIN model coupled with ARPEGE (experiment c), having as main difference the horizontal resolution, underestimate the observed precipitation amount (Fig. 13).The ALADIN model additionally generates two local maximums but they are positioned too southerly.The ECMWF model (experiment b) forecasts a more realistic precipitation amount.However, ALADIN coupled with the ECMWF model could not

Fig. 11 .
Fig. 11.Aladin model horizontal wind vertical cross section, along the line marked on Fig. 1, for 3 January 2008, 00:00 UTC.On the section, Bucharest is in the middle of the line A-B; wind speed [m/s] -color palette and isolines, wind vector -arrows.

Fig. 16 .
Fig. 16.Error growth function dependency on the forecast range.