J D´ıaz De Argandoña

A Ezcurra

J Sáenz

B Campistron

G Ibarra-Berastegi

F Saïd

Javier Díaz De Argandoña

Atmospheric tides over the Pyrenees: Observational study and numerical simulation

Barometric tides around the Pyrenees mountain range are analyzed by means of synoptic surface-station data recorded during one year, surface data from the Pyrenees Experiment (PYREX) and the CRA/LA Very High Frequency (VHF) wind profiler installed in the north of the range. Tides are decomposed into their diurnal and semi-diurnal components. Diurnal tides show a strong non-migrating component and are very dependent on local conditions. Semi-diurnal tides are more homogeneous and present a north-south asymmetry, also noted in the Alps. This cross-range asymmetry could be related to some interference effect caused by the mountain range in the migrating semi-diurnal tide wave. The asymmetry of the diurnal component presents a very strong seasonal variation, probably related to local diabatic effects. A three-month long simulation has been carried out with the National Center for Atmospheric Research's Weather Research and Forecasting (WRF) limited-area model to try to reproduce the tide structure. The validation of the results with wind-profiler data shows reasonable agreement with the observed diurnal tide and poorer results for the semi-diurnal component. At surface level, however, the model reproduces some of the features of the observed semi-diurnal tide, and especially the cross-range asymmetry.

Introduction

Atmospheric tides (also called atmospheric oscillations) are defined in analogy to oceanic tides as an atmospheric motion of the scale of the Earth, the periods of which are integral fractions of a day and affect pressure and wind. The diurnal, semi-diurnal or terdiurnal variations of pressure are usually masked in the midlatitudes due to the greater variations caused by transient synoptic weather systems [START_REF] Green | The tidal motion of the atmosphere[END_REF].

The main forcing agent for these oscillations, as opposed to ocean tides, is not the solar or lunar gravity pull but the daily variation in insolation and the thermal effects derived from it. The main components of the solar atmospheric tide are the semi-diurnal, with a 12 hour period, and the diurnal or 24 hour period component. Global-scale tides are usually referred to as migrating tides, and are the result of a gravity wave that travels westwards with the apparent motion of the Sun. Nevertheless, a significant part of the tide can be related to local characteristics, and this part is considered as the non-migrating component of the tide [START_REF] Lindzen | Atmospheric tides[END_REF].

Both semi-diurnal and diurnal components have been extensively studied in the last decades. Series of measurements at a great number of stations all around the world have been used to infer global patterns for tides, mainly for migrating tides [START_REF] Haurwitz | The diurnal and semi-diurnal barometric oscillations, global distribution and annual variation[END_REF][START_REF] Hamilton | Topographic effects on the solar semi-diurnal surface tide simulated in a very fine resolution general circulation model[END_REF][START_REF] Dai | Diurnal and semi-diurnal tides in global surface pressure fields[END_REF]. Also, series of data have been compiled for more restricted areas [START_REF] Willson | Atmospheric tidal motions over Australia below 20 kilometers[END_REF][START_REF] Mori | Semidiurnal tidal oscillation at Marcus island[END_REF][START_REF] Mass | Diurnal surface-pressure variations over the continental United States and the influence of sea level reduction[END_REF][START_REF] Bartzokas | Temporal variations of atmospheric tides over Athens, Greece[END_REF][START_REF] Riggin | Jicamarca radar observations of the diurnal and semi-diurnal tide in the troposphere and lower stratosphere[END_REF].

According to the classical theory of atmospheric tides [START_REF] Chapman | Atmospheric tides[END_REF], the semi-diurnal tide is mainly forced by heating in the ozone layer due to UV absorption, but also by absorption of infrared radiation by tropospheric water vapour. It presents little zonal variation, so its main component is the migrating tide.

Diabatic heating arising from the absorption of infrared radiation by tropospheric water vapour and from the ground is the main forcing of the diurnal tide. It is strongly influenced by land masses, topography and, in general, the characteristics of the land surface. As a consequence, the diurnal tide is generally much more irregularly distributed than the semi-diurnal one. It presents a more pronounced annual cycle, and is generally greater over land than over sea [START_REF] Mass | Diurnal surface-pressure variations over the continental United States and the influence of sea level reduction[END_REF][START_REF] Bartzokas | Temporal variations of atmospheric tides over Athens, Greece[END_REF][START_REF] Dai | Diurnal and semi-diurnal tides in global surface pressure fields[END_REF].

The classical theory has been complemented by globalscale numerical models that allow us to overcome some limitations of the classical theory, or by the checking of different forcing patterns [START_REF] Hagan | On modelling migrating solar tides[END_REF][START_REF] Hagan | GSWM-98: Results for migrating solar tides[END_REF][START_REF] Hagan | Migrating and nonmigrating semi-diurnal tides in the upper atmosphere excited by tropospheric latent heat release[END_REF]. These global models do not generally account for small-scale influences on the tides, such as ground topography or small-scale differences in surface sensible heat flux or atmospheric water-vapour content [START_REF] Kato | Atmospheric tides -a review[END_REF].

In a recent work, [START_REF] Hamilton | Topographic effects on the solar semi-diurnal surface tide simulated in a very fine resolution general circulation model[END_REF] carried out a simulation using a general circulation model with a very fine resolution of roughly 10 km in the horizontal. The authors analyze the influence of topography in the South American Andes region and in the Big Island of Hawaii. Their main finding is a decrease in the semi-diurnal tide amplitude from the eastern to the western side of the mountains, attributed to a shadowing effect of the mountain over the westward-travelling tidal wave.

A different topographic influence was described by [START_REF] Frei | Anomaly in the Alpine diurnal pressure signal: observations and theory[END_REF]. They observed an asymmetry in the tide amplitude at both sides of the Alps but along the longitudinal direction. By using a simple theoretical model, these authors showed that the mountain could create an interference pattern around the mountain in the migrating tide wave, causing this effect.

It is worth noting that in the above-mentioned simulation from [START_REF] Hamilton | Topographic effects on the solar semi-diurnal surface tide simulated in a very fine resolution general circulation model[END_REF] a similar north-south asymmetry can also be observed in the semi-diurnal tide pattern around Hawaii. From figure 12 of [START_REF] Hamilton | Topographic effects on the solar semi-diurnal surface tide simulated in a very fine resolution general circulation model[END_REF], the latitudinal tidal amplitude gradient far from Hawaii is about 0.05 hPa per degree, while over the island it is 0.1 hPa per degree.

The availability of extensive reanalysis databases has also attracted the interest of researchers wishing to find out about global-scale tide patterns. The usual archived time resolution of these databases is 6 h, the Nyquist frequency for the semi-diurnal tide. Nevertheless, the interpolation method developed by [START_REF] Van Den Dool | A temporal interpolation method to obtain hourly atmospheric surface pressure tides in reanalysis 1979-1995[END_REF] allows one to extract both diurnal and semi-diurnal tidal signals. [START_REF] Hsu | Tidal fluctuations as seen in ECMWF data[END_REF] found good agreement between experimental data and the semi-diurnal tide deduced from a European Centre for Medium-Range Weather Forecasts (ECMWF) analysis. However [START_REF] Ray | Comparisons of global analyses and station observations of the s 2 barometric tide[END_REF] and [START_REF] Ray | Barometric tides from ECMWF operational analyses[END_REF] found a bias in the phase of the semi-diurnal tide deduced from GEOS-1, the National Centers for Environmental Prediction/National Center for Atmospheric Research (NCEP/NCAR) reanalysis and the ECMWF operational analysis. In the latter the bias is about 20 min. This systematic error has not been explained yet.

These databases also present the opportunity of making use of them to force a limited-area model in order eventually to obtain a more detailed picture of the tide around topographic features, something not covered by classical theory or coarse-resolution simulations. This is a computationally expensive task, because tide calculation often implies a long integration period. However, due to the limited-area domain, it is much more affordable than simulations using global circulation models of comparable resolution.

The objective of the present work is to analyze the atmospheric tides around the Pyrenees mountain range, to find out whether a pattern similar to that found by [START_REF] Frei | Anomaly in the Alpine diurnal pressure signal: observations and theory[END_REF] in the Alps also exists in the Pyrenees and to check whether a limited-area model can reproduce the observed tides.

The next section presents the study area, the data used and the methodology employed for the tide calculation. Section 3 presents the observed tides, their seasonal variation and the characteristics of their amplitudes and phases. Section 4 describes the model simulation performed, the validation of its results against Very High Frequency (VHF) windprofiler data and the tides found from the simulation. The final section summarizes the main findings.

Study area, data and methodology.

The Pyrenees is a simple mountain range of elliptic shape lying on the Spanish-French border with its main axis in the zonal direction. Its length is about 400 km and its width 80 km in its central part, with a maximum height of about 3000 m in the same area. In the south of the Pyrenees, the Iberian range delimits the Ebro Valley, an almost flat river valley with a triangular shape. In the north of the Pyrenees the terrain is generally flat, except for the Massif Central in the northeast of the Pyrenees.

Regular pressure data from several synoptic stations around the mountain have been recorded from the internet server of the the Department of Atmospheric Science of Wyoming University (http://weather.uwyo.edu) every day during 2007. There are no regular meteorological stations inside the mountain range. For this reason we have also used some data from the Pyrenees Experiment (PYREX: [START_REF] Bougeault | Momentum budget over the Pyrénées: The PYREX experiment[END_REF]. The PYREX measurements were collected during 1990 and include surface measurements in locations not available in the synoptic network, inside the mountain range. From the PYREX experiment data, only the specially devised microbarograph transect, designed to compute the mountain drag, has been used in this study [START_REF] Bessemoulin | Mountain pressure drag during PYREX[END_REF]. Additionally, a three-month data record from the VHF wind profiler, installed at Lannemezan in the north-central part of the mountain chain, has also been used to assess the vertical structure of the tide and to validate the model results. In Table I we summarize the data used. Figure 1 shows the geographical domain and the location of the stations used in this work, together with the 1000 m height contour.

We have used the standard method to obtain the diurnal and semi-diurnal components of surface barometric tides. It consists of averaging the hourly pressure values during a given time span, creating a sort of composite day. The difference between each hourly average pressure and the total average gives the mean pressure perturbation at each hour of the day during the period considered. No reduction to sea level has been made. This calculation is frequently made using periods of time such as one month, three months or one year [START_REF] Chapman | Atmospheric tides[END_REF]. We have chosen a three-month period for the synoptic and profiler data and the whole PYREX observational period (two months).

Once we have obtained the diurnal pressure perturbation cycle, the diurnal and semi-diurnal components are worked out by Fourier fitting of the signal. The resulting components are characterized by their amplitude and phase, using p = S 1 sin(t + φ 1 ) + S 2 sin(2t + φ 2 ), Table I. Databases used in this study and their main characteristics.

No. of stations Period Resolution where t is the local mean solar time expressed in radians for the station referred to, S 1 and S 2 are the diurnal and semidiurnal solar tide amplitudes and φ 1 and φ 2 the phases. The residuals from this fitting are in general very small, as diurnal and semi-diurnal harmonics are the prevalent modes. An example is presented in Figure 2, where the fit and spectrum of the different harmonics are plotted. The errors of both amplitudes and phases have been estimated using the method proposed by [START_REF] Bartels | Über die atmosphärischen Gezeiten[END_REF] and outlined by [START_REF] Chapman | Atmospheric tides[END_REF]. This procedure is based on using several data sets to calculate different tidal parameters in the same location. The distribution of these parameters is then used to obtain the tide mean values and their uncertainty. We have obtained different data sets by randomly generating subsets of the available days for each period, e.g. randomly choosing 45 days from the 90 days in each of the three-month periods. From each subset, the tidal parameters have been calculated and the dispersion of these values has been used to estimate the error. We have tried different sizes for the subsets, and we can conclude that if the size of the sample is large enough (about 15 days) then the uncertainty obtained is virtually independent of the sample size.

After characterizing the tide using these experimental data, a three-month long integration has been carried out using NCAR's Weather Research and Forecasting (WRF) limited-area model. A first validation of the model wind results for the whole integration time has been made by comparing them with the wind observed by the profiler. We have used the correlation coefficient, the slope of the fitting line and the normalized root-mean-square (r.m.s.) error.

Next we have checked the ability of the model to reproduce the wind tide. A direct comparison of the local amplitudes and phases of the diurnal and semi-diurnal tides of the two components of the wind at different levels can obviously be made. Instead of using that procedure, we have first compared the cycles of the observed and modelled wind and then the diurnal and semi-diurnal components. This method allows us to obtain a better picture of the vertical diurnal cycle structure as a whole.

Finally, we have compared the observed and modelled semi-diurnal tides for the ground stations. It must be pointed out, however, that the observation and model periods come from different years. Hence we have only compared the semidiurnal tide, because it generally shows little interannual variability [START_REF] Cooper | Inferring solar UV variability from the atmospheric tide[END_REF][START_REF] Bartzokas | Temporal variations of atmospheric tides over Athens, Greece[END_REF]

Experimental results

Seasonal and spatial variation in tides

We have used three-month periods, centred on the reference month (e.g. June values correspond to the May, June and July period), to obtain the seasonal variation of the tides from the 2007 regular synoptic-station data. Figures 3 and4 show the amplitudes and phases of the resulting tides throughout the year. Following the criterion proposed by [START_REF] Cooper | Errors in atmospheric tidal determination from surface pressure observations[END_REF], some station data have been rejected, and Tarbes is not shown for compactness. The mean and standard deviation for the diurnal and semi-diurnal tidal amplitudes and phases for each station along the year are summarized in Table II.

For nearly all the stations, the semi-diurnal component is predominant over the diurnal one. A notable exception is the station of Zaragoza, for which the diurnal component in summer is double the semi-diurnal amplitude. In general the semi-diurnal component is much more uniform both through the year and between stations than the diurnal one.

The error bars in Figures 3 and4 must be interpreted with care when the errors in amplitude are greater than the mean values (e.g. in some diurnal components, especially in Biarritz). Notice that the errors are calculated from the distribution of the tidal vectors (amplitude and phase) and not separately from the amplitude and phase distributions. When the amplitude error is larger than the amplitude mean value, it does not of course mean that the amplitude could reach negative values, but usually reflects a large variability in the phases of the samples. As a consequence, when the error amplitude is greater than its mean value we have an almost complete indeterminacy in the phase. In these cases, no phase error bars have been drawn in Figures 3 and4. In s p i t eo ft h el a r g ee r r o rb a r s ,t h ec a l c u l a t e da m p l i t u d ec a n be significant in some cases, in the sense that the dispersion of the amplitude values will be in general much smaller than that indicated by the error bar, and the mean of the amplitudes much larger than the amplitude of the mean vector. In Figure 5 the samples from three stations have been represented for the June-July-August (JJA) period. It can be observed that the lower amplitude in Biarritz is mainly a consequence of the dispersion in the phase values.

Diurnal tide

The amplitude of the diurnal component shows some seasonal variation for all stations (Figure 3). Related to this fact, the uncertainties are in general quite large, indicating that the diurnal tide is far from being a constant phenomenon but can present a very large variation from month to month. The patterns in Zaragoza and Perpignan are similar, with a maximum in summer and a minimum in winter. This variation is especially intense in Zaragoza, for which the amplitude in summer is just four times the amplitude in winter. In Perpignan the variation is much less pronounced but nevertheless significant.

In contrast Biarritz, a coastal station, shows a completely different cycle, with a nearly zero diurnal tide from May-September and two maxima in April and October. Seasonal cycles in Pau and Tarbes (the latter not shown) are similar to Biarritz's, with larger summer amplitudes. Carcassonne shows a minimum in winter, and the diurnal tide in Barcelona varies from 0.2-0.4 hPa with minima in spring and autumn.

The diurnal tide amplitude is generally greater in continental stations than in some coastal ones (Barcelona and Biarritz), as shown by [START_REF] Dai | Diurnal and semi-diurnal tides in global surface pressure fields[END_REF]. The increase in the thermal inertia caused by the vicinity of the ocean is probably the cause of this weakening of the tide amplitude. The sea-breeze circulation, closely related to this effect, will also attenuate the pressure signal compared with a continental station. Perpignan, however, shows a bigger diurnal tide amplitude in summer than the French continental stations.

As mentioned above, Figure 5 shows that the low values of the summer diurnal amplitudes in Biarritz are related to the dispersion of the phases rather than to a weakening of the pressure variation amplitude. In this figure it can be seen that the amplitudes are comparable to those in Barcelona but the phases are mainly grouped around 0 • and 180 • .The scattering of the phases shows that the diurnal variation of pressure in Biarritz is not a well-organized phenomenon in the sense that although some diurnal oscillation exists it does not follow the same pattern from day to day. This diurnal oscillation cannot be considered as a proper diurnal tide, because it is not a periodical, repetitive oscillation. In contrast, the amplitude of the pressure cycle in Barcelona shows some variation, but the phase is nearly constant (i.e. the points in Figure 5 are grouped around the 0 • phase line). This feature suggests that the diurnal pressure oscillation in Barcelona is mainly driven by the diurnal heating cycle, developing a diurnal tide as in a continental station but not in Biarritz. This could be, at least partially, explained by the different topography of the two sites. The terrain is hilly in the surroundings of Barcelona, and this could contribute to the strengthening of the sensible heat flux diurnal cycle. In Biarritz, however, the terrain is quite flat, so the stable marine layer could more easily prevent the development of a heating-related diurnal pressure cycle. The analysis of the numerical simulation carried out for the summer period (see section 4) shows a well-developed sea breeze in Barcelona and a very weak one in Biarritz, sustaining this hypothesis. The difference in the climatology of the two stations, with Biarritz more exposed to synoptic perturbations from the Atlantic, could also contribute to the difference.

In general both the values and the seasonal pattern of the diurnal tide are very dependent on location. The large diurnal tide amplitude in Zaragoza, the most continental station, together with its great seasonal variability and summer maximum, very likely reflects the influence of the ground sensible heat flux in the non-migrating tide [START_REF] Dai | Diurnal and semi-diurnal tides in global surface pressure fields[END_REF][START_REF] Mass | Diurnal surface-pressure variations over the continental United States and the influence of sea level reduction[END_REF], probably enlarged by the thermal effect of the valley. D ´ıaz [START_REF] De Argandoña | Surface pressure disturbance in the Ebro Valley (Spain) produced by the Pyrenees mountains during PYREX[END_REF] showed that the mean topography-induced pressure perturbation at the south of the Pyrenees reaches maximum values of 7 hPa in Zaragoza, and, as shown above, the diurnal modulation in the pressure in Zaragoza could be as great as 1.4 hPa. The orographic pressure dipole is closely related to the wind circulations around the Pyrenees (Bènech, et al., 1998), so there must be a significant diurnal modulation in the momentum budget and the associated circulations around the Pyrenees.

All these features indicate that the diurnal tide in the surroundings of the Pyrenees is strongly influenced by the local conditions of the stations, so we can conclude that the main part of this tide is of non-migrating nature. The phase shift noted some days in Biarritz in summer with a reversal cycle has to be explained, however.

Semi-diurnal tide

The semi-diurnal component also presents an annual modulation, but the amplitude of this variation is much less important than the amplitude of the diurnal tide (Figures 3 and4). The errors are also much smaller than in the diurnal tides, and no significant dispersion in phase, similar to that in the Biarritz diurnal tide, can be found. The annual pattern is very similar in all the stations, with two maxima in spring and autumn and two minima in late winter and summer (Figure 3). The amplitude also presents a much weaker variability between stations than the diurnal component.

The phase of the semi-diurnal component is more or less constant through the year and from one station to another, as can be seen in Table II and in Figure 4, although again some modulation does exist. The mean value of the phase of the semi-diurnal component for all stations coincides fairly well with the 159 • value proposed by [START_REF] Haurwitz | The diurnal and semi-diurnal barometric oscillations, global distribution and annual variation[END_REF].

The features of the semi-diurnal component, such as its homogeneous values, its almost identical cycle in all the stations and its lower error, all indicate that this component is not much dependent on local conditions. Its origin must therefore be mainly the migrating tide. Nevertheless the coastal stations present a somehow weaker semi-diurnal signal than the continental stations at comparable latitude (e.g. Biarritz versus Pau or Barcelona versus Zaragoza). This difference is about 0.04 hPa in the north and 0.1 hPa in the south.

Asymmetry of the tide around the mountain range

In a study carried out around the Alps, [START_REF] Frei | Anomaly in the Alpine diurnal pressure signal: observations and theory[END_REF] found that the tide amplitudes were greater in the south of the Alps than in the north. The asymmetry was present in both tide components, but the cited authors emphasized the diurnal one. From a simulation by [START_REF] Hamilton | Topographic effects on the solar semi-diurnal surface tide simulated in a very fine resolution general circulation model[END_REF], a similar phenomenon seems also to be present in Hawaii.

This asymmetry can also be found in our data, as seen in Table II. Nevertheless this conclusion is supported only by data from two stations in the south, so it must be made with care. The mean values of the semi-diurnal tide amplitude for stations in the south of the Pyrenees (Barcelona and Zaragoza) are systematically greater than the corresponding values in the north, especially if we compare the coastal stations (Barcelona versus Biarritz) and the continental stations (Zaragoza versus others).

The asymmetry is not so clear in the diurnal amplitude and there is a key difference: the diurnal amplitude asymmetry is heavily seasonally dependent. For example, between Zaragoza and Biarritz the difference is about 0.05 hPa in January and nearly 1.4 hPa in July. This annual variation is not observed in the semi-diurnal amplitude meridional gradient and did not appear in the Alps [START_REF] Frei | Anomaly in the Alpine diurnal pressure signal: observations and theory[END_REF], although in this case all the stations are continental.

In the empirical distributions proposed for barometric tides [START_REF] Haurwitz | The diurnal and semi-diurnal barometric oscillations, global distribution and annual variation[END_REF][START_REF] Hamilton | Topographic effects on the solar semi-diurnal surface tide simulated in a very fine resolution general circulation model[END_REF][START_REF] Dai | Diurnal and semi-diurnal tides in global surface pressure fields[END_REF], the global distribution of the tide amplitude depends on the latitude. This variation can only account for 0.02 hPa and 0.04 hPa for the diurnal and semi-diurnal tide amplitudes in our case. This is clearly not enough to explain the asymmetry observed (nearly 0.2 hPa between Zaragoza and continental north stations, and 0.1 hPa between the two coastal stations for semi-diurnal tides).

No significant shadowing effect can be found between Perpignan and Biarritz in the semi-diurnal tide, but this is probably due to the geometry of the Pyrenees, lying in the zonal direction.

It has not been possible to obtain reliable diurnal tides from the PYREX data, because of the great dispersion of the tidal parameters and their very large errors. Related to this fact, and in contrast to the data from the synoptic stations, the tide spectra show that some of the higher order harmonics are larger than the diurnal one. We can conclude that the two-month time span is probably not long enough to obtain a significant diurnal signal from the data. A check performed in the three-month data recorded from the wind CRA/LA VHF profiler, also used in this study, results in a similar conclusion.

Figure 6 shows the semi-diurnal tidal amplitudes and phases along the PYREX microbarograph transect, obtained from the whole observational period (October and November 1990). The amplitudes are clearly greater in the south than in the north of the mountain range and this supports the presence of an asymmetry along the range.

The amplitude of the semi-diurnal tide presents a marked reduction with height, of about 0.35 hPa over 2000 m. A similar reduction of the semi-diurnal tide amplitude with height has been documented in several studies (e.g. [START_REF] Chen | Diurnal variation of pressureheights: A vertical phase shift[END_REF][START_REF] Hamilton | Topographic effects on the solar semi-diurnal surface tide simulated in a very fine resolution general circulation model[END_REF] and it is compatible with the classical theory of the semi-diurnal tide.

The tide asymmetry found implies that both the drag and the topographic pressure dipole must present a daily variation. Sure enough, both the drag calculated during the PYREX experiment and the pressure dipole intensity during the same period show tidal cycles. The tidal results for these two parameters are shown in Table III. Again we find that the diurnal and semi-diurnal modulation in drag, of the order of 0.2 Pa, is significant, as the mean absolute [START_REF] Frei | Anomaly in the Alpine diurnal pressure signal: observations and theory[END_REF] to explain the origin of the asymmetry in the diurnal tide over the Alps. The first group, comprising processes such as differential solar heating by albedo variations or slope exposition, is related to the spatial variation of the local thermal forcing. Although this group could produce some effect on the diurnal tides, it is improbable that the semi-diurnal tide is significantly affected by this kind of local effect [START_REF] Chapman | Atmospheric tides[END_REF].

The dynamic effect of the mountain on the incident air flow is very dependent on the conditions of the incident air mass. Some of these conditions, in particular the vertical stratification, could probably present a systematic day-night variation and thus produce a variation in the diurnal cycle of the pressure difference across the ridge. Hence the sign of this variation would depend on the direction of the incident air stream, and so on the sign of the drag. We have compared the tidal amplitudes and phases with the mean value of the drag over a number of samples and found no correlation.

A simple theoretical model used by [START_REF] Frei | Anomaly in the Alpine diurnal pressure signal: observations and theory[END_REF] predicts the existence of a mesoscale trapped Kelvin wave that interferes destructively (constructively) with the tidal migrating wave to the north (south) of the mountain range. The result of this interference will be the observed asymmetry in the tidal amplitude.

In the case of the Pyrenees, the seasonal variation of the diurnal tide is very irregularly distributed between the stations (compare Biarritz with Zaragoza, for example) and this irregularity leads to a strong seasonal variation in the asymmetry, hard to explain with only an interaction of a planet-scale wave with the topography. Hence, we do not think this explanation alone could be applied to our domain. The seasonal variation of the asymmetry is also compatible with a difference in the ground thermal balance. In our case the surface heat flux is added to the interference mechanism by the above-cited authors.

On the other hand, the annual cycle of the semidiurnal tide shows a similar pattern for all the stations, which strongly suggests a common origin for all these perturbations. From this point of view, an explanation similar to that proposed by [START_REF] Frei | Anomaly in the Alpine diurnal pressure signal: observations and theory[END_REF] for the diurnal tide seems admissible for the spatial structure of the semi-diurnal tide observed in the present study.

Numerical simulation

In the previous section, we have shown that the barometric tides around the Pyrenees are far from homogeneous. It is very likely that the topography and the soil type play an important role in the tide distribution, by interfering with the migrating tidal wave, by altering the thermal budget near the ground or by additional effects such as the influence of the topography on the distribution of vertically integrated water content (e.g. [START_REF] Smith | The influence of mountains on the atmosphere[END_REF], which will impact diabatic heating by infrared radiation. In order to depict the main features of the tidal distribution, we have used NCAR's WRF limited-area model, forced by ECMWF's ERA40 reanalysis data. We have chosen the period from June-August 2002, as this is the period of the available CRA/LA VHF windprofiler data, so we can use it to validate the model results. Additionally, it corresponds to the season of the year that shows the most intense spatial variation in the diurnal tide.

Model setup

Integration by means of the non-hydrostatic WRF version 2.2.1 mesoscale model has been arranged as follows. There exist two domains: the larger one (D1) has a horizontal resolution of 75 km, whilst the smaller one (D2) has a grid size of 15 km. The results of D2 have been analyzed in subsequent sections of this article, after discarding a five-point strip at the borders. The nesting between both domains is two-way. The boundary conditions are updated every 6 h from the ERA40 reanalysis data at a horizontal 1.125 • grid resolution. The WRF model has 35 vertical η levels reaching 20 hPa in the top layer. Vertical interpolation is linear in logarithm of pressure and the time step in the largest domain is 300 s. The microphysics parametrization used is the WRF Single Moment 3-class simple ice scheme [START_REF] Hong | A revised approach to ice microphysical processes for the bulk parameterization of clouds and precipitation[END_REF]. Radiation is computed according to the CAM scheme, both for long-and short-wave radiation, updating radiation computations every 10 min. The landsurface model consists of a thermally diffusive five-layer soil with a surface-layer scheme given by MM5 Monin-Obukhov similarity theory and a YSU planetary boundary-layer (PBL) parametrization scheme. The cumulus parametrization scheme used is the Betts-Miller-Janjic one [START_REF] Skamarock | A description of the advance research WRF Version 2[END_REF]. Some of the tidal signals already present in the reanalysis will be forced to the model by the boundary conditions. As the six-hour frequency is just the Nyquist frequency for the semi-diurnal tide, the intensity of this forcing will depend on the exact location of the external domain borders. In our case, only about 10% of the reanalysis semi-diurnal signal will be forced to the model, so virtually all the semi-diurnal tides present in the results will be produced by the WRF model itself.

Validation of the model results

Model wind results from the simulation are compared with the observed wind measured at three different levels in Figure 7. The model seems to capture well the main tendency of the observed wind, especially in the intermediate 6000 m level. To assess the accuracy of the model representation, we have calculated the correlation coefficient, the slope of the fitting line and the normalized r.m.s. error between the observed and modelled wind velocities at 37 vertical levels. The results for the different levels and wind components are shown in Figure 8. The best general agreement is observed between 4000 and 14 000 m, simultaneously indicated by the three parameters considered. In this region the r.m.s. error, normalized against the mean value of each component intensity, takes a value of 0.4 for the zonal component and 0.7 for the meridional component. The linear regression slope, of the order of 0.9 in both cases, shows that the model generally underestimates the observed wind.

In the lower levels the inaccuracy of the model output could be related to the strongly smoothed-down topography used by the model or a deficient representation of the PBL processes. The poor agreement above 14 000 m can be explained by the progressively weaker signal-to-noise ratio of the profiler above 12 000 m, so the measurement above this level could be of insufficient quality.

We can conclude that the model simulation reproduces the profile of the observed wind at Lannemezan with reasonable accuracy, in a layer between 4000 and 14 000 m. Below this layer the results are progressively worse. Aloft the results are uncertain, as the measured wind is probably not of high quality.

Model and observed wind tides

Both diurnal and semi-diurnal tidal signals are clearly present in the observed wind profile and the model results.

Using a technique similar to that used in the analysis of data from the surface, the hourly diurnal variation of the wind has been obtained for both the model and the VHF profiler data. Also, the diurnal variation due to either diurnal or semi-diurnal components and the sum of the two has been computed. Finally, as a quality check of the tidal model results, the difference between all these values has been obtained. The results are shown in Figure 9 (we will refer to the panels in this figure as (i, j), meaning the ith row and jth column).

The main feature of the observed diurnal cycle is a strong modulation in the meridional component of the wind, with maxima at 4 h (southwards) and 14 h (northwards). There is also a southwestward anomaly between 5000 and 10 000 m at about 21 h. The first two features are qualitatively well represented in the model cycle (Figure 9(2,1)), although the first appears earlier (at 2 h) and the second shows a discontinuity below 6000 m, not observed by the profiler. The last feature is also reproduced, but at a lower level (4000 m). These differences can also be observed in Figure 9(3,1), from which we can conclude that the diurnal cycle is reasonably well represented above 5000-6000 m but not below that level.

Comparing the diurnal component of the tide (Figure 9, column 2) we again find generally good agreement above 5000 m, with poor results below this level.

The semi-diurnal component is weaker than the diurnal in both observed and modelled data. The patterns are different in both, and the modelled tides are generally much weaker than the observed ones. In particular, the strong modulation in the wind direction at 4000-8000 m is not present in the simulation. The strongest amplitudes in the simulated tide appear in two layers, one below 5000 m and the other at 9000 m. The lower layer presents a variation of the meridional component, and the other of the zonal component.

Taking into account the comparisons carried out, the wind tides obtained from the model are far from perfect, especially for the semi-diurnal component. Nevertheless, some features of the tides are reasonably reproduced, especially considering that the diurnal component is predominant during the summer. As the tidal signal in the ground pressure can be considered as a consequence of the vertically integrated atmospheric variation, we can expect a reasonable representation of the barometric tides at ground level by the model.

Surface-pressure tides

The ground-pressure field obtained from the numerical simulation clearly presents diurnal and semi-diurnal tides, again with a negligible residual. In Figure 10, the amplitude and phase of both tides are represented for grid model points using a vector convention: the amplitude is represented by the vector's length and the phase by the direction of the vector, from east counterclockwise.

In the south of the Pyrenees, the semi-diurnal tide amplitude presents a very clear north-south gradient with minimum tides in the north of the mountain range. The position of this minimum, displaced to the north of the mountain axis, would in fact contribute to the tide asymmetry observed from our measuring stations. Note that the stations in the north are closer to the mountain range that those in the south.

The diurnal amplitude shows a more complex pattern and greater spatial variability. In the south there is a clear maximum very close to Zaragoza, which extends to the south through the Ebro valley. A symmetrical maximum also exists in the north, but away from our observation domain. The amplitude is clearly modulated by the terrain elevation, with a strong decrease in the mountains, for which the amplitude is about ten times smaller than in the valley. There is a clear phase shift over the mountains, similar to but smaller than that observed by [START_REF] Chen | Diurnal variation of pressureheights: A vertical phase shift[END_REF] over Taiwan.

Out of our measurement domain, both tides also seem to be affected by the relief of the Massif Central, in the northeastern corner of the domain represented in Figure 10.

Although the simulation period and the surface experimental data come from different years, the semidiurnal tide presents little interannual variation, so that both results can be compared. The comparison is made between synoptic stations and the closest grid point, both tides calculated for the summer period (JJA) and for the semi-diurnal tide. The result is shown in Figure 11. The model underestimates all semi-diurnal amplitudes, and this discrepancy is greater for continental stations (Carcassonne, Pau and Zaragoza). This can be explained by the relatively low upper model limit of 20 hPa, which will miss part of the ozone layer [START_REF] Zwiers | Simulation of solar tides in the Canadian Climate Centre general circulation model[END_REF]). As noted above, the north-south asymmetry is reproduced by the model. Both the simulated and experimental phases of the semi-diurnal tides are notably constant from one station to another, but there is a phase difference between experimental and model tides of about one hour.

In spite of the serious limitations pointed out, the results from the numerical model can capture some key features of the tidal distribution around the Pyrenees. Further research, with possible longer integration times and changes in the model parametrizations or the maximum height achieved by the simulation, is needed to refine these results and to obtain some clues about the causes of the distribution of the tide.

Conclusions

Atmospheric tides around Pyrenees range have been analyzed using a one-year record of synoptic stations, two months worth of data from the PYREX database and a three-month long record from the Lannemezan VHF profiler. Clear semi-diurnal and diurnal signals have been obtained from all the data. Both geographical and seasonal variations are observed in the distribution of barometric tides around the mountains. There is a particularly strong variation in the diurnal component and a much weaker one in the semi-diurnal tide. Generally, the semidiurnal component is predominant, but a stronger diurnal amplitude, twice as great as the semi-diurnal tide amplitude, has been found in Zaragoza, the most continental station, in summer. A clear asymmetry is also found, especially in the semi-diurnal component of the tide, between the north and the south of the mountains. The semi-diurnal asymmetry is almost constant through the year, but the diurnal one presents a very intense seasonal cycle. This asymmetry in the semi-diurnal component has been found in PYREX data.

The features observed in the tide cannot be explained by classical tide theory, which predicts only migrating, zonally constant planetary-scale tides and does not take into account any effects related to the terrain nature or topography. Several effects, summarized by [START_REF] Frei | Anomaly in the Alpine diurnal pressure signal: observations and theory[END_REF], can be invoked to explain this variability in the tide. According to the above-cited authors, interference of the migrating wave tide with the topography is the main mechanism active in the Alps. From our observations, this can hardly be the only explanation in the case of the Pyrenees, because it does not explain the great seasonal variability observed. It could probably be applied, however, to explain the semi-diurnal tide asymmetry.

A three-month simulation using NCAR's WRF limitedarea model has been carried out to check the ability of the model to reproduce the tides, and in particular the spatial pattern and the asymmetry observed. The model results have been validated using the Lannemezan VHF wind-profiler record. The model reproduction of the wind is generally good above 4000 m. The model is also able to reproduce the main features in the diurnal cycle of the wind, although the results are far from being quantitatively correct. The semi-diurnal barometric tide observed at the surface is in general well reproduced, and a north-south asymmetry is also obtained, although not as strong as in the observational data. There is a systematic bias in the phase. The diurnal tide is harder to compare, because the simulation and observational periods do not coincide and the temporal variability of this component makes any comparison uncertain. Nevertheless, a dipolar pattern is obtained, which leads to results qualitatively similar to those for the observed ground tides.

Figure 1 .

 1 Figure1. Geographical domain and location of the stations used in this study. Regular synoptic stations are marked with a square, the microbarograph network from PYREX with circles. The CRA/LA VHF wind-profiler situation is indicated with a star. The shadowed regions correspond to terrain higher than 1000 m.

Figure 2 .

 2 Figure 2. (a) An example of the fit obtained using the leading two harmonics (24 and 12 h) against the observed diurnal cycle. (b) Power spectrum of the harmonic analysis, showing the prevalence of the first two components. This example corresponds to the February-April period in Zaragoza.

Figure 3 .

 3 Figure 3. Seasonal variation of the diurnal and semi-diurnal tidal amplitude for synoptic stations. The error bars represent the probable error of the amplitudes, see text for details.

Figure 4 .

 4 Figure 4. Seasonal variation of the diurnal and semi-diurnal tidal phase for synoptic stations. The error bars are only drawn if the tidal amplitude is greater than its uncertainty.

Figure 5 .

 5 Figure5. Distribution of the diurnal tides obtained from 100 (Zaragoza and Barcelona) and 200 (Biarritz) random samples of 45 day size for the JJA period. The tide amplitude calculated for each sample is represented by the distance from each point to the origin and the phase by the angle from the horizontal counterclockwise. Notice that the great variability in the phase of the diurnal tide in Biarritz will lead to a very weak mean amplitude. Observe also the different patterns in the two coastal stations, Barcelona and Biarritz. The units of the label values are hPa.

Figure 6 .

 6 Figure 6. (a) Amplitude and (b) phase for semi-diurnal tides along the PYREX transect. The terrain height has also been represented.

Figure 7 .

 7 Figure 7. Comparison of the measured wind component with the modelled wind at three different levels for the Lannemezan station. The continuous line represents the measured wind and the dotted line the model output.

Figure 8 .

 8 Figure 8. Comparison between observed and modelled wind components at different levels for the Lannemezan station. The curves show the slope of the linear fit, the correlation coefficient and the normalized r.m.s. error.

Figure 9 .

 9 Figure 9. Perturbation of the wind velocity as a function of time and height. The first row represents measured data, the second row modelled data and the third row the difference between them. The first column represents the raw perturbation of the wind, the second column the diurnal component, the third column the semi-diurnal component and the last column the sum of both diurnal and semi-diurnal components. The vectors represent the direction of the horizontal wind perturbation (north upwards, west leftwards).

Figure 10 .Figure 11 .

 1011 Figure10. Amplitude and phase for the diurnal and semi-diurnal components of the pressure tide for the JJA period for the grid points of the simulation carried out. The vector's length is proportional to the amplitude and the phase is represented as the angle from the zonal eastwards direction, counterclockwise. Note that the amplitude scale is different between plots. The locations of synoptic stations are represented by squares.

Table II .

 II Values of the mean and standard deviation (in brackets) of seasonal amplitude and phase distributions of the tides for each of the synoptic stations.

		Diurnal		Semi-diurnal
	Site	S 1 (hPa)	φ 1 (rad)	S 2 (hPa)	φ 2 (rad)
	Biarritz	0.18 (0.1)	1.5 (1.3)	0.54(0.08)	2.8(0.2)
	Tarbes	0.19 (0.1)	1.2 (0.8)	0.55(0.07)	2.7(0.08)
	Pau	0.30 (0.1)	0.8 (0.5)	0.58(0.05)	2.9(0.2)
	Carcassonne	0.44 (0.1)	0.7 (0.2)	0.57(0.08)	2.9(0.2)
	Perpignan	0.39 (0.1)	0.6 (0.2)	0.58(0.04)	2.8(0.2)
	Barcelona	0.30 (0.1)	0.4 (0.4)	0.64(0.04)	2.8(0.1)
	Zaragoza	0.86 (0.4)	0.3 (0.1)	0.75(0.05)	2.7(0.08)

Table III .

 III Amplitude and phase of the diurnal and semidiurnal tides of the drag measured during PYREX and the pressure difference across the Pyrenees.

		Diurnal	Semi-diurnal
		S 1 (Pa)	φ 1 (rad)	S 2 (Pa)	φ 2 (rad)
	Drag	0.1380	-2.358	0.1914	-0.508
	p	65.08	-2.931	24.18	-0.925
	value measured during PYREX was 2.2 Pa with a standard
	deviation of 1.7 Pa. Unfortunately, the data from the PYREX
	experiment do not allow us to extend the analysis to the
	seasonal cycle of a whole year.
	Several physical processes have been proposed by

Acknowledgements

The authors are grateful for financial support by the ETORTEK Strategic Research Programme (Department of Industry, Trade and Tourism and Department of Transport and Civil Works of the Basque Government, Basque Meteorological Service-Euskalmet) through the EKLIMA21 project (ETORTEK IE08-217 and IE09-264). Funding was also received from the National R+D+i Plan, Spanish Ministry of Science and Innovation (CGL2008-03321/CLI).

The authors thank the ECMWF for granting access to reanalysis and operational data through the MARS archive system by means of the special project SPESIPRA. Comments by two anonymous reviewers have been very helpful and have led to a better final version of the paper.


F Youcef Ettoumi

H Sauvageot

A.-E.-H Adane

Statistical bivariate modelling of wind using first-order Markov chain and Weibull distribution

Keywords: Wind modelling, Stochastic Processes, Markov chain, Weibull distribution, Long-term probability

This paper studies the statistical features of the wind at Oran (Algeria). The data used are the wind speed and wind direction measurements collected every 3 h at the meteorological station of Es Senia (Oran), during the 1982/92 period. The eight directions of the compass card have been considered to build the frequency distribution of the wind speed for each month of the year and each direction. The three-hourly wind data have been modelled by means of Markov chains. First-order nine-state Markov chains are found to fit well the wind direction data, whereas the related wind speed data are well fitted by first-order three-state Markov chains. The Weibull probability distribution function has also been considered and found to fit the monthly frequency distributions of wind speed measurements. Two methods of wind data retrieval are thus made available. In fact, two models of chronological bi-series are obtained describing wind speed and wind direction.

1

Introduction

Wind contributes significantly to the thermal exchange processes between earth surface and atmosphere. It favours water evaporation at ground surface and carries various particles in suspension in the air. It also represents an appreciable source of energy. To account for the effect of wind erosion at the ground, knowledge of wind speed and direction is important. The wind can be considered as a short-time scale, random, non-stationary process. It is thus of great interest to assess its climato-statistical features. These characteristics can be studied in space and time by using the wind speed and direction measurements collected near the ground in various meteorological stations.

To describe the statistical features of the wind speed at a given location, the usual method is to fit the experimental wind data with probability distributions [START_REF] Weibull | A statistical distribution function of wide applicability[END_REF][START_REF] Khogali | Wind, solar energy potential in Makkah (Saudi Arabia). Comparison with Red Sea coastal sites[END_REF][START_REF] Hasan | Wind energy in West Bank and Gaza Strip[END_REF][START_REF] Rehman | Statistical characteristics of wind in Saudi Arabia[END_REF][START_REF] Morgan | Statistical distributions of wind parameters at Sydney, Australia[END_REF][START_REF] Garcia | Fitting wind speed distributions: a case study[END_REF].I n previous studies [START_REF] Ettoumi | Modelling and simulation of hourly wind in Algeria[END_REF][START_REF] Ettoumi | Optimum classification of wind energy potential in Algeria[END_REF], wind speed and direction data have been mapped for Algeria on a yearly basis. For a given location, the Weibull probability distribution function was shown to fit the annual frequency distribution of the three-hourly wind speed. The Weibull distributions obtained for the various meteorological stations of Algeria were classified by considering either the Ward method or the K-means algorithm [START_REF] Ward | Hierarchical grouping to optimize an objective functions[END_REF][START_REF] Kalkstein | An evaluation of three clustering procedures for use in synoptic climatological classification[END_REF][START_REF] Dixon | BMDP Statistical software[END_REF]. In more than 50% of the meteorological stations, mean wind speed was found to be lower than 3 m/s. For about 30% of the stations, the corresponding wind speed values ranged between 3 and 5 m/s. In the rest of the stations, wind speed was higher than 5 m/s. Thanks to this classification, the wind energy potential in Algeria can be estimated with a good accuracy.

Another way of estimating the wind statistical features consists in using autoregressive processes or Markov chains [START_REF] Castini | Stochastic modelling of wind velocities time series[END_REF]. This kind of modelling has the advantage of accounting for the time variations of wind speed and direction. Since the latter can easily be expressed in term of discrete states, it can be described efficiently by Markov chains or autoregressive models.

The goal of this paper is to present the results of the modelling, on a monthly basis, of the three-hourly wind speed and direction data recorded at the site of Es Senia (Oran) during the 1982/92 period using both Markov chains and Weibull probability distributions.

Data

The station of Es Senia is situated in Algeria at 0°37' W, 35°38' N and 90 m of altitude. At this station, the wind is measured by an anemometer and weathercock label, placed at 10 m above the ground. Wind direction angles are measured with respect to the North with an accuracy of 0.1°.

Wind speed and direction were averaged over three-hour periods. These measurements are available from 1982-1992. The resulting sequences of three-hourly wind data consist of fast random fluctuations superimposed on trend variations. The random fluctuations and trend variations are caused by atmospheric turbulence and geostrophical circulation respectively. Since the trend varies slowly, a month is a suf-ficiently short period to represent the three-hour wind variations by a stationary random variable [START_REF] Nfaoui | Stochastic simulation of hourly average wind speed sequences in Tangers (Morocco)[END_REF]. Hence, the obtained database has been divided into 12 subsets, corresponding to the 12 months of the year. A given monthly subset brings together all the 1982/92 three-hour wind data collected during the same months of the different years.

In meteorological applications dealing with wind measurements, the circle described by the weathercock in the horizontal plane is usually divided into eight sectors having the same spread of 45°. Every wind direction is assigned to one of the eight sectors by computing the sector corresponding to the observed direction.

Markov modelling

Physical aspects

The wind can blow in the same direction for a relatively long time. Then it suddenly changes from one direction to another, continues to blow in the new direction and so on. The resulting sequence of wind direction data is thus made of time series of stable states and abrupt changes between them. These states depend on the wind strength and can take all the weathercock directions. Wind speed and direction are both erratic and persistent. This means that wind variations in time can be described by first-order Markov chains [START_REF] Daniel | Stochastic simulation and forecasting of hourly average wind speed sequences in Jamaica[END_REF].

According to the division into eight sectors presented in section 2, wind direction measurements, expressed in degree, were classified following nine different states. These states were labelled by numbers 0-8 corresponding to no wind, north (N), north-east (NE), east (E), south-east (SE), south (S), south-west (SW), west (W) and north-west wind (NW) respectively. The nine states thus defined are those to be used for modelling wind direction data by first-order, nine-state Markov chains.

When analysing the experimental data, three states of wind speed can be distinguished, namely weak, mean or strong. The numbers 0, 1 and 2 labels these states respectively. They can be defined by considering the relative frequency distributions of wind speed obtained in a given site, for each month of the year and each of the eight sectors. Three typical states can be identified from the histograms of threehourly wind speed, which respectively correspond to: 0 Յ V (i) Ͻ 3 m / s (weak wind), 3 m / s Յ V (i) Ͻ 8 m / s (mean wind) and V (i) Ն 8 m / s (strong wind). Once these states are known, wind speed data can be modelled by first-order, three-state Markov chains.

First-order Markov model

A series of wind data collected every three hours, is considered, the number of intervals of 3 h being k at time t. X k is a random variable that may take one of the values 0, 1,…,n , …, N-1 at time t [START_REF] Ross | Stochastic processes[END_REF][START_REF] Stewart | Introduction to the numerical analysis of Markov chains[END_REF]. These values related to the N possible states of the phenomenon under study, are respectively denoted ⌽ 0 ,…, ⌽ n , …, ⌽ N -1 . Since the probability of observing ⌽ n state after 3 k h, is P (X k = n), a probability vector can be defined. It is written as:

Q ϭ [P(X k ϭ 0). …P(X k ϭ n).…P(X k ϭ NϪ1)],
(1)

with n P(X k ϭ n) ϭ 1. (2)
For the first-order Markov process, the probability of having a transition between the states occurring for two successive observations, is:

P ij ϭ P(X k+1 ϭ j / X k ϭ i), (3)
with i or j = 0, …, n, .... , N-1.

Considering all the possible transitions between the ⌽ n states, the Markov transition matrix:

M ϭ [P ij ], (4)
is obtained with

j P ij ϭ 1. (5)
If n ij is the number of transitions from ⌽ i to ⌽ j observed in the sequence of wind data, the maximum likelihood estimate of the transition probabilities is:

p ˆij ϭ n ij j n ij . ( 6
)
If Q k and Q k +1 , are the probability vectors observed after 3 k and 3 (k + 1) h respectively, the transition from Q k to Q k +1 is expressed as:

Q k+1 ϭ Q k M. ( 7
)
According to the ergodic properties of Markov chains, the probability vector Q, obtained after a long period of time, becomes independent from the initial state, and relation ( 7) can be written:

QM ϭ Q. ( 8
)
The solutions to this equation are the long-term Markov probabilities P (X k = n) which form the probability vector (3). These probabilities can also be estimated by calculating the corresponding a priori probabilities P n = n n / N M , where n n is the number of ⌽ n states observed over a month, and N M the whole number of observations in the month. To appraise the persistence of the phenomenon under study, the correlation between two three-hourly observations has to be calculated.

Thus, the autocorrelation coefficient corresponding to a lag of 3 q h, is given by:

C q ϭ NϪ1 i ϭ 0 NϪ1 j ϭ 0 ijp q ij P(X j ϭ i) Ϫ ͫ NϪ1 i ϭ 0 iP(X j ϭ i) ͬ 2 NϪ1 i ϭ 0 i 2 P(X j ϭ i) Ϫ ͫ NϪ1 i ϭ 0 iP (X j ϭ i) ͬ 2 (9)
where P q ij is the element (ij) of the matrix M q . The discrete model described above, is equivalent to a coding system which produces messages by combining N different symbols. The relative entropy of this system is, by definition:

H ϭ 1 log 2 (1 / N) N Ϫ 1 i ϭ 0 P i log 2 P i ( 10
)
The loss of information caused by the classification of wind data into N stable states, is then given by:

L ϭ 1ϪH. ( 11
)

Wind direction modelling

The wind direction data collected at Oran during the 1982/92 period, have been modelled by first-order, nine-state Markov chains, following two different approaches. The first one is based on the use of homogeneous Markov chains to describe the three-hourly wind direction fluctuations for each month of the year. This means that all the Markov parameters defined above are constant with time over the month.

To illustrate the resulting model, the Markov transition probabilities of wind direction data obtained in January, April, July and October are plotted against the various directions in Fig. 1. As shown in Fig. 2, the long-term Markov probability of having one of the nine states of the wind direction, calculated for each of the 12 months, is found to be almost equal to the corresponding a priori probability. More generally, homogeneous Markov models are seen to yield a good representation of the seasonal variations.

The second approach consists in setting up Markov models with time-varying parameters. The use of inhomogeneous Markov processes leads to improved models, which account for the daily variations. To obtain this kind of model, each of the 12 monthly subsets of three-hourly wind direction data has been divided into four samples, rather than considering only one subset per month. For each month of the year, these samples denoted S 1 ,S 2 ,S 3 and S 4, are made of all the measurements of the 1982/92 period performed in the intervals running from 0-6h ,6 -12 h, 12-18 h, and 18-24 h respectively. Then, each of the four samples has been modelled by first-order, nine-state Markov chains. The corresponding transition matrices are denoted M 1 ,M 2 ,M 3 and M 4 .

Fig. 3 gives, as an example, the results of the calculation of the Markov transition probabilities characterising the four January samples of wind direction data. The Fig. 2. Representation of long-term and a priori probability for hourly wind direction data. P 1 is the a priori probability corresponding to the sector 0. P 2 is the a priori probability corresponding to the sector 1. P(X = 1) is the long-term probability corresponding to the sector 0. P(X = 2) is the long-term probability corresponding to the sector 1.

product of these matrices has been performed, and the resulting matrix (M 1 M 2 M 3 M 4 ) has been compared with the transition matrix (M) obtained from the approach based on homogeneous Markov chains. For each of the 12 months, these matrices have been tested for homogeneity using the c 2 distribution with eight degrees of freedom. The level of confidence has been fixed at 1%. In this case, all the P ij values of these matrices are found to be falling within a 0.99 confidence range. This proves that the (M 1 M 2 M 3 M 4 ) and (M) matrices are identical.

The autocorrelation coefficient has been calculated using [START_REF] Ward | Hierarchical grouping to optimize an objective functions[END_REF]. The values of C q estimated for each month of the year and three different lags, namely 3, 6 and 9 h, are plotted in Fig. 4. These plots show that the autocorrelation coefficient decreases quickly when the lag (q) increases. When computing C q for various lags, the autocor- relation function has been found to be approximately equal to C 1 , where C 1 is the one-lag autocorrelation coefficient. This result means that the change in wind direction strongly depends on the immediate past. Hence, it justifies the limitation of the Markov chains used in this study to the first order.

The relative entropy was also computed for each of the 12 months. The results show that the loss of information due to the use of homogeneous first-order ninestate Markov model, only represents a few percentage points. The yearly variations of entropy for the inhomogeneous first-order nine-state Markov models are illustrated in Fig. 5. The mean monthly value of information loss is found not to exceed 11%. Hence, it may be concluded that such a model describes the wind direction variations for Oran with a reasonably good accuracy.

Wind speed modelling

First-order three-state Markov chains have been used to model the wind speed data collected every 3 h at Oran during the 1982/92 period. It has to be reminded that the three states are "weak wind" (state ⌽ 1 ), "mean wind" (state ⌽ 2 ) and "strong wind" (state ⌽ 3 ).

To illustrate the calculation of the Markov transition probabilities for wind speed, the plots of P 11 ,P 22 , and P 33 variations over the year are presented in Fig. 6. These diagrams show that the wind is weak at Oran almost all the time, and that strong winds are rare. This means that the wind energy potential is significant, but not very important in this region.

For each month of the year, the long-term Markov probability of having one of the three states of wind speed has been compared with the corresponding a priori probability. The results are plotted in Fig. 7 where it can be seen that both types of probability are identical for each of the three states. This result proves that the firstorder three-state Markov chains correctly fit the three-hourly wind speed variations.

The autocorrelation coefficient has been calculated for q = 1, 2 and 3 (i.e., lag equalling 3, 6, and 9 h respectively). The diagrams representing the variations of C 1 ,C 2 , and C 3 over the year are given in Fig. 8. The mean value of C 1 is close to 0.5, showing the strong correlation between two successive states. C 2 varies from 0.17-0.33, and C 3 is smaller than C 1 and C 2 . However, it sometimes reaches 0.22. More generally, the autocorrelation function is found to vary slowly. This means that the wind state observed at hour (h) is mainly influenced by the wind speed measured three hours before, i.e. h -3, but the measurements performed at the preceding hours, h -6, h -9,…, also bring non negligible information.

Such a result is confirmed by the values of relative entropy obtained over the year. Fig. 9 indicates that the yearly variations of the relative entropy are weak. Entropy is minimum in January, when it equals 75%, and maximum in May, when it reaches 89%. Hence, the resulting loss of information varies from 21-25%.

In spite of the equality of the a priori and long-term Markov probabilities mentioned above, the strong values of information loss imply that the first-order Markov model cannot account accurately for the wind speed fluctuations. To improve this kind of modelling, second and third-order transitions must be considered, since they are as important as the first-order ones.

Weibull distribution-based model

The implementation of second or third-order Markov models requires tedious computation.Then, instead of such approaches, we have considered the Weibull distribution to model the 3 h wind speed measurements performed at Oran [START_REF] Hennesey | Some aspects of wind power statistics[END_REF][START_REF] Tackle | Note on the use of Weibull statistics to characterise wind-speed data[END_REF][START_REF] Justus | Methods for estimating wind speed frequency distributions[END_REF][START_REF] Darwish | Wind energy potential in Irak[END_REF][START_REF] Løvseth | The renewable energy potential of Norway and strategies for its development[END_REF]. The Weibull probability density function is given by: where, v is the wind speed in m s -1 , k is a dimensionless shape factor, and A is a scale factor with the dimension of speed. The Weibull cumulative distribution function is:

f(v) ϭ k A ͩ v A ͪ kϪ1 exp ͩ Ϫ ͩ v A ͪ k ͪ (12)
F(V) ϭ ͵ V 0 f(v)dv ϭ 1 Ϫ exp ͩ Ϫ ͩ V A ͪ k ͪ ( 13
)
where F (V) is the probability for the wind speed to be equal or lower than V. This function can be transformed into a linear relation of the form y = ax + b by calculating:

y ϭ Ln{ϪLn[1 Ϫ F(V)]} ϭ kLnVϪkLnA. ( 14
)
The double logarithmic transformation has been applied to each of the monthly cumulative frequency distributions of the three-hourly wind speed data. The experimental data thus transformed have been smoothed using a least square method. The a and b parameters of the regression equation have been computed. Taking into account that a = k and b = -kLnA, the shape and scale factors of the corresponding Weibull distribution are obtained.

The plots of Figs 10 and 11 give the k and A values obtained at Oran for each of the eight directions and each of the 12 months. Fig. 10 shows that wind direction from west, north and south-west are dominant, with peak values in February, July, and June. For all sectors taken together, the peak is in April.

Fig. 11 also shows that the shape factor is always higher than 1. This means that the shape of the distribution is not exponential. However, Fig. 10 shows that, for the monthly curves, k is sometimes lower than 1.

For all the regression equations obtained in this study, the correlation coefficient is higher than 0.9 (Figs 12 and13). For the south, south-west and west directions, the correlation coefficient is higher than 0.97. For all the sectors together, it is higher than 0.96.

For a given month, each of the eight frequency distributions of three-hourly wind speed measurements is smoothed by a Weibull probability distribution function. All the results obtained with this kind of modelling, are in good agreement with those arising from experimentation. The results of the implementation of the usual hypothesis tests confirm that the Weibull probability distribution very well describes the monthly frequency distributions of hourly wind speed data.

Conclusion

This study confirms that wind is a Markov process, in which the main contribution arises from the first-order transitions. Two bivariate models are obtained. The first one is based on the first-order nine-state Markov chains for wind direction, and firstorder three-state Markov chains for wind speed. The second one is based on the first-order nine-state Markov chains for wind direction, and Weibull distribution, for wind speed.

The memory effect that characterises the Markov processes is well mirrored by the autocorrelation obtained for the two used variables.

The bivariate model consisting in the combination of only two Markov chains, correctly describes the three-hourly wind data. It generates two time series corresponding to two independent variables. However, this model leads to non negligible loss of information.

The three-hourly wind data are better described with the Weibull probability distributions associated to Markov chains. The tight link observed between experimental and retrieved series (95%) confirms this result. Moreover, this kind of bivariate model improves the usual Weibull-based wind modelling.

Using the Markov model of wind direction obtained in section 3.3, the Weibull model of wind speed generates joint times series. The addition of the wind direction model is useful since it enables to retrieve bivariate time series of the wind. This can be used for the reconstruction, simulation, and prediction of the wind.

As mentioned in section 3.4, the wind energy potential in the region of Oran is moderate, but sufficient to be used as an alternative energy source. However, the results of section 4 show that for some periods of the year and for some directions, the wind energy can be important.

Recently, both bivariate models have also been applied to hourly wind measurements collected in the meteorological station of Chasseral (Switzerland) from 1990-2000. This station is situated in the Jura mountains at 1600 m of altitude. The measurements were performed at 40 m above ground level. The obtained results are The wind data generated by the improved bivariate model could be useful to size the wind power generators to be installed in any location. These data are also usable for the evaluation of the performance of such generators.

Fig. 1 .

 1 Fig. 1. Transition probability P ii for the 3 h wind direction data collected at Oran and modelled by homogeneous first-order nine-state Markov chains.

Fig. 3 .

 3 Fig. 3. Transition probability P ii for 3 h wind direction data collected at Oran in January and modelled by inhomogeneous first-order nine-state Markov chains.

Fig. 4 .

 4 Fig. 4. Yearly variations of autocorrelation coefficients C 1 ,C 2 and C 3 for the wind direction, C 1 ,C 2 and C 3 are the correlation coefficients when two data are separated by three, six or nine h, respectively.

Fig. 5 .

 5 Fig. 5. Yearly variations of relative entropy for 3 h wind direction data.

Fig. 6 .

 6 Fig. 6. P 11 ,P 22 and P 33 transition probabilities for the 3 h monthly wind speed data collected in Oran and modelled by first-order three-state Markov chains.

Fig. 7 .

 7 Fig.7. A priori and Long-term probabilities for the 3 h monthly wind speed, noted P i and P(X = i) respectively with i = 1, 2 or 3. i represents the sector.

Fig. 8 .

 8 Fig. 8. Yearly variations of autocorrelation coefficients for the 3 h monthly wind speed. C 1 ,C 2 and C 3 are the correlation coefficients, when two data are separated by three, six or nine hours respectively.

Fig. 9 .

 9 Fig. 9. Yearly variations of relative entropy for the 3 h monthly wind speed.

Fig. 10 .

 10 Fig. 10. Shape (k) and scale (A) factors of Weibull distribution estimated at Oran for each of the eight directions and for months running from January to December.

Fig. 11 .

 11 Fig. 11. Shape (k) and scale (A) factors of monthly Weibull distribution estimated at Oran for all sectors taken together.

Fig. 12 .

 12 Fig. 12. Correlation coefficient between the experimental wind speed distribution and the corresponding values obtained from fitting the Weibull distribution for each sector.

Fig. 13 .

 13 Fig. 13. Correlation coefficient between the experimental wind speed distribution and the corresponding values obtained from fitting the Weibull distribution for all sectors together.

Acknowledgements

The authors wish to thank the National Office of Meteorology of Algeria for providing the data used in this study.