## Abstract

A case study of mountain-wave-induced turbulence observed during the Terrain-Induced Rotor Experiment (T-REX) in Owens Valley, California, is presented. During this case study, large spatial and temporal variability in aerosol backscatter associated with mountain-wave activity was observed in the valley atmosphere by an aerosol lidar. The corresponding along- and cross-valley turbulence structure was investigated using data collected by three 30-m flux towers equipped with six levels of ultrasonic anemometers. Time series of turbulent kinetic energy (TKE) show higher levels of TKE on the sloping western part of the valley when compared with the valley center. The magnitude of the TKE is highly dependent on the averaging time on the western slope, however, indicating that mesoscale transport associated with mountain-wave activity is important here. Analysis of the TKE budget shows that in the central parts of the valley mechanical production of turbulence dominates and is balanced by turbulent dissipation, whereas advective effects appear to play a dominant role over the western slope. In agreement with the aerosol backscatter observations, spatial variability of a turbulent-length-scale parameter suggests the presence of larger turbulent eddies over the western slope than along the valley center. The data and findings from this case study can be used to evaluate the performance of turbulence parameterization schemes in mountainous terrain.

## 1. Introduction

The structure of turbulence over flat, homogeneous terrain and under various atmospheric boundary layer (ABL) conditions has been investigated extensively (e.g., Wyngaard and Cote 1971; Kaimal et al. 1972, 1976; Oncley et al. 1996; Albertson et al. 1997; Piper and Lundquist 2004). These measurements and analyses have led to similarity theories that describe the flow over flat, homogeneous terrain very well. Although the need for the adaptation of the existing theory or the development of a new theoretical framework applicable to complex terrain has been recognized, the structure of turbulence in complex terrain remains relatively unexplored. Among a few papers devoted to this subject are those of Rotach et al. (2004), Weigel and Rotach (2004), and Rotach and Zardi (2007). Rotach et al. (2004) point out that, because turbulence parameterization schemes in atmospheric numerical models are based on measurements conducted over flat and horizontally homogeneous terrain, numerical models may have problems when dealing with complex terrain. Weigel and Rotach (2004) indicate that the turbulence structure is greatly influenced by the configuration of the underlying terrain. Rotach and Zardi (2007) argue that an extension and/or modification of the classical scaling for turbulence variables is required for which measurements in complex terrain are needed.

A correct parameterization of boundary layer turbulence is important for a variety of applications, such as air pollution dispersion and wind energy. In these applications, turbulence parameters are used as an input to atmospheric dispersion models and wind energy efficiency. Turbulent kinetic energy (TKE) is one of those parameters that describe the turbulence intensity in the ABL. In high-resolution numerical modeling, turbulent mixing is often parameterized using TKE. To advance the development of turbulence parameterization schemes, it is therefore important to investigate and document the temporal evolution and spatial structure of TKE in complex terrain.

In the Cartesian coordinate system, the mean TKE per unit mass is defined as the sum of variances of all three wind speed components (*u _{i}*,

*i*= 1, 2, 3):

The TKE budget equation is expressed as (e.g., Stull 1988)

where *U _{i}* and

*u*′

*are the mean and perturbation parts of velocity in the*

_{i}*x*direction, respectively, Θ

_{i}*and*

_{υ}*θ*′

*are the mean and perturbation part of the virtual potential temperature,*

_{υ}*p*′ is the perturbation part of the atmospheric pressure,

*ρ*is the mean air density, and

*ν*is the kinematic molecular viscosity. Repeated indices indicate summation, and bars denote time averaging. Our approach is to determine as many terms as possible. We unfortunately cannot determine the terms that describe horizontal advection and horizontal shear. We also cannot calculate the pressure correlation , which requires high-frequency pressure measurements that were not made during the Terrain-Induced Rotor Experiment (T-REX). The extent to which the TKE balance is closed can then provide some indication of the importance of horizontal heterogeneity and pressure redistribution. After rotating the coordinate system in the mean wind direction, the one-dimensional TKE budget equation can now be written as

where *U* and *V* are the mean values of the horizontal wind speed components (streamwise and transverse, respectively) and *W* is the mean vertical wind speed component. Term I represents a local storage of TKE, which is equal to the increase or decrease of TKE in time at a given location due to all of the TKE production and destruction terms. These terms include advection of TKE from the layers above or below by the mean vertical wind (term II); the buoyant production/consumption (term III), which depends on the sign of the heat flux ; the mechanical (shear) production (term IV), which is typically positive in the boundary layer because of opposite signs of the horizontal momentum fluxes and vertical wind shear; the vertical transport of TKE by turbulent eddies (term V); the dissipation *ɛ* of TKE into heat by molecular viscosity; and the residual Rs containing terms that describe horizontal advection and pressure redistribution of TKE and horizontal wind shear. Note that in most studies and models term II and the third part of term IV (which includes vertical convergence/divergence) are neglected. Our dataset offers the opportunity to calculate these terms, which might be significant because of the presence of upslope and downslope flows (see section 3 for details).

The overall evolution of TKE depends on the TKE production and destruction processes described by various terms in Eq. (3). Calculating these terms is nontrivial. For example, extremely high data sampling frequency is required for direct calculation of *ɛ* (e.g., Piper and Lundquist 2004). Nevertheless, there are some alternative methods for the evaluation of *ɛ* that are not as demanding, one of which will be described later in the text. In turbulence parameterizations used in atmospheric numerical models, *ɛ* is typically parameterized using TKE (e.g., Mellor and Yamada 1974). By quantifying the terms in the TKE budget Eq. (3) that determine the magnitude of TKE, turbulence parameterization schemes could therefore be evaluated. Moreover, *ɛ* provides insight into the local turbulence strength. This finds a variety of applications such as turbulence nowcasting at airports related to the presence of aircraft wake vortices (e.g., Frech 2007).

The objectives of this paper are to investigate the spatial and temporal variability of TKE, to investigate the applicability of a one-dimensional model of the TKE budget described by Eq. (3) during a mountain-wave event over complex terrain, and to obtain a dataset that can be used in the evaluation of turbulence parameterization schemes. The event took place during the first intensive observing period (IOP 1) of T-REX. This event was well documented by a scanning aerosol backscatter lidar (De Wekker and Mayor 2009) and a dense network of weather stations (Grubišić and Xiao 2006). In section 2, we briefly describe the T-REX field campaign and the data used to study the IOP 1 mountain-wave event described in section 3. This is followed by a detailed investigation of turbulence characteristics in section 4 and spatial and temporal characteristics of the TKE budget during this event in section 5. Section 6 summarizes the conclusions of this study.

## 2. Data

This study uses data collected during T-REX, which took place in Owens Valley, California, in March and April of 2006 in an effort to explore the structure and evolution of atmospheric rotors (Grubišić et al. 2008). Owens Valley lies east of the southern Sierra Nevada, which is the tallest, steepest, quasi-two-dimensional topographic barrier in the contiguous United States. Mountain waves and rotors are known to reach particularly striking amplitude and strength there. Climatological studies have shown that in the Sierra Nevada the months of March and April have the highest frequency of mountain-wave events (Grubišić and Billings 2007). Ground-based and airborne in situ and remote sensing measurements were conducted in Owens Valley during the 2-month period [for details, see Grubišić et al. (2008)].

In this paper, we primarily use data from the T-REX network of 16 automatic weather stations (TAWS) installed across Owens Valley by the Desert Research Institute (Grubišić and Xiao 2006), from three 30-m towers from the National Center for Atmospheric Research (NCAR) Integrated Surface Flux Facility (ISFF), and from the scanning Raman-Shifted Eye-Safe Aerosol Lidar (REAL; De Wekker and Mayor 2009). The location of the TAWS, the ISFF towers, and the REAL are shown in Fig. 1. The location and height of the three ISFF towers and REAL are listed in Table 1. The western tower (WT) was located on the western side of the valley, on an alluvial fan in the lee of the Sierra Nevada. The central tower (CT) and the southern tower (ST) were located along the valley’s central axis extending from north-northwest to south-southeast. The WT and CT were aligned across the valley. Between these two towers, REAL was positioned on the western slope, providing an undisrupted view in all directions within its range.

Each ISFF tower was instrumented with Campbell Scientific, Inc., model CSAT3 three-dimensional ultrasonic anemometers collecting data with a sampling rate of 60 Hz at heights of 5, 10, 15, 20, 25, and 30 m. Consecutive blocks of 6 data points were averaged to create a 10-Hz database. This averaging reduces the effects of aliasing in the calculation of the turbulence spectra (section 4).

Besides anemometers, pressure sensors were installed at a height of 1 m and temperature/relative humidity sensors were installed at heights of 5, 15, and 30 m. There were some additional instruments on the towers such as radiometers, but the data from these instruments are not used in this study. A detailed description of the ISFF towers in T-REX can be found in Grubišić et al. (2008) and online (http://www.eol.ucar.edu).

## 3. Mountain-wave event

During T-REX, 15 IOPs and 5 enhanced observing periods (EOPs) were performed (Grubišić et al. 2008). Our period of interest is IOP 1, which took place from 0000 UTC 2 March to 1500 UTC 3 March 2006 (local standard time is UTC − 8 h). IOP 1 was characterized by significant mountain-wave activity with strong surface winds. Satellite images indicate the presence of trapped lee waves over the study area (not shown). Before the onset of the strong surface winds, a convective boundary layer (CBL) developed in the valley (De Wekker and Mayor 2009) with thermally driven upvalley and upslope flows (Fig. 2a). A transition from a quiescent valley atmosphere to a valley atmosphere affected by mountain-wave activity occurred during midafternoon between 2300 and 0000 UTC (Fig. 2b). While southeasterly flow prevailed in the valley center the flow on the western slope became more irregular with downslope winds replacing the upslope winds and colliding with the upvalley flow (Fig. 2c). Southeasterly winds up to 15 m s^{−1} in a layer between 400 and 1500 m in the valley center and a sudden transition to northwesterly winds were also observed by a network of radar wind profilers in Owens Valley (not shown).

The transition from a quiescent period with a developing CBL to a disturbed ABL was accompanied by large spatial inhomogeneities in the aerosol backscatter detected by the REAL. Downslope-directed flows to the west of the REAL carried clean air with low backscatter intensity to central regions of the valley, where it collided with aerosol-laden air transported by the channeled upvalley flows. The resulting aerosol structures suggest highly turbulent conditions around the location of the REAL (Fig. 3; see also Fig. 6d in De Wekker and Mayor 2009). These turbulent conditions are confirmed by the raw 10-m data of east (*u _{x}*) and north (

*u*) wind speed components from the three ISFF towers, which show large fluctuations around the 5-min means (Fig. 4) and will be further explored in section 4b. Although at the CT and ST a mostly steady southwesterly (upvalley) flow is observed, the wind at the WT appears more irregular with rapidly changing wind direction from upslope to downslope flows around 2300 UTC. The change in wind direction at the WT is also apparent from the aerosol backscatter data in Fig. 3. A detailed analysis of the turbulence characteristics using data in Fig. 4 is presented in section 4.

_{y}The turbulence structure near the surface during and after transition periods from a convective to a stable boundary layer is poorly known and is potentially important for the modification of mountain-wave activity and the generation of rotors and subrotors on the lee side of mountains. To investigate the spatial and temporal characteristics of TKE and the terms in the TKE budget equation [Eq. (3)] surrounding the event, we investigate in detail a 6-h period from 2100 UTC 2 March to 0300 UTC 3 March. At the beginning of this period, approximately 29 min of data are missing at the WT and approximately 14 min of data are missing at the ST at all six levels. At the CT the time series are complete, but anemometers at 15 and 30 m did not operate properly, and data from these two levels could not be used. Turbulence characteristics during the transition period and the mountain-wave event will be discussed in detail in the next sections.

## 4. Estimation of turbulence from ISSF tower data

### a. Determination of the averaging time

For the analysis of near-surface turbulence, it is necessary to separate turbulence from the mean flow. This is a nontrivial task, especially in complex terrain and for a flow situation such as presented here. A separation between the turbulent scales and synoptic scales is typically defined by the presence of a spectral energy gap at the mesoscale (e.g., Metzger and Holmes 2008; Vickers and Mahrt 2003). This gap is not always obvious, however. In our case, spectra of the 6-h streamwise component *u* at all three towers and all levels do not show the gap (not shown here), which is not surprising because this is a single realization over complex terrain during a mountain-wave event. Therefore, we have to use another method that is more suitable for determination of the turbulent scales. Our task is to find a proper averaging interval for our dataset, which will be used for the definition of turbulent perturbations that are necessary for estimation of turbulent properties.

Figure 5a shows a log–linear representation of the streamwise *u* spectrum at all three towers at 10 m for the 6-h period of interest. For periods of longer than 15 min, there is much more energy at the WT than at the CT and ST. One possible explanation for the larger energy at the WT is the direct influence of the mesoscale transport associated with the mountain-wave-induced circulations (as seen in the REAL data; Fig. 3). Because our major interest is in the turbulent transport, we would like to keep the mesoscale influence to a minimum and therefore decided to consider only periods that are shorter than 15 min. We further substantiate this by using ogive plots of turbulence fluxes (e.g., Oncley et al. 1996; Berger et al. 2001; Lee et al. 2004). An ogive is defined as a cumulative integral by

where *p* and *q* represent any two variables, Co* _{pq}*(

*f*) is their cospectrum, and the integration goes from higher frequencies toward lower frequencies. If an ogive converges starting from a certain frequency, this is an indication that there is no more flux beyond this frequency or—equivalent—this period; thus, this period may be taken as the averaging period for the determination of turbulent perturbations. Figures 5b–d show ogive plots of the kinematic heat flux for each hour of a 6-h period for all three towers at three different levels (except for the first hour at the WT and ST because of the lack of 29 and 14 min of data, respectively, at the beginning of the period). Although at the CT and ST there is more or less clear convergence of ogives at periods that are shorter than or close to 15 min, at the WT there is at least a good tendency toward convergence. Similar behavior is evident in the ogive plots of momentum flux (not shown here). In conclusion, an averaging interval of 15 min seems appropriate for our study. In addition, in section 4e we discuss the sensitivity of the results to the length of the averaging interval.

To get detailed insight into the temporal variability of TKE within the 6-h period, we create 15-min intervals shifted by 1 min. A maximum of 346 intervals is available for each sonic anemometer. For each of these intervals, the coordinate system is rotated into the mean wind direction. Turbulent perturbations of the *u*, *υ*, and *w* velocity components are calculated as the residuals between the instantaneous velocity component and its mean value during each 15-min interval.

### b. Nature of turbulence during IOP 1

To gain insight into the nature of turbulence before processing large amounts of high-frequency turbulence data, it is helpful to assess the relative contributions of mechanical production [term IV in Eq. (3)] and buoyant production/consumption [term III in Eq. (3)] to the total TKE. This assessment can be made by evaluating the gradient Richardson number *R _{i}* given by (e.g., Stull 1988):

where the buoyancy production/consumption is represented by the squared Brunt–Väisälä frequency *N*^{2} in the numerator. Large positive values (*R _{i}* ≫ 0) indicate stable environments with weak vertical shear. Values close to zero (

*R*≈ 0) indicate weak stability and/or strong vertical shear. Large negative values (

_{i}*R*≪ 0) indicate statically unstable environments with weak vertical shear.

_{i}If the partial derivatives in Eq. (5) are approximated with finite differences using observations made at discrete height intervals, we can evaluate the bulk Richardson number *R _{B}* as follows (e.g., Stull 1988):

which is calculated using the mean wind speed and temperature data (15-min means), between 5 and 30 m for the WT and ST and between 5 and 25 m for the CT in the period from 1200 UTC 2 March to 1200 UTC 3 March (Fig. 6). During the 6-h period of interest (indicated by the vertical dashed lines in Fig. 6), the transition from unstable (*R _{B}* ≪ 0) to stable conditions (

*R*≫ 0) is clearly present at the WT (Fig. 6a). In the same period,

_{B}*R*≈ 0 at the CT (Fig. 6b) and ST (Fig. 6c), mostly because of strong wind shear, indicating strong local turbulent mixing along the valley. The wind shear is approximately 15 times as strong as

_{B}*N*

^{2}at the CT and ST, whereas at the WT the wind shear term and

*N*

^{2}appear to be of the same order of magnitude. Prior to this period, the situation is very unstable at the WT and ST while at the CT

*R*oscillates around zero. The

_{B}*N*

^{2}at the WT and ST is very small and negative, and wind shear is even smaller but positive. Therefore, the ratio of these two terms results in significant negative values of

*R*. At the CT, wind shear is much stronger than

_{B}*N*

^{2}in this period, leading to

*R*≈ 0. These observations are consistent with the turbulent structure revealed by the REAL (section 3).

_{B}The terms in Eq. (3) can be directly calculated with high-frequency measurements from the ISFF towers. Terms II, IV, and V contain vertical derivatives that are approximated using the measured values at two consecutive levels. The center of these levels, or midlevel, is the location at which all turbulent quantities calculated in the next sections will be specified. We start our investigation of the TKE budget with the calculation of the eddy dissipation rate *ɛ*. Estimation of this term relies on spectral analysis.

### c. Evaluation of eddy dissipation rate and spectral analysis

For evaluation of the eddy dissipation *ɛ*, we used the inertial dissipation method (IDM) based on Kolmogorov’s universal equilibrium hypotheses (e.g., Tennekes and Lumley 1972). For measurements that are performed in time at a single location, this method requires validity of Taylor’s hypothesis (TH) on frozen turbulence to transform from space (wavenumber) to time (frequency) domains. One can apply TH if the ratio between the standard deviation *σ _{M}* and the mean value

*M*of the horizontal wind speed in the observed time interval is less than 0.5 (e.g., Stull 1988):

It has been shown by Kolmogorov (e.g., Champagne 1978) that the energy spectrum of the velocity component that carries maximum variance follows the “ law” in the inertial subrange. In the log–log representation this can be expressed as follows (e.g., Večenaj et al. 2010):

where *k* is the wavenumber, *S _{u}*(

*k*) is the energy spectrum, and

*α*is the Kolmogorov constant (e.g., Stull 1988) for the streamwise velocity component. If the spectrum of the streamwise velocity component shows the behavior,

*ɛ*can be evaluated using TH from (e.g., Champagne et al. 1977):

where *U* is the mean streamwise velocity component. As mentioned before, this formulation of *ɛ* requires existence of the inertial subrange in which local turbulent isotropy is present. The behavior of spectra does not necessarily imply the presence of an inertial subrange, however. Champagne (1978) gives an overview of experiments in which the behavior even extended far into the low-frequency region where the flows were anisotropic. An indicator for the existence of the inertial subrange can be derived from Kolomogorov’s hypotheses (e.g., Tennekes and Lumley 1972). In theory, a ratio of the spectra of transverse (*υ*) to streamwise (*u*) velocity components and of the spectra of vertical (*w*) to streamwise (*u*) velocity components is anticipated in the inertial subrange (e.g., Batchelor 1959; Champagne 1978; Biltoft 2001).

We will now apply a basic spectral analysis to the 10-Hz wind speed data at all three towers. For this purpose a 1-h interval is chosen in which the time series of horizontal wind speed components at the WT (Figs. 4a,b) are most stationary (between 2221 and 2321 UTC 2 March). It is shown in Table 2 that the TH criterion is satisfied for all three towers in this interval. During this particular hour, a transition from unstable to stable conditions occurs at the WT while strong shear accompanied by weak stratification prevails at the CT and ST (Fig. 6). After rotating the coordinate system into the mean wind for each tower, the streamwise *u* and transverse *υ* velocity components are defined, and their corresponding spectra *S _{u}* and

*S*as well as vertical

_{υ}*S*are calculated as a function of frequency. Figure 7 shows the

_{w}*S*/

_{υ}*S*and

_{u}*S*/

_{w}*S*ratios over the chosen 1-h interval for all three towers at heights of 10 (Fig. 7a) and 25 (Fig. 7b) m. The

_{u}*S*/

_{υ}*S*ratios at the CT and ST approach near

_{u}*f*= 0.4 Hz and

*f*= 0.2 Hz at 10 and 25 m levels, respectively. At the WT, the ratio is smaller but still almost at the 25-m level. The

*S*/

_{w}*S*ratios converge to 1 at all three towers. Note that the ratio is anticipated theoretically and that, as pointed out by Biltoft (2001), there is no convincing experimental support for this theoretical ratio. Biltoft (2001) provides a number of possible explanations, including the presence of considerable scatter in lateral-to-longitudinal ratios of sonic-derived spectra caused by the data processing. A detailed investigation of the

_{u}*S*/

_{υ}*S*and

_{u}*S*/

_{w}*S*ratios for the T-REX dataset would require the analysis of a larger time interval than is considered here and is outside the scope of the current paper. We will continue our analysis of the eddy dissipation while assuming the presence of an inertial subrange. This assumption is supported by plotting compensated spectra

_{u}*f*

^{5/3}

*S*(

_{u}*f*), as shown in Figs. 7c and 7d. Compensated spectra show the inertial subrange as a horizontal line over a certain frequency band (e.g., Piper and Lundquist 2004). The compensated spectra of the

*u*component at all three towers at 10 and 25 m, shown in Figs. 7c and 7d, respectively, stabilize to a certain value at a frequency of approximately 0.4 Hz. Some more support of the assumption is provided by the behavior over frequencies of higher than 0.4 Hz in the power spectra in log–log representation (Fig. 8). For these reasons, we use Eq. (9) for estimation of

*ɛ*at all three towers.

The test of the TH criterion [Eq. (7)], which is required for the estimation of *ɛ*, is applied to all 346 of the 15-min intervals. This test resulted in the removal of a significant number of intervals from further analysis at the WT and removal of some intervals at the CT. More detailed information about the number of excluded 15-min intervals is given in Table 3. Only the perturbations of streamwise velocity component are used for the calculation of *ɛ*, because this component is least contaminated by self-induced distortion by the anemometer (e.g., Piper and Lundquist 2004).

Next, a specific procedure was followed to calculate *ɛ* using IDM. As shown in Fig. 8, the slope starts to appear at a frequency of 0.4 Hz. By expanding a frequency band from 0.4 Hz toward higher frequencies, successively adding one discrete point in the frequency domain produces a whole ensemble of frequency bands. The spectra at these bands are then linearly fitted, and the band at which the slope is closest to the is used for the calculation of *ɛ*. It appears that for most 15-min intervals the best fit is obtained in the band between approximately 0.4 and 4 Hz (Fig. 8). The mean values of the term *f*^{5/3}*S _{u}*(

*f*) are then calculated for each 15-min interval. Following Eq. (9),

*ɛ*is calculated, with the Kolmogorov constant

*α*taken here to be 0.53 (e.g., Oncley et al. 1996; Piper and Lundquist 2004; Večenaj et al. 2010). The final result is a time series of 15-min

*ɛ*values derived at all levels during the 6-h period of interest.

To evaluate the robustness of these *ɛ* estimates, we have used an alternative method to calculate *ɛ* that relies on the third-order structure function (e.g., Albertson et al. 1997). Results of the comparison of *ɛ* obtained by these two methods are given in appendix A, where a decent agreement between the methods is shown. The TKE from Eq. (1) as well as terms I–V in Eq. (3) are calculated in the next section.

### d. Evaluation of TKE and terms I–V

The TKE values are calculated using Eq. (1) for each 15-min interval. As an example, the time evolution of TKE for the CT at the 10-m level is shown in Fig. 9a. In the 6-h period, TKE is high between 2200 UTC 2 March and 0100 UTC 3 March, which approximately corresponds to the period of the coherent aerosol structures observed by the REAL (De Wekker and Mayor 2009, their Fig. 6).

The local storage term I from Eq. (3) is calculated for each level using the central-differencing method in time. It is found that this term is at least an order of magnitude smaller than any other calculated term and therefore can be considered to be negligible. According to Kaimal and Finnigan (1994), a balance between turbulence production and destruction (or local equilibrium) may occur in turbulent flows over complex terrain. Our estimations of term I indicate that this situation occurs for the period of interest.

The buoyancy term III is calculated as defined in Eq. (3) for each level. The calculation of the vertical advection term II, the mechanical term IV and the vertical transport term V from Eq. (3) requires two levels of data. These three terms are therefore calculated in the layer between the two consecutive levels using a second-order central-differencing method in height. The value of the terms in this layer is assigned to a height between the two levels or “midlevel.”

As an example, the results of the TKE budget analysis at the 7.5-m midlevel for the CT are shown in Fig. 9b. The mechanical term IV (black solid curve) and *ɛ* (gray solid curve) are largest and show a similar temporal variability. The two terms are highly correlated during the entire 6-h period, with a correlation coefficient squared *r*^{2} equal to 0.97. Therefore, in this particular case, the mechanical shear and molecular viscosity are the dominant physical mechanisms that determine the TKE. The mechanical term is governed by the vertical shear of the horizontal wind speed. The vertical wind convergence/divergence (third part of term IV) turns out to be negligible (in all of our cases presented later). Vertical advection, buoyancy, and the vertical transport term are denoted in Fig. 9b by the gray dashed, black dotted, and black dashed curves, respectively. Although vertical advection and buoyancy can be considered to be negligible during the observed period, the transport term fluctuates around zero, alternatively contributing to the TKE production and destruction. This example is consistent with the *R _{B}* analysis in section 4b, which indicates shear-dominated turbulence at the CT.

Analysis of the uncertainty of the estimated TKE and different TKE budget terms shows that estimations of term V, the vertical transport term, show the largest uncertainty. This uncertainty analysis is presented in appendix B using a representation of the statistical characteristics of each budget term. The results from the TKE budget analysis presented above will be discussed for all different heights and towers in detail in section 5.

### e. Sensitivity of the results to the choice of averaging interval

The optimal averaging time for the calculation of turbulence terms over complex terrain has been discussed in several papers but remains an open question (e.g., Sakai et al. 2001; Finnigan et al. 2003; Lee et al. 2004). To investigate the sensitivity of the TKE budget terms to the averaging interval, we compare the terms calculated using the 15-min averaging interval with the terms calculated using 1- and 30-min averaging intervals for the 6-h period of interest. The most sensitive variables appear to be TKE and the mechanical term (IV), and the least sensitive is *ɛ*. Figure 10 shows results of this analysis for TKE and *ɛ* at the 10-m level and for the mechanical production term IV at the 7.5-m midlevel. In Figs. 10d–f, term IV is positive and is placed in the first quadrant while *ɛ* is presented as negative and therefore is placed in the fourth quadrant. TKE values increase by using the 30-min averaging interval, especially at the WT and CT (Figs. 10a and 10b, respectively). This dependence of the TKE on the averaging interval is most likely due to the influence of mesoscale transport, which will be discussed in more detail in the next section. The mechanical term IV (first quadrant of Figs. 10d–f) is much smaller for the 1-min averaging interval than for the 15-min averaging interval, whereas the differences in the terms between the 30- and 15-min averaging intervals appear negligible. Dissipation rate *ɛ* (fourth quadrant of Figs. 10d–f) is least sensitive to the averaging time, with reduced scatter when larger averaging intervals are used.

## 5. Turbulence structure in Owens Valley

### a. Time evolution of the vertical turbulence structure

Time series of TKE from different levels within the selected 6-h period at all three towers are shown in Fig. 11. TKE at the WT (Fig. 11a) suddenly increases around 2315 UTC just around the time when downslope winds replace the former upslope winds (Figs. 2a,b). TKE increases with height at the WT with values from 9 m^{2} s^{−2} at 5 m to values of 16 m^{2} s^{−2} at 30 m. Along the valley axis [ST (Fig. 11c) and CT (Fig. 11b)], TKE is weaker and does not show a strong dependence on height, with a slight increase with height at the CT (from its maximum value of 6 m^{2} s^{−2} at 5 m to 7 m^{2} s^{−2} at 25 m) and a slight decrease with height at the ST (from its maximum value of 6 m^{2} s^{−2} at 5 m to 5 m^{2} s^{−2} at 30 m). A possible interpretation of these observations is that the TKE structure at the western alluvial slope is affected at the higher levels by mountain-wave activity whose influence is also felt by the higher levels at the CT (but not at the ST). Data obtained by the REAL show that mountain-wave-induced downslope advection of stable stratified air starts at the location of the WT around 2200 UTC 2 March (De Wekker and Mayor 2009, their Fig. 6b), with increased mixing aloft (above ~1 km). Time evolution of the aerosol structure indicates that this increased mixing is caused by the interaction of the downslope flow and the southeasterly flow along the valley, resulting in generation of coherent structures with a spatial scale of several hundred meters (De Wekker and Mayor 2009).

The time series of the TKE budget terms at different midlevels (Fig. 12) show that at the CT and ST the mechanical term and the dissipation rate *ɛ* are dominant mechanisms in TKE production/destruction at all levels. The mechanical term is considerably stronger than *ɛ*, however, especially at the lower levels (Figs. 12b,c). At the WT (Fig. 12a), the situation is considerably different. In general, the values of the budget terms are approximately one-half of those at the CT and ST, which contrasts with the stronger TKE at the WT when compared with the TKE at the CT and ST. Also, the vertical advection term (II) and the vertical transport term (V) play a significant role in the TKE budget at the WT. These observations confirm the different nature of turbulence at the western slope than along the valley axis.

Term I is at least an order of magnitude smaller than any of the other calculated terms (a maximum value of 0.004 m^{2} s^{−3} is obtained at the 25-m level at the WT) and therefore might be considered to be negligible. The question that arises now is how well the simplified TKE budget in Eq. (3) is closed by the calculated terms at the different towers and the different midlevels. In Fig. 13, we have summarized all of the positive and negative contributions of the terms in Eq. (3) to the TKE budget for each 15-min interval. The difference between the TKE production and destruction is also shown in Fig. 13 and is significant at all three towers at midlevels. At the CT and ST (Figs. 13b and 13c, respectively), the difference is mainly positive, which means that we measure more production than dissipation of TKE, whereas at the WT the difference changes sign. Overall, it can be concluded that the simplified TKE budget in Eq. (3) is not closed. The pressure correlation term, the horizontal TKE advection terms, and the horizontal wind shear terms (terms that we cannot calculate with our dataset and that are summarized in the residual term Rs) obviously play important roles in the TKE budget during this event.

The picture that emerges from these observations is that—in particular, at the WT, where TKE is strongest and the budget terms are smallest—TKE is not only produced locally but is also redistributed from layers aloft (where the turbulence activity was large) and advected from the area upslope of the WT toward the valley central axis. Snapshots and animations of the aerosol structure derived from the scanning backscatter lidar (Fig. 3 and De Wekker and Mayor 2009) along with the TKE behavior with height at the WT and CT support this picture.

### b. Turbulent-length-scale parameter

Many mesoscale models that parameterize the effects of turbulence use a local closure technique (e.g., Stensrud 2007). If a prognostic equation for TKE is included, a parameterization for *ɛ* is required (e.g., Mellor and Yamada 1974):

where is the mean value of TKE and Λ is an empirical length-scale parameter. There is unfortunately no unique rule as to how to determine this parameter. It is often chosen by nudging the simulated flow to the observed flow (e.g., Stull 1988). Using the data obtained during T-REX, we can determine TKE and *ɛ* for a range of turbulence conditions in complex terrain (see section 4) and use those values to estimate Λ from Eq. (10).

The validity of Eq. (10) is tested on the 6-h interval using the mean values of *ɛ* and TKE estimated on 15-min intervals. The test is performed at all midlevels at all three towers. A power law of the form

with the a priori power coefficient 1.5, is fitted to the scatter diagrams of *ɛ* versus TKE (e.g., Večenaj et al. 2010). The best fit between TKE and *ɛ* determines the value of *a*, the inverse of which is the length scale parameter Λ (Fig. 14a). The spatial variability of Λ suggests that the turbulent eddies along the valley axis (CT and ST) were much smaller than those at the western alluvial slope (WT), where Λ is on the order of a few hundred meters. This is another indication that mesoscale transport associated with the downslope winds and mountain-wave activity has a great influence on the turbulence structure. Furthermore, Λ increases with height—a result that is consistent with the concept that turbulent eddy sizes are limited by the proximity to the ground. To assess the robustness of the Λ calculations, the goodness of the fits that were used to retrieve *a* in Eq. (11) is also shown in Fig. 14a using *r*^{2}. The goodness of the fit deteriorates with height and is best at the ST, where *r*^{2} is larger than 0.5 at all levels. At the WT and especially at the CT, considerably smaller values of *r*^{2} are obtained. This indicates that the determination of a single length-scale parameter to relate the eddy dissipation to the TKE is questionable in our case.

Other approaches exist to estimate the parameter Λ. For example, following Kaimal and Finnigan (1994), Λ can be defined as

where *U* is the mean streamwise velocity component and *τ* is the integral time scale defined as an integral of the *u* wind speed component autocorrelation function. To be consistent with our 15-min averaging interval, we filtered the 6-h *u* component using 15-min moving averages and used the filtered time series for the autocorrelation function. Values of Λ calculated using Eq. (12) coincide at the 5-m level at all three towers (Fig. 14b), adding some credibility to the values obtained using Eq. (10). Toward higher levels, differences between the Λ values determined by the two methods become larger with much smaller values at the WT using Eq. (12). One possible explanation is that the mean horizontal wind direction during this 6-h period is poorly defined at the WT, resulting in very low values of *U* in Eq. (12).

The results presented above demonstrate the weakness in using a single turbulent-length-scale parameter to describe turbulence characteristics in complex terrain where nonlocal effects and mesoscale transport can become important.

## 6. Conclusions

In this study, data from T-REX IOP 1 are analyzed to investigate properties of the near-surface turbulence during a mountain-wave event. Temporal and spatial variability of the terms in the TKE budget are obtained from high-frequency data collected at multiple levels on three 30-m NCAR ISFF towers. A transition from a quiescent to a disturbed boundary layer is accompanied by large variability of TKE in space and time. Observations and analyses show that mesoscale transport associated with mountain-wave activity greatly enhances TKE on the western side of the valley. The contribution of the mesoscale transport becomes obvious when increasing the averaging interval for the turbulence calculations. In the central part of the valley, TKE is also strong but is mostly determined by local processes associated with the along-valley flow. Nevertheless, some evidence of the influence of scales that are larger than the local scale is still present at the highest levels on the tower that is located in the same cross-valley cross section as the western-slope tower.

We have shown that the TKE budget is poorly closed by the evaluated terms at all three locations. Possible reasons for this imbalance are the redistribution of TKE by the pressure fluctuations, the horizontal advection of TKE, and the horizontal wind shear—processes that we cannot estimate from our data.

Turbulence parameterizations in atmospheric numerical models often use a turbulence-length-scale parameter whose value is not well known. Using our calculations of the dissipation rate and the turbulent kinetic energy, we obtained values of this length-scale parameter. The spatial variability of the length-scale parameter suggests the presence of larger turbulent eddies at the western alluvial slope than along the valley axis. This finding is in agreement with the observation of aerosol structures from the scanning lidar. Turbulence parameterizations need to take into account the spatial variability of turbulent length scales in situations in which mountain-wave activity is apparent. Because boundary layer effects are important in the forecasting of mountain waves (e.g., Doyle and Durran 2007), an appropriate turbulence parameterization in these situations is highly needed. The data and results from our case study can be used to develop and evaluate turbulence parameterizations for complex terrain.

Among the many observations made during T-REX IOP 1, scanning aerosol lidar data proved to be particularly helpful in interpreting the results of the TKE budget analysis and in providing a more complete understanding of temporal and spatial variability of turbulence in mountainous terrain.

## Acknowledgments

This work was initiated while the first author, ŽV, was a visiting student at the University of Virginia (UVA). UVA is acknowledged for providing funding to ŽV to work on the T-REX dataset. The primary sponsor of T-REX is the U.S. National Science Foundation (NSF). This work was partially supported by the Croatian Ministry of Science, Education and Sports project BORA 119-1193086-1311, and by NSF Grant ATM-0524891 to the Desert Research Institute, which supported early stages of this work by VG. We thank Steve Oncley from NCAR for providing the raw turbulence data from the three ISFF towers, Tom Horst, Steve Oncley, and Branko Grisogono for helpful suggestions, and Ming Xiao from DRI for processing the TAWS data and generating Fig. 2. We also thank the three anonymous reviewers whose comments improved the quality of the manuscript.

### APPENDIX A

#### Comparison of *ɛ* Estimated Using IDM and 3SF

An alternative method for the evaluation of *ɛ*, applicable only to the inertial subrange in a strict sense and requiring isotropy, is based on the third-order streamwise velocity structure function (3SF). The method uses Kolmogorov’s four-fifths law (e.g., Piper and Lundquist 2004):

where *r* represents the spatial distance between the two measurements and the overbars denote averaging. Taylor’s hypothesis is used to convert Eq. (A1) to the time domain. If , where *τ* represents the time lag between the two measurements, then Eq. (A1) can be rearranged to

We used Eq. (A2) to estimate *ɛ* for each 15-min interval. The smallest period *τ* from our data for using 3SF is 0.1 s. According to Kolmogorov’s ⅘ law given in Eq. (A2), linear dependence of 3SF on *τ* is expected within the inertial subrange. We found that for different levels and different towers this linear dependence occurs when *τ* is between 5 and 30 s. Figure A1 shows the results of *ɛ* obtained using both IDM and 3SF for the CT at the 10-m level, indicating close agreement between the two methods (Fig. A1a) but with a slight underestimation of *ɛ* obtained by 3SF at higher values (Fig. A1b). Similar results are obtained at all levels at all three towers.

### APPENDIX B

#### Variability of TKE and TKE Budget Equation Terms

The variability of the magnitude of each individual term in the TKE budget equation can provide some insight to the uncertainty in the calculation of the terms. To this end, we created normalized histograms of certain term samples and then fitted a suitable probability density function (PDF) on these normalized histograms. Values of all terms are calculated for 346 of the 15-min intervals, each of which contains 9000 data points (10-Hz data). For each interval, TKE and TKE budget equation terms from II to V are calculated and then averaged. For the eddy dissipation *ɛ*, samples are much shorter (not exceeding 257 points) as a consequence of the limited frequency band that exhibits the behavior (see section 4c). As an example, Fig. B1 shows histograms of TKE and TKE budget terms. Also shown are the most suitable fitted distributions for one specific 15-min interval in the middle of the 6-h period of interest (between 2352 UTC 2 March and 0007 UTC 3 March) at the 10-m level and the 7.5-m midlevel at the CT. The most suitable distribution for TKE is the Weibull distribution (Fig. B1a); for terms II, III, and V, it is the normal distribution (Figs. B1b,c,e, respectively); for term IV, it is the exponential distribution (Fig. B1d); for *ɛ*, it is the lognormal distribution (Fig. B1f). In the background of a given panel in Fig. B1, PDFs of all 345 other 15-min intervals are shown. A lognormal distribution for *ɛ* has also been shown in other studies (e.g., Shaw and Oncley 2001; Goodman et al. 2006; Siebert et al. 2006). The vertical distribution of term V shows the largest scatter (Fig. B1e), and its calculation is therefore most uncertain.

## REFERENCES

## Footnotes

This article is included in the Terrain-Induced Rotor Experiment (T-Rex) special collection.