Impact of coherent eddies on airborne measurements of vertical turbulent fluxes

During the Hydrological-Atmospheric Pilot Experiment (HAPEX)-Sahel, which took place in Niger in the transitional period between the wet and dry seasons, two French aircraft probed the Sahelian boundary layer to measure sensible and latent heat fluxes. The measurements over the Niamey area often revealed organised structures of a few km scale that were associated with both thermals and dry intrusions. We study the impact of these coherent structures using a single day’s aircraft-measured fluxes and a numerical simulation of that day with a mesoscale model. The numerical simulation at high horizontal resolution (250 m) contains structures that evolve from streaks in the early morning to cells by noon. This simulation shows distribution, variance and skewness similar to the observations. In particular, the numerical simulation shows dry intrusions that can penetrate deeply into the atmospheric boundary layer (ABL), and even reach the surface in some cases, which is in accordance with the observed highly negatively skewed water vapour fluctuations. Dry intrusions and thermals organised at a few km scale give skewed flux statistics and can introduce large errors in measured fluxes. We use the numerical simulation to: (i) evaluate the contribution of the organised structures to the total flux, and (ii) estimate the impact of the organised structures on the systematic and random errors resulting from the 1D sampling of the aircraft as opposed to the 2D numerical simulation estimate. We find a significant contribution by the organised structures to the total resolved fluxes. When rolls occur, and for a leg length of about 30 times the ABL depth, the 1D sampled flux is shown to be sometimes 20% lower than the corresponding 2D flux when the 1D sampling direction is the same as the main axis of the M. Lothon (B) Centre de Recherches Atmosphériques, 8 route de Lannemezan, Campistrous 65300, France e-mail: lotm@aero.obs-mip.fr M. Lothon · F. Saïd Laboratoire d’Aérologie, UMR 5560 CNRS, Université Paul Sabatier Toulouse III, Toulouse, France F. Couvreux · S. Donier · F. Guichard · P. Lacarrère · J. Noilhan Météo-France/CNRM, Toulouse, France D. H. Lenschow National Center for Atmospheric Research, Boulder, CO, USA


Introduction
A fully turbulent, horizontally homogeneous and stationary medium is the basic framework assumed to describe the atmospheric boundary layer (ABL) measured from fixed ground stations, tethered balloons, or aircraft using a statistical approach. These observational platforms supply high resolution time series from which statistical moments are estimated. However, several problems have been raised about the eddy-correlation (covariance) method used to estimate the vertical fluxes in the ABL, due to sampling limitations in time and space.
For example, the surface fluxes deduced from the aircraft measurements made within the ABL by extending the vertical profiles of the fluxes down to the surface are often smaller than the estimates from ground measurements. This is commonly observed and only partially explained by filtering issues (Grunwald et al. 1998). One reason is that both aircraft and ground-based flux measurements are subject to sampling issues: the observations are limited in time and space to a range of sampling scales that determine the representativeness of the measurements. Lenschow et al. (1994) discussed the relationship between the length of aircraft flight legs and the associated systematic and random errors.
Coherent structures, which are organised structures of meteorological variables with a spatio-temporal periodicity, are at the core of this issue: their characteristic length scalesseveral times the depth of the ABL-belong to the upper size range of the spectrum of scales that contribute to the vertical transport of energy within the ABL. Consequently, they are difficult to sample by an aircraft leg that is a few tens of kilometres long. They are also fundamentally 3D structures that cannot be fully described with the 1D approach of an aircraft. Lohou et al. (1998) and Bernard-Trottolo et al. (2003) observed and characterised these structures and they showed with a numerical simulation (Lohou et al. 2000) how the anisotropy associated with these structures can make accurate flux estimation difficult.
During the Hydrological-Atmospheric Pilot Experiment (HAPEX)-Sahel experiment that took place in 1992 in Niger near Niamey (Goutorbe et al. 1997), two French aircraft flew at low altitudes to measure vertical fluxes in the lower troposphere. Organised structures were often detected, along with a skewed water vapour mixing ratio. Connections were found between the structures and boundary-layer entrainment-intrusion of dry air from above the boundary-layer top-as well as the presence of strong shear. Up to 40% relative error could be observed between the surface fluxes deduced from the aircraft measurements and the estimates made from the ground measurements, with higher probability of observing an underestimation of the surface flux with aircraft measurements.
Here we make use of a case observed during the HAPEX-Sahel experiment and simulated with a mesoscale model at a horizontal resolution as small as 250 m to provide further evidence of the impact of the coherent structures on aircraft measurements, and to stress the difficulty of measuring fluxes within an organised boundary layer.

A case study of the Sahelian boundary layer during the transition period between the wet and dry seasons
The HAPEX-Sahel experiment combined a long-term hydrological monitoring period from 1991 to 1993 and an 8-week intensive observation period that covered the end of the wet season and the beginning of the dry season in 1992 (Goutorbe et al. 1997). The daytime boundary layer observed on 8 October 1992 was numerically simulated in order to study the ABL during the transition between the wet and dry seasons. That day was part of an intensive observing period of the HAPEX-Sahel and occurred 23 days after the last rain in the one-degree square around Niamey (Wai et al. 1997), after a premature withdrawal of the intertropical convergence zone in Niger (Halpert and Ropelewski 1992). It was previously simulated with a mesoscale model at 10 km horizontal resolution by Taylor et al. (1997), who showed that such a simulation can represent relatively well surface heterogeneity and its effect on the energy budget of the atmosphere.

Numerical simulation set-up
The simulation was performed with the non-hydrostatic mesoscale Meso-NH model (Lafore et al. 1998). Three nested domains were defined, over which the horizontal resolution decreased from 5 km for the largest domain (D1) to 500 m for the medium domain (D2) and to 250 m for the smallest one (D3) (Fig. 1 (Masson et al. 2003), except for albedo, vegetation fraction and normalised difference vegetation index, which were obtained from the POLarisation and Directionality of the Earths Reflectances (POLDER) data analysis presented by Lacaze et al. (2003). The continental surface scheme is described by Noilhan and Planton (1989). The same turbulence scheme is used for the three domains (Cuxart et al. 2000), based on a prognostic equation for the turbulent kinetic energy. However, the 1D version of the scheme is used for D1, with the mixing-length parameterisation proposed by Bougeault and Lacarrère (1989), while the 3D scheme (Redelsperger and Sommeria 1986) is used for D2 and D3.

Aircraft measurements
On 8 October 1992, two aircraft flew at low altitudes to measure fluxes in the ABL between 0930 and 1200 UTC: the Merlin IV aircraft about 60 m above the ground and the ARAT aircraft about 300 m above the ground. The ARAT made two soundings at the beginning and end of the flight, west to Niamey at 13.5 • N, 2.2 • E. All legs were divided into 30 km segments to estimate the statistical moments of the wind components (u, v, w), the potential temperature (θ ), and water vapour mixing ratio (q) from the highest sample rate data (25 s −1 for the Merlin IV and 16 s −1 for the ARAT). A high-pass filter at 5 km wavelength was used to remove the mesoscale trend. The 24 segments flown by the ARAT aircraft are indicated in Fig. 2b. Hereafter, we use the letters 'a' and 'm' to denote segments flown by the ARAT and Merlin IV aircraft, respectively. For example, segment a10 is the 10th segment flown by the ARAT aircraft (see Fig. 2b). The Merlin IV aircraft flew a similar horizontal track to that of the ARAT and made a sounding at 14.5 • N, 2 • E. It also flew slightly different legs at the end of the flight, with a special final segment (m24) 30 m above Niger river. Segments m10 to m19 were flown approximately at the same time as segments a10 to a19, from 1050 to 1145 UTC. Figure 2a, b shows θ and q measured along the ARAT track at a 1 s sampling rate and corrected for a linear time evolution. A strong meridional gradient can be seen for both variables. The boundary layer is increasingly warmer and drier towards the north, in agreement with the climatology of the region. Using a least-squares fit along the track and assuming constant gradients along latitude, longitude and over time, we obtained a spatial variation of around 0.02 K km −1 and −0.03 g kg −1 km −1 with increasing latitude as well as a time variation of 1.8 K and −0.5 g kg −1 per hour for θ and q respectively.

Validation
Before using the numerical simulation to show the impact of thermals and dry intrusions, we compared the surface energy budget, the profiles of the mean meteorological variables, their horizontal gradients, and the distributions of their fluctuations. We found a good fit between the numerical simulation and the aircraft observations.
For instance, the model represents correctly the observed meridional gradient of increasing θ and decreasing q, as well as the diurnal warming and drying, at the time and level of the flight (300 m agl): meridional gradients of 0.01 K km −1 and −0.03 g kg −1 km −1 were found for θ and q, respectively, and time variations of 1.5 K and −0.6 g kg −1 per hour for the same variables. Figure 3 shows a good fit between the numerical simulation and the aircraft measurements for the profiles of θ and q at the two times when the soundings were made above  Niamey. The growth of the ABL is the same in both the observations and the numerical simulation. According to the latter, the ABL grew from a 400 m height at 0900 UTC to 2400 m at 1700 UTC, with significant drying ( 6 g kg −1 ) and warming ( 6 K). At the flight level (300 m a.g.l), however, the discrepancy between observed and modelled q in Fig. 3 is not negligible ( 1.2 g kg −1 ) and could be seen generally in the horizontal plane, despite a good fit for gradients. We also compared the fluctuation distributions of the meteorological variables for scales ranging from 250 m to 5 km. Figure 4 shows the distributions of w and q at these scales observed along the twenty-four 30 km-long legs of the ARAT aircraft and those found in the modelled 2D field at the same height averaged over the duration of the flight. In this figure, a low-pass 5 km Butterworth filter is applied to the aircraft measurements before comparison with the model. The observed distributions fit well with the averaged distributions simulated over the duration of the flight. A similar good match between model-and observation-derived distributions was already found by Couvreux et al. (2005) over the American Great Plains in late spring during IHOP (International H2O Project, Weckwerth et al. 2004). The modelled and observed variances of w and q are close to each other, with a large variability in q that is typical for this season (Said et al. 1997). The striking negative skewness of q observed with aircraft measurement (−0.9) is also found in the model (−0.6), but to a lesser extent. This negative skewness is very likely due to dry-air descent from the free troposphere into the ABL through entrainment, as shown in Fig. 10. During HAPEX, it was observed that the more negative the skewness, the larger the difference between aircraft-measured fluxes and surface-measured fluxes. (See Lenschow et al. 1994, who quantified the effect of skewness on both random and systematic errors.) The shoulder around 1.25 m s −1 in the modelled distribution of w may be explained by the inability of the simulation at this resolution to fully resolve the cores of the updrafts so Distribution of (a) q and (b) w observed with the aircraft (bars) and found in the modelled 2D field at the same height averaged over the flight duration (solid line). For aircraft data, all 30-km segments are considered. Only scales ranging between 250 m and 5 km are considered for both types of data that maximum velocities artificially accumulate near smaller values. It therefore leads to a slightly underestimated skewness. Couvreux et al. (2005) did not find this problem at higher resolution (100 m).

Error analysis
A dozen cases were investigated where vertical flux profiles were obtained from several legs flown at varying heights within the ABL. The linear profiles obtained for the sensible heat flux were used to estimate surface sensible heat flux by extrapolation to the ground (F. Saïd 2005, personal communication). Latent heat flux profiles did not allow this method because of large random errors, so surface flux was estimated from the lowest leg. Large differences (up to 40%) between aircraft estimates and ground station estimates were observed, with more probability of observing smaller flux with aircraft. Grunwald et al. (1998) reviewed the possible sources of errors when estimating surface fluxes using aircraft measurements: (i) filtering (usually an underestimate), (ii) flux divergence. During HAPEX, the organised structures were suspected to also increase errors in flux measurements, along with dry intrusions (linked to entrainment) penetrating deeply enough to make the ABL less homogeneous than expected.
The case shown here appropriately exemplifies this issue, since organised structures were observed and well-modelled. For this case, latent heat flux estimates were very uncertain due to large scatter (50%) and did not allow the use of the flux divergence method. Both heat flux estimates were approximately 20% smaller than those observed at the ground. Error due to filtering was not larger than 10%.
In the following, we define two different scales: the integral scales that allow us to estimate the systematic and random errors of flux estimates, and the characteristic scale length of the coherent structures when they occur.

Integral scales
The integral scale l w of w, which is a measure of the length over which w is relatively well correlated with itself, is defined from the autocorrelation function ρ w (r ) as: where r is the displacement and R w (r ) is the autocovariance funtion of w: Here angle braces denote an average over the length of the leg. Note that the autocorrelation function ρ w (r ) is the autocovariance function normalised by the variance σ 2 w = R w (0). A good estimate of l w can be obtained from the maximum of the running integral of (1) (Lenschow and Stankov 1986): which is reached at the first zero crossing of ρ w (r ). The integral scales of θ and q are similarly estimated with the autocorrelation functions ρ θ (r ) and ρ q (r ), respectively. We use the approximation found by Lenschow et al. (1994) to estimate the Eulerian integral length scale for covariance between w and the scalar s (here s = θ or q) because it is difficult to estimate this length scale when the covariance is small. Figure 5 shows the vertical profiles of l w , l θ , l q , l wθ and l wq found with the measurements of both aircraft. In this figure, we consider only the 10 collocated segments that were flown at the same time by the two aircraft (a10 to a19 and m10 to m19). For each segment, the flight height and the integral scale are normalised by the ABL depth z i . z i was previously interpolated at each segment using the temporal linear trend deduced from the two soundings of the ARAT aircraft and the meridional gradient of z i found in the numerical simulation. Profiles of measured l w and l θ fit well with the empirical law found by Lenschow and Stankov (1986): and where z * = z/z i . However, we find much larger l q than expected by (6), especially close to the surface. Dry intrusions penetrating down to the surface may explain this discrepancy (Mahrt 1991), which is also consistent with the negative skewness mentioned previously. The Eulerian integral scales l wθ and l wq in particular are larger (20-50%) than the empirical law found by Lenschow and Stankov (1986), especially close to the surface: 3.2 Systematic and random flux errors Due to the finite sample length in time and space, there is both a systematic and a random error in the measurements. Both depend on the ratio of integral scales to the sample length L.
According to Lenschow et al. (1994), the (relative) systematic error on the turbulent vertical flux of the scalar s is: where F is the true flux and F(L) is the ensemble average flux measured from samples of length L. In our case, we find the systematic flux error smaller than 1% for the lower height measurements (Merlin IV) and smaller than 1.5% for the ARAT measurements. The random error is usually larger than the systematic error. Lenschow et al. (1994) estimate this error as: where r ws is the correlation coefficient of w and a scalar s: We find a 15% random error for sensible heat flux for the Merlin IV aircraft measurements and 25% for the ARAT measurements. Thus we find about the same amount of random error as the discrepancy between aircraft and ground estimates. Consistent with the larger integral scales l q and l wq , we find twice as large random errors for latent heat fluxes. Lenschow et al. (1994) found that, in contrast to the systematic errors, the random errors are significantly affected by non-zero skewness. According to their study, the same relative random error in the second-order moment requires an average twice as long for a time series with skewness = −0.8 compared to zero skewness using their model for introducing skewness into a Gaussian process.

Coherent structure wavelength
An organised ABL shows periodic oscillations in w, θ and q. These oscillations are usually at scales several times larger than the integral scales, so that they are obvious in the autocorrelation function. We consider a wave w 2 = a 2 sin(2πr/λ), where r is the distance and λ the wavelength of the wave superimposed onto a turbulent sample of w (e.g. Fig. 6a).
As an example, Fig. 6b displays the autocorrelation function Rw(r )/Rw(0) that results from the addition of a wave w 2 of amplitude a 2 = 0.5 m s −1 and wavelength λ = 2.3 km to the fluctuations w measured along segment m1 (σ w = 0.75 m s −1 ). In this case, and for λ L, the autocovariance function of w 2 is periodic: and ρw(r ) has a slower decrease in autocorrelation at small separation distance than ρ w (r ) and a periodicity imposed by the wave. This shows that a periodic signal of a few km scale can be detected by a periodicity in the autocorrelation function although it is not so obvious in the time series (see Fig. 6a). It also shows that the estimates of integral scales based on the rate of the decrease in the autocorrelation at small separation distances, are increased by coherent structures of scales several times the integral scale of the series without coherent structures. In this example, lw 2 l w . Another example is given in Fig. 6c, where we consider a measured sample of q for which R q does not show any obvious periodicity (segment m6) and another sample with an obvious oscillation (segment m16). We find that we can simply model R q for segment m16 by superimposing on segment m6 a wave with λ = 2.8 km and a 2 σ q (σ q = 0.19 g kg −1 , σq = 0.21 g kg −1 ). The integral scale for segment m16 is 50% larger than for segment m6.

Evidence
The numerical simulations at 500 and 250 m resolutions showed very coherent structures in the morning at the flight height of the ARAT aircraft, organised in lines parallel to the flow direction (so-called convective rolls) at 1000 UTC, evolving to symmetric cells around noon. These "cloud-free streets" or "blue streets" observed at 1000 UTC are shown for q in Fig. 7 for R500 m and R250 m. For both, the banded structures are about 2 km apart, but for R250 m, they become more complex. Figure 8 displays the fluctuations in w, θ , q for R250 m at 1000 and 1200 UTC, after removing a 2D plane, at 1000 and 1200 UTC. The rolls at 1000 UTC are aligned with the mean wind, as previously observed by e.g. Weckwerth (1995). At 1200 UTC, the structures changed to cells (see Fig. 8), with a smaller mean ABL wind (3 m s −1 at 1200 UTC versus 6 m s −1 at 1000 UTC). Weckwerth (1997) found that −z i /L o (L o is the Obukhov length) determines the type of structures observed (also see Weckwerth 1999). She observed rolls with −z i /L o ≤ 25. In our study, −z i /L o increases throughout the day and becomes larger than 25 between 1100 and 1200 UTC, which is consistent with observed transformation to cellular structure.  Figure 9a displays the averaged modelled R w (r ) along east-west (black) and north-south (grey) axis over the entire R250 m domain at 1000 UTC. Rolls are responsible for the second maximum in both and for the difference in wavelength. The wavelength in the east-west direction is 4 km and 2.4 km in the north-south direction, as suggested by the orientation of the rolls in Fig. 8a. At 1200 UTC, these two functions merge and give the same characteristic wavelength in both directions, as suggested by Fig. 8b.
Even though the ARAT R w , R θ and R q are not as regular as the modelled R w , R θ and R q at 1000 UTC, they are consistent with organised structures at flight height. We find that about half of the 24 segments flown by the ARAT aircraft show significant organisation of w, θ and q, including six segments with an auto-correlation coefficient at the second maximum λ m 2 greater than 0.2 (which we use as a criterion for the existence of rolls). λ m 2 ranged from 1 km to 4 km along both flight directions. The structures were observed more often in the second part of the flight, so that we cannot prove that the rolls were already present around 1000 UTC as found in the numerical simulation. The segments flown within the small domain of R250 m around 0930 UTC do not show different λ m 2 in the two directions and they are too few to be representative. However, we found evidence of asymmetry in the structures later in the flight. Figure 9b shows the averaged ρ w over three north-south segments (a13-a15) and three east-west segments (a16-a18), which were flown around 1130 UTC. They show different wavelengths along the two axes (2 and 2.6 km, respectively), but the difference is smaller than at 1000 UTC in R250 m, consistent with a time closer to 1200 UTC, when the modelled structures have evolved to symmetric cells. Averaging ρ w , ρ θ and ρ q over all the segments removes the asymmetry, indicating that for most of the flight, the structures were mainly cellular, with a length scale similar to the model length scale (2-3 km). Possible explanations for the discrepancy between model and observation are (i) a different direction of the rolls in the real world, (ii) a different time of occurrence or (iii) a more complex field of structures.
We also observe a general asymmetry of the integral scales of w and θ , as shown in Fig. 5, where squares and circles distinguish one direction from another. On average, l w and l θ are smaller in the north-south direction, but l q does not show any asymmetry. This suggests that observed thermals and coherent structures (correlated fluctuations in w and θ ) may be elongated and that q behaves differently. Intrusions of drier air from the free troposphere down into the ABL likely complicates the organisation of q and of its covariance with w more than θ and w.  Figure 10a clearly shows dry intrusions (descent) from the top of the ABL. These intrusions, which entrain free-tropospheric air into the ABL, are characterised by much drier and warmer air and are responsible for the negative skewness of q. They are also clearly seen in the numerical simulation and penetrate deeply into the ABL, as shown in Fig. 10b. Figure 11 shows the vertical profiles of the mean horizontal wind measured with the aircraft around 1000 and 1200 UTC. The direction of the wind did not change significantly between 1000 and 1200, but the wind speed decreased by half within the monsoon flow, in broad agreement with the typical diurnal cycle of the ABL wind speed over the area (Parker et al. 2005). The (mostly directional) shear zone between the moist south-westerly monsoon flow and dry northeasterly Harmattan flow remains at the same height, while the ABL grows from 600 to 1000 m in height. Thus, the ABL is developing within the monsoon flux. Its top usually reaches the height of the shear zone between Harmattan and monsoon flows at midday or at the end of the morning, thus favouring dry intrusions and a positive feedback to the ABL growth during the afternoon.
In the numerical study, intrusions are more obvious around 1200 UTC than around 1000, consistent with the growth of the ABL between 1000 and 1200 UTC when the ABL is deep enough to reach drier layers and the shear zone (see Fig. 11), as observed by Wai et al. (1997), so that the entrainment of free-tropospheric dry air into the ABL is favoured. It is also consistent with the results of LeMone et al. (1999) who found a relationship between shear and rapid growth of the ABL.

Role of structures in fluxes and impact on measurements
Organised structures of scale λ 3 km can make the 1D aircraft measurements suspect for at least two reasons: (i) the scale of the coherent structures may be too large for them to be well sampled by a 30 km leg, (ii) the sample obtained by the aircraft might be biased depending on the flight direction and location. First, we estimate the contribution of the structures and dry intrusions to determine which part of the flux is missed if the structures are not well sampled. Second, we estimate the effect of a 1D approach relative to the two-dimensional organisation of the fluxes. Figure 12 displays the model covariance field obtained from the fluctuations shown in Fig.  8. It shows that the geometry of the structures affects the correlation field (similar line-typed patterns are found), and consequently the resulting fluxes.  Figure 13a shows the cross-correlation function normalised by the covariances between w and θ and between w and q averaged over the entire domain of R250 m at 1100. Figure  13b shows the two same functions obtained from the average of the 24 segments flown by the ARAT aircraft. In both the model and the observation, we observe a second maximum at a similar length of about 3 km. We also find that w and θ are less affected by the coherent structures than w and q, since the first minimum and second maximum of ρ wq are more emphasised. Note that the steeper slope of the observed ρ wθ and ρ wq is due to the fact that scales smaller than 250 m are not taken into account by the model. Figure 14 shows the effect of filtering scales larger than 1, 2, 4 and 5 km respectively on the flux estimate. The calculation of the sensible and latent heat flux profiles shown here was made from the 2D fields of fluctuations at each level of the simulation, after removing the 2D running mean. It demonstrates that when the scales >2.5 km are not included in the flux computation, the total flux is not retrieved, whereas filtering at cut-off scales larger than that does not make any change. Figure 14 also confirms the larger scale of water vapour transport. In order to estimate the relative contribution of the dry intrusions to the heat fluxes, we discriminated the dry intrusions (relative negative w and negative q) from the thermals (relative positive w and positive virtual potential temperature θ v ) from the remaining contributions to the total flux, based on criteria used by e.g. Berg and Stull (2004) or Williams and Hacker (1992). This conditional analysis defines three classes. Class 1 is for thermals (w > 0, θ v > 0), class 2 is for dry intrusions (w < 0, q < 0) and class 3 for the rest. Figure 15a displays the probability of each class (number of datapoints relative to the total number), while Fig. 15b, c display the contribution of each class to the total covariances w q and w θ respectively. For each class and diagram, we compare the ARAT aircraft measurements (all segments are considered) to what we find in R250 m at 300 m a.g.l, at 1000 and 1200 UTC. Note that the aircraft data are not low-pass filtered at 250 m wavelength. For both model and observation, we find a similar probability of observing dry intrusions or thermals of about 30% for both. We find a very large contribution of thermals in w θ , and a similar contribution of dry intrusions in w q . Between 1000 UTC and 1200 UTC in R250 m, there is a stronger contribution of fluctuations that are neither thermals nor dry intrusions. This may be due to a less organised ABL in late morning, as rolls evolve to cells. Comparing observation and model, we find a relatively similar distribution over the three classes. Figure 15a shows that class 3 is larger in the observations, which decreases the resulting fluxes. One reason for this is that the subgrid flux is not taken into account in the model estimates, while part of the observed datapoints are randomly distributed over the eight (= 2 3 ) classes that result from the combinations of the sign of the fluctuations of the three variables. In our simulation, only 5% of the total flux is not resolved in R250 m at 300 m. At 60 m where the Merlin IV flew, the contribution of the subgrid flux is about 20%. We also generally find a larger contribution of subgrid sensible heat flux than subgrid latent heat flux at the height of the two aircraft, which is consistent with larger l wq than l wθ as shown in Fig. 5. Inappropriate parameterisation of subgrid fluxes may lead to overestimated fluxes. There is also a possibility of underestimating the turbulent fluxes in R250 m because of unresolved maxima in fluctuations (see, for instance, the shoulder of the distribution of w in Fig. 4). Another reason for the difference between aircraft and model in Fig. 15 is linked with sampling errors.
As observed previously by Mahrt (1991) and Couvreux et al. (2007), we find negative skewness associated with the dry intrusions that contribute to the latent heat flux. Figure 16a displays the ratio of the sensible heat flux to the latent heat flux B * = H/L E as a function of the skewness of q, using the ARAT and Merlin IV aircraft measurements and R250 m. Note that B * is not rigorously a Bowen ratio, since the considered fluxes are not surface fluxes. We find larger negative skewnesses for smaller B * .
To quantify the link between dry intrusions and ABL growth, we relate the relative contribution of dry intrusions to entrainment velocity w e and convective velocity w * . and where subscript s indicates the surface and subscript i indicates the ABL top. wθ v is the buoyancy flux at the surface, θ v i the jump in θ v across the ABL top, and θ v m the mean θ v within the ABL. For both model and observations, the buoyancy fluxes at the surface and at the ABL top are deduced by extrapolating the profile linearly. For these observations, we have only one estimate of w e and w * using the soundings and all segments flown by the two aircraft. For R250 m, we can estimate w e and w * at several times, in particular later in the day. Figure 16b displays the contribution of dry intrusions to the latent heat flux divided by the surface latent heat flux as a function of the entrainment velocity w e . We find increasing contribution of the dry intrusions with increasing entrainment velocity, which is consistent with the process of entrainment associated with intrusions of dry tongues within the ABL (Couvreux et al. 2007). It is also consistent with a growing and warming ABL, confirmed by the increasing contribution of the dry intrusions with increasing w * (not shown). The vertical flux measured from 1D series of fluctuations may depend on the direction of the leg, because the reference mean removed from the signal can be higher or lower than the 2D field average, and the resulting fluctuations smaller. This may happen in particular if the leg is parallel to the longitudinal axis of the structures. In order to estimate the systematic and random errors at flight level due to the 1D approach relative to a 2D approach at 1000 UTC, we consider smaller domains in R250 m, that are 30 km squared (120 points × 120 points), centred every 2.5 km of the grid and from which we deduce a 2D covariance and 1D covariances of 30 km length along two directions. Figure 17 displays 2D covariances w θ and w q as functions of 1D covariances taken in two different directions at 1000 UTC: 60 • and 150 • clockwise from north, that are respectively parallel and transverse to the simulated rolls. One striking result is a random error that can be three times larger with the 1D approach than with the 2D calculation. Also, this random error varies with the direction: it is larger for a 1D computation along the roll axis. According to Wyngaard (1983), the ratio of the random error in 2D to that for 1D is 2 l ws L where L is the length of the side of the sampled square or of the sampled line respectively (L = 30 km here). The integral scales l ws found in the model are larger than in the observations, because of the unresolved scales. For both heat fluxes, we find integral scales of about 1100 m along the roll axis and 500 m transverse to it. This results in a scatter in 1D that is 11 times larger than the scatter in 2D in the parallel case and seven times larger in the transverse case, which is about what we observe in Fig. 17.
We find respectively about a 10 and 25% underestestimate of the 1D sensible and latent heat fluxes at 60 • , that is, along the roll axis. In contrast, the systematic error in heat fluxes between 1D and 2D estimates is <2% for a line transverse to the rolls. According to the theory (Eq. 8), this systematic error should be 3% for a line across the roll axis and 7% for a line along it. The latter is significantly smaller than what we actually find, which is likely due to the fact that the theoretical expression is valid for random turbulence only and in reality the flow contains non-random structures.
Finally, we extend this systematic flux error calculation to the entire depth of the ABL, for the same time (1000 UTC), as displayed in Fig. 18 Fig. 17 (a, b) 2D w θ covariance and (c, d) 2D w q covariance versus their corresponding 1D covariance along 30 km lines oriented with (a, c) 60 • heading (b, d) 150 • heading, from R250 m at 1000 UTC. 60 • is the direction of the rolls, while 150 • is perpendicular to them. For each cloud of datapoints, the average is represented with a large and thick circle, along with the standard deviations (solid line) for both x-axis and y-axis the 1D line oriented at 60 • , in the direction of the rolls (about 10% underestimate for the sensible heat flux and 25% underestimate for the latent heat flux). In the transverse direction, the latent heat flux is slightly underestimated (<5%) and the error for sensible heat flux ranges from −5% to +5%. Later in the day at 1300, when rolls change to cells, the relative underestimation of sensible heat flux along the 60 • direction relative to the 2D estimates decreases (not shown) and the differences between 60 • and 150 • are smaller: the sensible heat flux is underestimated by 3% along 150 • line and 6% along 60 • . The latent heat flux is underestimated by 10% along the 150 • line and 16% along the 60 • line. The fields of fluctuations are quite heterogeneous at this time, with cells of the same characteristic size only in 2/3 of the domain and no clear organisation elsewhere. The analysis is thus more difficult at these subsequent times.
In the case studied here, it is likely that the structures were not as simple as in the model and that the legs were not flown precisely along or transverse to the rolls. Therefore, the systematic error in measured fluxes due to this effect is probably smaller than 20% and closer to 10%.

Concluding remarks
During the HAPEX-Sahel experiment, the observations made with research aircraft over the Niamey area often revealed an organised ABL, with a few kilometre scale structures (that is a few times the ABL depth) associated with both thermals and dry intrusions. A numerical simulation at high horizontal resolution that was able to resolve these structures, evolving from rolls in the early morning to cells by noon. The simulation showed similar fluctuations of meteorological variables as those observed by the aircraft. In particular, rolls, thermals and dry intrusions were found in both. Dry intrusions were found to penetrate deeply into the ABL, sometimes even reaching the surface, which is manifested in the observed highly skewed q fluctuations. Dry intrusions, rolls and thermals organised at a few kilometre scale add to the heterogeneity observed by the aircraft with detrimental effects on the flux measurement accuracy. In particular, these structures can explain the larger random errors observed in that case. We also found from both observation and numerical simulation that the structures in q have more effects on random errors and contribute more to the latent heat flux than to sensible heat flux.
Even though spatial probing of an aircraft in an organised field may lead to smaller error than fixed point measurements from a ground station (Kanda et al. 2004), making accurate airborne measurements of fluxes in such a medium remains a challenge. We showed that the random error with a 1D approach is large relative to the 2D-computed flux, and that the estimates of fluxes from aircraft measurements with the presence of rolls in the ABL are considerably improved when the aircraft flies in the transverse direction compared to legs parallel to the roll axis. An extension of the study by Kanda et al. (2004) with the high horizontal resolution of LES should tell us more about the bias and random errors as a function of the scale, type, and direction of the structures, and allow us to optimise aircraft flight tracks.
This study suggests that stacked legs longer than 50 times the ABL depth are necessary to obtain accurate measurements of fluxes in the ABL. However, this is difficult to obtain with a single aircraft in a reasonable time. A possible complement to this would be a clear-air radar or a lidar capable of remotely detecting these coherent ABL structures.