Wind and Turbulence Statistics in the Urban Boundary Layer over a Mountain–Valley System in Granada, Spain

Abstract

Urban boundary layer characterization is currently a challenging and relevant issue, because of its role in weather and air quality modelling and forecast. In many cities, the effect of complex topography at local scale makes this modelling even more complicated. This is the case of mid-latitude urban areas located in typical basin topographies, which usually present low winds and high turbulence within the atmospheric boundary layer (ABL). This study focuses on the analysis of the first ever measurements of wind with high temporal and vertical resolution throughout the ABL over a medium-sized city surrounded by mountains in southern Spain. These measurements have been gathered with a scanning Doppler lidar system and analyzed using the Halo lidar toolbox processing chain developed at the Finnish Meteorological Institute. We have used the horizontal wind product and the ABL turbulence classification product to carry out a statistical study using a two-year database. The data availability in terms of maximum analyzed altitudes for statistically significant results was limited to around 1000–1500 m above ground level (a.g.l.) due to the decreasing signal intensity with height that also depends on aerosol load. We have analyzed the differences and similarities in the diurnal evolution of the horizontal wind profiles for different seasons and their modelling with Weibull and von Mises probability distributions, finding a general trend of mean daytime wind from the NW with mean speeds around 3–4 m/s at low altitudes and 6–10 m/s at higher altitudes, and weaker mean nocturnal wind from the SE with similar height dependence. The highest speeds were observed during spring, and the lowest during winter. Finally, we studied the turbulent sources at the ABL with temporal (for each hour of the day) and height resolution. The results show a clear convective activity during daytime at altitudes increasing with time, and a significant wind-shear-driven turbulence during night-time.

1. Introduction

In recent decades, there have been increasing concerns about climate change and the air quality, and how natural and anthropogenic processes affect them. According to a special report by the Intergovernmental Panel on Climate Change [1], human activities are estimated to have caused approximately 1.0 °C of global warming above pre-industrial levels, and this warming is likely to reach 1.5 °C between 2030 and 2052 if it continues to increase at the current rate. Moreover, wind field is another climate impact-driver affected by climate change, with an observed weakening of mean wind speed close to the surface which is expected to continue in the coming decades, together with an increase in severe wind storms and tropical cyclone peak winds. However, there are several atmospheric components and mechanisms that are not sufficiently understood, measured or whose effect still presents high uncertainties according to the Fifth and Sixth IPCC Assessment Report [2,3] and some related studies [4].

This is the case of the behavior of the lowermost region of the atmosphere, which has a crucial role in modelling and understanding climate and air quality. This layer, known as the atmospheric boundary layer (ABL), is the place where the emission of pollutants occurs, is directly responsible for the dispersion processes, and its correct modelling is essential for numerical weather prediction and climate models [5,6,7]. It is characterized by a turbulent behavior with significant temporal and spatial variations, which makes the accurate measurement and modelling of its internal mixing and its interactions with the surface and the rest of the atmosphere a challenging task.

Wind field is one of the key variables in understanding the complex processes in the ABL. Therefore, accurate measurements of wind speed and direction profiles are required by high-resolution numerical weather prediction models [8]. Moreover, those wind profiles measured over a particular site provide valuable information about local transport or to validate homogeneity assumptions for studies comparing close sites [9,10,11,12].

On the other hand, it is important to characterize and understand turbulence mechanisms and sources due to their complex interactions with other meteorological variables. Turbulent mixing is responsible for the redistribution of momentum, mass, temperature and humidity within the ABL [13]. The sources of turbulent mixing include buoyancy (that produces upwards convective mixing), wind shear (mechanical mixing) or radiative cooling in stratocumulus clouds producing top-down convective mixing [14]. Turbulence also has a role in new aerosol particle formation [15] and cloud microphysics [16].

An important application of social interest is the study and prediction of urban pollution events. In particular, urban areas with complex topography, such as basins, usually present pollution episodes due to the combined effect with weather patterns [17]. Bossioli et al. [18] showed the importance of correct parameterization of the ABL in a complex basin terrain for different seasons, with the simulation of pollution events under different conditions. These singular conditions are also the responsibility of special ABL features in Granada, a medium-sized city in southern Spain [19,20,21,22], causing serious air pollution events [23,24] especially during winter, with long periods with low wind, dry air and stable conditions. Therefore, thorough characterization of the ABL in this kind of locations is crucial for understanding the air quality.

In this context, lidar techniques represent a powerful tool to retrieve profiles of several ABL properties. Doppler lidars, those measuring the Doppler shift due to the aerosol particles motion by the wind, have been developed and improved in recent decades [25,26,27]. The strong potential of the technique prompted the European Space Agency (ESA) to design and launch the ADM-Aeolus satellite mission in 2018, with the first ever spaceborne Doppler wind lidar system contributing to improving the accuracy of numerical weather prediction. However, ground-based Doppler lidar systems are of great importance to validate this spaceborne measurements [28,29,30] and also to measure the 3D wind field inside the ABL and to retrieve turbulent properties with high temporal and vertical resolution, with particular applications in many fields [31].

Retrieval of the horizontal wind field has been widely used for wind energy, aviation and meteorology because of the possibility of studying effects such as wind shear [32,33], low-level jets [34,35,36] or wind gusts [37,38]. The high temporal resolution for the wind field has allowed many studies on ABL turbulence [39,40,41,42,43], often in combination with other remote sensing instruments [20,21,44,45]. Lane et al. [46] demonstrated the suitability of Doppler lidar to wind field profiles in urban environments, by comparing them with point measurements from anemometers. From that study, the importance of Doppler lidar in urban studies has been shown several times [47,48,49,50,51].

The different products retrieved from Doppler lidar measurements can also be combined to classify turbulence based on its source [52]. Manninen et al. (2018) developed a classifying method that represents a powerful tool to study and deeply analyze certain ABL scenarios, or to create long-term databases of vertically resolved ABL classification. They showed the performance of the ABL classification product when applied to two particular stations, namely, Jülich (Germany) and Hyytiälä (Finland). They presented an example of a clear-sky day in Jülich, where the convective ABL was properly developed, and a cloud-topped example in Hyytiälä, where the clouds avoided any possible convection. A statistical study was also performed for those stations, finding clear differences in the main sources and their seasonal and diurnal cycles due to the different locations.

This type of systematic algorithms is crucial for the development of networks and infrastructures such as Cloudnet (https://cloudnet.fmi.fi/, accessed on 14 December 2021) [53], which is part of the European ACTRIS (Aerosol, Clouds and Trace Gases Research Infrastructure, https://www.actris.eu/, accessed on 14 December 2021). There is currently high interest in creating and validating standardized and coordinated tools and datasets in order to exploit state-of-the art technology in an efficient and homogeneous way, as is the objective of initiatives recently funded by the European Cooperation in Science and Technology (COST), TOPROF [54] and PROBE [55].

In this study, we analyzed the first ever Doppler lidar measurements of wind field with high temporal and vertical resolution throughout the ABL over a southern European ACTRIS-Cloudnet station at Granada, Spain. The objective of this study was characterization of the ABL wind field and turbulent structures over this urban site and investigating their particular features. To this end, a statistical study was carried out to characterize the ABL in terms of mean values and probability distributions of horizontal wind and diurnal frequency of observed turbulence sources. We used a software toolbox developed by [56] to systematically derive different ABL properties for a 2-year database measured with Doppler lidar deployed at Granada.

The experimental site and the Doppler lidar system used for the study are presented in Section 2. In Section 3, the software toolbox to systematically derive different wind vector fields and ABL turbulent properties is described, as well as the data availability after the measurement pre-processing. The results of the wind and turbulence analysis are shown in Section 4, and the summary and conclusions are presented in Section 5.

2. Experimental Site and Instrumentation

This work was developed at the experimental observatory AGORA (Andalusian Global Observatory of the Atmosphere). In particular, the measurements for this study were performed at the UGR station, at the Andalusian Institute for Earth System Research (IISTA-CEAMA) in Granada (37.16°N, 3.61°W, 680 m above sea level, a.s.l.). This station is part of several instrumental and research networks as well as of ACTRIS.

Granada is a medium-sized city in south-eastern Spain located in a natural basin delimited by the Sierra Nevada mountain range, which reaches more than 3000 m a.s.l (Figure 1). This city is climatically characterized by large seasonal temperature differences, with winter average minima of 2 °C and summer average maxima of 33 °C (according to Spanish Meteorological Agency, AEMET, for the period 1981–2010). It is a dry area, in terms of low relative humidity (year average of 57%) and scarce precipitation (year average around 350 mm).

Mean surface winds are also light (less than 2 m s⁻¹ on average with more than 50% of calms), coming predominantly from west and northwest according to historical records [57,58]. An important factor characterizing the diurnal wind regime in Granada is the local mountain–valley thermal flow, which forms katabatic winds from Sierra Nevada mainly in the early night [59].

The vertical profiles of wind and of ABL turbulent properties in this study were obtained using the measurements of the Doppler lidar Stream Line (Halo Photonics), that is part of ACTRIS-Cloudnet [53].

Table 1. Technical details of the main optical elements of the Doppler lidar system.

Emission
Wavelength	1500 nm
Pulse energy	100 µJ
Pulse duration	200 ns
Pulse repetition rate	15 kHz
Reception optics
Telescope	monostatic optic-fiber coupled
Physical lens diameter	75 mm
Effective beam diameter *	(17.5 ± 1) mm
Lens divergence	33 µrad
Focal length *	(535 ± 35) m
Detection
Detection type	Heterodyne
Range resolution	30 m
Points per range bin	10
Sampling frequency	50 MHz
Velocity resolution	0.0382 m s−1
Nyquist velocity	20 m s−1

* Retrieved with the method of Pentikäinen et al. (2020).

gives a description of the instrument main features. The system consists of a solid-state pulsed laser emitting at 1.5 µm and a heterodyne detector using fiber-optic technology. The emission is done with low pulse energy (100 µJ) and high pulse repetition rate (15 kHz), what makes the instrument eye safe. The signal acquisition is performed continuous and autonomously in vertical stare mode with a temporal resolution around 2 s, and it also has full hemispheric scanning capability. For the regular measurements, conical scans with constant elevation of 75° and 12 equidistant azimuth points have been performed every 10 min.

For this study, the instrument was configured for operation with 30 m vertical resolution and an effective range from 90 m to 6000–9000 m. The focus of the optical system was experimentally found at (535 ± 35) m, using the method described in [60], increasing the instrument sensitivity at this height, but reducing it above 2 km. With the same method, the effective beam diameter was also determined to be (17.5 ± 1.0) mm. Note that the physical lens diameter in Table 1 refers to the physical size of the lens, whereas the calculated effective (1/e2) beam diameter describes how much is it illuminated by the laser beam. A more detailed description of the instrument can be found in [27].

3. Methodology

3.1. Halo Lidar Toolbox

The vertical profiles of wind field and turbulent properties within the ABL have been retrieved from vertical and scanning Doppler lidar measurements using several linked algorithms. The Halo lidar toolbox, which was developed at the Finnish Meteorological Institute, was used as the software processing chain. The aim of this toolbox is to produce harmonized Doppler lidar retrievals applied to measurements from different sites using robust methods presented in peer-reviewed articles [41,52,60,61,62,63,64,65,66,67,68] and providing consistent uncertainty estimates. Such harmonization on calculated quantities and used methods is essential to create a broad catalogue of ABL datasets in order to address the gap between the understanding of ABL physics and its representation in high-resolution climate and numerical weather prediction models [69].

In this study, only the horizontal wind product and the ABL classification product have been systematically analyzed. Nevertheless, the rest of the quantities calculated in the chain are required for the ABL classification; thus, they are briefly described below:

• Background correction has been applied to the raw data with methods published by [62,63]. A correction for the telescope focus of the instrument was also applied to the signal, as explained in [60]. The instrumental precision of radial velocities was estimated with the method given by [27,70], and attenuated backscatter coefficient (${\beta }_{att}$) with uncertainties were also calculated [64];
• From scanning measurements, three-dimensional wind vector profiles were calculated using the least squares method with singular value decomposition by assuming a stationary and horizontally homogeneous wind field, and radial velocity uncertainties were propagated to wind components as described by [65]. The wind component uncertainties were estimated with the method described by [66];
• The vertical velocity statistical momenta, i.e., variance, skewness, and kurtosis, were calculated from vertically pointing measurements at 3, 30, and 60 min resolutions. Statistics unbiased by random noise and sample size were calculated as given by [67] and standard errors were estimated with a bootstrap method described by [68];
• Wind shear vector, which can also be a source of turbulent mixing, is also calculated. This vector is calculated from the changes in $u$ and $v$ wind components with height [61], as:

  $$\overrightarrow{sh}=\frac{\Delta \overrightarrow{U_H}}{\Delta z}=\left(\frac{\Delta u}{\Delta z},\frac{\Delta v}{\Delta z}\right),$$

  (1)

  therefore, its module is:

  $$sh=\frac{\sqrt{\Delta u^2+\Delta v^2}}{\Delta z},$$

  (2)
• The dissipation rate of the turbulent kinetic energy, $\epsilon$, was calculated from vertically pointing measurements using the method presented by [41]. This quantity was defined as the rate at which the turbulence energy is absorbed by breaking the eddies down into smaller eddies until they are ultimately converted into heat by viscous forces [71], following the Kolmogorov hypothesis [72]. The method applies Taylor’s frozen turbulence hypothesis that eddies travel with the mean wind while maintaining their characteristics [73]. This quantity is then used as an indicator of turbulent mixing, instead of the combination of vertical skewness and variance [14]. The method used also provides an uncertainty estimate for $\epsilon$ [41];
• Finally, all the previously calculated quantities were combined following a decision tree to create a bitfield-based classification mask. This method was created by [64] following the profile-based Doppler lidar method introduced by [52] with the aim of objectively assigning a dominant source for turbulent mixing. The analyzed regions of the profiles were selected from calibrated ${\beta }_{att}$ to height ranges with sufficient atmospheric signal and no clouds (a threshold ${\beta }_{att}>$10 Mm⁻¹ sr⁻¹ was used based on the literature [14,52,74]). The presence of turbulence was obtained from $\epsilon$ with a threshold $\epsilon$ > 10⁻⁵ m² s⁻³ or $\epsilon$ > 10⁻⁴ m² s⁻³, depending whether the classified heights were below cloud or connected to the surface [64]. All range gates with surface-connected turbulent behavior during daytime were classified as dominated by convective mixing. During night-time, when ABL is assumed to be neutral or stably stratified [71], wind-shear derived turbulence is searched with a threshold $sh$ > 0.03 s⁻¹ [75]. Finally, range gates that are classified as turbulent but are unconnected to surface or clouds during daytime, and not related to wind shear during night-time, are labelled as ‘intermittent’ since turbulence is assumed to arise from other intermittent sources [76].

3.2. Data Processing

A 2-year database was gathered following the measurement procedure described in Section 2 with the Doppler lidar system from 3 May 2016 to 2 May 2018. To investigate seasonal changes in diurnal wind and turbulence cycle, the database was divided into meteorological seasons: winter (December to February, DJF), spring (March to May, MAM), summer (June to August, JJA) and autumn (September to November, SON). In order to ensure statistical representativeness of the study, an additional criterion on data availability was used. Wind vector retrieval was not available for certain range gates where the quality of the retrieval was not sufficient, where precipitation or fog was detected or if there were no measurements. With this in mind, we calculated the fractions of the data from the total analyzed period that were available for each hour of the day and range gate, and we selected those with more than 60% availability to perform the statistical analysis. The results of this calculation are shown in Figure 2 for each season. It can be observed that the highest altitudes with this criterion (red lines) were reached during summer and the lowest during winter. This fact was expected (if no technical issues are taken into account), because in Granada the ABL height, and consequently, the height with enough aerosol load to provide enough signal intensity, is higher in summer than in winter [19,20,22,77].

For each season, we averaged wind profiles within 1 h time intervals, but no range averaging (i.e., with the original range resolution of 30 m). The time interval of 1 h was selected following the averaging time used in most meteorological and air quality models [78]. With this approach, we aimed for a mean diurnal evolution of the horizontal wind speed and direction at different altitudes. The wind averaging was carefully tackled because of its vector nature. Depending on the application, the wind measurements may be vector- or scalar-averaged [79,80]. Scalar-averaging of wind direction may lead to incorrect results, due to the circular nature of this quantity. Meanwhile, vector-averaging may provide incorrect wind speeds when wind direction variance is large. Therefore, the wind averages we took in this study were scalar-averaged wind speed ($\overline{U_H}$, hereinafter simply referred to as $U_H$) and vector-averaged wind direction (${\overline{d}}_{RV}={\tan }^{-1}\left(\frac{\overline{u}}{\overline{v}}\right)+k$, hereinafter referred to as $d$, where $k$ is a correction to keep $d$ between 0° and 360° and depends on the signs of $\overline{u}$ and $\overline{v}$). However, they must not be interpreted as the polar components of any mean wind vector but have to be discussed as separated quantities.

Wind speeds were then handled as linear data, and means and standard deviations were calculated because they are the most used estimators in the literature, modelling, validation, etc. However, they are strictly only good estimators for normally distributed measurements. Indeed, wind speeds are experimentally observed to follow a Weibull distribution at any location [81]. This distribution is commonly used in the wind power industry to characterize wind regimes or to predict the production of a wind turbine, and was also chosen in this study to model the observed wind speeds at different heights. The probability density function described by the Weibull distribution with respect to wind speed depends on a scale parameter (β), proportional to mean wind speed, and a shape parameter (α) such that smaller values are given by the wider spread of winds around the mean wind speed.

$$f\left(U_H\right)=\frac{\alpha }{\beta }{\left(\frac{U_H}{\beta }\right)}^{\alpha -1}{e}^{-{\left(U_H/\beta \right)}^{\alpha }}$$

(3)

Wind direction data cannot be analyzed used linear statistics, but circular statistics was applied starting with calculating only circular mean as central estimator. Circular variance can be computed as a dispersion estimator [82], but it is not a robust estimator. Firstly, because it is highly affected by possible outliers, which was avoided in this study by filtering out measurements with wind speed less than 1 m/s, to avoid high error propagation when calculating ${\tan }^{-1}\left(\frac{\overline{u}}{\overline{v}}\right)$. Secondly, because in a strongly multi-modal circular distribution, a variance value indicating high dispersion could be obtained for a distribution with well-defined modes (with low variance each) in opposite directions, and in that case the circular mean might point in a direction of low sample frequency. Therefore, wind direction data were fitted to monomodal or bimodal von Mises distributions, also called Circular Normal, a function of the mean direction μ and the dispersion parameter κ that measures the concentration of unimodal circular data around the mean. The density function for such distribution is generally defined for a mixture of n modes as

$$f\left(d\right)={\displaystyle \sum }_{i=1}^n{p}_i\frac{e^{{\kappa }_i\cos \left(d-{\mu }_i\right)}}{2\pi I_0\left({\kappa }_i\right)}$$

(4)

where ${\mu }_i$ stands for the preferred direction of each mode, ${\kappa }_i$ represents the concentration parameter that indicates how directions cluster around ${\mu }_i$, and $p_i$ is the proportion of each component in the mixture. $I_0\left({\kappa }_i\right)$ is the Bessel function of the first kind and order 0. Wind direction data were fitted using Vector_Stats, a software tool calculating descriptive statistics and inference on circular variables [83].

In order to derive a statistical view of the complete wind vector, wind rose analysis of the database was also performed, without any seasonal division. To this end, we divided the hourly averaged profiles into three evenly spaced height ranges: 100–340 m a.g.l., 340–580 m a.g.l., and 580–820 m a.g.l. These ranges were selected so that there were at least 60% data available for the whole dataset. We also distinguished between two time intervals: ‘Daytime’, from 07:30 to 17:00 h UTC, and ‘Night-time’, from 21:30 to 03:00 h UTC. The reason for using these intervals is that they ensure that we are including only hours when it is daytime or night-time during the whole year, excluding twilight (defined when the sun has an elevation angle over −18° [84]) and transition times.

Finally, for the analysis of turbulent sources, the frequency of each source of mixing from the ABL-classification product was calculated for each time of the day (with 3 min resolution, as it was the maximum available by the software), for each season and range gate. This kind of analysis allows for characterizing the diurnal cycle of the ABL in terms of prevailing mixing sources.

4. Results

4.1. Horizontal Wind Field Characterization

We applied the seasonal analysis to the whole 2-year database, from 3 May 2016 to 2 May 2018. The results showed some seasonal differences during daytime, with calm winter and more windy spring and summer for all altitudes. During night-time, the winds were low (less than 1 m/s in average) for all seasons.

Figure 3a,d,f,i,k,n,p,s shows the time evolution of the mean hourly wind speed and direction profiles, from 100 m a.g.l. to the maximum available altitude gate (according to Figure 2). Two time intervals have been distinguished by directly identifying strong direction and/or speed gradients in the height–time average plots. The selected intervals are included in

Table 2. Selected hours for ‘Day’ and ‘Night’ time intervals for the averaged profiles at each season.

	DJF	MAM	JJA	SON
‘Day’ interval (h UTC)	10–18	10–22	10–21	9–20
‘Night’ interval (h UTC)	20–9	23–7	1–8	22–7

and depicted with orange (night) and blue (day) boxes in Figure 3a,d,f,i,k,n,p,s. The calculated mean profiles within such intervals are shown in Figure 3b,e,g,j,l,o,q,t), and the averaged profiles of wind speed standard deviations within all available measurements, $\Delta U_H$, are shown in Figure 3c,h,m,r. Wind direction standard deviations have not been computed because of their lower statistical meaning, as described in Section 3.2.

From wind speed plots in Figure 3a–c,f–h,k–m,p–r, we can observe that average $U_H$ increased with height during night-time for all seasons (from 2–3 m/s up to 5–6 m/s in winter, summer, and autumn and up to 9 m/s in spring). Meanwhile, daytime average $U_H$ started from 3–4 m/s at the lowest analyzed heights, then presented a slight increase around 200 m a.g.l. (especially in spring and summer) and remained constant before a final increase at high altitudes, reaching 6 m/s in winter, summer, and autumn and 10 m/s in spring. Wind speeds at the highest analyzed altitudes did not present diurnal differences for any season. The altitudes where this constant pattern was reached were around 440 m a.g.l, 940 m a.g.l., 1230 m a.g.l., and 1080 m a.g.l. for winter, spring, summer, and autumn, respectively.

Standard deviations also increased, with height, from around 1.5 m/s up to 5 m/s in all seasons with the exception of spring, when $\Delta U_H$ reached 7 m/s. This means that, during the analyzed period, wind speeds presented more diverse values at higher altitudes, especially in spring. However, no significant differences were observed between daytime and night-time standard deviations.

Direction plots in Figure 3d,e,i,j,n,o,s,t also show marked diurnal and seasonal differences, as well as contrast between the first 800 m and altitudes above. In winter, the mean daytime wind came from the NW, close to the surface, and from N at higher altitudes, whereas the mean nocturnal wind came from the SE, close to surface, and from the E at higher altitudes. In spring, most of the time there were winds coming from the W and NW, with the exception of the nocturnal wind close to the surface, which came from the S. Summer and autumn presented similar direction patterns, with mean daytime winds coming from the W and SW, and nocturnal winds coming from the E below 500 m a.g.l. and from the SW above this altitude.

After this analysis of mean and standard deviation profiles, the probability distributions described in Section 3.2 were fitted to the measured wind data. Wind speeds at each height were fitted to Weibull functions, obtaining the vertical profiles of scale ($\beta$) and shape ($\alpha$) parameters plotted in Figure 4. The distributions were calculated for two intervals (‘Daytime’, 07:30–17:00 h UTC and ‘Night-time’, 21:30–03:00 h UTC), and also using the complete database (‘All data’). It is observed how $\beta$, which is proportional to mean speed, is increasing with height, with higher daytime than night-time values. The profiles for $\alpha$, however, show a decrease with height in two different regimes: a faster decrease below a certain height (around 250 m a.g.l.) and then an almost constant behavior for higher altitudes. Moreover, higher daytime $\alpha$ values are observed at lower altitudes, whereas the behavior at upper altitudes is opposite. The study was also performed with a seasonal discrimination (see Appendix A). For all seasons, the scale parameter increased with height and was higher during daytime than during night-time up to a certain height, when there is no variability, or it is inverted. Moreover, this height increased in warmer seasons. Concerning the shape parameter, a decrease with height is observed for all seasons, with two differentiated height ranges in terms of decreasing rate.

As the next part of the study, Figure 5 includes all the wind roses from 1 h averaged wind data at three different height ranges and two common time intervals (‘Daytime’, 07:30–17:00 h UTC and ‘Night-time’, 21:30–03:00 h UTC), using the 2-year database. They were calculated with 22.5° angle intervals and the wind speed intervals: 0–0.5 m/s, 0.5–1 m/s, 1–1.5 m/s, 1.5–2 m/s, 2–5 m/s, 5–10 m/s, and ≥10 m/s. Wind directions at each height and time interval were fitted to von Mises distributions, finding that bimodal mixtures were able to model all wind roses with sufficient accuracy. The obtained distributions are plotted in Figure 5 as solid black lines, and the corresponding parameters are included in

Table 3. Parameters corresponding to bimodal von Mises distributions fitted to daytime and night-time wind directions divided into three height ranges.

	Daytime			Night-Time
Altitude a.g.l. (m)	580–820 (Figure 5a)	340–580 (Figure 5c)	100–340 (Figure 5e)	580–820 (Figure 5b)	340–580 (Figure 5d)	100–340 (Figure 5f)
${\mu }_1 (°)$	289.7	293.4	299.1	61.6	85.0	106.4
${\kappa }_1$	1.0	1.4	2.0	2.1	6.3	5.0
$p_1$	0.95	0.92	0.87	0.48	0.25	0.14
${\mu }_2 (°)$	213.8	205.8	205.3	233.3	0.0	91.5
${\kappa }_2$	15.0	10.0	5.0	0.9	0.03	0.5
$p_2$	0.05	0.08	0.13	0.52	0.75	0.86

.

At the lowest height interval (Figure 5e,f), there is a clear prevalence of NW daytime winds (main mode at 299.1°) with speeds mostly >2 m/s and a strong change to ESE (main mode at 106.4°) and weaker nocturnal winds. Analyzing those wind directions in the context of the local orography (see Figure 1), it is noticeable that this pattern is consistent with the katabatic winds expected at that height range due to slope effects [59], although a deeper analysis would be needed to detect such winds.

This pattern was not exactly preserved at higher altitudes (Figure 5a–d), although there are similarities. The daytime wind was also dominated by the NW direction, although the frequency distributions spread and there was a higher contribution of W and SW winds. The peak of the nocturnal winds moves towards the E and NE at higher altitudes (main modes at 85.0° and 61.6°) and a slight SW contribution (mode at 233.3°) also appears in the highest interval.

Those results at different heights represent an addition to the historical records mentioned in Section 2, and also to surface winds measured with a collocated anemometer (not included here for brevity). Such surface measurements were consistent with Doppler lidar results at lower heights with the predominant winds from the west and northwest during day-time and prevailing nocturnal winds from the east and southeast.

Finally, as a common feature for all heights and interval times, strong winds (more than 10 m/s, represented with brown color) are not frequent, but they come mostly from the SW when present. This is parallel to the Sierra Nevada mountain range axis, meaning that the wind coming from this direction is more likely to have higher speeds without orographic limitation.

4.2. ABL Turbulent Sources Characterization

Figure 6 displays the relative frequency of three of the turbulent source types, namely ‘convective mixing’ (Figure 6a,d,g,j), ‘cloud driven’ (Figure 6b,e,h,k) and ‘wind shear’ (Figure 6c,f,i,l), with respect to time of the day and altitude. Subplots from the different rows correspond to the distinguished seasons.

The convective mixing frequency plots show a clear diurnal evolution of the range gates with the most frequency (red colored, meaning more than 70%). This high frequency of convective mechanism usually starts with sunrise in the lowest heights and reaches higher altitudes throughout the day up to a maximum around afternoon, when the altitudes with high convective frequency start decreasing until sunset. Therefore, this statistical behavior is directly related with the process of the ABL growth through convection, but also with the processes that could shift or weaken such growth, such as temperature inversions or cloud presence. In spring (Figure 6d) and summer (Figure 6g), this maximum height seems to remain constant at 775 m a.g.l. from 14:30 h to 15:30 h UTC (for spring) and at 955 m a.g.l. from 12:30 h UTC to 16:30 h UTC (for summer), but it is likely due to lack of data over a certain height. Despite the expected presence of cases when the convection was strong enough to make the convective ABL to rise over those heights, the system might not receive signal with enough SNR to perform the complete analysis for the ABL classification above a certain altitude. This effect did not occur in winter (Figure 6a) and autumn (Figure 6j), when the maximum height with more than 70% convection was registered at 475 m a.g.l. at 14:30 h UTC and at 685 m at 14:45 h UTC.

Wind shear plots (Figure 6c,f,i,l) only presented frequency values during the night-time, because the classification algorithm assumes that convective mixing always dominates the surface-driven turbulence when it is present. The season with highest frequency of wind shear-driven turbulence was found to be autumn (Figure 6l), when this mechanism was responsible for more than 20% of the detected turbulence below 500 m a.g.l. around 20 h UTC. For spring and summer (Figure 6f,i), wind shear was responsible for around 10–15% of the detected turbulence in the same time interval (around 20 h UTC) at low altitudes, and in a separated layer around 500 m a.g.l. after 00 h UTC. It is likely that those high and stable wind shear frequencies between midnight and early morning in a decoupled elevated layer are related to the presence of nocturnal low-level jets [36,85], although a deeper investigation of this phenomenon would be needed to confirm it.

Cloud driven turbulence was more significant in winter and spring, with frequencies always less than 10% during the whole day, and mainly concentrated between 500 and 1000 m a.gl. In spring, nocturnal cloud-driven turbulence was also significantly detected close to the surface. During summer, this turbulence source was essentially only detected during night-time or in connection with convective clouds formed during ABL growth in the early morning (5–10 h UTC in Figure 6h).

Figure 7 investigates the diurnal and seasonal behavior of the assigned mixing sources, but taking into account all the possible given labels. For this analysis, we have taken five different range intervals and calculated the accumulated frequency of each source type with a 3 min time resolution. The amount of missing data at different heights presents diurnal dependence was directly related to ABL height, as expected because the SNR is usually much higher inside the ABL.

Figure 7 clearly shows the strong diurnal variation in convective mixing as well as its seasonal dependence, typical features for ABL in the mid-latitudes [20,22,86]. The frequency of this turbulent mechanism is almost 100% in the central hours for altitudes up to 700 m a.g.l. in summer, and a bit lower (but still high) in spring. This feature was not present in the analysis shown by [64] for Germany and Finland, where the surface-connected convective mixing was dominant during the central hours, but there was no negligible contribution of non-turbulent and cloud-driven turbulence cases.

Clouds frequency is very low for almost all subplots in Figure 7, and almost absent in summer. However, it must be noticed that the full cloud cover is not necessarily captured with this analysis, because there might be ABL-associated clouds above the maximum selected range of 1600 m a.g.l. In-cloud and cloud-driven cases were more frequent during winter and spring for altitudes above 1000 m a.g.l., with constant frequencies of 5–10%, as shown in Figure 6.

Wind shear was detected for all seasons during night-time at the lowest height interval, with a mean frequency (±standard deviation) of 8 ± 3%. For all seasons, mean wind shear frequency decreased to 3.5 ± 1.4% in the second height interval, and this source was only found to be relevant in the third height interval (700–1000 m a.g.l.) before sunrise in spring and summer, with 3.4 ± 1.3% frequency.

The ‘intermittent’ category, given to range gates with turbulence that was not related to the rest of sources, was the most frequently detected during night-time. This feature was in contrast to the high latitude scenario analyzed by [64], where the contribution of intermittent sources was mainly during midday hours in spring and summer. In contrast, for the mid-latitude site (Germany) analyzed in that study, the results were similar to our analysis for Granada, although the frequency of intermittent sources was much less, in favor of non-turbulent cases. A possible reason of the high frequencies might be the remaining presence of turbulence forming the residual layer, but we do not discard the need for finding a dynamic threshold for the turbulence presence, which could be too low for certain atmospheric situations. This category cannot thoroughly be investigated with the methodology used here; therefore, we are not able to discuss the reasons of this in depth. However, we found it important to highlight its presence and the need for further research.

5. Conclusions

In this study, we applied Doppler lidar technique to characterize wind field and ABL turbulent properties with high temporal and vertical resolution over an urban area with complex topography. The analyzed ABL properties were obtained with the Halo lidar toolbox, a software package developed in the framework of ACTRIS-Cloudnet that allows for continuous, automatic and standardized retrievals. With this methodology, we were able to gather a 2-year database of wind field and turbulence measurements for the Granada urban environment.

The statistical analysis of the horizontal wind over Granada revealed diurnal and seasonal patterns of the mean profiles from 100 to around 1000 m a.g.l. From this study, an important conclusion to be highlighted is the variation in observed mean winds with height, especially for wind directions. This feature supports the relevance of wind profiling within the ABL, because these variations could not have been inferred from surface measurements and are not probably captured by satellite measurements or global models. Moreover, we found that mean and standard deviation (or variance) were not sufficient to describe wind data, but Weibull and von Mises distributions were needed to model wind speeds and directions, respectively. The patterns showed two different wind regimes at altitudes close to the surface versus higher altitudes, in terms of Weibull shape parameter. Direction analysis revealed mean and modal daytime winds from the NW close to surface, with highest velocities in spring and summer, and weaker nocturnal winds from the SE. This clear pattern in wind direction is consistent with the expected katabatic wind at the analyzed altitudes due to the interaction between the city of Granada and the close mountain range of Sierra Nevada. Moreover, the NW–SE pattern was found to rotate counterclockwise with increasing altitude, as revealed by the wind roses and von Mises models.

Statistical analysis of the turbulent sources was performed for the same database. We could characterize the seasonal and diurnal behavior of the convective mixing, rising up to maximum heights around 14:30 h UTC for winter and autumn, and in a longer interval centered at the same time for spring and summer. Wind shear was found to be of similar importance (around 8%) at heights below 400 m a.g.l. for all seasons, with a significant presence up to 1000 m a.g.l. only for hours before sunrise in spring and summer. Cloud-driven turbulence was present mostly in winter and spring, with frequencies less than 10%.

Some conclusions can be also drawn from the combined analysis of wind and turbulence. In particular, we can observe that, in Granada, the nocturnal wind profiles are more stratified in terms of wind direction, whereas the diurnal profiles are more homogeneous in correlation with the hours of maximum convective mixing frequency. Moreover, regions where wind shear is frequently detected coincide with high average winds with direction gradients.

Finally, we believe that there are clear potential applications of this kind of information in air quality and pollutant dispersion studies. Poor air quality events are closely related to the capability of the ABL to mix and disperse the emitted pollutants, a complex issue that Doppler lidar measurements can clearly help to understand. Therefore, this study represents a solid basis to start specific and deeper analyses of the different patterns and features observed, e.g., finding different patterns of horizontal winds that are responsible for the analyzed mean behavior, or investigating the causes of the high nocturnal turbulence measured in Granada.