Interactive comment on “ The Chemistry of Atmosphere-Forest Exchange ( CAFE ) Model – Part 2 : Application to BEARPEX-2007 observations ”

The paper shows a comprehensive detailed evaluation of a multi-layer canopy exchange model, coupled to a representation of the lower boundary layer, by comparison with observations of fluxes and concentrations of a range of reactive compounds during the BEARPEX-2007 campaign in the Sierra Nevada, US. The paper is well written, addresses the application of an interesting dataset to not only conduct an quantitative but also a more qualitative analysis to further identify some of the main challenges involved in understanding canopy chemical exchange processes. It is clear from the presented detailed evaluation that in particular the dry deposition process including the


Introduction
Forest-atmosphere exchange of hydrocarbons, ozone, oxidized nitrogen and other reactive species impacts both atmospheric composition and ecosystem productivity, with broad implications for air quality and climate Isaksen et al., 2009;Fowler et al., 2009;Erisman et al., 1998). Quantifying deposition and emission to/from the forest, however, continues to present a significant experimental and theoretical challenge. Recent work has indicated that the air within and just above the canopy is highly oxidizing during the daytime (Farmer and Cohen, 2008;Holzinger et al., 2005;Lelieveld et al., 2008). This oxidative photochemistry affects the net biosphere-atmosphere exchange of many species. For example, the "escape efficiency" of highly reactive terpenoids is likely much less than unity (Ciccioli et al., 1999;Stroud et al., 2005;Bouvier-Brown et al., 2009a;Forkel et al., 2006), with consequences for scaling up leaflevel emissions for use in regional and global models. As a substantial in-canopy sink for oxidants like ozone (O 3 ), this chemistry could also contribute to downward O 3 fluxes (Goldstein et al., 2004;Kurpius and Goldstein, 2003;Fares et al., 2010a).
Chemistry can also influence surface fluxes of reactive nitrogen compounds, including NO x (≡NO + NO 2 ), acyl peroxy nitrates (APNs), alkyl nitrates (ANs) and nitric acid (HNO 3 ). Several measurement and modeling studies have demonstrated the influence of in-canopy gradients in radiation, O 3 and turbulent transport on fluxes of NO x (Gao et al., 1991;Dorsey et al., 2004;Duyzer et al., 2004) and the need to resolve canopy-scale processes in regional and global models (Ganzeveld et al., 2002a, b). One set of observations showing upward HNO 3 and APN fluxes over a young Ponderosa pine plantation suggests that, under certain conditions, intra-canopy chemistry may even alter the sign of fluxes traditionally assumed to be controlled by deposition (Farmer and Cohen, 2008). More recently,  observed downward APN fluxes at this same forest, but determined that the magnitude of the flux was sensitive to multiple in-canopy processes, including deposition, thermal decomposition and photochemical production.
Numerical modeling is an ideal tool for examining the interplay of physical and chemical processes contributing to net reactive gas exchange. Here, we apply the Chemistry of Atmosphere-Forest Exchange (CAFE) model to the comprehensive dataset from the Biosphere Effects on Aerosols and Photochemistry Experiment (BEARPEX) 2007 field campaign to investigate forest-atmosphere exchange at a young Ponderosa pine plantation. After a brief description of BEARPEX-2007 and a review of the key features of CAFE, we proceed with a detailed evaluation of observations from BEARPEX-2007. Our analysis focuses on the mechanisms controlling concentrations and fluxes of VOCs, hydrogen oxides, ozone, and reactive nitrogen.

Campaign and site description
BEARPEX-2007 was a multi-institutional collaborative research effort aimed at understanding the impact of forestatmosphere interactions on atmospheric composition. During the intensive measurement period of 15 August to 10 October 2007, a wide suite of chemical and meteorological observations were obtained within and immediately above a 17-yr-old Ponderosa pine plantation managed by Sierra Pacific Industries. The site is adjacent to the University of California's Blodgett Forest Research Station (BFRS), located in the western foothills of the Sierra Nevada Mountains, CA (38 • 58 42.9 N, 120 • 57 57.9 W, elevation 1315 m), and has been described in detail elsewhere . The BFRS overstory is primarily Ponderosa pine, with a few interspersed White fir, Douglas fir, Incense cedar, Black oak and Sugar pine, while the understory consists of Manzanita and Ceanothus shrubs.
For the current study, we simulate mean noontime (11:30-12:30 PST) observations from two sub-periods, designated "hot" (28 August-3 September, or day of year 240-246) and "cool" (13-18 September, or day of year 256-261). These windows were chosen because day-to-day meteorology (particularly temperature) is fairly uniform throughout each period, and because they contain the most overlap among chemical observations. Figure 1 illustrates near-surface temperature profiles for each period; the average temperature difference between the two periods is ∼10 • C. The hot and cool periods are representative of the general meteorological trend observed during BEARPEX-2007, characterized by a hot and dry August followed by a sharp transition to cooler, more humid conditions in September (Bouvier-Brown et al., 2009a;Wolfe et al., 2009); however, neither period is representative of the extreme conditions sampled during the campaign. Both periods are largely cloud free and remain under drought conditions, as the selected cool period precedes the first rain. Chemical observations from these periods are discussed in Sect. 3 and are summarized in Table S1 of the Supplement.

Model description
CAFE is a 1-D chemical transport model that resolves deposition, emission, chemistry and vertical diffusion throughout the canopy and mixed layer. The CAFE model and the details of its setup for BFRS are described in a companion paper (Wolfe and Thornton, 2011), and we will only briefly review the key aspects and modifications here. Table 1 lists important model parameters. The model domain consists of 86 layers in the vertical ranging from 0.01 m to 800 m, with noneven layer spacing that results in a fine-resolution grid of 36 layers within the forest canopy and 50 within the atmospheric boundary layer (ABL). High resolution within the canopy provides the greatest detail where vertical gradients in all  processes are largest and minimizes numerical artifacts from operator-splitting of chemistry and diffusion. Within each layer, the 1-D time-dependent continuity equation is solved to determine the rate of change for all chemical species: Terms on the right respectively represent rates of chemical production, chemical loss, emission, deposition, advection (horizontal mixing) and vertical turbulent flux divergence. In its current form, CAFE is designed to calculate steady-state concentration and flux profiles, i.e. it is not meant to perform as a prognostic model.
The resolved canopy includes an overstory with a height of 10 m, a one-sided leaf area index (LAI) of 3.2 m 2 m −2 and a leaf area dry mass (d) of 219 g m −2 , as well as an understory with a height of 2 m, an LAI of 1.9 m 2 m −2 and a d of 377 g m −2 . The leaf area density function (LADF), which describes how leaf area is distributed vertically, mimics observed vegetation structures (L. Misson, personal communication, 2008). Meteorological constraints, which are held constant throughout a model run, are taken as the mean noontime observations from the hot and cool periods (Table 2) with further parameterizations as outlined in the companion paper. Of particular importance for the current study is the imposed canopy microclimate, as this plays a large role in both chemistry and vertical mixing. Temperature is interpolated via a spline fit between four measurements ranging from 12.5 to 3 m and extrapolated to the ground (Fig. 1). This treatment is consistent with observations at 1.5 m recorded during BEARPEX-2009(R. Weber, personal communication, 2010. Friction velocity (u * ) and radiation were only measured above the canopy, thus their in-canopy profiles are estimated using standard parameterizations that assume an exponential decay as a function of cumulative leaf area. For BFRS, this results in a factor of 10 decrease in both u * and radiation between the top of the canopy and the ground. Schade et al. (2000) note that the top-down radiation attenuation algorithm is not optimized for coniferous forests as it does not explicitly account for "clumping" of needles. The flexibility of our simple exponential parameterization should minimize this problem via a judicious choice for the radiation extinction coefficient, k rad , as detailed in the companion paper.
Turbulent diffusion is represented using a first-order fluxgradient approach, also known as K-theory: Above 12.5 m, the eddy diffusion coefficient, K(z), is based on values used by Gao et al. (1993), scaled to an ABL height of 800 m. Below 12.5 m, K(z) is a function of friction velocity and canopy height and includes a correction factor to account for "near-field" effects of canopy elements on eddy diffusion (Makar et al., 1999;Raupach, 1989), though the latter is close to unity for the current study. As detailed further in Wolfe and Thornton (2011), diffusion coefficients are constrained by several metrics, including comparison with above-canopy flux-gradient relationships of nonreactive scalars and with previous estimates of the canopy mixing timescale. The resultant canopy residence time is ∼2 min for our conditions. Emissions of biogenic VOC (BVOC), including 2-methyl-3-buten-2-ol (MBO), isoprene (C 5 H 8 ), methyl chavicol (MCHAV, also known as estragole), and a suite of speciated monoterpenes (MT) and sesquiterpenes (SQT), are modeled in each canopy layer as a function of leaf density, light, temperature and vegetation density. For each emitted compound and in each layer, the emission rate is calculated in units of molecules cm −3 s −1 as E b is the basal emission rate in molecules per gram of leaf per second, C L (z) and C T (z) are dimensionless correction factors for light and temperature (Guenther et al., 1995), and the rightmost terms collectively represent the vertically-distributed leaf dry mass in grams of leaf per cubic centimeter. Basal emission rates are prescribed separately for the overstory and understory within the range of values reported for this forest (Bouvier-Brown et al., 2009b, c;Harley et al., 1998;Schade et al., 2000) and are adjusted to optimize model-measurement agreement during the hot period. Temperature and light corrections are taken from the literature (Bouvier-Brown et al., 2009c;Guenther et al., 1995;Harley et al., 1998) and are calculated as a function of the imposed canopy microclimate in each layer. Speciated MT emissions include α-pinene, β-pinene, limonene, 3-carene, myrcene, camphene, terpinolene, α-terpinene and γ -terpinene. SQT include α-bergamotene (ABERG), βcaryophyllene (BCARY), α-farnesene (AFARN) and unspeciated SQT (USQT). USQT are a proxy for the non-speciated SQT observations reported by Bouvier-Brown et al. (2009a, c). Soil NO emissions are a function of temperature assuming dry soil (Yienger and Levy, 1995;Williams et al., 1992) with a basal NO emission factor of 3 ngN m −2 s −1 . This gives temperature-corrected NO emission fluxes of 3.0 and 2.4 ngN m −2 s −1 for the hot and cool periods, respectively. Direct observations of soil NO fluxes are not available for BEARPEX-2007; however, modeled values are consistent with preliminary results from BEARPEX-2009measurements (E. Browne, personal communication, 2010. Deposition is calculated for 35 species using a standard resistance parameterization (Wesely, 1989;Zhang et al., 2003;Wesely and Hicks, 2000) and includes transfer across the laminar sublayer, stomatal and cuticular uptake and ground deposition. The stomatal resistance calculation includes environmental corrections for light extinction, temperature and vapor pressure deficit (Zhang et al., 2003) and is optimized to agree with observationally-constrained, "top-down" calculations of stomatal resistance during BEARPEX-2007. Cuticular resistances are based on standard values used in other models; however, as these are not well-constrained by observations, we will note in the discussion when uncertainties in non-stomatal deposition might influence interpretation of modeled fluxes. Deposition resistances (R dep (z)) for each species are calculated separately for the overstory and understory in each layer and scaled by LADF to give a first-order loss rate constant within each vertical layer: Multiplication of k dep (z) by a concentration yields the firstorder depositional loss rates. Chemistry in CAFE is based on a subset of the Master Chemical Mechanism (MCM) version 3.1 (http://mcm.leeds. ac.uk/MCM/) that includes all reactions stemming from oxidation of MBO, isoprene, α-pinene, β-pinene, propanal (C 2 H 5 CHO) and methane. MCM names and structures for key species mentioned in this study are listed in Appendix A. Our mechanism also includes a number of modifications and additions to the base MCM, most of which are described in the companion paper. Notably, CAFE incorporates a suite of 36 additional reactions for the initial oxidation of monoterpenes (excluding α-pinene and β-pinene), sesquiterpenes and MCHAV by OH, O 3 and NO 3 . Products of these reactions include small oxidized VOC with yields as reported by laboratory oxidation studies (Atkinson and Arey, 2003;Lee et al., 2006a, b), hydroxyl (OH) radicals from ozonolysis reactions, also with literature-reported yields (Atkinson and Arey, 2003;Lee et al., 2006a), and the generic peroxy radicals MTO 2 and SQTO 2 . The latter react with NO, HO 2 and RO 2 to form the species MTOX and SQTOX, which represent first-generation oxidation products of MT and SQT. Since these products are likely semi-volatile and their detailed chemistry is presently unknown, MTOX and SQTOX are given a deposition velocity equal to that of nitric acid (near the aerodynamic limit) and do not undergo further reactions. CAFE also incorporates isoprene dihydroxyepoxide chemistry (Paulot et al., 2009c) and assumes that the epoxide (IEPOX) deposits at the aerodynamic limit.
For the current investigation, we implement one significant modification to the default mechanism described in Wolfe and Thornton (2011). When VOC emissions are high (i.e. during the hot period), an additional "enhanced OH recycling" mechanism is required to bring modeled OH values into agreement with measurements. We employ a mechanism of the type where α is a stoichiometric constant. These reactions, listed in Table 3, are implemented only for first-generation MBO and isoprene-derived peroxy radicals (RO 2 = MBOAO 2 , MBOBO 2 , ISOPAO2, ISOPBO 2 , ISOPCO 2 and ISOPDO 2 ). The reaction "products" are those of the decomposition of the corresponding RO radicals, which are explicitly tracked. Values for α and k rec are tuned to optimize model-measurement agreement for OH and HO 2 ; for the current study, we choose α = 2.6 and k rec = 4.5×10 −11 cm 3 molec −1 s −1 . We discuss and evaluate the consequences of this mechanistic change further in Sect. 3.2.
Advection is treated as a simple mixing process in each model layer with a mixing rate constant (k mix ) of 0.3 h −1 (Dillon et al., 2002;Perez et al., 2009): Advection concentrations (C a ) are set constant throughout the model domain but are different for the hot and cool periods (Table 4). This scheme maintains reasonable concentrations for species that would otherwise build up to unreasonable values or decay below measured values during integration. Advection thus allows us to better constrain CAFE to BEARPEX-2007 observations. We will note when this term influences the interpretation of results, though it generally does not influence modeled exchange velocities as the timescale for advection is relatively long (several hours). For each period, meteorological observations (Table 2) are used to initialize diffusion parameters, emission rates, deposition velocities and chemical rate constants, which are held constant throughout a model run. Chemical concentrations are initialized with the same values used in the advection scheme (Table 4) and are thus constant throughout the vertical; these values are chosen to optimize model agreement with observations. Integration is accomplished via operator splitting using a Crank-Nicolson scheme to solve the diffusion operator and a forward Euler scheme for the chemical operator (Jacobson, 2005). Soil NO emission and ground deposition are incorporated into the diffusion operator, while canopy emissions, deposition and advection are represented in the chemistry operator. The model is run for two hours, which is sufficient time for relaxation of exchange velocity profiles. Fluxes and exchange velocities are calculated from concentration profiles at the end of a model run via

Results and discussion
In what follows, we compare CAFE model output with observations from the BEARPEX-2007 field campaign. Table S1 in the Supplement lists averaged chemical observations for each period. All quoted measurement values represent the means and standard deviations of 30-min averaged data. The hot period is typified by relatively high concentrations of BVOC, HO x , O 3 and oxygenated hydrocarbons and lower levels of NO 2 and acyl peroxy nitrates (APN = PAN + PPN + MPAN + . . . ); cold period data demonstrate the opposite trends. Differences in local atmospheric composition between the hot and cool periods are largely driven by temperature (as opposed to wind direction, for example), which controls emission rates and subsequent  photochemistry. Model results are examined with a particular focus on BVOC, RO x (= OH + HO 2 + RO 2 ), hydrogen peroxides, O 3 and reactive nitrogen (NO y ). We evaluate modifications implemented in the CAFE model (e.g. OH recycling) and provide an assessment of the relative contributions of deposition, emission and chemistry to above-canopy chemical fluxes for key species. Unless otherwise specified, model results discussed below are extracted from two "base" model runs, one each for the hot and cool periods. The base run for the hot period is carried out with the OH-recycling mechanism, while the base cold period run does not include OH-recycling. The reasons for this choice are detailed in Sect. 3.2. Reproducing observed concentrations is important for examining chemical contributions to fluxes, but we caution that the model is not strictly tailored towards reproducing all aspects of the chemistry (e.g. diurnal cycles) or, more importantly, horizontal transport. Modeled mixing ratios are, in a sense, constrained to observations via the advection term and the initial/advection concentrations. We run CAFE in this fashion because our primary goal is to understand the observed fluxes, though we also point out other interesting features in the model-measurement comparison when they arise. Modeled concentration and fluxes should not be taken as representative of daily or seasonal "average" conditions, but rather as mid-day "snapshots" from the two periods. An extended comparison table of modeled and measured concentrations can be found in Table S2 of the Supplement.

VOC
Within and immediately above the forest, concentrations of primary BVOC (MBO, isoprene, MCHAV, MT and SQT) are controlled by relative rates of emission and oxidation. Calculation of "bulk canopy" emission rates provides a means for validation of vertically-resolved emissions. Taking MBO as an example: integration of the hot-period MBO emission rate over the canopy height gives a bulk emission rate of 5.2 × 10 11 molec cm −2 s −1 (1.9 mgC m −2 h −1 ), or a boundary-layer average of 6.6 × 10 6 molec cm −3 s −1 . These rates are within the range of previous MBO flux measurements at BFRS (Baker et al., 1999;Schade et al., 2000) and of values employed by other models (Perez et al., 2009;Steiner et al., 2007). Though our emissions estimates generally agree with other literature values, the standard emission parameterization does have limitations. Previous work at BFRS has shown that basal emission rates can vary with temperature history and other factors (Gray et al., 2003(Gray et al., , 2006 and that tree-to-tree variability in emission rates can be substantial (Bouvier-Brown et al., 2009c).
Isoprene is not emitted in significant quantities from Ponderosa pine, Manzanita or Ceanothus (N. Bouvier-Brown, personal communication, 2009), but it can originate from less abundant vegetation within the forest stand and upwind, particularly Black Oak. Although direct measurements of above-canopy isoprene fluxes have not been conducted at BFRS, early isoprene gradient measurements and relaxed eddy accumulation observations in 1998 and 1999 indicated no significant emissions from the BFRS fetch (Dreyfus et al., 2002;Goldstein et al., 2001). Analysis of mixing ratio diurnal profiles at this site have determined that isoprene is primarily advected from a band of Oak located 30-40 km upwind (Dreyfus et al., 2002). The current construction of CAFE is unable to simultaneously reproduce the concentrations of isoprene and its main oxidation products, methyl vinyl ketone (MVK) and methacrolein (MACR), solely through our advection scheme. Thus, in addition to advecting isoprene at a rate of 1 ppbv h −1 , we invoke a substantial emission rate of isoprene (∼40% of the MBO emissions). While local (e.g. <500 m upwind) isoprene emissions are probably smaller than this, our isoprene emission rate is nearly identical to that used in the 4 km × 4 km grid cell of a three-dimensional model that contains BFRS (Steiner et al., 2007). The vertical profile of isoprene, and potentially its oxidation products, will depend somewhat on the nature of its sources (i.e. emission vs. advection). A small set of in-canopy isoprene gradients measured near the end of the BEARPEX-2007 campaign (after our cool period) suggest that in-canopy isoprene mixing ratios can exceed abovecanopy values by as much as a factor of 2 (J. Gilman, personal communication, 2010), but it remains unclear if this gradient can be attributed to local emissions. We will note when this issue affects our conclusions.
Since modeled mixing ratios of locally-emitted BVOC are primarily a function of the rates of emission and chemical loss (assuming roughly homogenous upwind emissions and a chemical lifetime of a few hours), it is worthwhile to compare modeled and measured concentrations of these (Table 5). MBO mixing ratios are reproduced to within 6% during both hot and cool periods, suggesting that the radiation and temperature adjustments are accurate for MBO emissions. MCHAV and total MT are predicted to within 10% during the hot period but are under-predicted during the cool period by 60-70%, while total SQT are over-predicted by 150% during the cool period. These errors could stem from the temperature corrections for emission rates, which become increasingly important at lower temperatures, or from deficiencies in modeled vertical mixing (note that the inlet for MT and SQT observations was moved from 1.5 m to 9.2 m between the hot and cool periods).
Even though total MT concentrations are well reproduced during the hot period, modeled terpene speciation differs from observations. The model generally under-predicts βpinene and 3-carene and over-predicts myrcene, camphene, terpinolene, α-terpinene and γ -terpinene (Table S2). Such discrepancies may arise from inaccurate estimates of emission speciation. Though we use the best estimates from leaflevel measurements (Bouvier-Brown et al., 2009c), terpene emissions are subject to plant physiological and environmental conditions that are not easily modeled. The relative terpene speciation has little impact on our conclusions regarding the chemical contribution to trace gas fluxes.
The terpene oxidation tracers MTOX and SQTOX show roughly the same seasonal trend as their VOC precursors (Table 5). Despite a fast deposition velocity, canopy-top concentrations of MTOX and SQTOX build up to 101 and 41 pptv, respectively, during the hot period. Many of the compounds represented by these tracers will contain alkenyl moieties and thus may still participate in oxidative chemistry. For both periods, near-surface MTOX and SQTOX gradients (not shown) match previous CAFE model results (Wolfe and Thornton, 2011).
In addition to the speciated single-height measurements, the BEARPEX-2007 dataset also includes vertical profiles of several classes of VOC acquired via proton-transfer mass spectrometry (PTR-MS). Details regarding instrumentation and measurement setup can be found elsewhere (Bouvier-Brown et al., 2009b;Holzinger et al., 2005). Figure 2 compares modeled BVOC profiles to four sets of PTR-MS measurements: total monoterpenes ( MT), MCHAV, the sum of MBO and isoprene and the sum of MVK and MACR. For clarity, observations and model results are presented from the hot period only and have been normalized to their canopy-top values; modeled and measured profiles exhibit similar normalized gradients for both the hot and cool periods. Modelmeasurement agreement is generally quite good, though the model under-predicts gradients of MT and MCHAV in the lower canopy. Potential explanations include an unidentified emission source near the ground such as decaying pine needles, as suggested by Stroud et al. (2005), or inefficient turbulent mixing in the lower canopy, which could lead to a buildup of BVOC emitted from the understory. The existence of a "weakly coupled" layer near the ground -or at least slower vertical mixing than modeled in CAFE -would also be consistent with model-measurement comparisons of NO 2 and PAN gradients (see Sect. 3.5). A lack of wind penetration into the deep canopy would seem the most likely cause (as opposed to a thermal inversion). To test this hypothesis, we conducted a simple sensitivity test where we uniformly reduced the eddy diffusivities (K-values) below z/ h = 0.5 (5 m). Results from this test (not shown) reveal that reducing K-values by a factor of three greatly improves agreement between modeled and measured MT gradients but increases the canopy residence time to 5 min.
The modeled MBO + isoprene gradient agrees fairly well with observations, increasing ∼20% between canopy top and the forest floor. This profile would be less steep if isoprene were primarily advected in CAFE, as the isoprene profile would then be more vertical and observed isoprene mixing ratios are 25-33% of MBO + isoprene at noon (Table 5). The mean observed MBO + isoprene mixing ratio at z/ h = 0.15 (1.5 m) is consistently 15% lower than that at 6.1 m (Fig. 2c). This feature persists even in individual 30-min gradient observations. As MBO (and isoprene, in CAFE) are primarily emitted from the overstory, this feature would be consistent with a depositional sink of MBO and isoprene near the ground (Stroud et al., 2005).
Profiles of the sum of MVK and MACR, which are firstgeneration oxidation products of isoprene, are fairly vertical in both the model and measurement (Fig. 2d), though the modeled profile shows a slight enhancement in the canopy due to production. Previous studies have suggested that MVK and MACR should also deposit to the canopy/ground with a deposition velocity similar to that of ozone (Stroud et al., 2005;Zhang et al., 2003;Karl et al., 2010). We do not include deposition of these compounds in CAFE since the observed mean profiles do not suggest strong deposition of oxidized VOC in this forest. The hydroxyl radical (OH) is the primary daytime oxidant for most VOC in the troposphere. OH reactivity (τ −1 OH ), or inverse OH lifetime, is defined as the sum of all OH loss rates divided by the OH concentration: Here, k i is the second-order rate constant for reaction of OH with species i having concentration C i . OH reactivity was measured directly during BEARPEX-2007 following the approach described in Mao et al. (2009) and is useful for constraining both VOC inventories and steady-state calculations of oxidant concentrations. Figure 3 compares model calculations of τ −1 OH with observed values. During the hot period, modeled (12.3 s −1 ) and measured (12.4 ± 2.0 s −1 ) τ −1 OH are in excellent agreement. About 63% of the modeled τ −1 OH is attributed to primary BVOC, with another 22% due to reactions with HCHO, CO, CH 4 and the first-generation oxidation products of isoprene (MVK and MACR) and MBO (IBU-TALOH and HOCH 2 CHO). The remaining 15% ("other") includes ∼300 reactions, each of which comprise <1% of τ −1 OH . During the cool period, modeled τ −1 OH (3.6 s −1 ) is lower than observations (6.8 ± 1.2 s −1 ) by almost a factor of 2. Model underestimates of MT and MCHAV mixing ratios during this period are not sufficient to explain this discrepancy, and the nature of this missing OH reactivity remains unclear. These results, including the missing OH reactivity during the cold period, are consistent with observationally- constrained bottom-up estimates of τ −1 OH . The latter study also demonstrated that measured anthropogenic VOC are a negligibly small contribution to τ −1 OH at BFRS.
Some of the missing OH reactivity might be attributed to a missing source of formaldehyde (HCHO). During the cool period, CAFE predicts HCHO mixing ratios of ∼1.3 ppbv, while measurements indicate a noontime mean of 12.5 ± 4.0 ppbv (Table S2); HCHO observations were not available during the hot period. Increasing HCHO mixing ratios to match observations (by raising the initial/advection HCHO concentrations) brings the modeled OH reactivity to 6 s −1 , which is within the range of observations. Maintaining this level of HCHO in the model leads to a 50% overprediction of hydroperoxy radical (HO 2 ) and hydrogen peroxide (H 2 O 2 ); OH increases by 15%. The sources of the elevated HCHO mixing ratios observed during the cold period are presently unknown but may be linked to oxidation of yet-unidentified reactive BVOC inferred from previous observations at BFRS (Choi et al., 2010;Holzinger et al., 2005). As constraining HCHO to measured values does not noticeably perturb the exchange velocities of key species in the model and because of a lack of measurements during the hot period, we retain the CAFE-predicted HCHO values for consistency between the two periods. A detailed evaluation of BEARPEX-2007 HCHO observations may be found elsewhere (Choi et al., 2010).

RO x
Cycling of hydrogen oxide radicals is driven by VOC and nitric oxide (NO). The sequence begins with reaction of OH and VOC to produce an organic peroxy radical (RO 2 ). Subsequent reaction of RO 2 with NO produces NO 2 and an alkoxy radical (RO). Typically, the latter reacts with O 2 to yield a hydroperoxy radical (HO 2 ) and a closed-shell aldehyde or ketone. In the CAFE mechanism, the latter two processes are combined. OH is regenerated upon reaction of HO 2 with NO to form NO 2 .
As a result of this cycling, we define the chemical families HO x = OH + HO 2 , RO x = HO x + RO 2 and NO x = NO + NO 2 . Partitioning within the RO x and NO x families is thus coupled by VOC abundance and reactivity with OH. Moreover, cross-reactions between RO x and NO x produce longer-lived reactive nitrogen species, the forest-atmosphere exchange of which can be sensitive to vertical gradients in this chemistry. Figure 4 depicts modeled profiles of OH, HO 2 and RO 2 radicals. For each period, the model was run both with and without the enhanced OH recycling mechanism discussed in Sect. 2.2. Both periods show small positive (increasing with height) in-canopy HO 2 gradients of ∼5%. The OH mixing ratio increases by ∼10% between the ground and the top of the canopy in the hot period and by ∼40% in the cool period. The relative gradients are mostly unaffected by the enhanced OH-recycling mechanism, though OH does exhibit a slight bulge maximizing at z/ h = 1.4 during the hot period with enhanced OH recycling, and RO 2 is ∼20% higher within the canopy than above for the same scenario. Enhanced OH recycling is required for replicating observations during the hot period. Excluding this mechanism leads to under-prediction of noontime OH by a factor of 6 and of HO 2 by ∼25%; RO 2 was not measured. Preliminary data from BEARPEX-2009 suggests measured OH mixing ratios from the 2007 campaign may be overestimated, perhaps by as much as a factor of 2 . Even in the face of such a systematic error, modeled values would still be too low by a factor of three and we would still require enhanced OH recycling.
Model-measurement mismatch of OH is a recurrent issue in investigations of RO x chemistry under conditions where BVOC such as isoprene are a dominant source of RO 2 (Hofzumahaus et al., 2009;Lelieveld et al., 2008;Thornton et al., 2002;Ren et al., 2008;Martinez et al., 2003;. Many of these studies, and others, have proposed mechanisms to augment radical production and propagation, including reduction in the formation rate of isoprene-derived organic hydroperoxides and/or enhancement of their photolysis rates (Thornton et al., 2002), additional production of OH during reactions of isoprene-derived first-generation RO 2 with HO 2 (Lelieveld et al., 2008;Thornton et al., 2002), inclusion of an unknown species "X" that reacts with RO 2 and HO 2 with the same efficacy as NO (Hofzumahaus et al., 2009), and RO 2 isomerization and decomposition (Peeters et al., 2009;Da Silva et al., 2010). We tested each of these mechanisms separately in CAFE, but found that no single mechanism could adequately reproduce observed HO x partitioning and abundance simultaneously with other key indicators, such as oxidized VOC abundance and speciation. For example, incorporation of the isoprene hydroxyperoxy radical isomerization/decomposition mechanism -as implemented in Stavrakou et al. (2010) with an OH yield of 3 from the photolysis of hydroxyperoxy aldehyde productsleads to a 30% over-prediction of HO 2 but a factor of three under-prediction of OH in the hot period. Furthermore, the postulated isomerization requires an allylic radical, thus firstgeneration MBO-derived peroxy radicals will not undergo analogous reactions. Using measured OH reactivity and concentrations, and assuming OH is in steady state (i.e. production equals loss), we estimate an observationally-constrained gross OH production rate (P OH ) of ∼4 pptv s −1 for noontime conditions during the hot period. Without OH recycling, modeled P OH for the hot period is ∼0.7 pptv s −1 and is mainly driven by O 3 photolysis and reaction of HO 2 with NO. As the model accurately reproduces measured OH reactivity during the hot period (Fig. 3), we conclude that the under-prediction of OH stems from inefficient recycling and/or excessive termination by RO x cross-reactions. Inclusion of the tuned OH recycling mechanism (Table 3) brings modeled OH and HO 2 to within the range of observations and increases RO 2 by a factor of ∼2.5. HO 2 is mostly conserved in the enhanced OH-recycling mechanism, thus the increase in modeled HO 2 is primarily due to a larger source from RO 2 reactions with NO.
Another potentially important OH source is ozonolysis of highly-reactive VOC not included in our emission inventory (Goldstein et al., 2004;Holzinger et al., 2005;Faloona et al., 2001). Holzinger et al. (2005) estimated that an average incanopy O 3 reaction rate of 5.25 × 10 8 molecules cm −3 s −1 , or 25 pptv s −1 , would be required to sustain the chemical contribution to in-canopy ozone fluxes inferred by previous studies (Goldstein et al., 2004;Holzinger et al., 2005;Kurpius and Goldstein, 2003). Given that our missing P OH is ∼3.3 pptv s −1 , an average OH yield of 13% from these reactions would be sufficient to sustain measured OH levels in the canopy. The resulting OH concentration from such a source, however, would lead to model overestimates of HO 2 and likely RO 2 . That is, such a source would still imply an incomplete understanding of RO 2 /HO 2 chemistry.
Our enhanced OH-recycling mechanism is similar to a blending of those proposed by Lelieveld et al. (2008) and Peeters et al. (2009). Our mechanism ties OH recycling to RO 2 + HO 2 reactions, but it is an additional process in competition with the peroxide-forming channel. The mechanism also simultaneously converts the primary MBO and isoprenederived RO 2 radicals into the relevant oxidation products as if passing through the respective RO radicals. Essentially, it is an enhanced RO 2 decomposition that yields OH and oxidized VOC but has little net effect on HO 2 . Failure to incorporate RO 2 destruction in the enhanced OH recycling mechanism leads to unrealistic RO 2 concentrations (>300 pptv), which in turn results in overestimation of several oxidation products -such as glyoxal and acetone -and unreasonably low NO/NO 2 ratios. With our enhanced OH recycling mechanism, model results are consistent with RO 2 and NO/NO 2 values derived from observationally-constrained steady-state calculations for this site Day et al., 2008) and with observations of total peroxy radicals at other forested locations (Cantrell et al., 1992;Qi et al., 2005). Furthermore, small-chain BVOC oxidation products agree reasonably well with BEARPEX-2007 observations (Table S2), though mixing ratios of these are also influenced by advection. Previous studies at BFRS have provided evidence for a temperature-dependent HO x source (Day et al., 2008;Farmer and Cohen, 2008). Our enhanced recycling mechanism is also consistent with this observational evidence as the OH production rate decreases with decreasing temperature by virtue of its reliance on RO 2 formed from BVOC.
In contrast to the hot period, modeled OH agrees with observations during the cool period without the need for additional OH recycling, while HO 2 is somewhat underpredicted. As CAFE underestimates measured OH reactivity by a factor of 2 during the cool period, however, this agreement is likely artificial. Incorporating enhanced OH recycling during the cool period leads to overestimation of OH by a factor of ∼2. Thus, by constraining modeled OH reactivity to the measured value and assuming the missing reactivity is caused by a non-methane hydrocarbon that is not MBO or isoprene, model-measurement agreement of OH concentrations during the cool period can be achieved when employing the enhanced recycling mechanism. As this result ultimately depends on the nature of the missing reactivity, and as OH, HO 2 , and RO 2 abundances are reasonably predicted by CAFE during the cool period without the recycling mechanism, we leave this issue for future investigation.

Peroxides
In high-RO x and high-VOC environments, peroxide formation is considered a key radical termination step: In the case of isoprene oxidation, further reaction of firstgeneration ROOH with OH can generate dihydroxyepoxides (Paulot et al., 2009c) or carbonyl-containing oxidized VOC: Comparison of modeled and measured peroxides/epoxides thus provides an additional indirect check on RO x abundance and chemistry in CAFE. BEARPEX-2007 observations include both hydrogen peroxide (H 2 O 2 ) and the sum of firstgeneration isoprene hydroxyhydroperoxides (ISOPOOH, see Appendix A) and isoprene dihydroxyepoxides (IEPOX). As shown in Table 5, ISOPOOH+IEPOX is over-predicted by a factor of 2 during the hot period, suggesting that production is too fast and/or losses -which include reaction with OH, photolysis and deposition -are too slow. Production of ISOPOOH depends on mixing ratios of HO 2 and first-generation isoprene-derived RO 2 . Observational constraints are only available for the former, thus we cannot rule out a model overestimation of isoprene-RO 2 , though this seems unlikely. Modeled OH and isoprene (RO 2 precursors) agree with measurements, and the over-prediction persists even without our enhanced OH-recycling mechanism -when OH and RO 2 concentrations are substantially lower. Isomerization of isoprene-RO 2 (Peeters et al., 2009) could mitigate production of ISOPOOH; however, as we noted in Sect. 3.2, this mechanism also leads to over-prediction of HO 2 in our model.
It is more probable that modeled sinks of ISOPOOH and/or IEPOX are too slow. The lifetimes of ISOPOOH against OH via Reactions (15) and (16) are 0.7 h and 8.5 h, respectively. Thus, 92% of this reaction proceeds through the epoxide-forming channel, and IEPOX comprises 57% of the ISOPOOH+IEPOX family during the hot period. IEPOX reacts with OH with a lifetime of ∼3 h using the recommended rate constant from Paulot et al. (2009c). A factor of 10 increase in the rate constants for reaction of OH with either ISOPOOH (via Reaction 16) or IEPOX reduces the model over-estimate to 29% or 17%, respectively. Deposition of ISOPOOH is currently implemented with a deposition velocity of ∼1.6 cm s −1 (Hall and Claiborn, 1997), while IEPOX is forced to deposit at the aerodynamic limit. An increase in deposition of ISOPOOH would reduce the model overestimate but is not a sufficient explanation on its own. Faster photolysis of ISOPOOH is also not likely a viable solution, as the rate for this process is ten times slower than reaction with OH in the base scenario. One potentially important sink for IEPOX is uptake to particles, which is not currently included in CAFE. Optimal model-measurement agreement during the hot period -assuming rates for other sinks are reasonably estimated -would require 88% of IEPOX be lost to particles.
In contrast to the hot period, ISOPOOH + IEPOX is underpredicted by 51% during the cool period. This may indicate an under-prediction of isoprene-RO 2 for this scenario, though we note again that both isoprene and OH are within the range of observations. At present, it is difficult to reconcile the large differences in model-measurement agreement between the two scenarios. One potential explanation may lie in advection. These compounds have been assigned initial/advection concentrations of 0 to simplify comparison to observations, but it is very likely that their "advection concentrations" are higher than we ascribe due to the upwind isoprene source. Increasing advection would improve modelmeasurement agreement during the cool period but degrade it during the hot period; thus, an increased advection source must also be coupled with a temperature-dependent sink.
Modeled H 2 O 2 concentrations agree with observations to within 3% and 16% for the hot and cold period, respectively. As noted in the companion paper (Wolfe and Thornton, 2011), CAFE forces H 2 O 2 to deposit at the aerodynamic limit by setting the effective Henry's Law coefficient (H * ) used in the deposition parameterization to 1 × 10 14 M atm −1 . This yields above-canopy exchange velocities of −4 to −5 cm s −1 . If we instead use the literaturerecommended H* of 1 × 10 5 M atm −1 (Seinfeld and Pandis, 2006), H 2 O 2 concentrations are over-predicted by as much as 50%. Though H 2 O 2 mixing ratios are somewhat dependent on our choice of initial/advection concentrations (Table 4), this finding is consistent with previous field studies that have reported diffusion-limited H 2 O 2 deposition over forests, much faster than predicted by the Wesely (1989) parameterization (Ganzeveld et al., 2006;Hall and Claiborn, 1997). Recent measurements at BFRS also suggest transportlimited H 2 O 2 deposition (Paulot et al., 2009a), a finding supported by our model results. Our choice to increase H * to match the observed H 2 O 2 exchange velocity does not necessarily imply that the molecular mechanism controlling H 2 O 2 surface loss is related to solubility; rather, we view the cuticular resistance (which is controlled partly by H * ) as a tunable parameter that could represent any number of yetundiscovered chemical or physical uptake processes.

Ozone
Deposition of O 3 is a major concern due to its impact on plant tissues (Darrall, 1989), which can reduce carbon sequestration (Sitch et al., 2007), enhance emissions of oxidized VOC  and alter uptake of other gas-phase species (Karl et al., 2010). As a terminal sink, deposition also influences the lifetime of gas-phase O 3 near the surface. Furthermore, ozonolysis of reactive BVOC can produce OH and oxygenated VOC within and immediately above the canopy (Ciccioli et al., 1999;Holzinger et al., 2005), stimulating gas-phase oxidative chemistry and growth of secondary organic aerosol (Li et al., 2011;Paulson et al., 1998). This chemistry may be fast enough to alter the net forest-atmosphere flux of both O 3 and BVOC (Goldstein et al., 2004;Kurpius and Goldstein, 2003;Stroud et al., 2005). Thus, characterization of the underlying mechanisms of forest-atmosphere O 3 exchange is critical for assessing both O 3 -induced ecosystem damage and our understanding of emissions and chemistry in this environment. Figure 5 compares model and measured O 3 concentrations, fluxes and exchange velocities. Mixing ratios agree reasonably well, decreasing slightly with height in the canopy region. Modeled downward (negative) fluxes and exchange velocities fall within the variability of observations but tend to under-predict mean values by ∼20% for both periods. CAFE successfully predicts the observed 20% increase in O 3 fluxes between the hot and cool period, which is due to an increase in stomatal conductance accompanying the reduced vapor pressure deficit and temperature (Table 2). This behavior is counter to the generally-positive correlation between temperature and O 3 fluxes observed on longer timescales at BFRS but is within the variability of measured O 3 fluxes from a 6-yr dataset (Fares et al., 2010a).
Previous work at BFRS has provided evidence that both deposition (stomatal and non-stomatal uptake) and in-canopy reactions with biogenic emissions can influence O 3 fluxes (Goldstein et al., 2004;Kurpius and Goldstein, 2003;Fares et al., 2010a). It is thus prudent to examine all processes contributing to O 3 fluxes, defined by the various terms in the mass balance equation (Eq. 1). Figure 6a shows verticallyresolved instantaneous rates for all processes during the hot period. The model predicts that deposition is the dominant process within the canopy region. Other rates are small but consistent with expected chemical behavior, which derives primarily from the NO x -O 3 equilibrium: Gross chemical O 3 production tracks light attenuation in the canopy, as it is rate-limited by NO 2 photolysis. Ozone photolysis controls gross chemical loss in the top half of the canopy, while reaction with soil-emitted NO dominates near the ground. The net chemical tendency (P +L) changes sign halfway through the canopy. These results are consistent with those from other canopy models (Stroud et al., 2005). Contributions from advection (A) and storage (∂C ∂t) are small. Rearrangement and integration of Eq. (1) shows that the vertical flux at any height, F (z), is the sum of the ground-up integrals of the rate of each process: Here, we group contributions into surface and "chemical" processes (for O 3 , E(z) = 0); these groups could also be thought of as "heterogeneous" and "gas phase." The integral over height of ∂C ∂t (the last term) is, by definition, the "storage" term employed in interpretation of flux observations (Rummel et al., 2007;Wolfe et al., 2009). Calculation of fluxes by this method yields the same values as those computed via Eq. (6), and normalization of each term by the modeled mixing ratio at any height gives the corresponding component of the exchange velocity. CAFE predicts that surface deposition controls the vertical flux of ozone (Fig. 6b), with chemistry inducing a slight positive slope on the O 3 exchange velocity profile above the canopy; the model yields similar results for the cool period. At face value, this result appears to refute previous claims of a substantial chemical ozone flux at BFRS; however, as discussed below, modeled ozone fluxes in CAFE are subject to potentially significant uncertainties in non-stomatal deposition and reactive BVOC emissions.
Since deposition dominates the modeled O 3 exchange, it is prudent to examine the contribution of various parameterized deposition pathways as enumerated in Table 6. Stomatal uptake, which is constrained by independent calculations based on observed latent heat fluxes (Wolfe and Thornton, 2011), accounts for 46% and 59% of the modeled exchange during the hot and cool periods, respectively. Non-stomatal deposition, which includes losses to both leaf cuticles and the ground, comprises the remaining 54% and 41%. Though the inclusion of non-stomatal deposition brings modeled ozone exchange velocities to within 24% and 20% of measured values, this agreement may be artificial as the true magnitude of non-stomatal deposition is not known. Constraints for nonstomatal deposition in CAFE are taken from a "big leaf" resistance model (Zhang et al., 2002(Zhang et al., , 2003 that assumes similar values across a fairly wide swath of ecosystems, and unlike the stomatal component, there is no simple way to validate this parameterization against observations. The possible influence of intra-canopy BVOC chemistry on observed above-canopy O 3 fluxes was not considered during the development of these parameterizations, thus the resulting nonstomatal deposition rates may be artificially high. Leaf and branch-level enclosure experiments have observed negligible cuticular ozone deposition to Loblolly pine  and citrus trees (Fares et al., 2010b). If our parameterization of cuticular and ground deposition over-estimates the magnitude of these processes at BFRS, then we are likely still missing a considerable O 3 sink within the forest. In the extreme case where non-stomatal deposition is ignored, CAFE underpredicts observed O 3 fluxes by 65% and 53% for the hot and cool periods, respectively. Modeled fluxes are also sensitive to prescribed values for eddy diffusion coefficients (Eq. 6); however, sensitivity tests indicate that diffusion would need to be significantly faster than currently modeled to increase the O 3 flux (Wolfe and Thornton, 2011).
Within the current model framework, in-canopy O 3 + BVOC reactions are not of sufficient magnitude to influence O 3 fluxes, seemingly at odds with inferred non-depositional O 3 fluxes at BFRS (Kurpius and Goldstein, 2003;Goldstein et al., 2004). In these studies, the authors note that emissions of very reactive BVOC -which can drive chemical O 3 fluxes -may not be included in current emission inventories. CAFE already contains emissions of highly reactive SQT that account for some of the missing ozone reactivity inferred previously (Kurpius and Goldstein, 2003;Goldstein et al., 2004). Bringing our model results into agreement with these studies, however, would require a substantial increase in emissions of the highly reactive SQT species or other yet-unmeasured BVOC. Such missing emissions have also been postulated from observations of "missing" OH reactivity at a forest in northern Michigan (Di Carlo et al., 2004) and from anomalously high HCHO concentrations measured at BFRS (Choi et al., 2010). For BEARPEX-2007, CAFE reproduces the observed above-canopy OH reactivity during the hot period (Fig. 3), when BVOC emissions should be highest. Thus, to affect ozone fluxes, unidentified emissions must react preferentially with O 3 in the canopy, similar to the SQT species β-caryophyllene and α-humulene (Bouvier-Brown et al., 2009c). The effects of very reactive BVOC on O 3 chemistry should be localized to the canopy or near-leaf airspace, consistent with the observation of markedly different terpene speciation between branch enclosure and ambient measurements (Bouvier-Brown et al., 2009a). Considerable non-stomatal ozone fluxes have been observed at several other forests (Hogg et al., 2007;Coe et al., 1995;Rondon et al., 1993), but whether such fluxes are driven by surface or gas-phase processes remains an open question. Targeted model sensitivity studies together with further in situ observations could shed additional light on such issues. Understanding the fate of O 3 in the forest must continue to be a priority, as the questions raised here are directly relevant to ecosystem health, aerosol formation and RO x chemistry in this environment.

Reactive nitrogen
The reactive nitrogen (NO y ) family encompasses a wide spectrum of atmospheric oxidized nitrogen compounds, including NO x , acyl peroxy nitrates (APNs), alkyl nitrates (ANs) and nitric acid (HNO 3 ), among others. Primary NO x sources in the troposphere include both anthropogenic (e.g. combustion and agriculture) and natural (e.g. soil and lightning) emissions (Jaeglé et al., 2005). The higher oxides of nitrogen are formed via reactions of NO x with RO x : We restrict our analysis to these four classes since they comprise the bulk of NO y at BFRS (Day et al., 2009;Ren et al., 2010) and their formation mechanisms are reasonably -if not yet quantitatively -understood. Table 7 compares modeled concentrations of NO y components to observations. Overall, the model is in decent agreement with measured NO y (= NO + NO 2 + PN + AN + HNO 3 ) and the temperaturedependence of NO y speciation, though PN (comprised primarily of APNs in the model) are slightly over-predicted and AN are slightly under-predicted. Speciation is discussed in more detail below. NO y measured during BEARPEX-2007 (1-1.3 ppbv) is about half that reported for previous years (∼2.5 ppbv), but the relative speciation is similar (Day et al., 2009;Murphy et al., 2006). This decrease is consistent with ∼10 ppbv lower average O 3 concentrations during BEARPEX-2007. Advection is important for maintaining reasonable concentrations of many NO y species, particularly NO 2 ; however, advection does not significantly affect modeled vertical fluxes.
Forest-atmosphere exchange of reactive nitrogen continues to be a significant uncertainty in assessing the influence of anthropogenic nitrogen emissions on forest productivity (Magnani et al., 2007;Thomas et al., 2010) and regional air quality (Steiner et al., 2006). Quantifying dry nitrogen deposition to forests remains a challenge for several reasons. Deposition velocities may vary by an order of magnitude for different classes of NO y (Farmer and Cohen, 2008;Turnipseed et al., 2006;Horii et al., 2006). As a result, deposition can alter the relative partitioning of remaining gas-phase NO y , which in turn affects NO y chemistry and deposition downwind. Furthermore, rapid in-canopy chemical transformations can alter the net forest-atmosphere exchange of NO y species Duyzer et al., 2004;Farmer and Cohen, 2008;Walton et al., 1997;Wolfe et al., 2009). Soil-emitted NO, often a primary NO x source in rural and remote regions , is rapidly converted to NO 2 by reaction with O 3 and peroxy radicals in the canopy (Gao et al., 1991), with implications for measuring the fluxes of NO x components. NO x partitioning within the canopy also affects the fate of APNs, which depends in part on the NO/NO 2 ratio. Oxidation of BVOC can enhance or alter the pathways for production of APNs and ANs, while temperature gradients can influence the decomposition of APNs to NO x , affecting fluxes of both of these components (Farmer and Cohen, 2008;Wolfe et al., 2009). In what follows, we examine the modeled vertical exchange for each class of NO y with a focus on the role chemistry plays in modifying the net abovecanopy flux. Figure 7a displays vertical profiles of NO 2 mixing ratios. NO 2 is lower during the hot period, likely because of faster conversion to HNO 3 and a decreased APN reservoir (Day et al., 2008). Concentrations increase near the ground due to fast conversion of soil-emitted NO via reactions with O 3 and peroxy radicals, as well as relatively enhanced thermolysis of APNs via the reverse of Reaction (R12). The measured NO 2 gradient is steeper than the model during the hot period, which may be symptomatic of a deficiency in our diffusion parameterization or of soil NO fluxes that are higher than our estimate. Gradients in the NO/NO 2 ratio (Fig. 7b) are driven by a balance between the soil NO emission rate, rapid establishment of the NO-NO 2 -O 3 equilibrium (Reactions R9-R11), and diffusion timescales. NO/NO 2 is lower during the hot period because of higher levels of RO 2 , which mainly convert NO to NO 2 via Reaction (R3) but also to ANs via Reaction (R13).

NO x
The mirrored shape of in-canopy flux profiles (Fig. 7cd) reflects the rapid inter-conversion of NO and NO 2 . Near the ground, low radiation results in net conversion of emitted NO to NO 2 , thereby increasing fluxes of NO 2 and decreasing those of NO. Flux profiles turn more vertical at z/ h ∼ = 0.6, when NO 2 photolysis becomes competitive with NO oxidation. The net NO x flux above the canopy is upward and nearly equal in magnitude to the soil NO emission flux, though it is mostly comprised of NO 2 . The canopy reduction factor, defined as the fraction of soil-emitted NO x that does not leave the canopy due to deposition and chemical transformations, is ∼11% for the current study. Because most of the soil-emitted NO x escapes the canopy and because the modeled NO x flux is almost entirely comprised of NO 2 , canopytop NO 2 fluxes likely provide an indirect check on our soil NO emission flux. Unfortunately these observations were not available during BEARPEX-2007, but we may draw a rough comparison to previous studies. During the hot period, CAFE predicts an above-canopy NO 2 exchange velocity of +3.5 cm s −1 , which is 50% higher than the +2.3 cm s −1 observed by Farmer and Cohen (2008) at BFRS in August 2004. Without soil NO flux measurements from either pe-riod, we are unable to conclude whether the higher NO 2 flux predicted for 2007 by CAFE is realistic.

APNs
APNs are a unique class of NO y in that their atmospheric residence time, determined partly by the chemical equilibrium (R12), is highly sensitive to temperature. Peroxyacetyl nitrate (PAN, CH 3 C(O)O 2 NO 2 ) comprises 70-90% of the observed APN budget during BEARPEX-2007 and evolves from a variety of anthropogenic and biogenic VOC precursors. Notable minor APNs include peroxypropionyl nitrate (PPN, C 2 H 5 C(O)O 2 NO 2 ) and peroxymethacryloyl nitrate (MPAN, CH 2 C(CH 3 )C(O)O 2 NO 2 ), which form during the oxidation of C 2 H 5 CHO and MACR, respectively. BEARPEX-2007 measurements (Table S1) include vertical concentration gradients and above-canopy fluxes of PAN, PPN and MPAN, as well as a separate measurement of total peroxy nitrates ( PN) that may contain contributions from other APNs and non-acyl species such as CH 3 O 2 NO 2 . Most of our analysis will focus on the speciated observations, since these include fluxes.
After accounting for both thermal decomposition via the reverse reaction in (R12) and subsequent loss of the RC(O)O 2 radical, canopy-top APN lifetimes for the current study range from ∼0.7 h during the hot period to ∼4.9 h during the cool period. Thus, both observed and modeled PAN concentrations increase by more than a factor of 2 between the hot and cool periods (Fig. 8a). PAN mixing ratios are over-predicted by ∼30% during the hot period but agree well with observations during the cool period. Model-measurement agreement here is partly coupled to our choice of initial/advection PAN concentrations (Table 4), though an assumption of negligible background concentration seems fair given the sustained high temperatures and steadily decaying daytime PAN maxima observed during this period ). This over-prediction is consistent with the steady-state analysis of , who suggested that overestimates of modeled PAN during warmer conditions may have resulted from underestimated sinks for CH 3 C(O)O 2 (PA) radicals. Additionally, our extensive chemical mechanism predicts several individually small sources of the PA radical, neglected by LaFranchi et al. (2009), that sum to ∼30% of the total PA production budget. In contrast, modeled PPN and MPAN concentrations agree with observations to within 5% for the hot period and within 17% for the cool period (Table S2).
Roughly 15% of modeled PNs can be attributed to compounds other than PAN, PPN and MPAN during the hot period. This result is consistent with the observational inter-comparison of Wooldridge et al. (2010), who showed that APNs other than PAN, PPN, and MPAN typically make up less than 10% of PNs. It is not consistent with many other models, however, including the base MCM, which would attribute a much larger fraction of PNs to species other than PAN, PPN, and MPAN (Perez et al., 2009). In our chemical mechanism, we neglect most anthropogenic VOC precursors (except acetaldehyde and propanal) and we invoke a rapid decomposition of βhydroxy acyl peroxy radicals (Wolfe and Thornton, 2011). These two characteristics explain the much lower predictions of PNs other than PAN, PPN, and MPAN. For example, if we were to exclude β-hydroxy acyl peroxy decomposition, modeled concentrations of PHAN (HOCH 2 C(O)O 2 NO 2 ) and C4PAN5 (HOC(CH 3 ) 2 C(O)O 2 NO 2 ), which are secondgeneration MBO oxidation products, lead to a 40% overprediction of PNs. The model overestimate of PNs during the hot period (Table 7) is mostly due the 30% overestimate of PAN.
The measured PAN gradient near the ground reveals that observations at 1.5 m are consistently lower than those at 5 m, with an average difference of 26 ± 14 pptv (mean ± 1σ ) between these two heights for the hot period. This difference constitutes a gradient of ∼17 ± 9% that is not captured in the modeled PAN profile. CAFE may not adequately represent deposition to the soil and ground litter, or diffusion near the ground, where chemical sinks are largest. We explore this in more detail below.
Modeled PAN fluxes and exchange velocities (Fig. 8bc) are under-predicted by 50-60% for both periods. Within CAFE, deposition of PAN occurs primarily through stomatal uptake (Table 6). It is possible that non-stomatal deposition is under-predicted, though laboratory measurements suggest that this term should be small compared to the stomatal component (Sparks et al., 2008). The 22% decrease in the observed PAN exchange velocities between the hot and cool periods suggests a temperature-dependent in-canopy loss process that is not represented in CAFE. Surface-facilitated thermal decomposition on sunlit canopy elements (which are warmer than the surrounding air) followed by loss of the PA radical seems a feasible mechanism, though the magnitude of this process would need to be larger than our total modeled deposition rate.
Alternatively, model-measurement disagreement may be related to the interplay of gas-phase chemistry and vertical mixing. Enhanced thermal decomposition due to the strong temperature gradient at the ground (Fig. 1) forces PAN out of chemical equilibrium (i.e. P < L), resulting in net chemical loss within the canopy and increasing its downward flux. This chemical perturbation, which we will call the chemical velocity (V c ) in analogy with the deposition velocity (V d ), amounts to −0.1 cm s −1 and −0.08 cm s −1 for the hot and cool periods, respectively (Table 7). Wolfe et al. (2009) estimated an average chemical velocity of −0.3 and −0.1 cm s −1 for two larger periods of BEARPEX-2007 that include our hot and cool periods. Thus, it is possible that CAFE underestimates this effect, especially during the hot period, which would also be consistent with the disagreement between near-ground modeled and measured PAN gradients as mentioned above. Decreasing diffusion in the lower canopy (as per the sensitivity test described in Sect. 3.1) does lower PAN mixing ratios by a few percent near the ground, but slower mixing ultimately leads to a decrease in the magnitude of the modeled above-cannopy PAN exchange velocity. Another potential cause for an under-estimate of the PAN chemical velocity is excess in-canopy production of PANor more specifically its precursor, the PA radical, which could be the case if modeled OH mixing ratios in the lower canopy are too high (recall that we have no constraint on these below z/ h = 0.9). This notion emphasizes the fact that APN exchange appears to be quite sensitive to the vertical profiles in chemical production and loss.
PPN and MPAN fluxes and exchange velocities (Table 8) are mostly within the large range of observed values, except for the PPN exchange velocity during the hot period. Measured PPN exchange velocities are quite fast (<−3 cm s −1 ) during the hot period, the possible implications of which have been discussed elsewhere . Within the CAFE model framework, an exchange velocity of this magnitude can only be obtained if PPN deposition rates are increased markedly by decreasing the cuticular resistance. The high variability of PPN and MPAN exchange velocities for the chosen observation windows precludes a more detailed model evaluation for these species.
Separating APN fluxes into chemical and surface (depositional) contributions can provide a more detailed look into the factors controlling forest-atmosphere APN exchange. Figure 9 compares the chemical velocities of PAN, PPN, MPAN and C4PAN5 for the hot period, as calculated according to Eq. (9). The latter species is the primary firstgeneration APN from MBO oxidation as predicted by the MCM and is an analog of MPAN, which derives from isoprene. Even though we do not have observations for comparison, we include C4PAN5 in this analysis for demonstrative purposes. Starting near the ground, the PPN chemical velocity diverges slightly towards less negative values than PAN, while the MPAN chemical velocity is even less negative. The chemical velocity of C4PAN5 shows the largest departure from PAN, becoming positive above z/ h = 0.3. These variations are not due to differences in deposition, as all APNs have the same deposition velocity in the model, which is also shown in Fig. 9. Deposition velocities may vary somewhat between APNs in reality, e.g. according to the functional groups on the carbon backbone. The diversity of modeled APN chemical velocity profiles is largely due to varying vertical distributions of their precursors. This is particularly evident in MPAN and C4PAN5 chemical velocities. Near the ground, the chemical velocity of MPAN and C4PAN5 is still controlled by thermal losses. Within and immediately above the canopy, however, oxidation of emitted isoprene and MBO leads to formation of methacrolein (MACR) and 2-hydroxy-2methylpropionaldehyde (IBUTALOH), respectively. In our MCM-based mechanism, these aldehydes are the sole precursors of MPAN and C4PAN5, respectively, and slightly enhanced levels of these precursors will enhance MPAN and C4PAN5 formation in the canopy. This production term continues to dominate above the canopy, whereas the thermochemical loss term becomes more constant as the temperature gradient is less pronounced here (note the near-vertical gradient of the PAN chemical velocity between z/ h = 1 and 2). Production of C4PAN5 is strong enough that CAFE predicts net emission of this compound from the forest, though this prediction is sensitive to the deposition term, which may be higher than we have modeled due to the hydroxyl functionality on C4PAN5. BVOC oxidation also produces PAN precursors, but these are not as specific as those of MPAN and C4PAN5. PAN production includes significant contributions from multi-generational oxidized VOC, such as acetaldehyde and methyl glyoxal, that are not as directly linked to BVOC emissions and thus are more evenly distributed in the vertical. The PPN chemical velocity is something of a special case, as the shape of its profile relative to that of PAN stems from an increased contribution from the storage term, ∂t dz. In this instance, the PPN precursor propanal (C 2 H 5 CHO) is evenly distributed in the vertical because its sole source in CAFE is advection. PPN produced aloft is transported into the canopy but cannot escape as easily, as decreased diffusion in the canopy can serve as a "trap" for gases with weak concentration gradients. This results in a slight buildup of PPN and thus a slight positive perturbation to the flux. Chemical flux contributions are slightly dampened during the cool period (not shown), though MPAN exchange velocities are still somewhat less negative than those of PAN and PPN (Table 8).
The strength of any APN production flux will depend on a number of factors, including BVOC emission rates, OH mixing ratios and canopy residence times. In particular, we noted earlier (Sect. 3.1) that a substantial isoprene emission rate is required to maintain agreement with isoprene observations, though previous studies have identified advection as the primary isoprene source at BFRS (Dreyfus et al., 2002). Replacing isoprene emission with advection would reduce the MPAN production flux and bring the modeled MPAN V ex closer to that of PAN or PPN, because the source of MACR would no longer be elevated in the canopy. Likewise, thermochemical APN loss fluxes depend on both the absolute temperature and the shape the temperature gradient in the canopy, which may change dramatically between a young and open forest like BFRS and a more mature forest. Moreover, APN production and loss are also subject to NO x concentrations and NO/NO 2 ratios. To expand the relevance of these findings to other ecosystems, future modeling work should probe the sensitivity of APN fluxes to such factors, particularly temperature gradients, BVOC emissions and soil NO emissions.
To conclude this section, we note that timescales for APN production and loss (30 min to several hours) are much longer than the modeled canopy residence time of ∼2 min; however, chemistry still obviously contributes markedly to the modeled fluxes of these compounds. It is not the absolute magnitudes of production and loss, but rather sustained gradients in these terms, that can perturb mixing ratio vertical profiles and thus drive diffusive fluxes. To the extent that our model captures the average mixing process, these results suggest that comparison of reactive and mixing timescales may not be an appropriate metric by which to assess the contribution of chemistry to forest-atmosphere exchange.

ANs
Alkyl nitrates (RONO 2 ) are formed as minor products during NO + RO 2 Reaction (R13), with typical branching ratios of 5-10% for AN formation (Atkinson and Arey, 2003). ANs are also produced during the oxidation of VOC by the nitrate radical (NO 3 ), though daytime NO 3 concentrations are generally too low to be important for the current study. Recent work suggests that chemical mechanisms may be incomplete with regard to AN chemistry, particularly concerning their ability to "recycle" NO x during oxidation by OH (Rollins et al., 2009;Perring et al., 2009a, b;Horowitz et al., 2007;Paulot et al., 2009b). Deposition of ANs is also poorly constrained by observations and may depend on the functional form of the R-group. BEARPEX-2007 measurements are limited to AN concentration profiles, though earlier observations at this forest have included AN fluxes. Figure 10a-c displays modeled mixing ratios, fluxes and exchange velocities for MBOANO 3 , which is a firstgeneration oxidation product of MBO that comprises ∼50% of the modeled AN budget. Concentrations are slightly higher during the hot period, consistent with the AN observations (Table 6) and with faster formation rates due to higher OH and BVOC concentrations. Vertical mixing ratio gradients for both periods are characteristic of strong deposition. Deposition velocities are tuned to match the value of 2.7 cm s −1 recommended by Farmer and Cohen (2008) by increasing the effective Henry's Law constant to 1 × 10 8 M atm −1 , effectively lowering the cuticular resistance. The slight increase in exchange velocities between the hot and cool periods results from the 8% rise in friction velocity, leading to a decrease in the modeled laminar sublayer resistance. Deposition thus dominates the flux and exchange velocity profiles of ANs with a rate that is mostly controlled by aerodynamic transfer to canopy surfaces. A small contribution from in-canopy production shifts MBOANO 3 fluxes towards less negative values by ∼10%. Chemistry-driven fluxes could become more important if deposition rates are lower than modeled -a distinct possibility considering the limited observational constraints on this process -or if AN yields are higher. Other primary ANs derived from local BVOC (i.e. MBOBNO 3 and the 4 isoprene-derived ANs, see Appendix A) exhibit the same vertical and seasonal patterns as MBOANO 3 .

HNO 3
Dry deposition of gas-phase nitric acid is a primary pathway for atmosphere-to-ecosystem nitrogen transfer. HNO 3 adsorbs readily to most surfaces, thus deposition is assumed to proceed at the aerodynamic limit. This view is generally supported by inferential (e.g. flux-gradient) measurements of HNO 3 fluxes over forests (Horii et al., 2006;Pryor and Klemm, 2004;Sievering et al., 2001), which report deposition velocities ranging from 2 to 10 cm s −1 . Previous eddy covariance measurements at BFRS (prior to BEARPEX-2007) have reported HNO 3 deposition velocities of 3-4 cm s −1 during winter but have also offered evidence that fast intra-canopy chemistry can influence HNO 3 fluxes, even to the point of creating a net upward flux (out of the forest) during the summer (Farmer and Cohen, 2008). Figure 11a-c illustrate modeled profiles of HNO 3 mixing ratios, fluxes and exchange velocities. HNO 3 concentrations are ∼2 times higher during the hot period relative to the cool period, consistent with faster production via Reaction (R14) due to higher OH levels. Fluxes and exchange velocities are fast and essentially driven by aerodynamically-limited deposition. As in the cases of ANs, changes in the laminar sublayer resistance give rise to different exchange velocities between the hot and cool periods. Modeled HNO 3 fluxes do include a small (∼5%) positive contribution due to in-canopy transformation of soil-emitted NO.

Oxidized nitrogen deposition
Ecosystem-scale nitrogen deposition affects biosphere productivity and represents a major pathway by which anthropogenic emissions influence the environment. Dry deposition typically constitutes ∼50% of total atmospheric N deposition, with the other half due to wet deposition (i.e. precipitation) (Bytnerowicz and Fenn, 1996;Sparks et al., 2008). Though HNO 3 is likely the dominant dry-depositing species, several studies have inferred that a significant fraction of the downward NO y flux is comprised of species other than nitric acid (Horii et al., 2006;Sparks et al., 2008). As detailed by the above discussion, inferring gross N deposition rates from net NO y fluxes without considering in-canopy chemistry can lead to errors. Figure 12 summarizes modeled above-canopy (z/ h = 2, or 20 m) NO y fluxes. Gross NO y deposition amounts to 23 and 14 pptv m s −1 (11 and 7 ngN m −2 s −1 ) for the hot and cool periods, respectively. This is within the range of other estimates of N deposition to California forests (Bytnerowicz and Fenn, 1996;Herman et al., 2003). HNO 3 constitutes 83% of deposited NO y during the hot period but only 72% during the cool period owing to decreased HNO 3 and increased APN and NO 2 mixing ratios. For both periods, upward NO 2 fluxes (driven by soil NO emissions) decrease the net modeled NO y flux by ∼30% relative to the gross deposition flux. Differences between total and depositional fluxes for individual classes of NO y are consistent with our earlier discussion. For example, APN fluxes are only 60% depositional during the hot period, while total AN fluxes underestimate the depositional flux by ∼10%. Note that our analysis is focused on gaseous oxidized nitrogen and thus does not consider dry deposition of ammonia (NH 3 ) or particulate ammonium nitrate (NH 4 NO 3 ). A small set of NH 3 flux observations recorded at BFRS in 2006 suggests an average NH 3 flux of −7.4 ngN m −2 s −1 for this location (Fischer and Littlejohn, 2007). If all of this flux is depositional, NH 3 uptake would be competitive with our estimated dry oxidized N deposition flux.
The picture presented in Fig. 12 should be interpreted with care. Relative NO y mixing ratios and deposition rates can vary widely by location and season. Deposition velocities are still highly uncertain for both APNs and ANs, largely because the mechanisms for uptake or heterogeneous loss are not understood. For example, given that CAFE underestimates PAN exchange velocities during the hot period, it is possible that non-stomatal PAN deposition is faster than represented by the standard resistance parameterization. It is also likely that some deposited species may be re-emitted as NO 2 or nitrous acid (HONO) rather than taken up by vegetation. In an analysis of HONO concentrations measured during BEARPEX-2007 require an unidentified HONO source of 1.6 ppbv day −1 , or 0.02 pptv s −1 , to reconcile observations with a steady-state estimate. Heterogeneous HONO production is generally thought to proceed via surface reactions of NO x (Goodman et al., 1999) and nitrate photolysis (Zhou et al., 2003;He et al., 2006). Assuming this missing source is purely heterogeneous (i.e. production occurs on canopy surfaces) and integrating over the canopy height, we estimate a HONO production flux of 0.2 pptv m s −1 , which is ∼67% of the modeled NO x deposition flux during the hot period. Similarly, if the total APN deposition flux during the hot period was treated as an emission of NO 2 , NO x fluxes could increase by as much as 20%. Recent measurements have even suggested that NO y emitted from canopy surfaces could originate from photolysis of deposited HNO 3 (Raivonen et al., 2006). Future efforts to close the N deposition budget should include direct field observations of speciated NO y fluxes and gradients, controlled laboratory experiments on uptake by vegetation and other surfaces (e.g. soil and ground litter), and detailed modeling work.

Conclusions
We have used the CAFE model to simulate observations from the BEARPEX-2007 field campaign at Blodgett Forest Research Station in the Sierra Nevada, CA. Our model results highlight a number of interesting features in the extensive BEARPEX-2007 dataset.
1. CAFE under-predicts OH concentrations by a factor of 6 during the hot period, requiring implementation of an enhanced OH-recycling mechanism that effectively converts key organic peroxy radicals to OH. During the cool period, model-measurement agreement for HO x is obtained without this mechanism; however, this agreement may be artificial as OH reactivity is underpredicted for this scenario. Thus, OH production is likely under-estimated by our base MCM mechanism for both periods.
2. Comparison of model results with H 2 O 2 observations suggests that H 2 O 2 deposition occurs at the aerodynamic limit, much faster than predicted by standard resistance parameterizations but corroborating recent direct observations.
3. Modeled O 3 exchange velocities under-predict observations by ∼20%; this discrepancy will grow if the magnitude of non-stomatal deposition is currently overestimated, as we argue. With known chemistry, incanopy chemical losses do not contribute significantly to the total modeled O 3 flux. Reproducing the chemical flux inferred from previous measurements at BFRS will likely require significant increases in BVOC emissions with high reactivity towards ozone. On a larger scale, such changes carry potential ramifications for quantifying ozone-induced ecosystem stress, BVOC oxidation pathways and intra-canopy oxidant sources.
4. Nitrogen oxide (NO x = NO + NO 2 ) fluxes are driven by soil emissions of NO. Roughly 89% of soil-emitted NO x escapes the canopy mostly as NO 2 , suggesting that NO 2 fluxes might be a proxy for soil NO emissions at this and similar young forests. Upward NO 2 fluxes from soil emissions reduce the net-downward above-canopy NO y flux by 30%. In other words, nitrogen deposition estimates from a total NO y flux measurement would be biased low by 30%.
5. Deposition and chemical loss contribute equally to modeled PAN fluxes, but the total modeled PAN flux is ∼50% lower than the observations. Coupled with a 30% over-prediction of PAN mixing ratios during the hot period, this suggests that in-canopy APN sinks are either under-estimated or masked by too much production. The underlying cause of model-measurement discrepancies in APN fluxes remains unclear, partly due to uncertainties in the parameterization of vertical mixing and non-stomatal deposition. Chemical production also influences APN fluxes, especially when their formation is closely tied to the oxidation of primary BVOC emissions. In contrast, fluxes of alkyl nitrates and HNO 3 are driven by deposition under our model conditions. 6. HNO 3 dominates model-calculated dry N deposition (which excludes NH 3 and particulate N) during the hot period, but other classes of NO y become non-negligible (∼28%) during the cool period. Such effects will carry implications for N deposition estimates from routine monitoring networks, which typically only measure wet and dry deposition of NO − 3 /NH + 4 and HNO 3 (Sparks et al., 2008).
It is clear that significant uncertainties still limit our understanding of forest-atmosphere exchange. First, chemical mechanisms fail to reproduce observed HO x concentrations under high-BVOC conditions. A number of OH-recycling schemes have been postulated to close this gap, but the underlying mechanisms remain unidentified, which will impede the predictive capability of any model aiming to track carbon through the emission and oxidation process. Second, Ktheory is a rough approximation to the true structure of turbulent transport within mature canopies, yet it persists as the standard for this type of model. A computationally efficient alternative to K-theory that accurately captures the key features of intra-canopy turbulence -particularly the influence of sweep-ejections and similar events -would improve confidence in future modeling efforts. Third, a lack of detailed experimental constraints on the mechanisms and efficiency of depositional processes and on BVOC emission inventories continues to prevent accurate parsing of fluxes into emission, chemistry and deposition. In particular, uncertainties in the magnitude of non-stomatal deposition for O 3 , APNs and other reactive species limit our ability to assess potential chemical perturbations on forest-atmosphere exchange and vice versa. In many instances, parameterizations are tuned so that observed trace gas fluxes are reproduced in models as being purely depositional or as direct emissions from the canopy to the atmosphere. This simplification will bear consequences for accurately modeling ecosystem responses to chemical and climate stresses, such as future changes in temperature and background ozone levels. Finally, our analysis of APN fluxes demonstrates that reactive timescales do not need to be faster than canopy mixing for chemistry to influence fluxes. Future efforts to interpret observed reactive trace gas fluxes must focus on both the absolute magnitudes and the gradients of production and loss that develop in a complex canopy environment. Modeling efforts should aim towards a comprehensive theory that allows incorporation of these chemical effects into surface exchange parameterizations of regional and global models, with the ultimate goal of fully representing the interaction of the forest with the atmosphere.