Evaluating Two Optical Methods of Woody-to-Total Area Ratio with Destructive Measurements at Five Larix gmelinii Rupr . Forest Plots in China

Accurate in situ leaf area index (LAI) estimates of forest plots are required to validate currently-used LAI map products. Woody-to-total area ratio (α) is a crucial parameter in converting the plant area index estimates of forest plots obtained by optical methods into LAI. Although optical methods for estimating the α of forest canopy have been proposed, their performance has never been assessed. In this study, five Larix gmelinii Rupr. forest plots with contrasting plot characteristics (i.e., tree age, tree height, management activities, stand density, and site conditions) were selected. The performance of two commonly used optical methods, namely, multispectral canopy imager (MCI) and digital hemispherical photography (DHP), in estimating the α of L. gmelinii forest plots was evaluated by using the reference α of the selected forest plots. The reference α of forest plots was measured via destructive method by harvesting two or three representative trees in each plot. Large variations were observed amongst the reference α of the selected forest plots (ranging from 0% to 56%). These α were also highly correlated with the site conditions and management activities in these plots. The effective α (αe) or α estimated using the leaf-on and leaf-off periods MCI or DHP images with or without consideration of the clumping effects of canopy element and woody components were 1.57 to 4.63 times the reference α in the five plots. The overestimation of α or αe was mainly caused by the preferential shading of woody components by the shoots in the leaf-on canopy. Accurate α estimates for the L. gmelinii forest plots with errors of less than 20% can be obtained from MCI when the clumping effects of canopy element and woody components are considered in the estimation.


Introduction
Leaf area index (LAI) is typically used in many models to describe vegetation-atmosphere interactions as crucial variable controlling processes, such as photosynthesis, respiration, and rain interception [1].LAI is defined as half of the total green leaf area per unit of the flat ground area [2,3].LAI can be measured in situ and retrieved from remote sensing images.In situ LAI measurements usually provide LAI estimates at the plot scale (e.g., 20 × 20-100 × 100 m).By contrast, LAI maps provide LAI estimates at the regional, national and global scales.Over the past two decades, several commonly used LAI map products have been published, including CYCLOPES [4], GLOBCARBON [5], GLASS [6], MODIS [7], and GEOV1 [8].However, several studies highlighted huge differences amongst these products [9][10][11][12].Therefore, accurate and reliable in situ LAI measurements are needed to assess and improve the accuracy of these products and for them to match the accuracy requirements of several applications, such as the Global Climate Observing System (GCOS).The GCOS has specified that LAI map product values must be within 20% of in situ LAI measurements and be improved to within 5% for future applications [1,13].
Methods for estimating in situ LAI measurements can be categorized into direct and indirect methods.On the one hand, direct methods are time consuming, laborlabor-intensive and can sometimes be destructive to canopies [14].On the other hand, indirect methods are widely used for estimating in situ LAI of vegetation canopies due to their high efficiency, low cost, and nondestructiveness [15].Amongst these methods, the optical ones are most commonly used, including digital hemispherical photography (DHP), LAI-2000/LAI-2200 (Li-Cor, Lincoln, NE, USA), tracing radiation and architecture of canopies (TRAC) (3rd Wave Engineering, Winnipeg, Manitoba, Canada), HemiView (Dynamax, Houston, TX, USA), SunScan (Delta-T Devices, Cambridge, UK), multispectral canopy imager (MCI) [14], and multiband vegetation imager (MVI) [15].In situ LAI measurements can be estimated based on the gap fraction or radiation attenuation measurements collected by using optical methods.For forests, the gap fraction or radiation attenuation measurements collected by optical methods include the contribution of leaves or shoots and woody components (e.g., stems, branches, flowers, and fruits).Therefore, the estimates derived from indirect methods are plant area index (PAI) instead of LAI.PAI is the sum of LAI and woody area index (WAI).Optical methods require a parameter called woody-to-total area ratio (α = WAI/PAI) [16] to convert PAI into LAI.
The direct method estimates the α of forest canopies by directly measuring the α of several representative trees in the plots via destructive methods [16].However, this method is labor-intensive and time-consuming and has been rarely used in practice [14].The indirect methods for estimating the α of forest canopies mainly include optical methods, such as LAI-2000/LAI-2200 [17,18], DHP [19,20], MCI [14], and MVI [21].Amongst these methods, LAI-2000/LAI-2200 and DHP estimate the α of forest canopies based on the gap fraction or radiation attenuation measurements collected at leaf-on and leaf-off periods [17][18][19][20].Meanwhile, the effective WAI (WAI e ) obtained from LAI-2000/LAI-2200 or DHP at leaf-off period forest canopies is represented as the WAI e of leaf-on period forest canopies in estimating α [17][18][19].Other studies have estimated the α of forest canopies by calculating the ratio of pixel numbers between woody components and the sum of leaves or shoots and woody components in a group of single-channel DHP images [13,22].However, these traditional optical methods are unable to discriminate leaves or shoots, branches, stems, fruits and sky.Two multispectral canopy imaging devices (i.e., MVI and MCI) have been developed and applied to estimate α based on the special spectral characteristics of canopy element at visible (VIS) and near-infrared (NIR) bands [14,21].Previous studies have reported a wide range and relatively large values of α (ranging from 0.03 to 0.41) [14,17,18,23,24].Therefore, ignoring α can result in approximately 3% to 40% estimation errors in LAI estimates obtained by optical methods.Thus, α must be considered in the in situ LAI estimation of forest canopies.
Although several optical methods have been developed for estimating the α of forest canopies, no study has evaluated the estimation performance of those methods.Moreover, three aspects of those methods require further investigation.Firstly, a constant α value is often used in estimating the LAI of forest plots with the same tree species [16,25].However, the variations in the α of plots with the same tree species yet varying stand densities, tree ages, site conditions and management activities are often ignored and whether such variations can be accepted for obtaining accurate in situ LAI measurements is yet to be verified.Secondly, for most natural forest canopies, the spatial distribution of the canopy element and woody components in the space largely deviates from the random distribution assumed by the inversion models used by optical methods.Zou et al. reported that the most effective plant area index (PAI e ) and WAI e estimates derived from optical methods are 35% to 75% and 63% to 90% of the PAI and WAI of leaf-on and leaf-off forest scenes, respectively [15].Therefore, the clumping effects of canopy element and woody components should be considered in estimating the PAI and WAI of forest canopies.However, no study has evaluated how the clumping effects of canopy element and woody components influence the α estimation of forest plots.Thirdly, several optical methods have been developed for estimating the α of forest canopies.Some studies have attempted to cross-compare the α estimates obtained by different optical methods [19,20].However, to the best of our knowledge, no study has attempted to validate the estimation performance of those methods based on the α estimates obtained by direct methods, such as the destructive method.
In this study, five permanent Larix gmelinii Rupr.plots, with varying tree ages, tree heights, stand densities, and site conditions, were selected as elementary sampling units to validate LAI map products.The reference α of each plot was obtained by using the destructive method and the variations amongst the α the five selected forest plots were analyzed.Two classical optical methods (i.e., MCI and DHP) were used to simultaneously estimate the α of these five sites.Factors that affect the estimation performance of MCI and DHP to estimate the α forest plots were analyzed, including inversion model and the clumping effects of canopy element and woody components.The α estimates derived from DHP and MCI were compared with those obtained by the destructive method.The performance of these two optical methods was evaluated afterwards to determine the best method for estimating the α of L. gmelinii forest plots.

Effective PAI (PAI e ), Effective WAI (WAI e ), PAI, and WAI
The PAI of leaf-on forest canopies can be estimated by inverting Beer's law based on the gap fraction measurements [15,26,27] as follows PAI = −ln(p e (θ))cos(θ)γ e G e (θ)Ω e (θ) where p e (θ), G e (θ), and Ω e (θ) are the canopy element gap fraction, canopy element projection coefficient and canopy element clumping index (Ω e ) at zenith angle θ. γ e is the needle-to-shoot area ratio for quantifying the effect of needle clumping within a shoot.Previous studies have reported that the G e (θ) of leaf-on forest canopy is approximately intersected at the zenith angle near 57.3 • and can be assumed to be equal to 0.5 at such zenith angle [13,15,28].Therefore, many studies have attempted to estimate the LAI and PAI of forest canopies based on the Beer's law and gap fraction measurements collected at this zenith angle or the narrow zenith angle ranges near 57.3 • [13, 15,27,29,30].The PAI of leaf-on forest canopies can be calculated based on Equation (1) as follows (57.3) where p e (57) and Ω e (57) denote the canopy element gap fraction measurements and the Ω e estimates obtained at zenith angle 57 • or the narrow zenith angle ranges near 57 • , respectively.Aside from Beer' law, the PAI of leaf-on forest canopies can be estimated based on Miller theory [31] without prior knowledge of G e (θ) (Miller): Based on the gap fraction measurements collected from the five annuli of LAI-2000/LAI-2200, a modified Miller theory was used to estimate the PAI of leaf-on forest canopies [15,32] (LAI-2200): where θ i is the center zenith angle of the ith annulus and p e_i (θ i ), G e_i , and W i are the canopy element gap fraction, canopy element projection coefficient, and weight factor of the ith annulus, respectively.Ω e_i and G e_i were calculated by averaging the Ω e (θ) and G e (θ) estimates with the zenith angles covered by the ith annulus, respectively.The zenith angle ranges of the five annuli of LAI-2000/LAI-2200 were 0-13    , respectively [15,32].W i can be calculated as follows When normalised to 1.0, the values of W i in Equation ( 5) were computed as 0.041, 0.131, 0.201, 0.29, and 0.337, which correspond to the five annuli with center zenith angles of 7    , and 68 • , respectively [15].
The modified Miller theory also can be adopted by MCI to estimate the PAI of leaf-on forest canopies (MCI_0-85) as follows By using Equation ( 5), the W i for the nine annuli of MCI were computed as 0.0038, 0.0303, 0.0596, 0.0872, 0.1120, 0.1335, 0.1510, 0.1640, and 0.2750 whilst the zenith angle ranges covered by these annuli were 0-5     If the beyond-(Ω e ) and within-shoot clumpings (γ e ) are not considered in estimating the PAI of leaf-on forest canopies (where Ω e (θ), Ω e (57), Ω e_i , and γ e are assumed to be equal to 1, then the estimates derived from Equations (2), (3), (4), and (6) are PAI e_57.3 , PAI e_Miller , PAI e_LAI−2200 , and PAI e_MCI_0−85 , respectively.Given that the spatial distribution of the canopy element for most natural forest canopies usually deviates from the random distribution.Optical methods usually fail in detecting small gaps between needles within shoots [14,33].Therefore, both Ω e and γ e should be considered in estimating the PAI of leaf-on forest canopies.Several algorithms have been developed to estimate the Ω e and woody components clumping index (Ω w ) of leaf-on and leaf-off forest canopies, including logarithmic averaging (LX) [34], modified logarithmic averaging (LXW) [35], gap size distribution (CC) [36,37], modified gap size distribution (CMN) [38], combination of gap size and logarithmic averaging (CLX) [39], and Pielou's coefficient of spatial segregation (PCS) [40].Previous studies show that PCS tend to produce error Ω e and Ω w estimates [15,35,38].Moreover, high similarities are observed between CC and CMN as well as between LX and LXW [15].For conciseness, only the three algorithms of CC, LX and CLX recommended by Zou et al. [15] were used in this study.Table 1 presents the Ω e and γ e estimation formulae used in this study.

Table 1.
The Ω e and γ e estimation formulae used in this study.F m (0, θ) is the measured total gap fraction of the canopy element at zenith angle θ (i.e., accumulated gap fraction from the largest to the smallest gaps), F mr (0, θ) is the total gap fraction after removing those large gaps resulting from the nonrandom distribution of the canopy element [14,36,41], p e (θ) is the mean canopy element gap fraction of all segments at θ, ln[p e (θ)] is the mean logarithmic canopy element gap fraction for all segments at θ [34,41], n is the number of segments, p e_k (θ) is the canopy element gap fraction of segment k, Ω e_CC_k (θ) is the Ω e of segment k at θ estimated by CC [35,39], and A n is half of total needle area in a shoot.A p (0 • , 0 • ), A p (45 • , 0 • ), and A p (90 • , 0 • ) are the shoot projection areas measured by projecting the shoot at zenith angles 0 • , 45 • , and 90 • and azimuth angle 0 • , respectively [16].γ is the effective needle-to-shoot area ratio derived based on the A p (0 • , 0 • ), A p (45 • , 0 • ), A p (90 • , 0 • ), and A n of each shoot sample.γ e is the needle-to-shoot area ratio that corrects the overestimation of γ for woody components in estimating the PAI of leaf-on forest canopies [14].

2.2.
Effective Woody-to-Total Area Ratio (α e ) and Woody-to-Total Area Ratio (α) If Ω e and Ω w are not considered in estimating the α of leaf-on forest canopies, then the estimates derived from optical methods are computed as α e (effective woody-to-total area ratio) [14,19,42]: Sea et al. [22] proposed another method to estimate the α e of leaf-on forest canopies based on the canopy and woody component cover fractions ( f c and f w ) as follows (percentage) where f c is the ratio of the total pixel numbers of leaves or shoots and woody components to the pixel number of the classified DHP or MCI images, whereas f w is the ratio of the pixel number of woody components to the pixel number of classified DHP or MCI images [13,22].
If Ω e and Ω w are considered in estimating the α of leaf-on forest canopies, then α can be estimated as follows α = WAI /PAI (14) where PAI and WAI are were obtained from DHP or MCI.For the destructive method, harvesting all trees in each forest plot is impossible, especially for those plots that have been established as long-term observation plots for validating LAI map products.Some studies have attempted to estimate the α of forest plots based on the α of several harvested trees Forests 2018, 9, 746 6 of 26 located in or nearby these plots [16,23,43].This study also applied such method to estimate the α of forest plots based on the α of destructive trees and the basal area of the trees in these plots: where B d , B c , and B s are the basal area of the dominant, codominant and suppressed trees in the plots, respectively; and α d , α c and α s are the estimated α of the representative dominant, codominant, and suppressed trees as measured by the destructive method, respectively.

Plots Description
The five forest plots selected in this study are located in the Saihanba National Forest Park in Hebei Province, China (Figure 1).These plots are monocultural with a flat terrain and a size of 25 m × 25 m.Forest inventories were carried out during the LAI measurement campaign (Table 2).For plots 1 and 2, the majority of the branches below the live canopies were harvested during the selective thinning performed in the winter of 2014.For plots 3 and 4, the dead branches with heights of less than about 2.5 m were low-and medium-level harvested before the LAI measurement campaign, respectively.For plot 5, all branches of forest canopies were not harvested before the LAI measurement campaign.
Forests.2018, FOR PEER REVIEW 6 of 4 of several harvested trees located in or nearby these plots [16,23,43].This study also applied such method to estimate the  of forest plots based on the  of destructive trees and the basal area of the trees in these plots: where   ,   , and   are the basal area of the dominant, codominant and suppressed trees in the plots, respectively; and   ,   and   are the estimated  of the representative dominant, codominant, and suppressed trees as measured by the destructive method, respectively.

Plots Description
The five forest plots selected in this study are located in the Saihanba National Forest Park in Hebei Province, China (Figure 1).These plots are monocultural with a flat terrain and a size of 25 m × 25 m.Forest inventories were carried out during the LAI measurement campaign (Table 2).For plots 1 and 2, the majority of the branches below the live canopies were harvested during the selective thinning performed in the winter of 2014.For plots 3 and 4, the dead branches with heights of less than about 2.5 m were low-and medium-level harvested before the LAI measurement campaign, respectively.For plot 5, all branches of forest canopies were not harvested before the LAI measurement campaign.The typical shoots used for measuring the mean element width and needle-to-shoot area ratio were collected by using a ladder at the three height classes of forest canopies (i.e., top, middle and bottom).Two to four shoot samples were randomly obtained at each height class of the five plots.The mean element width of each shoot sample was calculated by using the method described in Leblanc et al. [44].The mean element width of each plot was derived by averaging the mean element width estimates of all typical shoots in the plot.
The shoot projection areas were recorded by using a Canon 6D camera equipped with a Canon 24-70 mm lens and a flat, levelled white panel with two rulers laid on its top surface.Three shoot projection images were obtained for each shoot by rotating the shoot main axis at zenith angles 0 • , 45 • , and 90 • and azimuth angle 0 , 0 • ), and A p (90 • , 0 • ) were obtained for each shoot based on the three shoot projection images.The A n of each shoot was measured by applying the volume displacement method described in Chen et al. [33].Afterwards, the γ of each typical shoot was calculated by using Equation (10).The γ of each plot was estimated by averaging the γ of all sampled typical shoots.The γ e of each plot was calculated by using Equation (11) based on the γ and reference α of the plot.

DHP Images Acquisition and Processing
A cross-pattern sampling scheme with nine sampling points (with a 5 m distance interval between each sampling point, except for the central sampling point) was used to collect DHP images (Figure 2).The sampling points were marked with wooden stakes to facilitate the two periods collection of leaf-on and leaf-off DHP images at each sampling point.These DHP images were collected by using a Canon 6D camera with a height of about 1.2 m and equipped with a Sigma 8 mm fisheye lens.The DHP images had a resolution of 5472 × 3678 pixels.The exposure parameter was set manually to avoid the overexposure or underexposure of DHP images.All images were collected before sunrise, after sunset, or under an overcast sky.Both leaf-on and leaf-off periods DHP images were collected for each plot (Figure 3a,b, respectively).The leaf-on period DHP images were collected between 11 August and 2 September 2017 at the maximum PAI (maximum LAI) period of the five plots.The leaf-off period DHP images were collected between 28 October and 15 November 2017 at the minimum PAI (equals to WAI and the LAI equals to 0) of the five plots.A total of 90 DHP images were obtained from these plots.
The DHP images were preprocessed by using the same procedure described in Gonsamo and Pellikka [45].The image preprocessing procedure includes cropping the circular areas of the DHP images, selecting the blue channel, applying a correction for gamma, and thresholding to classify the canopy element and sky.The in-house software called Measurement Tools of Vegetation Structural Parameters (MTVSP, version 2018, Fuzhou, China), which was developed using the MATLAB programming language, was used to calculate the p e (θ), p w (θ), Ω e , Ω w , WAI e , PAI e , WAI, and PAI of the leaf-on and leaf-off forest plots based on the classified DHP images and MCI image pairs.For the Forests 2018, 9, 746 8 of 26 DHP method, the calculation procedures for p e (θ), p w (θ), p e_i (θ i ), p w_i (θ i ), p e (57), and p w (57) were the same as those described in Zou et al. [15].An assumption of spherical projection function was made for the canopy element (G e ) and woody components (G w ) projection functions in the PAI and WAI estimation because no G e and G w measurements were collected in the field for the five plots.
The Ω e (θ) and Ω w (θ) estimates at 10 zenith angles ranging from of 0 • to 9 • with a 1 • interval were assumed to be equal to Ω e (10) and Ω w (10), respectively [15].This assumption was adopted to obtain 91 Ω e (θ) and Ω w (θ) estimates that match the p e (θ) and p w (θ) measurements at the same zenith angle range of 0 • to 90 • with a 1 • interval.The Ω e_i or Ω w_i in Equation ( 4) was calculated by averaging the Ω e (θ) or Ω w (θ) estimates with the zenith angles covered by the ith annulus.The Ω e (57) or Ω w (57) in Equation ( 2) was calculated by averaging the Ω e (θ) or Ω w (θ) with zenith angles ranging from 52 • to 62 • with a 1 • interval.Three groups of PAI e and WAI e estimates of PAI e_57.3 and WAI e_57.3 , PAI e_Miller and WAI e_Miller , as well as PAI e_LAI−2200 and WAI e_LAI−2200 were calculated for each plot by using Equations ( 2) to (4) based on the p e (θ), p w (θ), p e_i (θ i ), p w_i (θ i ), p e (57), and p w (57) measurements.Afterwards, three α e estimates were obtained for each plot by using Equation ( 12) based on the three groups of PAI e and WAI e estimates.A total of 21 groups of PAI and WAI estimates were calculated for each plot by using Equations ( 2 Vegetation Structural Parameters (MTVSP, version 2018, Fuzhou, China), which was developed using the MATLAB programming language, was used to calculate the   () ,   (),   ,   ,   ,   , WAI, and PAI of the leaf-on and leaf-off forest plots based on the classified DHP images and MCI image pairs.For the DHP method, the calculation procedures for   (),   (),  _ (  ),  _ (  ),   (57), and   (57) were the same as those described in Zou et al. [15].An assumption of spherical projection function was made for the canopy element (  ) and woody components (  ) projection functions in the PAI and WAI estimation because no   and   measurements were collected in the field for the five plots.
For the DHP method, six segment sizes were used by LX (5°, 15°, and 30°) and CLX (15°, 30°, and 45°) to estimate the   and   of leaf-on and leaf-off forest canopies.The reasons for selecting the six segment sizes are explained in detail in Zou et al. [15].Seven   () or   () estimates were obtained by using CC, LX, and CLX at each zenith angle ranging from 10° to 90° with a 1° interval at each plot (one estimate for  _ () or  _ () and three estimates each for  _ () and  _ () or  _ () and  _ () ).The   () and   () estimates at 10 zenith angles ranging from of 0° to 9° with a 1° interval were assumed to be equal to   (10) and   (10), respectively [15].This assumption was adopted to obtain 91   () and   () estimates that match the   () and   () measurements at the same zenith angle range of 0° to 90° with a 1° interval.The  _ or  _ in Equation (4) was calculated by averaging the   () or   () estimates with the zenith angles covered by the i th annulus.The   (57) or   (57) in Equation ( 2) was calculated by averaging the   () or   () with zenith angles ranging from 52° to 62° with a 1° interval.Three groups of   and   estimates of  _57.3 and  _57.3 ,  _ and  _ , as well as  _−2200 and  _−2200 were calculated for each plot by using Equations ( 2) to (4) based on the   (),   (),  _ (  ),  _ (  ),   (57), and   (57) measurements.Afterwards, three   estimates were obtained for each plot by using Equation ( 12) based on the three groups of   and   estimates.A total of 21 groups of PAI and WAI estimates were calculated for each plot by using Equations ( 2

MCI Images Acquisition and Processing
Two representative locations were selected for the collection of MCI images in each plot.The same procedures for collecting MCI image pairs adopted in previous works were used in

MCI Images Acquisition and Processing
Two representative locations were selected for the collection of MCI images in each plot.The same procedures for collecting MCI image pairs adopted in previous works were used in this study with two key differences [14,41].Firstly, the MCI image pairs were collected in zenith directions ranging from 0 • to 80 • with a 10 • interval (i.e., nine zenithal measurements at each azimuthal direction).Secondly, the manual exposure model was used to capture MCI images to avoid overexposure and underexposure.The MCI instrument height was about 1.2 m and MCI image resolution was 3488 × 2616 pixels.The MCI measurements were collected before sunrise, after sunset or under an overcast sky.To minimize the effect of forest canopy movement on the collection of MCI image pairs, the MCI measurements were collected under windless or breezy conditions [14,41].Both leaf-on and leaf-off periods MCI measurements were collected in each plot (Figure 3c-f).The leaf-on period MCI measurements were collected from 20 August to 12 September 2017 at the maximum PAI (maximum LAI) period of the plots.The leaf-off period MCI measurements were collected from October 31 to 16 November 2017 and 4 May 2018 at the minimum PAI (equals to WAI and the LAI equals to 0) of the plots.A total of 1080 MCI image pairs were obtained from the five plots.
The procedures used for processing the MCI image pairs were the same as those described in Zou et al. [14], except that the images were subdivided into two sub-images to reduce the amount of work needed.The MCI image pairs were eventually classified into sky, woody components, and shoots [14].For the MCI method, p e (57) and p w (57) were calculated by using the classified MCI images of the seventh annuli collected with a center zenith angle of 60 • .
The full size of each pixel of the MCI images was calculated based on the field of view of MCI, image resolution, mean tree height of the plot and the zenith angle where these images were collected [14,21].The pixel number of the segment required as the key parameter for LX and CLX in the Ω e and Ω w estimation of leaf-on and leaf-off forest plots was calculated based on the full size of each pixel of the MCI images.A total of 27 Ω e (θ) or Ω w (θ) estimates were obtained for each annulus of the MCI measurements at the nine zenith angles covered by the annulus with a 1 • interval (three estimates each for Ω e_CC (θ), Ω e_LX (θ), and Ω e_CLX (θ) at nine zenith angles).The center zenith angle of the nine zenith angles was defined as the zenith angle where the MCI images were collected.The Ω e_i or Ω w_i in Equation ( 6) was calculated by averaging the nine Ω e (θ) or Ω w (θ) estimates of the i th annulus.The Ω e (57) or Ω w (57) in Equation ( 2) was treated as the Ω e_i or Ω w_i estimate of the seventh annulus.Two groups of PAI e and WAI e estimates of PAI e_57.3 and WAI e_57.3 as well as PAI e_MCI_0−85 and WAI e_MCI_0−85 were calculated for each sampling point by using Equations ( 2) and ( 6) based on the p e_i (θ i ), p w_i (θ i ), p e (57), and p w (57) measurements.Two α e estimates were calculated for each sampling point by using Equation (12) based on the two groups of PAI e and WAI e estimates.The α e of each plot was calculated by averaging the α e of the two sampling points.Six groups of PAI and WAI estimates were calculated for each plot by using Equations ( 2) and ( 6) based on the two groups of PAI e and WAI e estimates (PAI e_57.3 and WAI e_57.3 as well as PAI e_MCI_0−85 and WAI e_MCI_0−85 ) and the three Ω e and Ω w estimates derived from CC, LX, and CLX (three estimates each for group PAI 57.3 and WAI 57.3 as well as PAI MCI_0−85 and WAI MCI_0−85 ).Six α estimates were obtained for each plot by using Equation ( 14) based on the six groups of PAI and WAI estimates.
Classifying the leaf-on period DHP images into sky, shoot, and woody components by using the method described in Sea et al. [22] entails a large uncertainty.Therefore, the classified leaf-on period MCI images were used to estimate the α e of leaf-on forest plots if the percentage method was adopted.The f c or f w was calculated by taking the ratio of the sum of the canopy element or woody components pixels to the sum of the image pixel number of the classified leaf-on period MCI images.Afterwards, the α e of each sampling position was calculated by using Equation (15) based on the obtained f c and f w estimates.The α e of each plot was then calculated by averaging the two α e estimates of the two sampling points.

The α Measurements of The Representative Trees and Plots
The trees in each plot were categorized into dominant, codominant and suppressed classes based on DBH.The representative trees harvested for each plot were selected according to the three DBH classes.Two or three representative trees with a distance of less than 100 m between the trees and plots were harvested to measure the α of each tree.One tree for each DBH class was harvested for each plot except for plot 2. For plot 2, only one dominant and one codominant tree were harvested because all suppressed trees were felled in the winter of 2014.Those trees located within plots were not harvested because these plots were established as long-term observation sites.The canopy structure characteristics of the harvested trees were guaranteed to be similar to those trees inside each plot.All representative trees were felled in the middle of August and the beginning of September in 2017.
The branches (live or dead) and shoots of each harvested tree were categorized according to their height from 1.2 m above ground (equals to instrument height) to the tree top with a 1 m interval.The branches and shoots located below 1.2 m were classified as individual height class.Both the stem and branch sections were assumed to be circular truncated cone.The top and base diameters of the stem section of each height class were measured by using a tape or digital caliper in the field.The branches and shoots of each height class were bagged and labeled in the field before being immediately shipped to the laboratory.After clipping the shoots from the live branches, the branch areas of both dead and live branches were measured by using a digital caliper.In this study, the branches attached directly to the stem were defined as level 1 branches whilst those attached to level 1 branches were defined as level 2 branches.In this way, the branches of each height class of the harvested tree were classified into six levels.The levels 1 and 2 branches were divided into branch sections with a branch length of 15 cm to facilitate the measurement of branch area.Meanwhile, the other four levels branches were divided into branch sections with a branch length of 10 cm.The bottom and top diameters of each branch section were measured by using a digital caliper.The branch base and tip diameters as well as the branch length were measured for those branch sections with a branch tip.Aside from stems and branches, the fruit area of each height class for the harvested trees was measured under the assumption that these fruits are spheroids.The long and short axes of these fruits were measured by using a digital caliper.
For each harvested tree, approximately 15 typical shoots were randomly selected from each 1 to 5 adjacent height classes to measure the specific leaf area (SLA).The adjacent height classes with similar structural characteristics of shoots and needles shared the same SLA.The number of adjacent height classes sharing the same SLA ranged from two to five for all harvested trees.Approximately 300 to 350 typical needles were randomly selected from the sampled typical shoots to measure the SLA in the laboratory.Twisted, incomplete, and abnormal needles were discarded during the SLA measurements.The hemisurface area of the needles was measured using the volume displacement method described in Chen et al. [33].Typical needle for measuring SLA, whilst the needles of each height class of the harvested trees were dried in an oven at 67 • for approximately 48 hours until the weight of the needles was almost unchanged.Afterwards, the dry mass of the typical needles and the needles of each height class were weighted by using high-precision electronic balances.
Two leaf and woody components area estimates were obtained for each harvested tree with measurement heights of 0 m and 1.2 m, respectively.A measurement height of 0 m or 1.2 m means that all needles and woody components or those needles and woody components with a height of above 1.2 m were considered in calculating the leaf and woody components areas of each harvested tree.The leaf area of each height class for the harvested trees was calculated by multiplying the dry mass and SLA of each height class.The leaf area of each harvested tree is equal to the sum of the leaf area of all height classes which height exceeded the measurement height.Meanwhile, the woody components area of each harvested tree is equal to the sum of the woody components areas of the branch and stem sections of those height classes exceeding the measurement height.Two α estimates were obtained for each harvested tree with the measurement heights of 0 m and 1.2 m.Afterwards, two α estimates with the measurement heights of 0 m and 1.2 m were obtained for each plot by using Equation (15) based on the α d , α c , and α s estimates of the harvested trees and the B d , B c , and B s estimates of the plot.

α Estimates Obtained from The Destructive Method
The proportions of branch and stem areas to the woody components area of each harvested tree range from 16% to 32% and from 63% to 84% in the five plots, respectively (Figure 4).The branch areas of the harvested trees mainly comprise of the branch areas of B1, B2 and B3 (Figure 4).Meanwhile, the proportions of B6 and sum of B4, B5, and B6 branch areas to the woody components area of each harvested tree are smaller than 0.1% and 3% in the five plots, respectively (Figure 4).Interestingly, the proportions of the fruit area to the woody components area of each harvested tree vary obviously amongst the five plots (Figure 4).For example, the proportions of the fruit area to the woody components area of each harvested tree range from 7% to 14% in plot 1 yet range from 0% to 7% in the other plots (Figure 4).Meanwhile, the proportions of the sum of B3, B4, and B5 branch areas and fruit area to the woody components area of each harvested tree range from 5% to 28% in the five plots (Figure 4).
Table 3 shows that the proportions of the woody components area of those woody components with heights of <1.2 m to the woody components area of the harvested trees in each plot slightly vary amongst the five plots at a range of 2% to 6%.However, the proportions of the woody components area of the branches located below live canopies to the woody components area of the harvested trees in the five plots largely vary at a range of 1% to 21% (Table 3).The proportions of the woody components area of all branches located below live canopies to the woody components area of all harvested trees in each plot are the largest in plots 4 and 5 (0.18 and 0.21, respectively; Table 3) and smallest in plots 1 and 2 (0.01; Table 3).
Large variations can be observed amongst the α of the harvested trees in plots 3 to 5, especially plot 3 (Table 4).The variations in value and proportion of the α of harvested trees are 0.19 and 112% in plot 3 and <0.03 and 19% in plots 1 and 2 (Table 4).Large differences are also observed in the α of the five plots estimated with a measurement height of 0 m (Table 4).For example, the difference in proportion of the α in plots 1 and 5 reaches as large as 56% (Table 4).Compared with the α of plots 1 and 5, the difference in proportion of the α in plots 3 and 4 is relatively smaller and equivalent to 19% (the measurement height is 0 m).Meanwhile, only minor differences are observed between the α of plots 1 and 2 (2%; the measurement height is 0 m) (Table 4).The α of plots 4 and 5 are larger than those of the other three plots (Table 4).In sum, the variations in proportion of the α in the five plots estimated with two measurement heights of 0 m and 1.2 m range from 2% to 6% (Table 4).Hereinafter, only the plot α estimated using the destructive method with a measurement height of 1.2 m is regarded as the reference α of the five plots.Table 3. Proportions of the woody components area of those woody components with heights of <1.2 m (including stem and branches) or the branches located below live canopies to the woody components area of the harvested trees in each plot.Table 3. Proportions of the woody components area of those woody components with heights of <1.2 m (including stem and branches) or the branches located below live canopies to the woody components area of the harvested trees in each plot.

α e Estimates Obtained from DHP and MCI
Figure 5a shows that the proportions of the woody components pixel number to the classified leaf off period MCI image pixel number at zenith angles ranging from 0 • to 80 • with a 10 • interval are remarkably larger than the proportions of the woody components pixel number to the classified leaf on MCI image pixel number at the same zenith angle in the five plots.For example, the proportions of the woody components pixel number to the classified leaf-on and leaf off periods MCI image pixel number in plot 1 are 0.15 and 0.62 at zenith angles 0 • as well as 0.50 and 0.96 at zenith angles 80 • , respectively (Figure 5a).For leaf-on forest canopies, the proportions of woody components shaded by shoots at zenith angles ranging from 0 • to 80 • with a 10 • interval are large, especially at small zenith angles, and tend to decrease along with zenith angles in the five plots (Figure 5b).For example, the proportion of woody components shaded by shoots can reach as high as 0.70 at a 0 • zenith angle, and decreases to 0.35 at a 80 • zenith angle in plot 1 (Figure 5b).Meanwhile, the proportions of woody components shaded by shoots at a zenith angle ranging from 0 • to 80 • with a 10 • interval are larger at those plots with fewer branches located below live canopies (plots 1 and 2) compared with those plots with a large number of branches (plots 4 and 5) (Figure 5b).5a shows that the proportions of the woody components pixel number to the classified leaf off period MCI image pixel number at zenith angles ranging from 0° to 80° with a 10° interval are remarkably larger than the proportions of the woody components pixel number to the classified leaf on MCI image pixel number at the same zenith angle in the five plots.For example, the proportions of the woody components pixel number to the classified leaf-on and leaf off periods MCI image pixel number in plot 1 are 0.15 and 0.62 at zenith angles 0° as well as 0.50 and 0.96 at zenith angles 80°, respectively (Figure 5a).For leaf-on forest canopies, the proportions of woody components shaded by shoots at zenith angles ranging from 0° to 80° with a 10° interval are large, especially at small zenith angles, and tend to decrease along with zenith angles in the five plots (Figure 5b).For example, the proportion of woody components shaded by shoots can reach as high as 0.70 at a 0° zenith angle, and decreases to 0.35 at a 80° zenith angle in plot 1 (Figure 5b).Meanwhile, the proportions of woody components shaded by shoots at a zenith angle ranging from 0° to 80° with a 10° interval are larger at those plots with fewer branches located below live canopies (plots 1 and 2) compared with those plots with a large number of branches (plots 4 and 5) (Figure 5b).Table 5 shows that the α e estimated by using the leaf-on and leaf-off periods DHP or MCI images is 2.29 to 4.63 times the reference α of the five plots regardless of the optical method and inversion model used in the α e estimation.For example, the α e derived using the leaf-on and leaf-off periods DHP images at plot 2 was ~4.63 times the reference α of plot 2 if the LAI-2200 inversion model was adopted in the α e estimation (Tables 4 and 5).Similarly, the α e derived using the leaf-on and leaf-off periods MCI images at plot 2 was ~3.5 times the reference α of plot 2, if the MCI_0-85 inversion model was adopted in the α e estimation (Tables 4 and 5).The RMSE and MAE of the α e derived using the leaf-on and leaf-off periods MCI or DHP images were larger than 0.41 (205%) and 0.40 (203%) in the five plots, respectively (Table 6).Similarly, the α e derived from the percentage method by using the leaf-on period MCI images was 2.19 to 2.57 times the reference α of the five plots (Table 5).However, the α e obtained by using the leaf-on period MCI images based on 57.3 and MCI_0-85 was closer to the reference α of the five plots whilst its RMSE and MAE were <0.05 (26%) and 0.04 (22%), respectively (Tables 4-6).Small differences (<7%) were observed between the α e estimates obtained by using the leaf-on period MCI images based on 57.3 and MCI_0-85 in the five plots (Table 5).Similarly, the variations in proportion between the α e estimates obtained using the leaf-on and leaf-off periods MCI images based on 57.3 and MCI_0-85 were <13% in the five plots (Table 5).By contrast, large differences were observed between the α e estimates obtained by using the leaf-on and leaf-off periods DHP images based on 57.3 and Miller (3% to 27%), Miller and LAI-2200 (11% to 40%), and LAI-2200 and 57.3 (5% to 15%) in the five plots (Table 5).

α Estimates Obtained from DHP and MCI
The α e estimates derived without consideration of Ω e and Ω w by using the leaf-on and leaf-off periods DHP or MCI images were systematically larger than those derived with consideration of Ω e and Ω w by using the same inversion model in the five plots (with the exception of the α e estimates obtained by using the leaf-on and leaf-off periods DHP images based on the Miller inversion model in plots 2 and 3) (Tables 4, 5, 7 and 8).Specifically, those α estimates obtained with consideration of Ω e and Ω w by using the leaf-on and leaf-off periods DHP or MCI images were 1.57 to 4.13 times larger than the reference α of the five plots (Tables 4, 7 and 8).The RMSE and MAE of these α estimates derived using the leaf-on and leaf-off periods DHP or MCI images with consideration of Ω e and Ω w were >0.30 (153%) and 0.30 (150%) or 0.21 (107%) and 0.20 (103%) in the five plots, respectively, regardless of the Ω e and Ω w estimation algorithm and inversion model adopted in the estimation (Tables 9 and 10).Based on RMSE, MAE, and R 2 , the accuracies of the α estimates derived using the leaf-on period MCI images with consideration of Ω e and Ω w showed little to no improvement compared with those estimates derived from the same inversion model but without consideration of Ω e and Ω w in the five plots (Tables 6 and 10).The RMSE and MAE of the α estimated with consideration of Ω e and Ω w by using the leaf-on period MCI images were <33% and 27% in the five plots, respectively (Table 10).
For the DHP method, only slight differences (0% to 16%) were observed in the α estimated using the leaf-on and leaf-off periods DHP images based on the same inversion model but with different Ω e and Ω w estimation algorithms in the five plots (Table 7).No significant differences (0% to 22%) were also observed in the α estimated using the leaf-on and leaf-off periods DHP images based on the two inversion models of LAI-2200 and 57.3 but with the same Ω e and Ω w estimation algorithm in the five plots (Table 7).By contrast, large variations (0% to 30%) were observed in the α estimated using the leaf-on and leaf-off periods DHP images based on two groups of inversion models (Miller and 57.3 as well as Miller and LAI-2200) but with same Ω e and Ω w estimation algorithms in the five plots (Table 7).
For the MCI method, only small differences (2% to 18%) were observed in the α estimated using the leaf-on period or leaf-on and leaf-off periods MCI images based on the same inversion model but with different Ω e and Ω w estimation algorithms (except for CC) in the five plots (Table 8).However, large differences (0% to 82%) were observed in the α estimated using the leaf-on period or leaf-on and leaf-off periods MCI images based on the same inversion model but with different Ω e and Ω w estimation algorithms (include CC) in the five plots (Table 8).Meanwhile, only slight differences (0% to 9%) were observed in the α estimated using the leaf-on period or leaf-on and leaf-off periods MCI images based on the two inversion models with the same Ω e and Ω w estimation algorithms (excluding CC) in the five plots (Table 8).By contrast, obvious differences (0% to 45%) were observed in the α estimated using the leaf-on period or leaf-on and leaf-off MCI images based on the two inversion models with CC in the five plots (Table 8).

Factors That Affect the Accuracy of the Reference α Estimates
The destructive method requires considerable time and effort to measure the reference α of forest plots.To reduce workload, some studies have simplified the measurement procedures, such as by discarding the fruits and parts of branches of the harvested trees [16].For example, all branches (excluding B1, B2, and the end level branches) and fruits were not considered in estimating the reference α of forest plots in Chen [16].However, the sum of the proportions of the B3, B4, and B5 branch areas and fruit areas to the total woody components area of the harvested trees ranges from 5% to 28% in the five plots (Figure 4).Therefore, simplifying the measurement procedures of the destructive method as can be seen in Chen [16] can introduce large bias in the reference α estimation of L. gmelinii forest plots.Such simplification cannot be directly used to measure the reference α of forest plots with different plant functional types.At least, this simplified method [16] is not accurate enough to estimate the reference α of L. gmelinii forest plots.However, based on the results in Figure 4, B4, B5, and B6 can all be neglected in estimating the reference α of L. gmelinii forest plots due to their small contributions to the woody components areas of the harvested trees in the five plots.By contrast, fruits should be considered in estimating the reference α because they are widely distributed across the forest canopies (Figure 4).
The mismatch between the measurement heights of the destructive and optical methods has usually been ignored in previous studies [16,23,33].Specifically, these studies have estimated the reference α of forest plots with a measurement height of 0 m, but their optical instruments are usually operated and levelled above ground in the plots.Meanwhile, in this study, the DHP and MCI instruments were mounted on a tripod with a measurement height of about 1.2 m above ground.Therefore, the woody components and leaves or shoots of forest canopies with a height of < 1.2 m cannot be sampled by optical methods.Neglecting the mismatch in the measurement heights of harvesting and optical methods can introduce approximately 2% to 6% bias into the estimation of the woody components area and α of the harvested trees in the five plots (Tables 3 and 4).Therefore, such mismatch must be considered when evaluating the performance of optical methods in estimating the α of forest plots.

Impact of Tree Age, Stand Density, Site Conditions, and Management Activities on the Reference α of Forest Plots
Management activities, such as removing the branches below live canopies are usually performed by forest managers to promote tree growth in forest plots.The proportion of the woody components area of all woody components located below live canopies to the woody components area of the harvested trees in plot 5 is 0.21 (Table 3), which indicates that management activities affect the reference α of forest plots obviously.If the woody components of the trees located below live canopies were not harvested for plots 1 and 2 during the thinning campaign prior the experiment and if the proportions of the woody components area of the woody components located below live canopies to the woody components area of the trees in plots 1, 2, and 5 were assumed to be the same, then the reference α of plots 1 and 2 would become equal to 0.19 and 0.19, respectively.The variations in proportion between the recalculated α of plots 1 and 2 and the reference α of plot 5 were approximately 21%.Given that no obvious differences can be found in the site conditions of the three plots, the variations between the recalculated α of plots 1 and 2 and the reference α of plot 5 were mainly attributed to tree age and stand density.On the other hand, the large variations in the reference α of plots 1, 2, and 5 also be indicative of the larger role of management activities (relative to tree age or stand density) in increasing the large difference in the reference α of the five plots.
Plots 3 and 4 sharing similar characteristics, including tree age, mean tree height, and average DBH (Table 2).Furthermore, those branches with a height of ~≤ 2.5 m were low-and medium-level harvested in plots 3 and 4 prior the experiment, respectively.Plot 4 had poorer site conditions compared with plot 3 because this plot had many yellow needles distributed amongst its canopies and had a larger dry soil depth as compared with plot 3.However, the reference α of plot 4 was 20% larger than that of plot 3 thereby indicating that site conditions affect the reference α of forest plots obviously (Table 4).
In conclusion, tree age, stand density, site conditions, and management activities contribute to the large differences in the reference α of the five selected L. gmelinii plots (Tables 2 and 4).Therefore, these four factors must be considered when estimating the reference α of L. gmelinii plots.Amongst these factors, site conditions and management activities have the largest contributions to the variations in the reference α of the L. gmelinii plots.Meanwhile, the reference α of those forest plots with the same tree species and similar characteristics of stand density and tree age yet with different levels of branch removal and site conditions must be dissimilar.A specific α estimate should be obtained for those forest plots (Tables 2-4).The same can be said for those forest plots with the same tree species, branch removal levels and site conditions yet with different stand density and tree age (Tables 2-4).However, those forest plots with the same tree species and similar characteristics of stand density, site conditions, tree age, and management activities (e.g., plots 1 and 2) must show slight variations in their reference α (Tables 2 and 4).Chen [16] also reported small differences (13%) in the reference α of two Picea mariana plots with similar stand density, tree age, tree height, and management activities.This argument is consistent with the findings of this study.

Factors That Affect the α Estimation of the Two Optical Methods
Figure 5 and Tables 5-10 reveal that the performance of DHP and MCI in estimating the α of forest plots is influenced by four factors, namely, the inversion model, the clumping effects of canopy element and woody components, the Ω e and Ω w estimation algorithms, and the preferential shading of woody components by shoots.Amongst these factors, the inversion model only slightly influences the α e or α estimation of forest plots, if 57.3, LAI-2200 and MCI_0-85 are used in the estimation (Tables 5,  7 and 8).Meanwhile, the α e or α estimates obtained using the leaf-on and leaf-off periods DHP images based on Miller and 57.3 as well as Miller and LAI-2200 showed relatively larger differences compared with those estimated by 57.3 and LAI-2200.This finding is consistent with the conclusions of Zou et al. [15], who reported small differences (<10%) in the WAI e , PAI e , WAI, and PAI of leaf-on and leaf-off forest canopies estimated using LAI-2200 and 57.3.Relatively large differences (2% to 17%) were also observed in the WAI e , PAI e , WAI, and PAI that were estimated using Miller and 57.3 as well as Miller and LAI-2200 [15].Such variations can be attributed to the fact that the gap fraction measurements collected from DHP at zenith angles near the horizon tend to approach zero.When using Miller in estimating PAI and WAI, the null gap fraction measurements could introduce errors into the estimations even if the null gap fraction measurements processing solutions were applied [15].Zou et al. [15] pointed out that several inversion models, including LAI-2200 and 57.3, can avoid such effect.Given that accurate α estimates were obtained by using leaf-on period MCI images based on MCI_0-85 (similar to LAI-2200) and 57.3 in the five plots (Tables 8 and 10), both 57.3 and MCI_0-85 suggest to be used in the α estimation of forest plots.
The α e or α estimates obtained with or without consideration of Ω e and Ω w using the leaf-on and leaf-off periods DHP or MCI images were 1.57 to 4.63 times the reference α of the five plots (Tables 4, 5, 7 and 8).Given that the WAI e and WAI of leaf-on forest plots were derived by using leaf-off period DHP or MCI images, the overestimation of α e and α indicates that the PAI e and PAI of leaf-on forest plots are significantly underestimated when leaf-on period DHP or MCI images are used in the PAI e and PAI estimation.Moreover, the underestimation of the PAI e and PAI of leaf-on forest plots was mainly driven by the underestimation of WAI e and WAI because of the obvious differences in the proportions of woody components pixel number to the leaf-on and leaf-off periods MCI image pixel numbers at the zenith angles ranging from 0 • to 80 • with a 10 • interval in the five plots (Figure 5a).The underestimation of the WAI e and WAI of leaf-on forest plots can also be reflected in the large proportions of woody components shaded by shoots at zenith angles ranging from 0 • to 80 • with a 10 • interval (Figure 5b).The shading of woody components by shoots in canopies increases the woody components gap fraction of forest canopies as their WAI should be, thereby resulting in the underestimation of WAI e and WAI.However, the overestimation of α e or α can be largely reduced if only the leaf-on period MCI images are used in the α e or α estimation (Tables 4, 5 and 8) because the underestimation of the denominator (PAI e or PAI) and numerator (WAI e or WAI) cancels out the α e or α estimation errors to a large degree.The underestimation of the WAI e and WAI of leaf-on period forest plots as estimated using leaf-on period MCI or DHP images reveals the obvious overestimation of α e and α in the five plots, which were obtained using leaf-on and leaf-off periods DHP or MCI images, can be attributed to the fact that the woody components are preferentially shaded by shoots.In L. gmelinii forest canopies, the shoots are directly located on the branch and stem surfaces, thereby leading to preferential shading of branches and stems.
The woody components of forest canopies include stems, branches, and fruits.Given that these woody components are not randomly distributed in space similar to leaves or shoots, both the clumping effects of the canopy element and woody components should be considered in estimating the α of forest plots.However, minor differences were observed between the RMSE and MAE of the α e and α estimates obtained using the leaf-on period MCI images without or with consideration of Ω e and Ω w (Tables 6 and 10).Therefore, the clumping effects of the canopy element and woody components of leaf-on forest canopies are not effectively evaluated by the adopted Ω e and Ω w estimation algorithms regardless of the Ω e and Ω w estimation algorithm and segment size used in the α e and α estimation.Previous studies have also reported that the Ω e and Ω w estimation algorithms used in this study do not always perform well in estimating the Ω e and Ω w of all leaf-on and leaf-off forest canopies with different plant functional types, tree species composition, PAI, and WAI [15,27,46].For example, Zou et al. [15] observed small differences between the majority of the WAI estimates derived using 57.3 or LAI-2200 with consideration of Ω w (derived from CC, LX, and CLX) and the WAI of leaf-off birch forest plots.However, large differences were observed between the majority of the PAI estimates (derived by using the same inversion model and Ω e estimation algorithm) and the PAI of leaf-on birch forest plots [15].This report indicates that the Ω w of leaf-off birch forest plots was effectively estimated by using the three Ω e and Ω w estimation algorithms; by contrast, the Ω e of leaf-on birch forest plots was not effectively estimated by using these algorithms [15].Accurately estimating α depends on both accurate PAI and WAI estimates whilst accurately estimating PAI and WAI estimates, in turn, relies on both accurate estimates of Ω e and Ω w .Therefore, the mismatch between the performance of Ω e and Ω w estimation algorithm in estimating the Ω e and Ω w of leaf-on and leaf-off forest plots can introduce errors into the estimation of α.The fact of the preferential shading of woody components by shoots can also influence the performance of Ω e and Ω w estimation algorithms in evaluating the clumping effects of the canopy element and woody components of leaf-on forest canopies.

Determining Whether Accurate α Estimates Can Be Obtained from DHP and MCI
Tables 4-10 reveal that relatively accurate α estimates with estimation errors of lower than 20% can be obtained from MCI in the L. gmelinii plots.If no errors are exist in the LAI estimation of L. gmelinii forest plots except in estimating α, then MCI can be considered an effective method that meets the required maximum LAI estimation error (20%) for in situ LAI estimates specified by GCOS.Unlike DHP, MCI can discriminate leaves or shoots, woody components and sky based on the NIR and VIS band images of forest canopies.DHP has the disadvantage of the insufficient sampling of forest canopies at zenith angles close to the zenith due to the intrinsic defect of hemispherical photography.The proportions of woody components shaded by shoots at zenith angels ranging from 0 • to 85 • with a 10 • interval are large at those zenith angles close to the zenith in the five plots (Figure 5b).Therefore, for DHP, the estimation errors in the p e (θ), p w (θ), Ω e (θ), and Ω w (θ) measurements with zenith angles close to zenith can introduce bias into the α estimation of forest plots.By contrast, MCI can obtain images with sufficient sampling of forest canopies at all directions of the upper hemisphere and is thereby recommended for estimating the α of forest plots.
The methods adopted in previous studies [17][18][19][20] for estimating the α of forest plots by using leaf-on and leaf-off periods DHP images or radiation attenuation measurements can introduce significant errors into the estimations (Tables 4-10).Therefore, these methods would not suggested to be used in estimating the α of forest plots, at least for L. gmelinii forest plots.The percentage method also overestimated the α in the five plots (Table 5), thereby indicating that this method must not be used in estimating the α of L. gmelinii forest sites.This conclusion contradicts the report of Woodgate [13], who found a 3% difference in the α estimated by the percentage method and the reference α of eucalypt forest scenes.The contrasting canopy structures of the tree species covered in Woodgate and in this study may explain such contradiction.

Limitations and Perspectives
Given the huge amount of time and effort required by the destructive method to measure the α of forest sites, only two or three representative trees were harvested for each plot in this study.Considering the obvious variations between the α of the harvested trees in each plot (Table 4), more representative trees must be harvested to obtain highly reliable and accurate estimations of the reference α of forest plots.For example, three trees can be harvested for each class of dominant, codominant, and suppressed trees in each plot.Given that only five plots were covered in this study, the impact of tree age or stand density on estimating the reference α of forest plots was not evaluated.More forest plots (>5) must be covered in the future to obtain highly reliable conclusions regarding the influence of tree age, stand density, site conditions and management activities on estimating the reference α of L. gmelinii forest plots.More plant functional types must also be covered in the future to check whether MCI can effectively measure the α of forest sites with different tree species compositions.

Conclusions
The conclusions of this study are as follows.(1) The woody components area of harvested trees mainly comprises of the woody components areas of branches and fruits.Caution is needed when parts of woody components, especially fruits, are discarded in estimating the α of harvested trees.The mismatch in the measurement heights of destructive and optical methods should also be considered in the α estimation of forest plots.(2) Large variations (0% to 56%) are observed in the reference α of the five plots with the same tree species yet with different tree age, stand density, site conditions, and management activities.Moreover, the impact of site conditions and management activities on the reference α of forest plots is larger than that of tree age and stand density.(3) The performance of DHP and MCI in estimating the α of forest plots is obviously affected by the inversion model and the preferential shading of woody components by shoots in leaf-on canopies.The large overestimation of the α (i.e., 1.57 to 4.63 times larger than the reference value for this study) that was derived using leaf-on and leaf-off periods DHP or MCI images can be mainly attributed to the preferential shading of woody components by leaves or shoots.The impact of the adopted Ω e and Ω w estimation algorithms on estimating the α of forest plots is relatively small.(4) Relatively accurate α estimates (with estimation errors of <20%) can be obtained from MCI if the appropriate inversion model and Ω e and Ω w estimation algorithm are used in the α estimation of L. gmelinii forest plots.
Ω e -canopy element clumping index Ω e (θ)-canopy element clumping index at θ Ω e (57)-canopy element clumping index at zenith angle 57 • Ω e_CC (θ)-canopy element clumping index estimated using the gap size distribution algorithm at θ Ω e_CC_k (θ)-canopy element clumping index of segment k at θ estimated using CC Ω e_LX (θ)-canopy element clumping index estimated using the logarithmic averaging algorithm at θ Ω e_CLX (θ)-canopy element clumping index estimated using the combination of gap size and logarithmic averaging algorithm at θ Ω e_i -Ω e estimate of the ith annulus Ω w -woody component clumping index Ω w (θ)-woody component clumping index at θ Ω w (57)-woody component clumping index at zenith angle 57 • Ω w_CC (θ)-woody component clumping index estimated using the gap size distribution algorithm at θ Ω w_LX (θ)-woody component clumping index estimated using the logarithmic averaging algorithm at θ Ω w_CLX (θ)-woody component clumping index estimated using the combination of gap size and logarithmic averaging algorithm at θ Ω w_i -Ω w estimate of the ith annulus γ e -needle-to-shoot area ratio γ-effective needle-to-shoot area ratio α-woody-to-total area ratio α e -effective woody-to-total area ratio

Figure 1 .
Figure 1.Location of the Saihanba National Forest Park in Hebei Province, China.Figure 1. Location of the Saihanba National Forest Park in Hebei Province, China.

Figure 1 .
Figure 1.Location of the Saihanba National Forest Park in Hebei Province, China.Figure 1. Location of the Saihanba National Forest Park in Hebei Province, China.
) to (4) based on the three groups of PAI e and WAI e estimates (PAI e_57.3 and WAI e_57.3 , PAI e_Miller and WAI e_Miller , as well as PAI e LAI −2200 and WAI e_LAI−2200 ) and seven Ω e and Ω w estimates derived from the combinations of the three algorithms of CC, LX, and CLX and the four segment sizes of 5 • , 15 • , 30 • , and 45 • (seven estimates for each group of PAI and WAI estimates, including PAI 57.3 and WAI 57.3 , PAI Miller and WAI Miller , as well as PAI LAI−2200 and WAI LAI−2200 ).Afterwards, 21 α estimates were derived for each plot by using Equation (14) based on 21 groups of PAI and WAI estimates.
) to (4) based on the three groups of   and   estimates ( _57.3 and  _57.3 ,  _ and  _ , as well as    −2200 and  _−2200 ) and seven   and   estimates derived from the combinations of the three algorithms of CC, LX, and CLX and the four segment sizes of 5°, 15°, 30°, and 45° (seven estimates for each group of PAI and WAI estimates, including  57.3 and  57.3 ,   and   , as well as  −2200 and  −2200 ).Afterwards, 21  estimates were derived for each plot by using Equation (14) based on 21 groups of PAI and WAI estimates.

Figure 2 .
Figure 2. Sampling scheme of the digital hemispherical photography (DHP) method.

Figure 2 .
Figure 2. Sampling scheme of the digital hemispherical photography (DHP) method.

ForestsFigure 3 .
Figure 3.The DHP or multispectral canopy imager (MCI) images collected from the same sampling point or the same sampling point with the same zenith and horizontal angles at leafon and leaf-off periods in plot 1.(a) Leaf-on period DHP image, (b) leaf-off period DHP image, (c) leaf-on period visible (VIS) band MCI image, (d) leaf-on period near-infrared (NIR) band MCI image, (e) leaf-off period VIS band MCI image, and (f) leaf-off period NIR band MCI image.

Figure 3 .
Figure 3.The DHP or multispectral canopy imager (MCI) images collected from the same sampling point or the same sampling point with the same zenith and horizontal angles at leaf-on and leaf-off periods in plot 1.(a) Leaf-on period DHP image, (b) leaf-off period DHP image, (c) leaf-on period visible (VIS) band MCI image, (d) leaf-on period near-infrared (NIR) band MCI image, (e) leaf-off period VIS band MCI image, and (f) leaf-off period NIR band MCI image.

Figure 5 .
Figure 5.The proportion of the woody components pixel number to the classified MCI image pixel number (a) and the proportion of the woody components shaded by shoots, which is calculated by using the classified leaf-on and leaf-off periods MCI images (b) at zenith angles ranging from 0° to 80° with a 10° interval in the five plots.The proportion of the woody components shaded by shoots at each zenith angle is derived by taking the ratio of the shaded woody components pixels (calculated by subtracting the woody components pixel number of the classified leaf-on period MCI images from the woody components pixel number of the classified leaf-off period MCI images) to the total classified leaf-on period MCI image pixel number.

Figure 5 .
Figure 5.The proportion of the woody components pixel number to the classified MCI image pixel number (a) and the proportion of the woody components shaded by shoots, which is calculated by using the classified leaf-on and leaf-off periods MCI images (b) at zenith angles ranging from 0 • to 80 • with a 10 • interval in the five plots.The proportion of the woody components shaded by shoots at each zenith angle is derived by taking the ratio of the shaded woody components pixels (calculated by subtracting the woody components pixel number of the classified leaf-on period MCI images from the woody components pixel number of the classified leaf-off period MCI images) to the total classified leaf-on period MCI image pixel number.

Table 4 .
(15)α of the harvested trees estimated via the destructive method with a measurement height of 0 m in each plot.The α of each plot was calculated by using Equation(15)with measurement heights of 0 m and 1.2 m.

Table 5 .
α e estimates obtained from DHP and MCI in the five plots.

Table 6 .
Correlation statistics between the α e estimated from DHP and MCI and the reference α of the five plots.Statistics are given at the 95% confidence level from the two-tailed Student's t-test.

Table 8 .
α estimates obtained with consideration of Ω e and Ω w by using the leaf-on period MCI images in the five plots.

Table 9 .
Correlation statistics between the obtained α (estimated with consideration of Ω e and Ω w by using the leaf-on and leaf-off periods DHP images) and the reference α of the five plots.Statistics are given at the 95% confidence level from the two-tailed Student's t-test.

Table 10 .
Correlation statistics between the obtained α estimates (derived with consideration of Ω e and Ω w by using the leaf-on period MCI images) and the reference α of the five plots.Statistics are given at the 95% confidence level from the two-tailed Student's t-test.