Water vapor density and turbulent fluxes from three generations of infrared gas analyzers

Fast-response infrared gas analyzers (IRGAs) have been widely used over three decades in many ecosystems for 10 long-term monitoring of water vapor fluxes in the surface layer of the atmosphere. While some of the early IRGA sensors are still used in these national and/or regional eco-flux networks, optically-improved IRGA sensors are newly employed in the same networks. The purpose of this study was to evaluate the performance of water vapor density and flux data from three generations of IRGAs – LI-7500, LI-7500A, and LI-7500RS (LI-COR Bioscience, Inc., Nebraska, USA) – over the course of a growing season in Bushland, Texas, USA in an irrigated maize canopy for 90 days. The energy balance ratio, which is the 15 sum of turbulent fluxes divided by the sum of surface available energy, was used to assess systematic biases of the IRGA sensors for evapotranspiration (ET). Water vapor density measurements were in generally good agreement, but temporal drift occurred in different directions and magnitudes. Means exhibited mostly shift changes that did not impact the flux magnitudes, while variances of water vapor density fluctuations were occasionally in poor agreement, especially following rainfall events. LI-7500 variances were largest compared to recent LI-7500RS and LI-7500A results manifesting in widened cospectra, 20 especially under unstable and neutral static stability. Agreement among the sensors was best under the typical irrigation-cooled boundary layer, with a 14% interinstrument coefficient of variability under advective conditions. Generally, the smallest variances occurred with the LI-7500RS, and high-frequency spectral corrections were larger for these measurements resulting in similar fluxes between the LI-7500A and LI-7500RS. Fluxes from the LI-7500 were best representative of growing season ET based on a world-class lysimeter reference measurement but using the energy balance ratio as an estimate of systematic 25 bias corrected most of the differences among measured fluxes.


Introduction
The eddy covariance (EC) method is a standard way to monitor water vapor flux between the surface and atmosphere at most spatial scales and environments, including marine (Honkanen et al., 2018;Takahashi et al., 2005), forest (Novick et al., 2013;and Paw, 2011). In water-limited regions, the need to conserve a subsurface source, such as the U.S. Ogallala Aquifer, serves as motivation for agricultural producers to estimate the crop water use for daily irrigation scheduling (Xue et al., 2017). Current crop production involves innovative water saving measures, such as variable rate irrigation management, requiring high quality evapotranspiration (ET) data to supplement efforts to calculate the correct amount of water to apply to crops (O'Shaughnessy 35 absorption of water vapor from a dewpoint generator over a range of temperatures from 17 to 41°C. Based on the manufacturer calibration and re-calibration sheets (after a certain period the IRGA is returned to the manufacturer for re-calibration), the span drift is primarily a function of temperature, whereas the zero drift is chiefly influenced by the measurement range of water vapor density.
In addition to the IRGA, a sonic anemometer is necessary to determine water vapor flux. This pair of instruments 70 introduces systematic error due to their physical separation, which is a source of high frequency turbulent signal loss (Massman, 2000). The magnitude of flux attenuation is enhanced by lighter wind speed and a greater ratio of horizontal separation to sensing height (Horst and Lenschow, 2009). The expected cospectra, or eddy flux in the spectral domain, can be estimated analytically with a series of transfer functions (Massman, 2000;Moncrieff et al., 1997) that account for signal loss at low and high frequencies. A spectral correction factor can often be determined based on how this modeled cospectrum 75 departs from the measured cospectrum, indicating the degree of flux loss for a given observation period and EC system.
To address offset errors of water vapor density from an IRGA, data are typically compared to another type of sensor. In a comparison to the enclosed-path EC155 system (Campbell Scientific, Logan, UT, USA), errors in water vapor density were generally between -3 and 3 g m -3 (Novick et al., 2013). Such errors were largest in early to mid-morning hours coinciding with Water vapor density data among the three infrared open-path IRGAs were compared in a fashion similar to Mauder et al. (2006). The following characteristics of variance ( " ′ " ′ '''''''' ) and covariance ( ′ " ′ '''''''' ) were of interest: regression intercept (a), slope (b), and coefficient of determination (r 2 ); root mean square deviation (rmsd); and bias (d). Comparability between LI-7500RS and the other two models was found using rmsd, defined as: where 2,3 is the i th observation for the LI-7500/A and 56,3 is the i th observation for the LI-7500RS. Interinstrument variability was also determined, which is like rmsd except the mean of the EC systems is the reference value. For fluxes, interinstrument variability was expressed relative to flux magnitude using the coefficient of variation (CVI-I). Implausible values of 20 Hz data, defined as greater than 30 g m -3 or less than 2 g m -3 , were removed prior to taking half-hourly means. Additionally, while the LI-7500A and LI-7500 were calibrated in 2014 and 2015, a correction to these data was made based on a factory calibration 135 after data was collected. Otherwise, no additional conditioning was performed on the raw data. Given the interest in sensor sensitivity, comparisons were also made between collocated HMP155A and IRGA(s) at each tower, which were assumed to be sensing identical air parcels containing equal water vapor density.
To further ascertain the performance of IRGAs, (co)spectral density of " ( " ) measurements were calculated for each of three EC systems using Welch's periodogram method (Blanken et al., 2003). The distribution of power across frequencies, 140 particularly signal loss at high frequencies, can indicate differences in flux characteristics with an expectation that latent heat would be underestimated. Of particular interest are results from an advective environment in which high frequency variation is enhanced (Prueger et al., 2012). Data were conditioned by linear detrending on half-hour (36,000 points) segments (Zhang et al., 2010). Spectral density ( : ; ) was calculated across these segments with a Hamming window length of 360 and overlap of 180 observations. Then the spectra were averaged into 100 evenly spaced bins on the logarithmic scale. The same procedure 145 was repeated for the cospectra of vertical velocity and water vapor density, indicating the behavior of water vapor flux in the spectral domain. Finally, ogives were calculated to summarize differences in cospectra across wavelengths by integrating the cospectra from low-frequency energy to high-frequency energy on a scale from 0 to 1. The (co)spectra and ogives were multiplied by the frequency and normalized by mean (co)variance to make the data dimensionless.
After examining raw variances and covariances, water vapor fluxes (E) were processed using Eddypro (v6.2.0) software 150 (LI-COR Bioscience, Lincoln, Nebraska, USA) for half-hour averaging periods when availability of data exceeded 90% ( and " were recorded for at least 32,400 of 36,000 possible observations). Prior to computing fluxes, a statistical screening of time series data was implemented. Spikes were detected using the median absolute deviation for each half-hour (Mauder et al., 2013) and replaced with the half-hour mean of non-outlier observations. Then data was detrended by block average and corrections were made to account for sensor separation, tilt of the sonic anemometer via double rotations (Fratini and Mauder, 155 2014), and spectral energy loss in both low (Moncrieff et al., 2004) and high (Moncrieff et al., 1997) frequency ranges. Based on spectral losses and other corrections, E was calculated iteratively. The original water vapor flux was multiplied by the spectral correction factor of < " < '''''' before adding WPL density fluctuation terms (Kaimal and Finnigan, 1994). Sensible heat (H) was then corrected for humidity effects that arise from using sonic temperature in place of air temperature (Van Dijk et al., 2004). Finally, this corrected H was multiplied by its spectral correction factor, and the WPL term was added to the corrected water vapor flux to create a final E or λE. Approximately 13.5% of available data were removed through results of steady-state and fully developed turbulence tests (Mauder and Foken, 2004). The acquisition ratio of each half-hour was obtained by dividing the count of non-filtered fluxes by the maximum number of observations (Kim et al., 2015).
Intercomparison of λE and its systematic error (d) and random uncertainty (e) components was conducted on half-hourly and daily timescales. The measured λE is assumed to be the difference between the actual flux and these errors (Lasslop et al., 165 2008). Systematic error can be evaluated in the context of surface energy balance, such that d is zero when turbulent flux equals the available energy measured through solar radiation, ground heat flux, and heat storage during a given period (Mauder et al., 2013). The estimate of systematic error is then and 170 where the terms in the numerator are independent (H is sensible heat flux, and is latent heat flux) for each EC system and those in the denominator are shared among the EC systems. J was calculated as the sum of soil and photosynthesis heat storage since the other components of heat storage contribute negligibly to instantaneous energy balance in this ecosystem (Kutikoff et al., 2019). Random error associated with sampling was quantified with the method of Finkelstein and Sims (2001), which 175 calculates the variance of the covariance using the raw timeseries data for each averaging period. Together, error quantification can indicate if half-hour fluxes from the three EC systems statistically differ for half-hours in which turbulent flux measurements are reliable.
Water vapor flux was compared using the equivalent total water depth ET for daily totals. Gap filling, following Reichstein et al. (2005), was done for half-hours that were flagged for any of the three EC systems based on steady-state and developed 180 turbulence tests (Mauder and Foken, 2004), occurrence of precipitation, and high relative humidity (RH > 95%). Total gapfilled ET was close to the sum of the half-hour observations, with approximately a 3% greater flux for each EC system. Flux accuracy of the three EC systems was assessed in relation to a large weighing lysimeter, which has an accuracy of 0.05 mm hr -1 (Evett et al., 2012b). Located within 30 m of the EC system, lysimeter ET was computed using a soil water balance approach from a subsection of the same field. Briefly, the mass change of water measured by the weighing lysimeter was 185 calculated and converted into a flux based on the surface area of the lysimeter and density of water. Description of the lysimeter data can be found in Moorhead et al. (2017).

Results
The findings of the study are presented in three subsections, including water vapor density mean and fluctuations, spectra and cospectra, and fluxes. All were influenced by irrigation and precipitation events. Water added to the field included 498 mm 190 from 33 separate subsurface drip irrigations (SDI) (Evett et al., 2019) and 238 mm of precipitation (Evett et al., 2018), consistent with an average growing season (Gowda et al., 2009;Tolk et al., 2013). However, much of that rainfall (88%) occurred after 1 August, and combined with crop maturity, eliminated the need for irrigation after 18 August.
Data filtering also impacted all comparisons. After all threshold and precipitation screenings, 3,577 out of a possible 4,320 half-hour observations are available for analysis. The acquisition ratio was comparable to similar studies (Wu et al., 195 2015). Between 9:00 AM and 9:00 PM (LST), the ratio exceeded 92%, whereas EC system issues reduced availability in the predawn hours to as low as 61% for the half-hour ending at 7:00 AM (Fig. 1).

Water vapor density validation
The long-term zero drift of water vapor density for the three IRGAs was evaluated as the three-month change in bias D " . As the study period began, the reference value of water vapor density ",M ranged from 3 to 18 g m -3 . Accordingly, the measured 200 values ",$ for the LI-7500 and LI-7500RS were biased low and the LI-7500A was biased high. After applying the postcorrection to the LI-7500 and LI-7500A data, all ",$ were between 0.11 and 1.31 less than ",N (Fig. 2). At the end of the study period, all IRGAs clearly showed an increased bias relative to the HMP155. Interestingly, the LI-7500 and LI-7500A had moved towards larger values, whereas the LI-7500RS moved towards smaller values (Fig. 2). That resulted in the LI-7500 8 To investigate the unexpected large drift exclusive to the LI-7500RS on DOY 191, biometeorological data were assessed.
Light southerly winds and moderately humid conditions were observed when Δ ",56 increased from -0.96 to -2.45 between 8:30 and 9:30 PM LST. While nothing unusual occurred meteorologically, a 3°C drop in temperature and 10% increase in RH 225 accompanying the loss of daytime heating was noted. It was instructive to look at the variation in RH as estimated using vapor and ambient pressure from the IRGAs and sonic temperature from the CSAT3. While the magnitude of RH did vary slightly among the sensors, the increase in RH was similar for the LI-7500 and LI-7500A while being less than half for the LI-7500RS.
In the hours immediately prior and after, the slopes of Δ ",56 among the IRGAs and HMPs are nearly in lockstep. Unlike other deviations that exist on a subdaily timescale, this new offset continued until DOY 197.
Step changes are a dominant feature in 230 the linear regression between ",OP/2 and ",OP56 .
Differences between the means and fluctuations of " are summarized in Fig. 4 as a function of day of year. Since variance of the " time series reflects the mean of squared fluctuations " < 8 '''' , greater variance in the half-hourly data reflects larger fluctuations " < . While the LI-7500 tended to have consistently greater "

Spectra and cospectra
Since the three analyzers had the same specifications and were configured to measure turbulence in the same fashion, any deviations in spectral characteristics would be an indication of possible drift. Returning to the distinct LI-7500RS error on 255 DOY 191, spectra were examined during the interval from 8:00-9:30 PM (LST), which consisted of three spectra corresponding to consecutive flux averaging periods. Overall, as evident from Fig. 5a-c, the shapes of spectra were in close agreement during the daytime, whereas the nighttime peak frequency was shifted to lower frequencies indicating the predominance of large eddies after sunset. At 8 PM, the three spectra were nearly identical and matched the predicted -2/3 slope (Fig. 5d). In the following hour, the spectra of the LI-7500A and LI-7500 remained nearly identical, whereas the LI-7500RS spectra were 260 greatly modified. Based on the 20 Hz timeseries, air humidity began to decrease suddenly at roughly 8:40 PM in concert with a doubling of fluctuation amplitude. As the other two IRGAs and HMPs continued to indicate increasing air humidity, ",STJOPUU56 steadily rose for nearly one hour until ",STJOPUU56 again agreed with the other instruments. Because only the averaging period between 9 and 9:30 PM is affected by increased variance water vapor, the spectrum corresponding to that half-hour is the period with a shift towards higher frequencies. 265 Cospectra were viewed through the lens of atmospheric stability because it predicts their shape according to Monin-Obukhov similarity theory (Kaimal and Finnigan, 1994). For all cospectra, the LI-7500 tends to have greater energy in the production and dissipation spectral regions while being nearly identical in the inertial subrange, and these differences translate into higher latent heat fluxes (Fig. 6). Lower frequency components of flux were clearly greater, especially in unstable and neutral conditions, as observed by the LI-7500 (the oldest version), compared to the LI-7500A and LI-7500RS. While the two 270 newer sensors exhibited similar behavior and relatively smaller fluxes than the LI-7500, under unstable conditions the LI-7500RS showed a difference in performance from the LI-7500A at high frequencies. For all three IRGA, co-spectra dipped at 2.5 Hz, which should not occur in any desired instruments (Kaimal and Finnigan, 1994). Strong turbulent motions were likely captured more by the LI-7500A within the surface layer. These cospectra were shifted towards lower frequency compared to those in neutral and stable conditions, favoring larger eddy sizes with a smaller percentage of energy accumulated in the inertial 275 subrange (Fig. 6b). This middle frequency range is where the IRGAs were most similar. Regardless of sensor, unstable conditions featured a flattened peak and more energy towards lower frequencies, as expected for various scalar fluxes measured with the same instrumentation (Wolf and Laca, 2007). However, in an irrigated cropland environment, the surface layer is prone to become stable more often than the surrounding area due to a temperature inversion forced by the relatively wetter, cooler canopy. A previous study demonstrated this effect by using simultaneous sensing over adjacent irrigated cotton and 280 non-irrigated winter wheat fields, where energy production as depicted by : ; was two orders of magnitude smaller for the irrigated field than the non-irrigated field (Prueger et al., 2012). Accordingly, in the present study, variability among cospectra was small under these conditions with relatively few large eddies (Fig. 6e). In contrast, under neutral and unstable conditions, the LI-7500 departed largely from the other two sensors with energy contribution from low frequency eddies. https://doi.org/10.5194/amt-2020-302 Preprint. Discussion started: 27 August 2020 c Author(s) 2020. CC BY 4.0 License.

Water vapor fluxes 285
For much of the study period, lE from the LI-7500RS and LI-7500A were similar with slightly larger magnitude than the LI-7500. Overall interinstrument variability CVI-I of lE was 20%, about that of the underlying water vapor variance, and errors on average were less during daytime hours than nighttime (Table 2). For an average diel cycle, the largest CVI-I occurred during the middle of the night, rapidly declined after sunrise, reached its smallest value of 10% at 4 PM, and then increased at a relatively slow rate after sunset. On a seasonal basis, there was a slight, nonlinear increase in CVI-I over time, with mean 290 values increasing from approximately 16% to 24%. Overall, the LI-7500 measured a 15% greater flux than the LI-7500RS both on average and during only daytime hours. Meanwhile, LI-7500A and LI-7500RS fluxes were nearly identical, with 0.5% less flux measured by the LI-7500A and an additional 0.2% difference during the daytime. While the daily bias was as equally positive as negative, the LI-7500A tended to underestimate flux through the first and last third of the study period although possible rainfall effects exist. Greater flux was observed on 27 of the 41 days from DOY 196 -226, which coincided with 295 greater accumulated ET (Fig. 7). Relative error varied little by time of day. An increase in variability during advective conditions was due to greater mean (co)variance. Under advective conditions, the coefficient of determination was particularly small (see Table 2), but this coincided with large turbulent fluxes including downward sensible heat that was also slightly biased towards increased magnitude.
The 90-day ET (Fig. 7) was in good agreement among the three IRGAs, with slightly greater seasonal flux from the LI-300 7500, consistent with the larger variance in the timeseries of " . Systematic underestimation of ET for all IRGAs is consistent with advective conditions, especially in the earlier part of the growing season where the gap in daily ET is particularly large for a similar magnitude of ET (Fig. 7). Even if all spectral loss is corrected for, based on the conservation of water vapor and eddy covariance theory, the measured EC flux should be less than the true flux under advective conditions. Approximately 16% of accumulated ET was underestimated from LI-7500A or LI-7500RS relative to the accumulated lysimeter ET at the end 305 of the growing season (Fig. 7). However, only less than 5% of accumulated ET was underestimated from the oldest LI-7000 analyzer (Fig. 7). Furthermore, the EC and lysimeter should differ more with increasing mean ET because the advective component of ET, not captured by EC systems, is more likely to be elevated (Alfieri et al., 2012). The greater flux from the LI-7500 occurs nearly symmetrically on a diel basis, with relative differences smallest during the day. The mean daytime error of measured flux λE between the LI-7500A and LI-7500RS systems was 4.5%, with the LI-310 7500A estimating greater ET than the LI-7500RS on approximately three out of every four days. Systematic error d averaged 0.08 mm for the LI-7500RS system, which is rather large considering the mean measuring flux of 0.2 mm. Larger systematic error is typically associated with greater flux underestimation due to failure to capture all low frequency signals, consistent with the observed cospectra (Vickers and Mahrt, 1997). In contrast, daily λE differed by 18.6% between LI-7500 and LI-7500RS systems and the magnitude from the LI-7500RS only exceeded that of the LI-7500 on a single day. Comparing daily 315 declining ET (Fig. 8). Random error e was overwhelmingly similar among the sensors, indicating that uncertainty due to sampling has little effect on differences in estimated ET.

Water vapor variance errors 320
Water vapor variance and flux were compared from three similar eddy covariance systems yielding similar results in rain-free This is encouraging despite demonstrated substantial errors in the water vapor density measurements. More important for flux is the (co)variance of the water vapor density. Despite screening the data for quality, several outliers were observed in the ρ W < 8 '''' 325 which contributed to deflated r 2 values and notable discrepancies in water vapor fluxes. Overestimation of water vapor variance could contribute to overestimated flux but is not necessarily the case (Mauder et al., 2006). At noon on DOY 190, LI-7500A overestimated ρ W < 8 '''' by 5.9 g 2 m -6 ; while corresponding values were only 0.63 and 0.11 for LI-7500 and LI-7500RS, respectively. '''' (3.5 g 2 m -6 and 0.84 g 2 m -6 greater relative to LI-7500A and LI-7500RS). However, in a vast majority of cases, large ρ W < 8 '''' was observed with both the LI-7500 and LI-7500A relative to the LI-7500RS and were associated with a recent rainfall event. For instance, a large discrepancy in ρ W < 8 '''' among the three IRGAs occurred an hour after light rain on DOY 211, which suggests that thick water droplets may have been still evaporating from the mirror surface. Antecedent conditions were dry and with the cessation of precipitation, a sudden increase in mean wind 335 speed from under 3 to 5 m s -1 and a wind shift from east to south enabled sensible heat advection as clouds began to dissipate.
Although air humidity decreased by the end of the half-hour for all IRGAs, the magnitude measured by the LI-7500 was much smaller at the start of the averaging period than at the end, in contrast to observations by the LI-7500A and LI-7500RS. Further, we observed that the LI-7500A air humidity began decreasing within the first 15 minutes, suddenly increased by approximately 5 g m -3 , and then began a rapid decrease. This pattern is different than what was observed by the LI-7500RS, which initially 340 increased and then quickly decreased at an earlier time than for the LI-7500A (not shown). Large variability of air humidity in time and space caused large errors of water vapor density. The LI-7500 " '''' decreased to 7.01 g m -3 while the LI-7500A " '''' increased to 17.92 g m -3 . These corresponded to Δ " of 8.83 g m -3 and 1.95 g m -3 , respectively. While the LI-7500RS performance during this time was markedly better than that of the other two sensors, the -1.13 g m -3 bias was still different from its long-term offset. Resulting water vapor fluxes were smallest for the LI-7500A and largest for the LI-7500, with 345 sampling by the LI-7500RS seeming to best reflect the variations in eddies during a period of substantial air mass change. A similar event occurred on DOY 196. However, for the half-hour of interest, a relatively small difference in ρ W < 8 and LI-7500A resulted in a larger flux difference, in which a large, likely overestimated flux was measured by the LI-7500.
Interestingly, 20 Hz fluctuations for all systems were dampened during roughly the first half of this averaging period, showing signs of low frequency atmospheric motion. Once turbulence became more typical of a well-mixed boundary layer, the 350 amplitude of " < then grew with a larger variance noted in the LI-7500A and LI-7500 compared to the LI-7500RS. This is exactly what was observed on DOY 211 during its relevant averaging period. The effect of rainfall may linger depending on its timing. All the sensors exhibited some degree of non-stationarity in the " timeseries from late night on DOY 224 into the early morning of DOY 225. However, only the LI-7500A continued to exhibit this behavior for several more hours while the other two sensors showed constant flux. Because " was so similar between the LI-7500 and LI-7500RS, it seems that the LI-355 7500A was uniquely sensitive to the intermittent turbulence during this calm period. Based on the combination of high relative humidity and light winds, these observations were subject to increased random error as expected.

Water vapor flux errors
In the context of ET measurement, total daily magnitude is of prime importance for practical applications. Therefore, flux errors during the daytime, roughly between 09:00 and 17:00 LST, contribute to the vast majority of ET variation. The similarity 360 between the LI-7500A and LI-7500RS fluxes is reflected by the lack of scatter in covariance data. As expected, errors were larger during advective periods than for other times, but overall correlation between ρ W Uncorrected fluxes were assessed to assure that the data processing steps did not appreciably affect our findings. Postprocessing of turbulent fluxes could increment fluxes while causing greater error (Irmak et al., 2014). The magnitudes of a, d, 365 and rmsd were slightly smaller for all comparisons, and b and r 2 were nearly identical, indicating that the corrections contributed little to measurement uncertainty. For instance, the rmsd decreased by 6.8% and 7.3% for daytime fluxes against the LI-7500 and LI-7500A, respectively. Among the corrections, sensor separation and frequency response were of most interest for the LI-7500RS and LI-7500A pair since they are newer different optical analyzers. This may be why among the three generations of IRGAs, the LI-7500RS consistently had a larger spectral correction factor by approximately 2 to 4%, but 370 again, this served to only slightly decrease flux error. Its midday mean value of 1.11, though slightly larger than for the LI-7500 and LI7500A, was still less than reported in a feedlot for an LI-7500 and CSAT-3 EC system (Prajapati and Santos, 2017). This suggests that high frequency attenuation was relatively minor when turbulent intensity was large, and any missing flux was more attributable to low frequency. While the LI-7500 high frequency energy compared more favorably to the LI-7500A than the LI-7500RS, a large departure from the LI-7500A and LI-7500RS pattern was clearly observed at low 375 frequencies (Fig. 6).
It has previously been shown that turbulent flux error partitions into primarily random error, with daytime systematic error only as large as 0.018 mm (30 min -1 ) (Alfieri et al., 2011). In contrast, Sect. 3.3 demonstrated that the magnitudes of systematic error were generally large in response to daytime energy balance residuals. The different findings are based on different The authors declare that they have no conflict of interest. Tables   Table 1: