Measuring N2O Emissions from Multiple Sources Using a Backward Lagrangian Stochastic Model

Nitrous oxide (N2O) emissions from agricultural soil are substantially influenced by nitrogen (N) and field management practices. While routinely soil chambers have been used to measure emissions from small plots, measuring field-scale emissions with micrometeorological methods has been limited. This study implemented a backward Lagrangian stochastic (bLS) technique to simultaneously and near-continuously measure N2O emissions from four adjacent fields of approximately 1 ha each. A scanning open-path Fourier-transform infrared spectrometer (OP-FTIR), edge-of-field gas sampling and measurement, locally measured turbulence, and bLS emissions modeling were integrated to measure N2O emissions from four adjacent fields of maize production using different management in 2015. The maize N management treatments consisted of 220 kg NH3-N ha−1 applied either as one application in the fall after harvest or spring before planting or split between fall after harvest and spring before planting. The field preparation treatments evaluated were no-till (NT) and chisel plow (ChP). This study showed that the OP-FTIR plus bLS method had a minimum detection limit (MDL) of ±1.2 μg m−2 s−1 (3σ) for multi-source flux measurements. The average N2O emission of the four treatments ranged from 0.1 to 2.3 μg m−2 s−1 over the study period of 01 May to 11 June after the spring fertilizer application. The management of the full-N rate applied in the fall led to higher N2O emissions than the split-N rates applied in the fall and spring. Based on the same N application, the ChP practice tended to increase N2O emissions compared with NT. Advection of N2O from adjacent fields influenced the estimated emissions; uncertainty (1σ) in emissions was 0.5 ± 0.3 μg m−2 s−1 if the field of interest received a clean measured upwind background air, but increased to 1.1 ± 0.5 μg m−2 s−1 if all upwind sources were advecting N2O over the field of interest. Moreover, higher short-period emission rates (e.g., half-hour) were observed in this study by a factor of 1.5~7 than other micrometeorological studies measuring N2O-N loss from the N-fertilized cereal cropping system. This increment was attributed to the increase in N fertilizer input and soil temperature during the measurement. We concluded that this method could make near-continuous “simultaneous” flux comparisons between treatments, but further studies are needed to address the discrepancies in the presented values with other comparable N2O flux studies.


Introduction
At present, more than 50% of non-CO 2 greenhouse gas (GHG) emissions (mainly N 2 O and CH 4 ) come from agricultural sources [1][2][3][4]. Nitrous oxide emitted from agricultural soils is mainly attributed to the application of nitrogen (N) fertilizers [5,6]. Nitrogen losses via N 2 O emissions not only give rise to a negative impact on the environment (strong GHG input into the atmosphere) but also reduce Figure 1. The experimental site was at the Purdue Agronomy Center for Research and Education (ACRE). The field experiment was operated by (a) a scanning OP-FTIR, and seven retroreflectors created seven line-sampling (LS) concentration sensors (LS-1-LS-7) to "scan" four treatments (T1-T4) managed by contrasting field and N practices shown in Table 1. Three point-sampling (PS) inlets were deployed at the edges of fields (W, N, and E) to measure the background N2O concentrations or provided additional downwind sensors for fields of interest. (b) The prevailing wind came from the SSW, and the mean wind velocity was 3.6 ± 1.8 m s −1 . Table 1. Four treatments (T1-4) represented four N2O sources. Field practices were no-till and chisel plow. With the total N rate of 220 kg NH3-N ha −1 , NH3 was applied as a full (220 kg N ha −1 ) or a split (total N rate was equally split and applied) application at two timings of the fall in 2014 and the spring before planting in 2015.

Treatment (ha) Tillage
Anhydrous Ammonia (Total = 220 kg NH3-N ha −1 ) 2014 Fall 2015 Spring T1 (1. 2) No-till 220 0 T2 (1.1) No-till 110 110 T3 (1.1) Chisel plow 110 110 T4 (1. 2) No-till 0 220  Table 1. Three point-sampling (PS) inlets were deployed at the edges of fields (W, N, and E) to measure the background N 2 O concentrations or provided additional downwind sensors for fields of interest. (b) The prevailing wind came from the SSW, and the mean wind velocity was 3.6 ± 1.8 m s −1 . Table 1. Four treatments (T1-4) represented four N 2 O sources. Field practices were no-till and chisel plow. With the total N rate of 220 kg NH 3 -N ha −1 , NH 3 was applied as a full (220 kg N ha −1 ) or a split (total N rate was equally split and applied) application at two timings of the fall in 2014 and the spring before planting in 2015.

Treatment (ha) Tillage
Anhydrous Ammonia (Total = 220 kg NH 3  OP-FTIR spectrometer (MIDAC Model2501-C, MIDAC Corporation, Irvine, CA, USA) bounded by a retroreflector array of 26 corner cubes. The LS lengths were between 100 and 150 m at the height of approximately 1.5 m above ground level (a.g.l.). The OP-FTIR was mounted on a scanner with horizontal and vertical rotaries (YUASA computer numerical control, CNC) to scan seven retroreflectors deployed in or on borders of four fields to create seven open path sampling lines (LS-1-LS-7), and it took nearly thirty minutes to complete one cycle of the scanning (i.e., from LS-1 to LS-7). The scanner and OP-FTIR were not synchronized. A dwell time of each path was three minutes-assuring two-to-three approximately 1 min FTIR spectra per retroreflector. Some LS were placed in the fields to increase the sensitivity of emission from the field (e.g., relatively high concentrations from an active emission source) [29]. Each 1 min spectrum was acquired by co-adding 64 interferograms (IFG) at a resolution of 0.5 cm −1 , and a co-added IFG was apodized using a triangular function and converted to a single beam (SB) spectrum by Fourier-transform using AutoQuant Pro4.0 software package (MIDAC Corporation, Irvine, CA, USA). Details of acquiring OP-FTIR spectra were described in Lin et al. (2019) [30].
A difference frequency generation (DFG) mid-IR laser-based N 2 O/H 2 O analyzer (IRIS 4600, Thermo Fisher Scientific Inc., Waltham, MA) sequentially measured N 2 O concentrations from air streams collected through a gas sampling system (GSS). The DFG N 2 O analyzer measured N 2 O concentration with the precision (1 standard deviation) of 0.15 ppbv for three-minute averages. The analyzer was calibration checked with a standard N 2 O gas with a concentration of 520 ppbv every four hours. The GSS system included a sampling pump and four DC power solenoids connected to four sampling lines to collect gas samples by drawing the air through lines. Each sampling line consisted of a 9.5 mm-diameter Teflon ® tube fitted with a 1.0 µm Teflon ® filter [39]. Then, lines were located on a border of the field (N, W, E) ( Figure 1). In addition to three PS lines, a 50 m synthetic open path sampling system (SOPS) was made up of ten tubes (i.e., ten inlets) and deployed along one of the OP-FTIR paths (LS-4 between T2 and T3) to measure the path-averaged N 2 O concentration. Since the SOPS can measure N 2 O concentrations with high precision (±0.15 ppbv) and accuracy (frequent certified gas calibration), the path-averaged concentration measured by the SOPS was used as a benchmark to calibrate the bias of the concentration measured by the OP-FTIR. Thus, four sampling lines included three PS lines on the edges of the field and one SOPS line between the T2 and T3 treatments (i.e., NT and ChP). Each line was used to collect gas samples at the flow rate of approximately seven liters every minute, and the sampling interval was five minutes, which was controlled by electrifying solenoids.
Wind information (i.e., direction, velocity, and turbulence statistics) was measured using a 3D sonic anemometer (Model 81000, RM Young Inc., Traverse City, MI, USA), which was installed at 2.5 m height on the meteorological mast and recorded at 16 Hz. Ambient temperature, humidity, and barometric pressure were measured using a temperature-humidity sensor (Model HMP45C, Vaisala Oyj, Helsinki, Finland) and a pressure sensor (278, Setra, Inc., Boxborough, MA, USA), which were installed at 1.5 m a.g.l. mast. These meteorological data were averaged every five minutes by a data logger (Model CR1000, Campbell Scientific, Logan, UT, USA).

Open-Path FTIR
Calculations: A synthetic SB background was generated using IMACC software (Industrial Monitoring and Control Corp., Round Rock, TX, USA) to convert sample SB spectra to absorbance spectra. The N 2 O concentrations were quantified by partial least square (PLS) algorithms based on the Beer's law quantitative spectral processing software (Thermo Fisher Scientific TQ Analyst Version 8.0). Sixty spectra of mixed-gas standards (water vapor + N 2 O) were generated using a laboratory FTIR spectrometry joined with a multipass gas cell (i.e., optical path = 33 m) and used to build the PLS quantitative model. The analytical region for N 2 O calculations was used the combination window of 2188.7-2204.1 cm −1 and 2215.8-2223.7 cm −1 .
Calibrations: Fluctuations of ambient humidity, temperature, and path length substantially influenced the accuracy and precision of N 2 O calculations. The bias of the N 2 O concentrations derived from the OP-FTIR was calibrated by the path-averaged N 2 O concentration measured from the SOPS. The relative error between the FTIR-and SOPS-measured N 2 O concentrations (i.e., bias = [(x FTIR /x SOPs ) − 1] × 100%) was calculated to infer the accuracy of N 2 O concentrations determined by OP-FTIR [30]. The variations in ambient humidity and temperature are regarded as leading error sources for OP-FTIR quantification [31]. We assumed that ambient humidity and temperature were relatively constant within thirty minutes; as a result, the offset of the path-averaged concentration between the SOPS and OP-FTIR (LS-4) was considered as a constant within each 30 min interval and applied to calibrate other OP-FTIR paths (LS-1 to LS-7).
Uncertainty: It was assumed that ambient humidity, temperature, and N 2 O concentrations remained stationary in a well-mixed surface layer within thirty-minute intervals. For each path length (i.e., optical path = 100, 200, and 300 m), ten OP-FTIR spectra were acquired in the 30 min interval and calculated for N 2 O concentrations. A standard deviation of N 2 O concentration from ten spectra represented the uncertainty in N 2 O concentrations measured by the OP-FTIR using a PLS algorithm in a specific path length [31].

Emission Measurements
Gas measurements began just after the NH 3 application from 01 May to 11 June 2015. Emission rates (Q) of N 2 O were estimated by a Lagrangian stochastic dispersion (LSD) model (Windtrax2.0: Thunder Beach Scientific, http://www.thunderbeachscientific.com) operating in backward mode while N 2 O concentrations were estimated by an LSD model operating in forward mode.

Assumptions for Field Measurements
Assumptions required for this experiment included: (1) Stationary turbulence and flux-the averaged quantities of variables were invariant over thirty minutes. (2) Homogeneous emissions and turbulence-the 30 min averaged quantities of variables are invariant under a spatial translation.
(3) The dominance of N source-N source for N 2 O production was predominately from the applied NH 3 within sixty days after applications and atmospheric N deposition, (a)biotic N fixation and (in)organic N residues (i.e., soil organic matter, or 2:1 clay minerals) were negligible.

Backward Lagrangian Stochastic (bLS) Model
The upwind source emissions can be predicted based on the measured downwind gas concentrations using the LSD model in "backward" mode (bLS). The bLS model used the measured turbulence statistics to simulate the ratio of the LS-or PS-measured gas concentration (C) to the surface-emission rate (Q) by "backward" tracing trajectories of gas particles from sensors (downwind) to source (upwind). Emission rates can be inferred based on the deviation of concentration between the downwind fields and background (C bg ) (i.e., ∆C downwind−C bg ), which is described as: where (C/Q) sim is the bLS-simulated ratio of concentration to emission rate, Q is the inferred surface emissions from the source of interest, and ∆C is the deviation of gas concentrations. For the PS sensors (Figure 1), the simulated ratio can be calculated as the following equation: Atmosphere 2020, 11, 1277 6 of 20 For the LS and the SOPS sensors (Figure 1), thirty equidistant points were defined along the path, and the same number of particles (N) released from each point. Consequently (C/Q) sim (Equation (1)) was determined for each source j and path i as: where N is the number of the released particles from a sensor, and w 0 is the vertical velocity of "touchdowns" based on the turbulence statistics. The concentration at each point or path can be expressed as a function of a background concentration and the emissions from the four adjacent sources as: where a ij represents the ratio of concentrations to emission rates simulated ((C/Q) sim ), from the contribution of the jth source (j = 1, 2, . . . , m) to the ith sensor (i = 1, 2, . . . , n). The four N 2 O sources (m = 4) and ten gas sensors (n = 10) resulted in an over-determined [29] system of simultaneous equations of the form of Equation (4): A standard singular value decomposition (SVD) algorithm was used to solve the set of equations to produce the best estimate of the background concentration and the emissions of each source.

Forward Lagrangian Stochastic (fLS) Model
The downwind gas concentrations can also be predicted based on the measured emission rates of the upwind source using the LSD model in "forward" mode (fLS). In the fLS model, the measured turbulence statistics and emission rates are used to simulate the gas concentrations along each path based on (C/Q) sim (Equation (2)).

Emissions Quality Assurance
Numerous filtering criteria have been published for quality assurance of the bLS-estimated emissions. Thresholds of low wind conditions (i.e., u* < 0.15 m s −1 ) and exclusion of extremely stable and unstable atmosphere (i.e., |L| ≤ 5 m) are likely universal criteria for using the bLS model. Other published criteria vary with different experimental configurations. Although a TDF > 0.1 was the most commonly used threshold in the literature, the threshold for TDF will be influenced by the configuration of sensors and emission sources (e.g., varied from 0.1 to 0.7) [20,34,35,38]. This study investigated the effect of changing the TDF threshold on the emissions estimations. In addition, the sensitivity of N 2 O concentrations determined by OP-FTIR (i.e., σ ppbv) was used to provide additional criteria for checking the model performance. Ideally, the minimum concentration of the measured N 2 O among sensors (i.e., 7 LS + 3 PS) should be identical to the modeled background (C bg ) (Equation (5)). The C bg was compared with the lowest N 2 O concentration within all of the gas sensors (including SOPS), and the quality of the modeled emissions was considered poor and excluded if the deviation between C bg and the lowest concentrations was out of the range of ±3σ ppbv (|N 2 O min -C bg | > 3σ ppbv). The likely cause of this was either emissions non-stationarity over the one-half hour or spatial heterogeneity of the emissions source.

Advective Interferences
In a multi-source configuration, air parcels (an air parcel refers to a body of air in which thermodynamic properties were approximately uniform) advected from adjacent fields can be detected by sensors deployed in a field of interest and interfere with gas emission estimations. The output TDF from the bLS emissions model and the known areas of emission fields (T1-T4) were used to evaluate the magnitude of air parcels from adjacent fields. For instance, T2 was interfered with by T1, T3, and T4 when the wind came from NW ( Figure 1). In order to evaluate the magnitude of interferences from T1, T3, and T4, only the LS sensors of LS-2, LS-3, and LS-4 deployed in the T2 plot were used to calculate TDFs. Table 2 shows the procedure for the estimation of advective interference. First of all, the T1-T4 treatments were assumed to have the same emission rates. Based on the same turbulent statistics, the TDFs were calculated by changing the sizes of areas stepwise (i.e., T2 only, T2 + T1, T2 + T3, and T2 + T4). Given the known TDF and the specified area, the released particles from LS-2, LS-3, and LS-4 sensors covered 0.36, 0.76, 0.93, and 1.12 hectares of the areas of T2 only (1.5 ha), T2 + T1 (2.7 ha), T2 + T3 (2.9 ha), and T2 + T4 (2.6 ha), respectively ( Table 2). The fraction of air parcels (FRAC air , %) with touchdowns in each combination of fields (field of interest and an adjacent field) during each measurement period was compared to that of the field of interest. For example, FRAC air of T1, T3, and T4 represented potential adjacent interferences to the T2 treatment for a given emissions measurement. In the example of Table 2, more interferences came from the T4 (36% of FRAC air ) than the T1 (19% of FRAC air ) treatment for the T2 treatment at NW wind. Table 2. An example of determining the fraction of air parcel (FRAC air ) that traveled from adjacent fields (i.e., T1, T3, and T4) to the field of interest (i.e., T2), representing the magnitude of interferences, under a particular wind condition (wind direction: NW (317 • ); mean wind velocity: 3.7 m s −1 ; friction velocity: 0.32 m s −1 ; Monin-Obukhov length: −37 m). TDF: touchdown fraction. TD cover † : touchdowns cover area for the specified field, calculating by multiplying the area (e.g., 2.9 ha) by TDF (e.g., 0.32) of the specified field (e.g., T2 + T3). § : TD cover at each plot was calculated based on the difference of the TD cover area between specified fields. For example, TD cover areas of T2 only and T2 + T3 were 0.36 and 0.93, respectively. Thus, the TD cover area of T3 only was 0.57 (i.e., 0.93 subtracted 0.36). FRAC air (%): a fraction of air parcels, the percentage of TD cover area at each plot among four plots (treatments).

Emissions Uncertainty Analysis
The uncertainty in N 2 O concentrations in the scanning system (Section 2.3) was introduced in the fLS emissions model to estimate the uncertainty in bLS emission measurements. The procedure for estimating the uncertainty in emission rates is shown in Figure 2.
First of all, the known emission rates of the T1-4 treatments were used to "forward" calculate N 2 O concentrations of seven LS and three PS sensors based on the measured turbulence statistic using the fLS model. The background (C bg ) of N 2 O was assumed 350 ppbv, and 50,000 particles (N) were considered in the bLS model (step 1). Second, the uncertainty in N 2 O concentrations (σ ppbv shown in Section 2.3) was introduced into the fLS-calculated concentrations to infer emission rates using the bLS model. It was assumed that the uncertainties in N 2 O measured from three PS sensors were negligible. Only the uncertainty from seven LS sensors was considered in the model. As a result of two possibilities (N 2 O ± σ ppbv) for each LS sensor, there were 128 possible permutations (2 7 )

Statistical Analysis
Statistical analysis with PROC MIXED procedure in SAS [40] was used to test the effect of the treatment on the bLS-measured N2O emissions using ANOVA. Outliers (values outside 1.5 times interquartile range [IQR] below lower quartile [Q1] and above upper quartile [Q3]) were removed, and the data were shifted and transformed using the square root transformation to meet the normality assumption prior to statistical comparisons. The Fisher's least significant difference (LSD) test was used for multiple comparisons among population means of N2O emissions measured over the study period (α = 0.05).

OP-FTIR Sensitivity
The LS N2O concentrations were derived from the measured OP-FTIR spectra. Both the precision and accuracy of N2O concentrations increased with increasing path lengths. Uncertainties in N2O concentrations (standard deviation) reduced from ±3.3 to ±1.4 ppbv by increasing optical path lengths from 100 m to 300 m ( Figure 3). Nitrous oxide concentrations were underestimated by 14%, relative to the background concentration of 349.0 ppbv, at the optical path length of 100 m. By increasing path lengths to 300 m, the quantitative accuracy increased by approximately 10% (bias reduced from 14% to 4%). The increased path length increased the path-integrated concentrations (i.e., concentration x path length) as well as improved gas quantification [31]. Since the optical path lengths were between 200 (i.e., LS-2 and LS-6) and 300 m (i.e., LS-1, -3, -4, -5, and LS-7) in this scanning system, the uncertainties determined from the optical paths of 200 m (±2.3 ppbv) and 300 m (±1.4 ppbv) were

Statistical Analysis
Statistical analysis with PROC MIXED procedure in SAS [40] was used to test the effect of the treatment on the bLS-measured N 2 O emissions using ANOVA. Outliers (values outside 1.5 times interquartile range [IQR] below lower quartile [Q1] and above upper quartile [Q3]) were removed, and the data were shifted and transformed using the square root transformation to meet the normality assumption prior to statistical comparisons. The Fisher's least significant difference (LSD) test was used for multiple comparisons among population means of N 2 O emissions measured over the study period (α = 0.05).

OP-FTIR Sensitivity
The LS N 2 O concentrations were derived from the measured OP-FTIR spectra. Both the precision and accuracy of N 2 O concentrations increased with increasing path lengths. Uncertainties in N 2 O concentrations (standard deviation) reduced from ±3.3 to ±1.4 ppbv by increasing optical path lengths from 100 m to 300 m ( Figure 3). Nitrous oxide concentrations were underestimated by 14%, relative to the background concentration of 349.0 ppbv, at the optical path length of 100 m. By increasing path lengths to 300 m, the quantitative accuracy increased by approximately 10% (bias reduced from 14% to 4%). The increased path length increased the path-integrated concentrations (i.e., concentration × path length) as well as improved gas quantification [31]. Since the optical path lengths were between 200 (i.e., LS-2 and LS-6) and 300 m (i.e., LS-1, -3, -4, -5, and LS-7) in this scanning system, the uncertainties determined from the optical paths of 200 m (±2.3 ppbv) and 300 m (±1.4 ppbv) were averaged to Atmosphere 2020, 11, 1277 9 of 20 calculate the overall uncertainty in N 2 O concentration of each LS sensor (i.e., ±1.8 ppbv). Concentration biases (~4%) were calibrated by the SOPS before flux calculations.
Atmosphere 2020, 11, x FOR PEER REVIEW 9 of 20 averaged to calculate the overall uncertainty in N2O concentration of each LS sensor (i.e., ±1.8 ppbv). Concentration biases (~4%) were calibrated by the SOPS before flux calculations.

Emission Measurements
The influence of the TDF threshold for a valid bLS emissions estimate was evaluated relative to the potential for advective interference in the emissions estimate. We investigated the relationship between advective interferences (denoted by FRACair) and TDFs ( Figure 4). The measured emissions were sorted into three categories based on the TDF threshold: high (>0.9), median (0.9 > TDF > 0.5), and low (0.5 > TDF > 0.1) to examine the relationships between emissions and TDFs ( Figure 5).

Emission Measurements
The influence of the TDF threshold for a valid bLS emissions estimate was evaluated relative to the potential for advective interference in the emissions estimate. We investigated the relationship between advective interferences (denoted by FRAC air ) and TDFs ( Figure 4). The measured emissions were sorted into three categories based on the TDF threshold: high (>0.9), median (0.9 > TDF > 0.5), and low (0.5 > TDF > 0.1) to examine the relationships between emissions and TDFs ( Figure 5).

High TDF Threshold (>0.9)
When the TDF values of the field of interest are nearly 1, the four treatments (T1-T4) alone had FRAC air values (mean ± SD) of 89 ± 12%, 87 ± 14%, 84 ± 16%, and 80 ± 19% respectively ( Figure 4). The targeted field received a "clean" upwind background and resulted in high TDFs (~1.0) when FRAC air~1 00%. For instance, T1 received a clean background when the wind came from SW (202.5-247.5 • ), as well as high TDFs (1.0) ( Figure 4). Likewise, the SE wind (112.5-157.5 • ), NE wind (22.5-67.5 • ), and NW wind (292.5-337.5 • ) led the TDFs of 1.0 for T2, T3, and T4 treatments, respectively. Increased interference (0-20% of the FRAC air of the adjacent upwind field) reduced TDFs from 1.0 to 0.9. For instance, T2 started interfering with T1 when the wind directions shifted from SW to S (~12% of FRAC air from T2); likewise, T4 began to interfere with T1 at the direction of 270.0-292.5 • (~6% of FRAC air from T4) (Figure 4a). The mean (± SD) emissions of the T1 through T4 treatments were 3.9 ± 3.2, 0.8 ± 1.1, 1.5 ± 1.1, and 1.5 ± 2.1 µg m −2 s −1 , respectively, over the 41 days of measurement ( Figure 5). The uncertainty in N 2 O concentrations (±1.8 ppbv shown in Figure 3) was used to calculate the uncertainty in N 2 O emission rates, and the corresponding uncertainties in the emission estimates (mean ± SD) of the treatments T1 through T4 were 0.6 ± 0.2, 0.3 ± 0.1, 0.1 ± 0.1, and 0.5 ± 0.3 µg m −2 s −1 , respectively ( Figure 6). The emission rates indicated that the T1 treatment (i.e., the combination of NT and full N rates applied in the fall) resulted in the highest N 2 O emissions in the spring.    3.2, 0.8 ± 1.1, 1.5 ± 1.1, and 1.5 ± 2.1 µ g m −2 s −1 , respectively, over the 41 days of measurement ( Figure  5). The uncertainty in N2O concentrations (±1.8 ppbv shown in Figure 3) was used to calculate the uncertainty in N2O emission rates, and the corresponding uncertainties in the emission estimates (mean ± SD) of the treatments T1 through T4 were 0.6 ± 0.2, 0.3 ± 0.1, 0.1 ± 0.1, and 0.5 ± 0.3 µ g m −2 s −1 , respectively ( Figure 6). The emission rates indicated that the T1 treatment (i.e., the combination of NT and full N rates applied in the fall) resulted in the highest N2O emissions in the spring. Negative emissions implied an N2O "sink". Soil uptake of N2O, however, was unlikely when the soil was enriched with available N-substrates after NH3 application. One of the assumptions of the bLS model required a substantial difference in emission rates between the source (hotspot) and background to make the "delta concentration" detectable (Equation (1)). Negative emissions measured from a clean background and minimum interferences (TDF > 0.9 shown in Figure 5) were attributed to emissions lower than the detection limit. The standard deviation of the "negative" emissions of T1-4 (TDF > 0.9) was 0.4 µ g m −2 s −1 , so the MDL of this scanning open path method was calculated as ±1.2 µ g m −2 s −1 (3σ) (gray bars shown in Figure 5).
Adjacent field interference not only led to increased emissions uncertainty but also increased bias. The emissions bias became substantial if the upwind treatment emitted more than the downwind field. For instance, high N 2 O emissions were measured from T1 at the S wind (180-225 • ) (Figure 5a) substantially interfered with the T4 emissions measurements (with 60% of the FRAC air from T1 shown in Figure 4d), resulting in "negative" emissions (less than the lower boundary of MDL of −1.2 ug m −2 s −1 ) from the T4 treatment (Figure 5d). Even when the TDF values were relatively high, the T1 treatment (upwind) emitted more N 2 O than T4 (downwind), so the air parcels from the T1 treatment likely carried higher N 2 O concentrations over the T4 plot (N 2 O advection) and introduced substantial biases in the downwind treatment emissions estimate. Emissions lower than −1.2 µg m −2 s −1 were considered as incorrect estimations and excluded.
At low TDF thresholds, the bias in adjacent field interference substantially influenced the field with low emission rates. For instance, emissions lower than −1.2 µg m −2 s −1 (lower bound of the MDL) were measured for T4 under an S wind because of the interferences from the T1 treatment emissions (Figures 4 and 5). Emissions from the T3 treatment were subjected to interferences from the emissions from all the other treatments when the wind direction was 180-270 • (prevailing SW wind in 2015 measurements). In contrast to the emissions from the T4 treatment, most of the emissions from the T3 treatment were still measurable and greater than the MDL (Figure 5c) even under the three-source interferences of the T1, T2, and T4 emissions.

Integrated Effects of Wind Direction and Speed
As expected, the adjacent field interference of the bLS-emission estimates was a function of the wind speed and direction. Increased wind speed resulted in greater bias (i.e., emissions < −1.2 µg m −2 s −1 shown in Figure 7a) and increased uncertainty (Figure 7b) in the emissions measurements. For instance, the increased wind speed amplified the interferences from T1 as well as emission biases for the T4 measurements at the S wind (Figure 7a). Advective interferences and the corresponding uncertainties in emissions were highly wind-speed-dependent. For instance, the NE wind (0-90) resulted in high TDFs (>0.9) and low uncertainties (0.1 ± 0.1 µg m −2 s −1 ) for T3 measurements. Interferences and emission uncertainties (0.5 ± 0.5 µg m −2 s −1 ) substantially increased at the wind from SW (180-270 • ) for T3 measurements. These uncertainties increased with increasing the wind speed: linear regression R 2 value of 0.20 for winds from 0-90 • ; R 2 = 0.68 for winds from 90-180 • ; R 2 = 0.90 for winds from 180-270 • ; R 2 = 0.42 for winds from 270-360 • (Figure 7b).
Although this study only investigated the effect of uncertainties in N 2 O concentrations on the propagation of uncertainty in the bLS-measured emission, it is worthwhile to mention that the geometry of the source and sensor locations is another source of the uncertainties in bLS measurements [27][28][29]. Flesch et al. (2009) [29] mentioned that the source-sensor geometry is critical, and the ideal sensor deployment would be: each sensor only "sees" one emission source. The "optimal" deployment suggested by Flesch et al. (2009) [29] was also implemented in our study; however, the downwind sensors in the field of interest still detected the blended plume from the adjacent field via advection under different wind conditions. Therefore, the quantified advective interferences (Figure 4) associated with the propagated uncertainty in emission estimations ( Figure 6) were used to assure data qualities in this study.
"optimal" deployment suggested by Flesch et al. (2009) [29] was also implemented in our study; however, the downwind sensors in the field of interest still detected the blended plume from the adjacent field via advection under different wind conditions. Therefore, the quantified advective interferences (Figure 4) associated with the propagated uncertainty in emission estimations ( Figure  6) were used to assure data qualities in this study. Figure 7. The effect of the mean wind speed on interferences from the upwind fields: (a) the T1 treatment interfered with the T4 treatment at the S wind, and these interferences increased with increasing the wind speed. Increased wind speed led to great biases (emissions < −1.2 µ g m −2 s −1 ) in the T4 measurements; (b) the increased wind speed increased uncertainties (standard deviation, SD) in the T3 measurements.

Data Filtering and Statistics
There was a total of 1968 possible half-hour bLS emission values (48 × 41 days) for each treatment over the period of 01 May to 11 June 2015. The OP-FTIR did not operate during intermittent rain events and field operations (e.g., herbicide spraying) due to the risk of FTIR mirror contamination, retroreflector contamination, or thermal IR signal absorption. The continuous operation was commonly impossible due to thermal IR signal absorption of condensed water on the spectrometer and/or retroreflector during the night or early morning. Strong wind-induced vibrations resulted in poor optical path alignment and low-quality OP-FTIR spectra resulting in poor N2O concentration estimation. Lastly, infrequent dusty wind-blown events also resulted in low spectra quality and N2O concentration determination. As a result, only 675 of the 1968 possible measurement periods (approximate 34%) were available to assess treatment effects on emissions in 2015.
Estimated emissions were excluded during ill-defined turbulence conditions of u* < 0.15 m s −1 , |L| < 5 m, when TDF < 0.1, and when |N2Omin− | > 3σ ppbv (exclusion criteria 1 shown in Figure  8). Nearly 51% of the total collected half-hour emissions (n = 345 out of 675) remained for each treatment. Most of the excluded emissions were nighttime (i.e., 20:00-06:00, local time, LT) values when the boundary layer was stable with low winds. Estimated emissions less than the MDL (−1.2 µ g m −2 s −1 ) were excluded due to the severe impact of the advection-induced bias (exclusion criteria 2 in Figure 8). An additional 62 and 101 half-hour periods of emissions were removed from analysis for the T2 and T4 treatments, respectively, mainly resulting from the predominant interference from the T1 plot under the prevailing SSW wind (Figure 1). For the T1 and T3 measurements, high emission plots, only 11 and 20 half-hour emissions were less than −1.2 μg m −2 s −1 . Although emission estimations of the T3 treatment were interfered with by emissions from three adjacent fields (T1, T2, T4), most of the emissions from T3 were measurable (>MDL). The range of the advection-induced uncertainties was from 0.1 to 2.0 µ g m −2 s −1 under different wind conditions and increased with increasing advective interferences ( Figure 6). Estimated emissions were excluded if the uncertainty

Data Filtering and Statistics
There was a total of 1968 possible half-hour bLS emission values (48 × 41 days) for each treatment over the period of 01 May to 11 June 2015. The OP-FTIR did not operate during intermittent rain events and field operations (e.g., herbicide spraying) due to the risk of FTIR mirror contamination, retroreflector contamination, or thermal IR signal absorption. The continuous operation was commonly impossible due to thermal IR signal absorption of condensed water on the spectrometer and/or retroreflector during the night or early morning. Strong wind-induced vibrations resulted in poor optical path alignment and low-quality OP-FTIR spectra resulting in poor N 2 O concentration estimation. Lastly, infrequent dusty wind-blown events also resulted in low spectra quality and N 2 O concentration determination. As a result, only 675 of the 1968 possible measurement periods (approximate 34%) were available to assess treatment effects on emissions in 2015.
Estimated emissions were excluded during ill-defined turbulence conditions of u* < 0.15 m s −1 , |L| < 5 m, when TDF < 0.1, and when |N 2 O min − C bg | > 3σ ppbv (exclusion criteria 1 shown in Figure 8). Nearly 51% of the total collected half-hour emissions (n = 345 out of 675) remained for each treatment. Most of the excluded emissions were nighttime (i.e., 20:00-06:00, local time, LT) values when the boundary layer was stable with low winds. Estimated emissions less than the MDL (−1.2 µg m −2 s −1 ) were excluded due to the severe impact of the advection-induced bias (exclusion criteria 2 in Figure 8). An additional 62 and 101 half-hour periods of emissions were removed from analysis for the T2 and T4 treatments, respectively, mainly resulting from the predominant interference from the T1 plot under the prevailing SSW wind (Figure 1). For the T1 and T3 measurements, high emission plots, only 11 and 20 half-hour emissions were less than −1.2 µg m −2 s −1 . Although emission estimations of the T3 treatment were interfered with by emissions from three adjacent fields (T1, T2, T4), most of the emissions from T3 were measurable (>MDL). The range of the advection-induced uncertainties was from 0.1 to 2.0 µg m −2 s −1 under different wind conditions and increased with increasing advective interferences ( Figure 6). Estimated emissions were excluded if the uncertainty was 50% greater than the emission rates (exclusion criteria 3 shown in Figure 8). Low emissions, however, were more vulnerable to these criteria and excluded (Figure 8).
Atmosphere 2020, 11, x FOR PEER REVIEW 14 of 20 was 50% greater than the emission rates (exclusion criteria 3 shown in Figure 8). Low emissions, however, were more vulnerable to these criteria and excluded (Figure 8). By excluding the data measured during the ill-defined turbulence conditions (exclusion 1), the mean N2O emissions of the four treatments over 41 day measurement period were 2.7 (T1), −0.2 (T2), 2.3 (T3), and −0.5 (T4) µ g m −2 s −1 , respectively (Figure 8). The mean "negative" emissions of the T2 and T4 treatments appeared to be largely due to the interference from T1 ( Figure 5). These mean emissions became 2.8 (T1), 0.4 (T2), 2.5 (T3), and 0.4 (T4) µ g m −2 s −1 by further removing the half-hour emissions less than −1.2 µ g m −2 s −1 (exclusion 1&2), showing a substantial improvement for emission estimations from the T2 and T4 treatments. Exclusion criteria 3 reduced the measured emissions to only 24% of the total collected samples (n = 675). These left-over and relatively high emissions increased mean N2O emissions (i.e., 3.9 [T1], 0.7 [T2], 3.7 [T3], and 0.7 [T4] µ g m −2 s −1 ) and likely overestimated the overall average emission rate over the measurement period. Therefore, we only implemented the exclusion criteria 1 and 2 to clean the data from each treatment. Moreover, the skewed distribution of the estimated emissions showed that the "instantaneous" high emissions substantially influenced the overall mean emission values. These high emissions were considered outliers (the value greater than the upper whiskers shown in Figure 8

Treatment Comparisons
The split-N rate application (T2) resulted in lower mean N 2 O emissions over the study period than the full-N rate application (T1) in the NT field (p < 0.05). This was consistent with the results of other studies [7,[41][42][43]. For split N applications, chisel plowing the field before planting corresponded with higher mean N 2 O emissions in the growing season than when NT was done (T2 vs. T3, p < 0.05). This tendency for lower emissions under the NT management is consistent with the results of all prior N 2 O emission studies using soil chamber measurements at the ACRE [12,13,44]. The effect of tillage practices on N 2 O emissions are variable and based on climate regime and soil types [7,42,45]. This study suggested that NT management can mitigate N 2 O emissions, consistent with some field studies [42,46], but in contrast with others [47,48]. This inconsistency of the NT effect on N 2 O emissions may have been due to the duration of NT adoption. A long-term (>10 years) NT adoption was reported to mitigate N 2 O-N loss and the overall net global warming potential [42,48].

Flux Comparison with Other Studies
In this study, the quality-assured half-hour N 2 O emissions were observed from −1. The bLS method has been extensively qualified in comparison with flux gradient [24], ZINST [53], mass balance [21], integrated horizontal flux method [54], and Eulerian inverse model [55], and these methods were in good agreement.
Compared with other micrometeorological studies measuring N 2 O emissions from cereal crop production systems [56][57][58], this study showed an increase in the estimated emission by a factor of 1.5~7.  [58] used the eddy covariance (EC) method to measure N 2 O-N emissions from the ammonium-nitrate-and urea-fertilized soils and showed the range of half-hour emission from 0 to 1.4 µg m −2 s −1 . Not only do different N management practices significantly influence N 2 O emissions (Section 3.6), but the reduced N rate likely decrease emission rate [56]. These studies [56][57][58] usually had a lower annual total-N input from 98 to 180 kg-N ha −1 input with the mean N rate of 130 ± 28 kg-N ha −1 than this study (a total of 220 kg-N ha −1 annual input). Furthermore, soil microbes are sensitive to soil temperature, and the decreased temperature tends to reduce soil microbial activities as well as N 2 O productions. These studies [56][57][58] had a lower soil temperature (usually between 0 and 20 • C) during measurements than this study, which soil temperature was from 8.6 • C to 30.6 • C with the mean temperature of 18.2 ± 3.9 • C at the 10 cm soil (data not shown).
The EC method is one of the most common micrometeorological techniques often used to monitor long-term N 2 O-N loss and compare with the static chamber method. Several studies showed that EC and chamber methods had a similar flux magnitude for N 2 O-N measurements, but the EC method usually observed higher N 2 O-N emissions than chamber measurements [15,59,60]. The measurement of N 2 O fluxes by static and dynamic chambers has been more widely used to compare gas emissions among treatments and investigate the management effects on N 2 O emissions. Previous chamber studies cited in Section 3.6 [7,12,13,41,[43][44][45][46][47][48] showed the mean emission values ranged from 0.01 to 0.32 µg m −2 s −1 . Increased emissions measured using micrometeorological approaches are likely due partially to differences in the resistance to gas transfer through the still air in the soil pore space, through in the laminar layer of air over the soil, and through the turbulent boundary layer of air over the soil. Different soil gas permeability would then create differing resistances to gas flow in the soil, and differences in atmospheric turbulence would create differing resistance to gas flow in the turbulent atmosphere. Several studies showed increased wind speeds potentially increased N 2 O emissions, which was also shown in Figure 7a [16,[61][62][63]. Poulsen et al. (2017) [62] and Pourbakhtiar et al. (2017) [63] indicated that the wind-induced emissions could be 20-100 times larger than the emissions that are only driven by the molecular diffusion (e.g., most of the chamber measurements).

Challenges and Potentials for Improvement
It is challenging to quantify the overall N 2 O-N loss over a long-term period (days to years) because of the high temporal variability in field N 2 O emissions. The OP-FTIR plus bLS system can measure gas fluxes with high temporal resolution and ideally have a better estimation for calculating the cumulative N 2 O-N loss. This method, however, could not measure emissions in low wind conditions commonly during nighttime (i.e., 20:00-06:00, LT, in this study) and had the potential for advective interferences under high wind conditions. For assessing "representative" cumulative N 2 O emissions, estimates of the emission during low and high wind conditions are required. Simply scaling up the mean emission using this method to estimate the overall N 2 O-N loss over the measurement period (e.g., conversion from µg-N m −2 s −1 to kg-N ha −1 over 41 days in this study) likely resulted in a substantial positive bias as well as misleading conclusion about the overall N 2 O-N loss. Although our studies showed that the scanning OP-FTIR plus bLS model could make "simultaneous" flux comparisons between treatments, the causes for the relatively high half-hourly emissions measured using this method need to be further studied.
For the method improvement, advective interference issues can be minimized by (1) adding measurements between treatments or reducing treatments, (2) measuring larger fields or by careful selection of adjacent treatments to prevent expected high and low emission treatments from interacting. Furthermore, the emission estimates under low wind conditions (e.g., largely nocturnal) must be included to evaluate long-term cumulative emissions, which can be achieved by incorporating different measurement techniques and appropriate "gap-filling" methods for missing data points [56][57][58]. For instance, chamber methods typically have a higher sensitivity to measure low emission rates than micrometeorological methods and can be operated in most weather conditions (e.g., rain events or low winds) [16]. Grant and Omonode (2018) [64] also showed that chamber measurements tend to be comparable to micrometeorological emission measurements during low winds and nocturnal conditions. The integrated uses of chambers and this or other micrometeorological methods would improve the temporal representativeness of seasonal to annual N 2 O-N loss estimations.

Summary and Conclusions
Measuring the above-canopy emissions in adjacent treatment plots was shown to require careful assessment of advection from adjacent fields. The commonly used threshold for sufficient sampling of an identified source (i.e., TDF > 0.1) in bLS modeling of emissions did not constrain the model conditions enough in this medium plot multi-source situation. In this study, TDF thresholds greater than 0.9 indicated that a field of interest received a clean upwind background and resulted in low emission uncertainties (0.5 ± 0.3 µg m −2 s −1 ) for N 2 O flux estimations. The TDF threshold values ranging from 0.5 and 0.9 typically had one predominant upwind field influencing the emission estimates of the downwind field of interest and increased emission uncertainties (i.e., 0.6 ± 0.4 µg m −2 s −1 ). The TDF thresholds less than 0.5 typically had all three fields contributing to upwind interferences and resulted in relatively high emission uncertainties (i.e., 1.1 ± 0.5 µg m −2 s −1 ). Under these conditions, the advection of N 2 O emitted from the upwind field likely negatively biased the downwind field emission estimations.
After minimizing this advective influence on the bLS-modeled emissions from the four adjacent fields, the effects of different N and tillage management practices on N 2 O emissions were generally consistent with prior studies based on soil-chamber emissions estimates. A single fall application of N led to higher growing season N 2 O emissions than the split-N rates applied in the fall and spring.
For the same split-N applications in the fall and spring, the ChP practice tended to increase N 2 O emissions compared with the NT treatment.
Implementing the bLS method in a multi-source area is still ongoing research and becoming interesting to many users. Although the discrepancies in the presented values with other micrometeorological studies need to be further investigated, this study provides the current or future users with the approaches to identify the issues, potential solutions, and uncertainty analyses in a multi-source situation.