Application of Physically-Based Slope Correction for Maximum Forest Canopy Height Estimation Using Waveform Lidar across Different Footprint Sizes and Locations : Tests on LVIS and GLAS

Forest canopy height is an important biophysical variable for quantifying carbon storage in terrestrial ecosystems. Active light detection and ranging (lidar) sensors with discrete-return or waveform lidar have produced reliable measures of forest canopy height. However, rigorous procedures are required for an accurate estimation, especially when using waveform lidar, since backscattered signals are likely distorted by topographic conditions within the footprint. Based on extracted waveform parameters, we explore how well a physical slope correction approach performs across different footprint sizes and study sites. The data are derived from airborne (Laser Vegetation Imaging Sensor; LVIS) and spaceborne (Geoscience Laser Altimeter System; GLAS) lidar campaigns. Comparisons against field measurements show that LVIS data can satisfactorily provide a proxy for maximum forest canopy heights (n = 705, RMSE = 4.99 m, and R = 0.78), and the simple slope correction grants slight accuracy advancement in the LVIS canopy height retrieval (RMSE of 0.39 m improved). In the same vein of the LVIS with relatively smaller footprint size (~20 m), substantial progress resulted from the physically-based correction OPEN ACCESS Remote Sens. 2014, 6 6567 for the GLAS (footprint size = ~50 m). When compared against reference LVIS data, RMSE and R for the GLAS metrics (n = 527) are improved from 12.74–7.83 m and from 0.54–0.63, respectively. RMSE of 5.32 m and R of 0.80 are finally achieved without 38 outliers (n = 489). From this study, we found that both LVIS and GLAS lidar campaigns could be benefited from the physical correction approach, and the magnitude of accuracy improvement was determined by footprint size and terrain slope.


Introduction
Forest ecosystems are a substantial carbon sink, biodiversity reservoir, and driver of microclimate and ecological processes [1][2][3][4].We are interested in quantifying forest physical characteristics (e.g., canopy heights and stand volumes) since such measures provide reasonable proxies of carbon storage and ecosystem dynamics (e.g., productivity, diversity, and mortality) [5].However, continuously and accurately monitoring forests over a large spatial domain represents a significant challenge.Improvements in remote sensing techniques have addressed this challenge by using active or passive sensors that are sensitive to the structural attributes of forests [1,5,6].The passive system in optical remote sensing has produced biophysical estimations for Leaf Area Index (LAI), biomass, gross primary productivity, and net primary productivity (e.g., [7][8][9][10]).However, these measures are only derived from theoretical conjugations between forest structures and indirect observations (i.e., surface reflection, absorption, and re-emission of solar radiation).In recent studies, active sensors (light detection and ranging (lidar) or radio detection and ranging (radar)) have become more attractive to the remote sensing community because they overcome some constraints of the passive system [11,12] (e.g., indirect measures, dependence on time of day and season, insufficient energy for certain wavelengths, and sensor saturation at highly dense forest [13][14][15]).
Lidar sensors are used to calculate surface heights by measuring the time taken between emission and return of laser pulses.Lidar sensors can use discrete or full-waveform recording systems [16].The discrete lidar detects multiple laser scattering events within a small-footprint: the first return represents the top of obstacles and the last return represents the ground surface (if the laser pulse can penetrate through all canopy components) [17,18].One or more returns could be located in between the first and last return.Forest canopy heights and vertical/horizontal characteristics can be retrieved from the 3-dimensional locations of these scattering events.On the other hand, the full waveform lidar system records a continuous vertical profile of a waveform within a given footprint [19,20].The vertical profile is related to the level of fluctuation in the magnitude, or brightness, of lidar returns as laser pulses pass through the forest canopy.
The Land, Vegetation, and Ice Sensor (a.k.a., the Laser Vegetation Imaging Sensor; LVIS) and the Geoscience Laser Altimeter System (GLAS) are among the leading edge full waveform lidar altimeter systems capable of generating vertical information of the regional and global land surface.
The LVIS instrument, with a relatively small-footprint, is suitable for precise regional-scale mapping of forest canopy height and biomass [21,22].The GLAS, with a larger-footprint and wide spatial coverage, has provided practical means for monitoring various global forest attributes [23][24][25].The performance of LVIS and GLAS for retrieving canopy vertical structure is generally stable over flat terrains where the disturbance by topographic features (e.g., slope and roughness) is minimized.However, retrieving forest canopy height from backscattered lidar signals over the steep slopes or uneven topography found in mountainous regions is not easy [26][27][28][29][30][31][32][33][34][35], because the vertical extent of waveforms collected by the sensor increases as a function of terrain slope and footprint size [27,32].Also, the combination of a laser pulse pointing at an off-nadir angle and terrain slope alters waveform extent in both directions (widening or narrowing) [32,33].Thus, quantifying and minimizing slope and off-nadir angle effects on lidar waveform are key challenges to overcome in the application of lidar remote sensing in accurate forest structure mensuration.
To understand complex waveform characterization and to obtain more accurate vegetation canopy height estimation, several radiative transfer models have been developed for demonstrating lidar waveform formulation over various sensing environments [31][32][33].Amongst models, the extended Geometric Optical and Radiative Transfer (GORT) model of Yang et al. [32] provides a useful mean to quantify the influence of topographic features, footprint size, and off-nadir angle on lidar waveform characterization (i.e., height level and magnitude of backscattered waveform).Their physically-based approach is able to correct the impact of surface topography and footprint size on vegetation height estimates with an assumption that vegetation canopy is uniformly distributed within each footprint.A subsequent study by Lee et al. [29] implements this correction approach to circumvent slope-related errors in the forest canopy height estimation using LVIS and GLAS data over a small homogeneous forested area.The simple physically-based correction has resulted in a significant accuracy improvement and quantified slope effects on height estimation, but their approach still needs to be tested over a larger area with diverse forest types and topographic conditions.
The objective of this study is therefore to practice the simple physically-based slope correction approach of Lee et al. [29] for the estimation of maximum forest canopy heights over additional study sites based on the two different footprints of LVIS and GLAS data.We examine the LVIS height metrics (with/without the slope correction) based on comparisons against field measurements.Then, we investigate the accuracy improvement in GLAS height retrievals using the proposed slope correction approach.Proving the efficacy of the slope correction method for maximum canopy height retrieval, we attempt to quantify terrain slope and footprint size effects on canopy height estimation by relating a gradient of slope to observed errors in both lidar sensors.Off-nadir pointing effect on height retrieval is not considered in this study, because the mean off-nadir angles of LVIS and GLAS data used are negligible (2.3° ± 1.4° for LVIS and 0.35° ± 0.05° for GLAS) and there is a lack of accurate information of sensor zenith and azimuth angles and slope orientation.Instead, we briefly discuss possible off-nadir angle effects as well as additional error sources in a separate section.
This paper is organized as follows: Section 2 lists data used in this study and Section 3 presents a thorough description of methods for retrieving maximum forest canopy heights from LVIS and GLAS altimetry.Section 4 compares results of maximum forest canopy height estimates from field, LVIS, and GLAS datasets and discusses a practical aspect of slope correction method, additional error sources (e.g., geolocation and off-nadir pointing), and limitations of our study.In Section 6, we summarize and conclude our study.

Field Measurements
We used field-measured maximum canopy heights from the North American Carbon Program, which is a multi-disciplinary research program designed to understand North America's carbon sources and sinks [36].These data were obtained from the Oak Ridge National laboratory Distributed Active Archive Center [37].Field measurements in the datasets are from five experimental forest stations of New England for the year 2009: Bartlett Experimental Forest (BF), Harvard Forest (HF), Howland Research Forest (HO), Hubbard Brook Experimental Forest (HB), and the Penobscot Experimental Forest (PE).Sierra National Forest (SN) for the year 2008 is an additional site of field measurements (Figure 1).Main plots (50 × 200 m) of the five New England sites consist of 16 subplots (25 × 25 m).The maximum forest canopy heights were derived from extracting the tallest individual tree in each subplot.For the SN study site, the main plot is 100 × 100 m, with nine subplots (33.3 × 33.3 m).In this study, 1010 field-measured subplots were used.In both New England and SN sites, Global Positioning System (GPS) locations are provided for four corners of each subplot, which were used to find co-located lidar footprints.Table 1 lists plot information and corresponding LVIS data acquisition years for each field site.

LVIS
The LVIS instrument is an airborne and full-waveform lidar sensor with a 1064 nm wavelength laser (footprint size = ~20 m), developed by the National Aeronautics and Space Administration (NASA) [38].This sensor has been used to monitor canopy vertical structure and surface topography over several forested sites since 2003.We used 13 LVIS datasets whose geographical locations are overlapped with field measurements and GLAS data (Figure 1).The LVIS datasets were separated into two groups for this study.The first was used to examine whether the physically-based correction method improves the accuracy of LVIS maximum forest canopy heights across different study sites, and to showcase whether LVIS derivations can be used as an appropriate proxy of field measurements (Table 1).The second was then used to evaluate the GLAS data with larger-footprints by following the approach of Lee et al. [29] (Table 2).All acquisition years of field and LVIS data are identical (Table 1), but those of LVIS and GLAS data are variant (up to 3 years of differences) (Table 2).Note that two LVIS groups are independent of each other and data acquisition dates could be different even though their sites are identical.Details of LVIS uses are provided in Section 3.

Table 2. Datasets for inter-comparisons between LVIS derived heights and GLAS height metrics.
track, with a maximum cross track separation of 15 km at the equator [39].A 1064-nm laser pulse is implemented in the GLAS instrument, but the returned laser energy of GLAS forms an ellipsoidal footprint.The National Snow and Ice Data Center provides GLAS data that are acquired from 2003-2009.In this study, the latest release of GLAS laser altimetry data (Release 33) was used.Among six altimetry products (GLA 05, 06, 12, 13, 14 and 15) of release 33 GLAS data, the level-2 altimetry product (GLA14) was employed for the maximum forest canopy height estimation.GLA14 is useful for obtaining land surface elevation, laser footprint geo-location, and waveform parameters such as signal beginning and echo energy peaks.This study selected GLAS data taken from May to October for the years 2005 and 2006, as this period best approximates the growing season in the Northern hemisphere [40].This selection scheme may avoid underestimations of forest canopy height due to the backscattered signal from leafless deciduous species.Moreover, the GLAS laser campaign (i.e., Laser-3) in 2005 and 2006 has the smallest ellipticity of its footprints (nearly circular shape [39]) and this stability in shape and size of GLAS footprints can reduce uncertainties in the surface height estimation [39].Based on the length of major axis and ellipticity of the GLAS dataset, this study assumed all the GLAS footprints were circular with a 50 m diameter [41].

Ancillary Datasets
Three ancillary datasets used in this study were elevation, slope, and land cover classes over the continental USA (CONUS).National Elevation Dataset (NED) of the United States Geological Survey provides seamless and nationally consistent Digital Elevation Model (DEM) data with various spatial resolutions.This study used a standard NED product at the one arc-second resolution (~30 m), which covers the CONUS and island territories except for Alaska.Overall absolute vertical accuracy of the NED DEM is 2.44 m root-mean-square-error (RMSE) [42].Spatially continuous slope data were generated using the ArcGIS 10.0 software.Holmes et al. [43] reported that the accuracy of NED-derived slope data is less than 3° at the one arc-second scale.DEM and slope of ground surface associated with LVIS and GLAS were calculated based on the average of all the grid cells intersecting with a LVIS or GLAS footprint.The National Land Cover Database (NLCD) 2006 land cover data over the CONUS at the 30 m spatial resolution is derived from the Landsat 7 Enhanced Thematic Mapper Plus (ETM+) and Landsat 5 Thematic Mapper (TM) archive [44].Among 16 available land cover classes in the NLCD, three dominant forest classes (deciduous, evergreen, and mixed forests) were the interests of our study.

Methods
We tested the slope correction method across different study sites and for two different lidar footprint sizes (i.e., LVIS and GLAS).Additionally, we quantified terrain slope and footprint size effects on canopy height retrieval.Our study was performed in two sequential steps (Figure 2).
The first step made comparisons of LVIS estimates (with/without slope corrections) against field-measured maximum forest canopy heights.Building upon the previous study of Lee et al. [29], we added more sample sites to explore the feasibility of the physically-based approach over different forest types and terrain conditions.In the second step of our research, the slope correction method was applied to the estimation of GLAS maximum forest canopy heights.It is important to note that the LVIS forest canopy heights were the only appropriate estimates for the evaluation of GLAS metrics.The reasons are three-fold; (1) a subplot of field measurement was too small (25 × 25 m or 33.3 × 33.3 m) to reflect all of the trees within a GLAS footprint (~50 m diameter).Even if several subplots were merged, its spatial extent did not matched a GLAS footprint); (2) statistically meaningful evaluation using field measurements was impossible since so few subplots in our field data were spatially intersected with GLAS footprints; and (3) a relatively large number of LVIS footprints overlapped each GLAS footprint and LVIS has been shown to match well with field-measured maximum forest canopy heights.Lee et al. [29] and Choi et al. [40] have also validated the GLAS metrics using LVIS estimates.From the above steps, we quantified slope related errors and analyzed how terrain slope and footprint size influence lidar canopy height retrieval.

GLAS Preprocessing
GLAS data quality is influenced by cloud contamination and atmospheric saturation.Cloud-and saturation-free GLAS waveform data were selected in this study based on three filters; (1) a minimum signal to noise ratio of 15; (2) an elevation difference between the GLAS and NED data less than 50 m; and (3) a signal saturation level of 0. The second filter is useful to remove GLAS data distorted by low-cloud contamination.An additional filter was derived from the NLCD data to identify GLAS footprints over deciduous, evergreen, and mixed forests.

Maximum Forest Canopy Height Retrieval from LVIS and GLAS
To retrieve maximum forest canopy height from LVIS and GLAS, this study implemented a direct height calculation method based on decomposed waveform metrics.This direct method relies on the accurate identification of elevations for both canopy top and canopy bottom of the tallest plant within a LVIS or GLAS footprint.Usually, signal start and ground peak elevations of a waveform are used to approximate canopy top and canopy bottom, respectively [34,35].Signal start elevation is defined as the first elevation at which the waveform energy exceeds a threshold.The threshold was set as 4.5 times the background noise level in GLAS [45] and three times the background noise in LVIS [38].Ground peak elevation is determined by identifying the last Gaussian peak from decomposed multiple Gaussian distribution curves.With direct methods, maximum canopy height over flat terrain can be expressed by relative height from ground surface to signal start elevation where 100% of cumulative lidar backscattered energy is recorded (a.k.a., RH100 waveform metric).
RH100 of waveform can be distorted by terrain slope and off-nadir pointing angle [32,33].However, we implemented the simpler slope correction approach of Lee et al. [29] with a nadir viewing assumption due to negligible off-nadir pointing angle of both LVIS and GLAS data used in this study.We expected that the higher slopes with larger footprint sizes would lead more interference on the maximum forest canopy height estimation when compared to the lower slopes with smaller footprint sizes.The terrain slope effects on RH100 with nadir-viewing lidar can be expressed by the following equation (see Figure 3): where, z s and z g are the elevation of signal start and ground peaks, respectively.H max is the actual maximum forest canopy height, θ represents the slope, and d refers to the diameter of the waveform footprint.From the above equation, H max can be retrieved from RH100, with given footprint size, d, and slope (θ), such that: ( We applied this equation to LVIS and GLAS waveform data in order to correct H max estimates.This study also provided uncorrected H max values for LVIS and GLAS data to show the level of slope effects and accuracy improvements using the slope correction approach.

Comparison of LVIS Metrics to Field Measurements
We explored the physically-based slope correction approach for LVIS forest canopy height estimation across different sites.LVIS height metrics (with/without slope corrections) were compared against the field-measured heights from six study sites.For the plot-level comparisons, valid LVIS footprints (i.e., co-located with field plots) should overlap with a subplot.This analysis selected LVIS data whose center points fell inside each subplot.On average, seven LVIS footprints were registered over each subplot.The maximum RH100 of valid LVIS footprints was calculated and defined as the uncorrected LVIS forest canopy heights (LVIS RH_UC ) of each subplot, whereas the corrected LVIS estimate (LVIS RH_C ) was derived from Equation (2).Note that d in Equation (2) refers to the diameter of a LVIS footprint (~20 m).We evaluated the agreement between LVIS and field estimates using four statistical metrics: bias, mean-absolute-errors (MAE), RMSE and coefficient of determination (R 2 ).Based on these comparisons, we attempted to quantify the slope effect on LVIS height retrieval by relating a gradient of slope to the difference between LVIS RH_UC and field measured maximum height.

GLAS Maximum Height Estimation and Evaluation
We calculated RH100 for GLAS data using Equation (1).In the GLA14 product, two standard altimetry variables represent signal begin and ground peak elevations: (a) signal begin range increment, SigBegOff and (b) centroid range increment for the last Gaussian Peak, gpCntRngOff 1. Theoretically, gpCntRngOff 1 are assumed to represent the ground level elevation within a GLAS field-of-view, while SigBegOff refers to the highest point of a surface.In practice, (SigBegOff − gpCntRngOff 1) can be converted into RH100.From RH100, we obtained both uncorrected (GLAS RH_UC ) and corrected GLAS (GLAS RH_C ) maximum canopy height estimates from Equation (2).Note that here d in Equation (2) refers to the diameter of a GLAS footprint (~50 m).We then performed footprint-level comparisons over seven study sites based on corrected LVIS maximum canopy heights spatially corresponding to each GLAS footprint.The three tallest LVIS estimates (slope corrected) for each GLAS footprint were then averaged and compared to the GLAS metrics.We also calculated bias, MAE, RMSE, and R 2 to evaluate this analysis.

Comparison between LVIS and Field Measurements
Among 1010 available subplots, 705 valid subplots (i.e., co-located with LVIS shots) were used for the comparison between field measurements and LVIS metrics (LVIS RH_UC and LVIS RH_C ).Table 3 lists statistics for valid subplots and field-measured maximum forest canopy heights.The maximum height of most field measurements (664 subplots) ranged from 5-45 m, but there were extremes in the height distribution (45-80 m) for 41 subplots, which were mostly from mature stands at SN.The quality of LVIS height estimates and the impact of slope on LVIS heights were evaluated by comparing corrected and uncorrected LVIS heights with field measured maximum heights.Both LVIS RH_UC and LVIS RH_C showed moderate linear relationships with the field-measured heights as depicted in Figure 4a-b.Overall, both LVIS RH_UC (bias = 2.42 m) and LVIS RH_C (bias = 1.53 m) were slightly overestimating maximum heights when compared to the field measurements.However, we found that the simple slope correction could slightly improve the accuracy of LVIS maximum forest canopy height estimation (MAE = 4.08, RMSE = 5.39, and R 2 = 0.778 for LVIS RH_UC ; MAE = 3.72, RMSE = 4.99, and R 2 = 0.782 for LVIS RH_C ).According to Table 4, the application of this physical correction across different field locations improves estimates of maximum canopy height, which was previously shown for a single study site [29].PE was the only site that did not show improvement.Among these sites, HB and SN show a larger improvement (RMSE of 1.11 m and 0.71 m improved for HB and SN, respectively) than other sites because both contain large regions of high-sloped terrain (average slope = 9.77° ± 5.42° for HB and 7.77° ± 4.59° for SN).In the case of PE, the slope correction degraded the height retrieval accuracy (RMSE of 0.23 m increased).A plausible reason for this degradation is that 183 subplots at PE are located in relatively flat terrain (overall, 2.44° ± 1.46°) compared to the other sites, which leads us to believe that non-slope related errors were more dominant in this study site.We also observed relatively larger RMSEs in BF, HB, and HO sites than in HF, PE, and SN sites.This might be explained by larger geolocation uncertainties of field campaigns (discussed in Section 4.4).Figure 4c represents the error distribution as a function of slope gradients for the LVIS height retrieval.Filled and unfilled markers refer to LVIS RH_UC and LVIS RH_C , respectively.Although LVIS footprint size (~20 m) is much smaller than GLAS' (~50m), LVIS height could also be affected by terrain slope.In this figure, height difference between uncorrected LVIS estimates and field measurements increases with slope indicating a clear influence of slope on LVIS RH100.Disparities between field-measured and LVIS RH_UC heights were significantly dependent on terrain slope (R 2 = 0.06 (p < 0.001); slope of fitted line = 0.25).On the other hand, we obtained less slope-dependent errors (R 2 = 0.01 (p = 0.064); slope of fitted line = 0.069) after the physical slope correction process.A F-test confirmed that the coefficient of the fitted line was significantly changed by the slope correction scheme (F = 11.88 (p < 0.001)).

Comparison between GLAS and LVIS
In the case of GLAS height retrieval, we identified 527 valid GLAS footprints that are spatially co-located with at least three LVIS estimates.Selected GLAS metric (GLAS RH_C ) ranged from 2.98-65.04m (mean = 32.96m; std.dev.=12.12 m).Both GLAS metrics (GLAS RH_UC and GLAS RH_C ) showed a tendency to overestimate when compared to the reference LVIS maximum canopy heights (Figure 5a,b).Biases for GLAS RH_UC and GLAS RH_C were 7.22 m and 1.98 m, respectively (Table 5).These overestimations of GLAS forest canopy heights also occurred in the previous studies [28][29][30]35,40].Similarly to the smaller lidar footprint (i.e., LVIS), the physical slope correction approach produced significantly better RMSE and R 2 for the GLAS height estimate: RMSE = 12.74 m and R 2 = 0.54 for GLAS RH_UC ; RMSE = 7.83 and R 2 = 0.63 for GLAS RH_C .The use of slope correction could explain an additional 9.1% of the variability in the reference LVIS height and with a better RMSE.This improvement of RMSE (4.91 m) and R 2 (0.091) for GLAS RH_C was much greater than of LVIS RH_C (RMSE = 0.39 m and R 2 = 0.004 improved).If outliers are removed based on Cook's distance threshold (>3 times of mean Cook's distance) [46], we obtain a RMSE = 5.32 m and R 2 = 0.80 for GLAS RH_C (n = 489) without 38 outliers.
The distribution of errors by terrain slope gradients generated similar results when compared to the analysis for field-measured and LVIS heights.As shown in Figure 5c, errors in the uncorrected GLAS metric (i.e., GLAS RH_UC ) were largely dependent on the topographic condition (R 2 = 0.50 (p < 0.001) and slope of fitted line = 0.71 for GLAS RH_UC ).Not surprisingly, the simple physical correction satisfactorily improved the bias and global mean of errors in the GLAS metric: changes in (a) bias of −5.24 m; (b) MAE of −3.27 m; and (c) RMSE of −4.91 m.The fitted line of errors also showed that GLAS RH_C is more independent from terrain slope conditions over the footprints (R 2 = 0.06 (p < 0.001) and slope of fitted line = 0.17).A F-test confirmed that the slope correction process significantly changed the coefficient of the fitted line (F = 153.47(p < 0.001)).

Physical Slope Correction and Lidar Footprint Size
The sloped terrain generally lengthens the full extent of lidar waveform and decreases the level of laser energy at the forest canopy and ground peaks [32].Figure 6a simulates the relationship between waveform lidar altimetry and theoretical effects of slope (0°-70°) and footprint size (10-100 m) using the equation (=0.5 × footprint size × tan(slope)) of Lee et al. [29].This figure clearly represents the increased slope effect on lidar waveform as a function of terrain slope and footprint size.As expected, the larger footprint size and greater slope tend to generate more errors in the retrieved lidar forest canopy heights [29,30,32,33,40].In actual slope effects observed from our study, systematic errors in the uncorrected LVIS (LVIS RH_UC ) and GLAS (GLAS RH_UC ) metrics agreed well with the simulated slope effects (Figure 6b,c), although some variations occur within bins (range = 5°).Figures 4c and 6b supported the use of this simple slope correction for the smaller lidar footprint (LVIS) and Figures 5c  and 6c displayed the effectiveness of slope correction for the larger lidar footprint (GLAS).From the physically-based slope correction approach, we achieved RMSE improvements of 0.34 m and 4.91 m from LVIS and GLAS height retrieval, respectively.In the case of GLAS, the use of slope correction explained an additional 9.1% of the variability in the reference LVIS height, while explanatory power of LVIS based height estimation showed minimal improvement (0.4%).This discrepancy in the magnitude of accuracy improvement of both lidar datasets can be explained by the magnitude of the slope effect, which is a function of lidar footprint size and terrain slope.Furthermore, this explanation supports stronger slope-dependency (R 2 = 0.50 (p < 0.001); slope of fitted line = 0.71) of GLAS than that of LVIS (R 2 = 0.06 (p < 0.001); slope of fitted line = 0.25).Despite the efficacy of the physical slope correction approach, errors remain in LVIS and GLAS height estimation.Possible additional error sources and limitations are discussed in the following section.

Limitations and Further Plans
With available in-situ, airborne and spaceborne lidar datasets, this study attempted to test the proposed physical slope correction scheme across different footprint sizes and study sites.We found that both smaller-and larger-footprint campaigns could benefit from the physical correction approach.However, other error sources beyond the scope of our study might influence the results.First, geolocation errors may induce larger uncertainties in both LVIS and GLAS practices.The horizontal geolocation error for GLAS L3a was 2.4 ± 7.3 m (mean ± 1 standard deviation) [26] and that of LVIS was less than 2.0 m [47].Additionally, average geolocation errors of in-situ datasets were 1.7 m, 5.0 m, 6.0 m and 6.4 m for PE, HO, BF, and HB sites, respectively (available only at these sites) [37].From a simple comparison shown in Figure 7, geolocation errors are closely related to the overall accuracy of LVIS height retrieval.Although we could not account for the uncertainty contribution of this error source, it is important to evaluate the geolocational accuracy for remotely estimated forest canopy height from LVIS and GLAS altimetry datasets.
Second, forest condition, which is closely related to the shape of lidar waveform [33], was not explicitly considered in this study.Our study was performed with the assumption that similar vegetation structure is evenly distributed within each LVIS (~20 m) and GLAS (~50 m) footprint.However, Hyde et al. [48] reported that unevenly distributed canopy structure might induce larger uncertainty when calculating canopy height using waveform lidar.For instance, if the taller trees are at the edge of a lidar footprint, they may not be detected because of the low laser intensity at the edge.This effect will cause errors in the height estimation, precluding comparisons between remotely sensed and field measurements (in addition to the geo-location problems mentioned above).This study only used the simple slope correction scheme due to the negligible off-nadir angles and the lack of accurate information for the sensor zenith/azimuth angle and slope orientation [29].However, to fully quantify lidar height distortion on slope, topography information (slope and aspect) and sensor zenith and azimuth angle are important [32,33].For instance, we assumed that LVIS and GLAS data are acquired under the nadir-viewing condition, but in practice, both instruments might have a slight off-nadir angle.According to Yang's model [32] and ranges of the off-nadir angle for LVIS (3°-5°) and GLAS (0°-2°) from previous studies [29,34], estimated height errors may range from 0.1-1.0m for a 2-degree off-nadir angle in GLAS and from 0.2-1.5 m for a 5-degree off-nadir pointing angle in LVIS.For our case, the off-nadir angle effect was less than previously reported, but this error does magnify the uncertainties to some degree.
Additionally, the asymmetric shape and size of individual footprints may introduce uncertainty in lidar height retrieval, especially in GLAS.In our study, we assumed that every GLAS shot has a 50-m circular footprint.However, the eccentricity (= 1 − b 2 /a 2 , a = major axis and b = minor axis) and the major axis of the GLAS footprint vary from 0.48-0.63 and 51.20-55.41m respectively, according to the summary of Laser Profile Array (LPA) [41].These disparities may induce slight geo-locational disagreement between the laser footprint and actual forest structure and it might propagate into overall uncertainty.Actual footprint sizes of LVIS and GLAS also vary by flight altitude and sensing elevation, which may contribute to uncertainty.
With respect to the use of this physically-based slope correction method in a global scale study, accurate DEM information plays a critical role in correcting slope effect on GLAS waveform data.In this study, we used the NED dataset over CONUS that is derived from a photogrammetry-based approach.NED data are probably more accurate than other possible global DEM data, especially the Shuttle Radar Topography Mission (SRTM) dataset.In fact, several previous studies exploited the difference between NED and SRTM to retrieve forest canopy height since SRTM gives the height of the scattering center of the forest [49,50].It is unclear whether the physically-based slope correction method using SRTM would get results as good as using NED.This will need to be determined before using our method in global scale studies due to the limited spatial coverage of photogrammetric DEMs.
Important future studies will require: understanding of forest conditions, considering sensor and topography conditions, and asymmetric footprint shape and size.With an accurate construction of the physically-based slope correction scheme, we can test the effectiveness of using SRTM as a precursor for a global scale study.Also, this study focused only on site-specific and plot-footprint level tests of waveform lidar, and did not attempt to investigate large-scale tests and mapping forest canopy heights.Spatially and temporally continuous mapping is critically important in forest monitoring, ecological modeling, and long-term forest management.Therefore, our subsequent study will pursue the generation of forest canopy height maps using appropriate relationships between geo-predictors (e.g., temperature and precipitation) and forest vertical structures [23][24][25]51] since both airborne and spaceborne lidar data are still not perfectly continuous in space and time.

Conclusions
We tested a physical slope correction approach for waveform lidar data to generate accurate maximum forest canopy heights.Building on the previous research [29] constrained on a single homogeneous forest stand, this study was performed over six additional study sites with various terrain conditions and forest types.Two types of lidar altimetry were derived from airborne (Land, Vegetation, and Ice Sensor; LVIS) and spaceborne (Geoscience Laser Altimeter System; GLAS) sensors.The LVIS height metrics (with/without slope corrections) showed a moderate linear relationship with field-measured maximum forest canopy heights (n = 705; R 2 = 0.778 and RMSE = 5.39 m for the uncorrected; R 2 = 0.782 and RMSE = 4.99 m for the slope-corrected LVIS).Like the results of Lee et al. [29], the simple physical approach provided slight benefits on the accuracy of LVIS heights whose footprints are relatively small (size = ~20 m).Apart from random data errors due to the spatial mismatches between LVIS footprints and field subplots, we obtained improved RMSE (changed by -0.34 m) using the physical slope corrections.Distribution of errors as a function of slope gradients showed that disparities between field-measured and uncorrected LVIS heights were likely dependent on terrain slope effects (R 2 = 0.06; slope of fitted line = 0.25), but this dependency was successfully removed by the slope correction approach (R 2 = 0.01 (p = 0.064); slope of fitted line = 0.069).Similarly to LVIS, the physical slope correction performed well with the larger-footprint lidar (i.e., GLAS with ~50 m diameter).The use of the slope correction method on GLAS altimetry produced better measures of maximum forest canopy heights (n = 527; R 2 = 0.54 and RMSE = 12.74 m for the uncorrected; R 2 = 0.63 and RMSE = 7.83 m for the slope-corrected GLAS).After excluding 38 outliers (n = 489), we could explain 80% of variability in the reference LVIS heights with RMSE = 5.32 m.Errors in the uncorrected GLAS metrics were largely dependent on the terrain slope (R 2 = 0.50 (p < 0.001) and slope of fitted line = 0.71).Despite the correction, GLAS metrics produced slight overestimates (R 2 = 0.06 (p < 0.001 and slope of fitted line = 0.17) as in previous studies [22,25,33].Systematic errors in the uncorrected LVIS and GLAS metrics agreed well with the simulated theoretical slope effects as a function of slope and footprint size.Our results showed that both smaller and larger footprint lidar estimates of maximum canopy height significantly improved the simple physical slope correction.

Figure 1 .
Figure 1.Geographical locations of field measurements and LVIS dataset.

Figure 2 .
Figure 2. Overall scheme of physically-based slope corrections for LVIS and GLAS maximum forest canopy height retrievals.

Figure 3 .
Figure 3. Schematic definitions of actual maximum forest canopy height and lidar-derived maximum canopy height for a waveform over the flat (solid blue line) and slopped (dash red line) terrains.Lower left and right panels demonstrate impact of terrain slope on lidar height retrieval.H max is actual maximum canopy height and H lidar (=RH100) is lidar-derived maximum canopy height.d is lidar footprint diameter, and θ is the terrain slope.

Figure 4 .
Figure 4. Comparisons of (a) uncorrected (LVIS RH_UC ) and (b) corrected (LVIS RH_C ) LVIS metrics against field-measured maximum forest canopy heights.Colors correspond to each field site; (c) Distribution of LVIS height errors by slope gradients.Filled markers with red fitted line represent errors for LVIS RH_UC , while unfilled markers with blue fitted line correspond to errors for LVIS RH_C .

Figure 5 .
Figure 5. Comparisons of (a) uncorrected (GLAS RH_UC ); and (b) corrected (GLAS RH_C ) GLAS metrics against LVIS derived maximum forest canopy heights.Outliers are determined by Cook's distance threshold (within three times the mean Cook's distance); (c) Distribution of GLAS height errors by slope gradients.Filled markers with red fitted line represent errors for GLAS RH_UC , while unfilled markers with blue fitted line correspond to errors for GLAS RH_C .

Figure 6 .
Figure 6.(a) Simulated theoretical terrain slope effects across different footprint sizes (10-100 m) and slope gradients (0°-70°).Errors in slope-uncorrected (b) LVIS and (c) GLAS metrics used in this study along with the simulated slope effects.Bar histograms represent the frequency (%) of LVIS and GLAS data given each slope range (numbers in plot represent number of lidar shots in each bin (not %)).Circles are the simulated slope effects of given footprint size and terrain slope (refer to Figure 6a) and box-whiskers plots correspond to errors in the LVIS and GLAS estimates.

Figure 7 .
Figure 7. Impact of geolocation errors on overall accuracy of LVIS height retrieval.This comparison includes Bartlett Experimental Forest (BF), Howland Research Forest (HO), Hubbard Brook Experimental Forest (HB), and Penobscot Experimental Forest (PE) sites only (where geolocation accuracy assessments were performed).

Field Measurement Site LVIS SiteTable 1 .
Datasets for inter-comparisons between field measured and LVIS waveform derived heights. a

Table 3 .
Statistics of valid subplots and field-measured heights corresponding to LVIS footprints.
a BF: Bartlett Experimental Forest; HF: Harvard Forest; HO: Howland Research Forest; HB: Hubbard Brook Experimental Forest; PE: the Penobscot Experimental Forest; SN: Sierra National Forest; b Min.: Minimum of maximum forest canopy heights for each field measurement sites; c Max.: Maximum of maximum forest canopy heights; d Ave.: Averages; e Std.: Standard deviations; * Total number of valid subplots (co-located with LVIS footprint) for six field measurement sites.
c RMSE: Root Mean Squared Errors; d R 2 : Coefficient of Determination.
a Without outliers: statistics without outliers based on Cook's distance.Number of outliers: GLASRH_UC and GLASRH_C = 39 and 38, respectively; b MAE: Mean Absolute Errors; c RMSE: Root Mean Squared Errors; d R 2 : Coefficient of Determination.