Edinburgh Research Explorer A systematic evaluation of multi-resolution ICESat-2 ATL08 terrain and canopy heights in boreal forests

The launch of NASA's Ice, Cloud, And Elevation Satellite-2 (ICESat-2) in September 2018 provides the scientific community an opportunity to observe high-resolution and three-dimensional surface elevations with global coverage. ICESat-2's Land and Vegetation Height (ATL08) data product focuses on the along-track terrain and canopy heights observations at a 100 m × 11 m spatial resolution. This work expands on past ATL08 validation studies to assess a higher spatial resolution (30 m × 11 m) version of ATL08's height product. This new dataset enables higher resolution mapping and fusion with Landsat data, but has not previously been validated across large geographic extents. In this paper, we examine the accuracy of multi-resolution ICESat-2 ATL08 across North America boreal forests using Land, Vegetation, and Ice Sensor (LVIS), an airborne laser ranging system as reference datasets. Overall, strong agreements of terrain elevation and canopy height were found between ATL08 and LVIS at both 100 m × 11 m (RMSE terrain = 2.35 m; bias terrain = (cid:0) 0.17 m; RMSE canopy = 4.17 m; bias canopy = 0.08 m) and 30 m × 11 m (RMSE terrain = 3.19 m; bias terrain = 0.49; RMSE canopy = 4.75 m; bias canopy = 0.88 m) spatial resolutions. We found the accuracy of high-resolution terrain and canopy height measurements were constrained by sensor and external conditions during the time of data acquisition with lower uncertainties observed from samples along high-intensity ground tracks and with low topography/slope variabilities. Through this work, we provide insight into the use of multi-resolution ICESat-2 ATL08 for terrain and canopy heights characterization in northern forests. The results found in our study serve as a benchmark for end users to select high-quality ATL08 for a variety of scientific applications.


Introduction
Boreal forests represent one of the largest carbon reservoirs in the terrestrial ecosystem, encompassing ~30% of the global forest area (Brandt et al., 2013), containing approximately 40% of the terrestrial carbon (Alcaraz-Segura et al., 2010;Gorham, 1991;Kasischke, 2000), and approximately 15% of global aboveground biomass (Santoro et al., 2021). Approximately 30% of the boreal forests are found in North America with a total area encompassing over 6 million km 2 (Wells et al., 2020). An accurate characterization of the extent and variability of North American boreal forests is critical to understand current and future global carbon cycle dynamics under the rapidly changing climate. Among all the forest attributes, canopy height has been identified as an important geophysical parameter with strong association with the aboveground biomass (Drake et al., 2002;Lefsky et al., 2002;Neuenschwander et al., 2020). Therefore, spatially-explicit observations of forest height structure not only provide scientists a lens to understand the spatial and temporal distribution of forest carbon stock but also help policy makers to effectively propose strategies for carbon emission mitigation.
Over the past several decades, there has been significant development in the use of laser altimeters, installed on ground- (Anderson et al., 2018;Hancock et al., 2014) air- (Duncanson and Dubayah, 2018;Duncanson et al., 2014;Zhao et al., 2018) and spaced-based Neuenschwander et al., 2020) platforms for multi-scale forest structure parameter estimation and mapping. The laser pulse emitted from a laser altimeter is capable of penetrating the dense canopy surface with the forest vertical profile being recorded from the returned laser energy. In addition to its promising capacity on vegetation structure observations, lidar sensors are less constrained by the solar angle than other optical remote sensing platforms (Patenaude et al., 2004;Zheng and Moskal, 2009).
In 2018, NASA launched two spaceborne laser altimeter missions: ICESat-2 and GEDI (Global Ecosystem Dynamics Investigation; Markus et al., 2017). Since GEDI is onboard the International Space Station with the spatial bounds between 51.7 • N and S, ICESat-2 is the only current spaceborne laser ranging system capable of capturing three-dimensional measurements globally. While the primary goal of ICESat-2 is to measure changes in the cryosphere (Neuenschwander and Pitts, 2019;Smith et al., 2020), one of the science objectives defined by Markus et al. (2017), is to measure the vegetation canopy height as a basis for large-scale biomass and biomass change estimation. Particularly, the unique data-collecting strategy adopted by ICESat-2 mission with intensive high-latitude samplings being recorded (Montesano et al., 2015;Patterson et al., 2019) also makes it a suitable platform for vegetation structure and biomass monitoring of boreal forests. In contrast to the full waveform data from GEDI and the former ICESat/GLAS mission, ICESat-2 carries a photon-counting laser ranging system named ATLAS (Advanced Topographic Laser Altimeter Systems) with detection sensitivities at a single photon level (Neuenschwander and Pitts, 2019;Neumann et al., 2019). The reduced laser power requirement, in addition to the high repetition rate during the normal operation, allows ATLAS to provide consistent and high-resolution observations of the surface in the along-track direction (Neuenschwander and Pitts, 2019).
Along with field inventory data, high-precision terrain and canopy height measurements from airborne lidar systems (ALS) serve as a gold standard for validating spaceborne lidar observations Duncanson et al., 2021). Since the pre-launch stage, the performance of photon-counting lidar for terrain and forest height estimation has been evaluated using both simulated and real ICESat-2 observational data Montesano et al., 2015;Mulverhill et al., 2022;Neuenschwander et al., 2020;Neuenschwander and Magruder, 2016). Neuenschwander and Pitts (2019) integrated the ATL08 algorithm and DRAGANN (Differential, Regressive, and Gaussian Adaptive Nearest Neighbor), a signal finding method, on simulated ATLAS datasets to test its potential on terrain and canopy height estimation from two different ecosystems, located in Sonoma County, CA and Alaska Tundra/Taiga respectively. Their results demonstrated a high accuracy (RMSE<1 m) of terrain height estimation in Alaska Tundra/Taiga ecotone which is mainly occupied by sparse vegetation. For Sonoma county, characterized by high topographic relief and a wide range of height and canopy covers, errors were slightly higher but still relatively low with <2 m RMSE. However, the canopy height estimated from the simulated ATLAS underestimated the true canopy height in both study sites, with mean Bias and RMSE of − 4.12 m/2.65 m (Alaska Tundra/Taiga) and − 3.96 m/3.19 m (Sonoma County, CA) respectively, primarily due to the vertical sampling error on heterogeneous canopy surface (Neuenschwander and Magruder, 2016). Uncertainties in Sonoma County also varied as a function of canopy structure, with higher uncertainties for tall trees and dense canopies . Gwenzi et al. (2016) generated the pseudo-ATLAS data from Multiple Altimeter Beam Experimental Lidar (MABEL), an airborne photon-counting lidar sensor to pre-validate the potential of ICESat-2 mission on canopy height estimation in savanna ecosystems. The results illustrated a modest correlation with RMSE ranging from 2.9 to 4.4 m with respect to segment size. After its launch, on-orbit ICESat-2 observations have also been tested against reference datasets across multiple ecoregions (Guerra-Hernández et al., 2022;Liu et al., 2021;Neuenschwander et al., 2020;Neuenschwander and Pitts, 2019). For example, the ATL08 canopy height accuracy was evaluated in the region of Extremadura across Central-West Spain with ALS as reference datasets (Guerra-Hernández et al., 2022). The absolute RMSE and Bias were found to vary across ecosystem types with range from 0.95 m to 2.24 m, and − 0.56 m to − 0.14 m respectively. ATL08 terrain and canopy height estimates were also assessed in boreal forests across southern Finland . The mean absolute errors of terrain height estimates were found to range from 0.49 to 0.59 m, while the mean bias of canopy height measurements was estimated to be 3.05 m. Although strong agreements were found among those studies on the promising capability of ICESat-2 for terrain and canopy height characterizations, the accuracy is observed to be limited by multiple conditions (e.g., day/night, beam intensity, snow coverage, canopy cover, etc.), explained by the number of photons captured by ATLAS from certain land surface types and solar background noise during daytime operations. Therefore, post-stratified assessment of ATL08 becomes a pressing need for end users to understand the potential biases of the datasets and properly select the high-quality ATL08 observations under certain conditions.
Although the standard ATL08 data product is reported at a fixed step size of 100-m, terrain and canopy labeled photons were remapped to a 30 m along-track to estimate higher resolution terrain and canopy heights in support of a boreal biomass mapping project (Duncanson et al., 2021). This higher spatial resolution version of ATL08 provides a direct training of Aboveground Biomass Density (AGBD) models with available field plot data in boreal systems (which are typically circular or square ~ 30 m plots) following the protocols used with the GEDI (Duncanson et al., 2022). This higher-resolution product also allows for a more direct comparison of airborne full-waveform lidar (e.g., LVIS) as reference datasets for the validation of the entire canopy height profile. In addition, the downscaling of ATL08 to the resolution equivalent to Landsat imagery enables the integration of multi-sourced remotely sensed observations for large-scale and wall-to-wall mapping of forest structure attributes. However, the downscaling of ATL08 could potentially cause increased uncertainties on ground and canopy height estimation due to reduced number of photons within each segment and comprehensive accuracy assessments of high-resolution ATL08, thereby could serve as a benchmark for the appropriate selection of high-resolution ATL08 datasets. A primary objective of this study was to validate this higher spatial resolution product which has had much less attention than the official 100 m × 11 m resolution product.
Although prior studies have examined the ICESat-2 ATL08 data product in several sites across the globe, our current knowledge of ICESat-2 performance in terms of its potential for terrain and canopy heights estimation is still quite limited, mainly due to 1) the lack of highquality spatially and temporally coincident reference data, and 2) the limited spatial coverage of airborne reference datasets across large topographic and vegetation structural gradients. Further, the higher spatial resolution product has yet to be validated. In order to address the above-mentioned issues, this study aims to examine the accuracy of ICESat-2 ATL08, measured at different spatial resolutions and across the entire North American boreal forest using all available LVIS observations acquired in the 2019 airborne campaign as reference datasets. The specific objectives of this work include: 1) examine the accuracy of ATL08 terrain and canopy heights at fixed step size of 30 and 100 m respectively, 2) assess the impacts of sensor specifications (e.g., beam intensity, snow coverage, day/night etc.) on ATL08 data quality in North America, and 3) explore accuracies across the topography/slope, canopy height and land-cover gradients.

Study area
Our study is in the North American boreal forest, and primarily in the NASA Arctic Boreal Vulnerability Experiment (ABoVE) domain (Miller et al., 2019). As part of the circumpolar boreal zone, the boreal forest biome is occupied by a large proportion of intact, primary forests with a total area over 6 million km 2 in North America (Wells et al., 2020). The forested lands are comprised of deciduous and coniferous species, dominated by white spruce (Picea glauca), black spruce (Picea mariana), balsam fir (Abies balsamea) and trembling aspen (Populus tremuloides) (Matasci et al., 2018). The main disturbance agent found in the area is fire (Kurz and Apps, 1999) which plays a critical role in vegetation succession, energy regulation and climate change. Fig. 1 demonstrates the geographic extent, Digital Elevation Model (DEM) and land cover types of our study region.

Airborne lidar data
LVIS Facility data were acquired between July-August 2019 with a total of ~78,000 km flown, as shown in Fig. 1 (a). LVIS Facility is a fullwaveform scanning laser altimeter capable of measuring the terrain height, vegetation height and sub-canopy structure (Hofton et al., 2000b) within the footprint diameter of ~10 m. The laser altimeter was onboard the NASA Gulfstream V plane with flight altitude of ~12.5 km and scanning angle of ±6 • around nadir (Hofton et al., 2000a). For our comparison to ICESat-2 data, we used the LVIS Facility level 2 products available from National Snow & Ice Data Center (NSIDC, https://nsidc. org/data/LVISF2/versions/1) with the footprint-level terrain elevation and vegetation RH metrics being incorporated. Specifically, the ZG and RH98 from level 2 products were used to represent the ground and canopy heights within each LVIS facility footprint.

ICESat-2 ATL08 product
The ATLAS instrument onboard the ICESat-2 is a photon-counting laser altimeter (Neuenschwander and Pitts, 2019). Although ATLAS also relies on the same laser ranging concept with the travel time of emitted laser pulse in tandem with the satellite position and pointing information for the geodetic measurements over the earth surface, the unique photon-counting system featured with low laser energy requirement, high repetition rate and single-photon level detection sensitivity makes great advantage for ICESat-2 in terms of its consistent and high-resolution observations in the along-track direction. Operated in the green wavelength (532 nm) of the electromagnetic spectrum with laser repetition rate of 10 kHz, each transmitted laser pulse from ATLAS is split into six individual beams parallel to the ascending and descending orbital path (92 ̊ inclination), arranged in three pairs with an energy ratio of 1:4 (Markus et al., 2017). Each emitted laser pulse can illuminate a surface area of ~10-11 m in diameter (Magruder et al., 2020), with adjacent footprints of 70 cm across the along-track direction.
The ATL08 product is generated from ATL03 geolocated photon datasets through a series of noise filtering and photon classification procedures, adopted at both ATL03 and ATL08 levels, to provide the along-track measurements of terrain and canopy heights (Malambo and Popescu, 2021;Neuenschwander and Pitts, 2019) at segment size of 100 m × 11 m. The key outputs of the ATL08 product include the terrain elevation, canopy height, and related RH metrics describing the structural profiles of the canopy across the vertical plane. All of the terrain and canopy parameters correspond to certain geolocated segments at fixed step size of 100 m. In addition, there are multiple flags linked to each segment describing the instantaneous conditions (e.g., day/night, beam intensity, snow/no-snow, etc.) of the data acquisition.
In addition to the official ATL08 Version 5 datasets, one of the primary objectives of this study is to validate a high-resolution ATL08 dataset having a 30 m × 11 m resolution across the North America boreal forests with spatially coincident LVIS Facility observations as reference datasets. The ATL08 data product at 30 m step size was created by mapping the photon labels by using the photon indices provided on the ATL08 data product back to the ATL03 geolocated photon data.
Contrary to the standard ATL08 product, the 30 m × 11 m product produced here includes all of the signal photons (i.e. canopy and ground) for the relative height (RH) metrics calculation, whereas only classified canopy photons are used to calculate RH metrics on the 100 m × 11 m ATL08 data product. Including the ground photons in the RH metric calculation allows more comparable RH metrics to those from waveform lidar (e.g. LVIS, GEDI) which calculate RH metrics as waveform percentiles, including the ground portion of waveforms. This study incorporated ATL08 datasets acquired from Oct, 2018 to Nov, 2020 (100 m × 11 m ATL08), and from Oct 2018 to Sep 2020 (30 m × 11 m), respectively. Fig. 2 illustrates the spatial patterns of ATL08 tracks in Central Alaska. The profiles of classified ATL03 photons (ground photon, canopy photon, top of canopy photon), along with the labeled ATL08 canopy top height at fixed step size of 100 m and 30 m are demonstrated in Fig. 3. For this project, we will only analyze the ATL08 measurements north of 51.6 • N.

Data processing
The primary data processing procedures involved in this work were carried out on the Multi-mission Algorithm and Analysis Platform (MAAP), a cloud-computing platform jointly developed and implemented by NASA and ESA. As a virtual open and collaborative environment, MAAP is equipped with high-performing computing capabilities, a set of open-source tools and algorithms, along with direct access to Earth science data acquired from Earth observing satellites, sub-orbital platforms and field campaigns (Albinet et al., 2019(Albinet et al., , 2020. The ATL08 and LVIS Facility datasets used in this study were directly retrieved from MAAP to the local workspace based on the given spatial bounds. The validation is conducted following five main steps: 1) ATL08 pre-processing using the basic set of filters, 2) spatial matching between the ATL08 and reference datasets, 3) exclusion of all the ATL08 segments with stand-replacing disturbance occurring between the date of acquisition for the ATL08 and reference datasets, 4) extraction of terrain and canopy metrics from ATL08 and LVIS, and 5) accuracy assessment. The methodological framework was demonstrated in Fig. 4.

ATL08 noise filtering
Prior to pairing the ATL08 segments and LVIS footprints given their central coordinates, we applied a basic set of filters to pre-process all the ATL08 granules to remove noisy estimates of vegetation height. Specifically, the multiple scattering warning flag (msw_flag) from ATL09 was utilized to exclude all the segments with cloud contamination during the time of acquisition. Once all the cloud-free observations were  extracted, the canopy height (h_canopy) and number of canopy photon (n_toc_ph + n_can_ph) parameters were analyzed and all the segments without top of canopy photons presented after ATL08 photon classification procedures were removed for canopy height validation. We also carried out a land-cover filtering process and removed all the highresolution ATL08 with canopy height above a given threshold associated with specific land-cover types.

Spatial matching and metrics extraction
The next step carried out was to spatially pair the LVIS Facility footprints with all available ATL08. After the quality-filtered ATL08 were extracted through the filter, we segmented all the datasets into 1̊ latitude by 1̊ longitude boxes with all the ATL08 segments falling within each latitudinal/longitudinal boundary grouped. The same operation was conducted on the LVIS Facility data based on the central coordinates of each footprint while the LVIS trajectories passing through multiple latitudinal/longitudinal boxes were further split into several individual tracks. This step is to improve processing efficiency of the geospatial operation given both ATL08 and LVIS Facility spanned across a broad geographical domain with huge data volumes.
Both datasets were then converted to the WGS84 system with UTM projection and each ATL08 segment was created based on the central coordinates, inclination angle and segment/footprint sizes. Finally, the LVIS Facility footprints with centers falling within the ATL08 segments were spatially joined (Fig. 5). The terrain height (ZG) of all the LVIS Facility footprints within each segment were averaged while we maximized the canopy height (RH98) from all the footprints to represent the true canopy height at a segment level.

Disturbance filtering
Although the temporal window between the acquisitions of reference datasets and ATL08 is relatively short (< 3 years), the land-use/landcover conditions over all the examined ATL08 granules were checked to exclude the potential biases introduced by natural and anthropogenic disturbances. To address this issue, we adopted the most up-to-date National Terrestrial Ecosystem Monitoring System (NTEMS) forest change map which includes pixel-level stand-replacing disturbances attributed by fire and harvest between 1985 and 2020 across Canada. All the ATL08 segments that spatially matched pixels with stand-replacing disturbances occurring between 2018 and 2020 were removed for further analysis. Due to the lack of up-to-date forest disturbance map in Alaska, the historical fire data between 2018 and 2021 were retrieved from the Alaska Large Fire Database and all ATL08 segments falling within the fire perimeters were excluded.   5. Spatial configuration of LVIS Facility footprints and ATL08 segments in spatial join operation.

Accuracy assessment
The accuracy of ATL08 was examined by comparing with the terrain and canopy height measurements from the reference datasets. The accuracy metrics we calculated for the validation include 1) Bias, 2) Root Mean Square Error (RMSE), 3) Relative Bias (Bias%), and 4) Relative Root Mean Square Error (RMSE%). Once the ΔH Terrain and ΔH Canopy per segment were obtained, we calculated the above mentioned statistic items through following equations: Where n referred to the number of ATL08 segments, H ATL08 represents the terrain and canopy height metrics from the ATL08 while H LVIS represents the terrain and canopy metrics extracted from the paired LVIS Facility footprints.
The accuracy of ATL08 was also checked across the environmental gradients which includes topography/slope, land-cover types and reference canopy height. To achieve this, we gridded the LVIS-Facility terrain (ZG) and canopy height (RH98) to 30 m spatial scale and calculated the terrain slope and reference canopy height within each ATL08 segments. The Copernicus Land Cover Datasets were also incorporated to assess the accuracy of ATL08 terrain and canopy heights across land-cover gradients. The terrain slope of all available ATL08 observations ranges from 0 ̊ to ~50 ̊ and we evenly split this into 10 ̊ intervals (e.g. 0 ̊ to 5 ̊; 5 ̊ to 10 ̊) . 19 intervals were evenly created across reference canopy height gradients (0 m to 38 m) as well to evaluate its impact on ATL08 canopy height accuracy. The intervals of each environmental variable were chosen given the criteria that 1) sufficient samples were included within each bin, and 2) the number of intervals were large enough to capture the impact of environmental variables on variability of the ATL08 data quality. Finally, the statistics describing the correspondence between ATL08 and reference datasets were calculated as a function of the environmental variables.

Terrain elevation accuracy
Through spatial matching and disturbance filtering processes, we derived 57,424 and 215,487 valid ATL08 samples at 100 m * 11 m and 30 m* 11 m spatial scales. Prior to an evaluation of ATL08 terrain height accuracy post-stratified by different external conditions, we summarized the proportion of standard ATL08 samples under combinations of pa- (Fig. 6). This procedure facilitates our understanding in how seasonality impacts the number of data acquisitions under certain conditions, and its interrelation to the reliability of post-stratified analysis. We found a majority of ATL08 samples acquired during winter are with presence of snow coverage (> 95%), among which ~75% of observations were collected during nighttime. In contrary, external conditions of ATL08 acquired during the growing season are opposite, with daytime and no-snow data accounting for the majority of ATL08 samples.
Strong agreements of terrain elevation between 100 m × 11 m ATL08 and LVIS were observed across the entire North America boreal forests, with overall RMSE and bias of 2.35 m and − 0.17 m, respectively. In comparison, the accuracy of terrain elevation from 30 m × 11 m ATL08 slightly decreased, with RMSE and bias of 3.19 m and − 0.49 m, respectively. The results above include all the ATL08 with LVIS footprints spatially overlapped regardless of ATLAS specifications during the time of acquisition. Fig. 7 (a, b) illustrates the histograms of terrain residuals post-stratified by beam intensity and we found the ATL08 acquired through strong beam consistently showed higher accuracy on terrain height estimation (RMSE 100m = 2.07 m; Bias 100m = − 0.13 m; RMSE 30m = 2.92 m; Bias 100m = − 0.27 m). Thereafter, we filtered out weak beam dataset and focused on accuracy of high-intensity ATL08 terrain height observations under sets of external parameters (Snow/Nosnow, Day/Night). Fig. 6. The ratio of monthly 100-m ATL08 at certain conditions to all available 100-m ATL08 acquired monthly. Fig. 7. Histograms of terrain height residuals at 100 m (7a; 7c; 7e) and 30 m (7b; 7d; 7f) step size stratified by beam Intensity (7a; 7b), day/light condition (7c; 7d) and snow/no-snow conditions (7e; 7f) with accuracy metrics (RMSE & Bias) displayed.
For the ATL08 at 100 m × 11 m spatial resolution, RMSE of terrain height for snow & strong beam conditions is 2.21 m (mean bias = 0.03 m), which is comparable to the accuracy of nighttime & strong beam observations (RMSE = 2.14 m, mean bias = 0.04 m). This can be attributed to the unique geographic characteristics of our study region where most nighttime observations in high-latitude North were acquired during snow season. Interestingly, the ATL08 acquired during nighttime showed higher uncertainties on terrain elevation estimation and an increase in absolute RMSE by 0.12 m. The high-resolution ATL08 data is also found to demonstrate similar response to sets of external conditions, with a lower RMSE and higher bias being observed during daytime (RMSE = 2.02 m; bias = 0.19 m). The long snow season in North American boreal, along with presence of deep snow on the ground is considered as the primary reason the daytime observations outperformed nighttime & snow data on terrain height characterizations.
In addition to the sensor and external conditions, the ATL08 terrain elevation estimates were further examined across slope gradients at both spatial scales to evaluate the potential uncertainties introduced by topographic/slope variability (Figs. 8 & 9). We summarized the number of valid 100 m ATL08 samples used for slope analysis, as shown in Table 1. Overall, larger uncertainties of terrain height estimates were found in regions with high topographic relief, as described by higher RMSE and substantial variabilities of terrain height residuals within slope intervals. Contrary, ATL08 terrain height measurements across flat terrain is more robust. The day/night condition and beam intensity displayed lower impact on terrain height accuracy across slope gradients in our study site, potentially due to 1) lack of nighttime samples in highlatitude North during the ATL08 acquisition window, and 2) long snow season with presence of snow cover dramatically increasing the surface reflectivity and sufficient number of ground photons can be captured by ICESat-2/ATLAS even along low-intensity tracks for terrain height characterizations.

Canopy height accuracy
The canopy height estimated from ATL08 at different spatial resolutions displayed a similar relationship to sensor and external parameters. In brevity, the strong beam demonstrated significantly lower RMSE at 100 m (RMSE = 3.85 m) and 30 m (RMSE = 4.22 m) step size ( Figure 10). Therefore, we filtered out all the weak beam datasets and focus on the high-intensity observations under sets of external conditions. We found the ATL08 acquired under snow condition had higher accuracy on canopy height characterization, at both 100 m (RMSE = 3.78 m; Bias = − 0.08 m) and 30 m (RMSE = 3.92 m; Bias = − 0.73 m) step size. This could be explained by an increase of canopy photons received by ATLAS through the snow season when the canopy surface is covered by high-reflectance snow. Compared with ATL08 acquired during daytime, the nighttime observations yielded lower absolute RMSE on canopy height measurements (RMSE 100 m = 3.74 m; RMSE 30 m = 3.80 m), which explains the significance of reduced background noise on high-quality canopy height measurements.
The ATL08 canopy height estimates was also examined across a slope gradient at both spatial scales to evaluate the potential errors introduced by topographic variabilities (Figs. 11 & 12). The increase of slope was found to introduce additional uncertainties on canopy height characterization at both spatial scales, with an RMSE of ~7.5 m (100 m × 11 m ATL08) and ~ 20 m (30 m × 11 m ATL08) in high-relief terrain (slope > 45 • ). Though the accuracy statistics of ATL08 canopy height estimates in high-relief terrain can be less representative due to limited sample size, a dramatic increase of uncertainties along slope gradient demonstrated strong impact of terrain variabilities on ATL08 canopy height measurement errors. We further post-stratified the uncertainties by sensor and external conditions and found the night & snow & strong beam observations displaying higher accuracy on canopy height estimates.
The accuracy of high-quality ATL08 canopy height data was further examined across reference height and land-cover gradients, as shown in Fig. 13. A relatively consistent RMSE% (~30%) and increasing absolute RMSE ranging from 2.5 m to 17.5 m was observed for ATL08 canopy height datasets with corresponding reference height over 10 m, indicating higher uncertainties of ATL08 canopy height measurements in forest stands occupied by tall vegetation. The bias is found to be strongly height-dependent with tall trees more prone to be underestimated (Bias <0 m), compared with a trend of overestimation (Bias >0 m) for those ATL08 segments occupied by short vegetations (< 10 m). The uncertainties of ATL08 on canopy height characterization were also evaluated across forested land-cover gradients. We found ATL08 underestimated the canopy height in both closed-canopy deciduous broadleaf and mixed forests, which are mainly dominated by tall vegetation with mean canopy height of 15.1 m and 14.9 m, respectively. In contrast, an overestimation (0.21 m < Bias <0.81 m) was observed in all other land-cover classes (e.g., closed-canopy evergreen needleleaf forest, open-canopy deciduous broadleaf forest, etc.) with mean canopy height ranging between 5.1 m and 10.8 m, which showed correspondence to our finding of the height-dependent ATL08 canopy height accuracy.
Finally, we analyzed the potential uncertainties from reference datasets for canopy height validation. Fig. 14 illustrated the accuracy of filtered ATL08 canopy height (Night & Strong Beam & snow) after corresponding low-quality reference datasets were excluded. Only ATL08 samples with at least 5 (100 m × 11 m ATL08) or 3 (30 m × 11 m ATL08) footprints contained were extracted for analysis. We observed significant improvements on canopy height agreements once all the ATL08 linked to low-quality reference datasets were filtered.

ICESat-2 for terrain elevation measurement
The ATL08 elevation datasets spanning across the North American boreal biome were examined to evaluate the potential of ICESat-2/ ATLAS on terrain height characterization. We found a strong agreement of mean terrain elevation derived from ATL08 classified ground photons and LVIS-Facility footprints at segment size of 100 m × 11 m and 30 m × 11 m, respectively. Notably, the downscaling of ATL08 to a higher spatial resolution resulted in a slight increase of uncertainty on terrain elevation estimates (RMSE All 30m ATL08 = 3.19 m), compared with original ATL08 dataset at fixed step size of 100 m (RMSE All 100m ATL08 = 2.35 m). This difference indicates the impact of scaling effect on characterizing the heterogeneity of terrain relief and the adequate representation of relief variability can serve as a critical basis for high-quality canopy height measurements.
The results found in this study enable us to draw some interesting and unexpected conclusions in terms of the performance of the ATL08 terrain product under different sensor and external conditions. We observed the accuracy of ATL08 terrain elevation is less constrained by day/night conditions at the time of data acquisition given the daytime observations can even outperform the nighttime data with higher degree of accuracy observed. Note, however, that due to the high northern latitude of the study area, the nighttime data were also associated with more snow cover during the winter months, and therefore there is some degree correlation between day/night acquisitions, seasonality, and snow coverage. Particularly, the presence of snow is expected to bias the ground surface height with an overestimation of terrain elevation. While ATL08 acquired during snow & nighttime conditions have higher signal to noise ratio (SNR) due to more ground photons and lower background noises captured by ICESat-2/ATLAS, the dynamic ground surface height throughout the long snow season can introduce additional biases on ground height estimation. Therefore these findings may not be applicable to lower latitudes. Meanwhile, some evidence has been found that the ATL08 classification algorithm is more prone to misclassify noise photons as surface with a relatively low SNR (Malambo and Popescu, 2020) and some unexpected results such as daytime observations Fig. 8. RMSE of ATL08 terrain height estimates at 100 m (8a; 8c; 8e) and 30 m (8b; 8d; 8f) step size stratified by day/light condition (8a; 8b), snow/no-snow condition (8c; 8d), and beam intensity (8e; 8f).
outperforming nighttime data can partially be attributed to this effect (Malambo and Popescu, 2021). In addition, the quality of ATL08 terrain product is found to be related to a combination of error factors associated with landscape structure. Some studies have been carried out to explore the impact of topographic variability on ATL08 terrain height accuracy at a 100 m × 11 m spatial scale (Liu et al., 2021;Neuenschwander et al., 2020). Our study suggests that topography/slope is one of the dominant factors regulating terrain elevation accuracy and highresolution ATL08 datasets is observed to be more susceptible to slope variability with larger uncertainties observed in high-relief terrain.

ICESat-2 for canopy height estimation
One of the primary objectives of this study was to evaluate the quality of multi-resolution ATL08 canopy height estimates in the boreal. Strong correlations were observed between ATL08 and reference canopy height datasets at both spatial resolutions. We found the ATL08 samples acquired with strong beams consistently demonstrating superior capacity on canopy height estimations at both spatial resolutions, which was expected as the weak beams were not designed for vegetation. Since each beam pair comprised of a strong and weak beam (energy ratio 4:1), strong beam indicated more photons reflected from surface, which makes ATL08 ground and canopy finding algorithm less susceptible to background noise photons (Neuenschwander and Pitts, 2019).
The ATL08 canopy height measurement errors found in this work are slightly higher than results reported by Guerra-Hernández et al. (2022) with a study region in Extremadura across Central-West Spain. Higher uncertainties observed in this study can be attributed to 1) the huge spatial coverage of our study areas with complex forest structure, and 2) increased spatial resolution of our 30-m ATL08 dataset validated. Our study areas spanned across huge geographic domains across North America boreal with forests characterized by substantial variabilities in terms of canopy cover, vertical structure and environmental conditions. Canopy shape of primary species in North America boreal regions are also more pointy than those found in Extremadura, while this unique canopy architecture, along with low sampling density of ICESat-2/ ATLAS can potentially result in lower possibilities that top of canopy being detected by returned photons (Neuenschwander and Magruder, 2016). In addition, downscaling of ATL08 to a higher spatial resolution can introduce additional uncertainties due to limited ground and canopy photons contained within individual segments for ATL08 ground and canopy finding. Results found in this work showed strong agreements with previous studies carried out in boreal forests, indicating that canopy height measurements acquired with strong beam / nighttime conditions yield higher accuracy Varvia et al., 2022). This is primarily due to an increase of terrain and canopy photons, as well as decrease of solar background noises captured by ICESat-2/ATLAS for the ATL08 ground and canopy finding. Our results regarding the influence of snow cover on canopy height estimation showed that snow cover generally increased the accuracy of canopy heights at both spatial resolutions. Since canopy height is reported as the 98th percentile with respect to the ground, presence of snow on the ground surface increases the ground elevation, thus shrinking the trees. However, shallow snow would not impact the canopy height substantially while providing a more accurate ground elevation, thus decreasing uncertainty in heights. This result is likely due to increased accuracy in terrain elevation estimation from the increased reflectivity of snow covered terrain, which subsequently increases the accuracy of height estimation. It should be noted, however, that much of our study area has snow on the ground for most of the year, and therefore the samples between snow and no-snow data are latitudinally dependent. In lower latitude boreal forests, canopy structures are often more complex, with relatively denser and taller forests. Therefore our assertion of increased surface albedo yielding more accurate terrain elevations may be convoluted by the relationship between snow cover, latitude, and forest structure. Further research to unpack these relative drivers of ATL08 canopy height accuracies is desirable.
The ATL08 canopy height estimates were found to be either overestimated or underestimated depending upon segment-level forest structure variability and plant function types characterized by different land-cover classes. Although previous studies suggested ATL08 is more prone to underestimate true canopy height due to the vertical sampling uncertainties, our findings indicated the canopy height can be overestimated by ATL08 in forested regions occupied by short vegetation. This is attributed to the misclassification of solar background photons (i. e. noise photons) or low lying fog or clouds as canopy photons.

Progresses, challenges and prospects
This study provides the first look at high-resolution ATL08 datasets with a focus on its capacity for large-scale terrain and canopy height characterization across North America boreal forests. We found that ATL08 resampled to 30 m can provide consistent and high-quality terrain and canopy height measurements with accuracies similar (but slightly lower) to the standard ATL08 at 100 m × 11 m spatial resolution. The performance of 30 m resolution ATL08 datasets acquired with different operational and external parameters were examined through stratified analysis and we found the quality of ATL08 is highly dependent upon sensor configuration and external condition (day/night & snow/no-snow & Beam Intensity), as well as terrain and canopy structure characteristics along tracks. Through this analysis, we obtained a combination of conditions which yielded high-quality terrain and canopy height estimates, which provides the scientific community a baseline to carry out strategic filtering processes for the extraction of highquality ATL08 observations. Although ICESat-2/ATLAS is not designed primarily for vegetation structure, its demonstrated utility for canopy height characterization can meet our ever-increasing demand for largescale forest structure and biomass observations, especially in highlatitude boreal ecosystem outside the orbital bounds of GEDI (>51.6 • N).
Slope effects have been found to be the primary source of uncertainties on ATL08 terrain and canopy height measurements at both spatial scales. We consider the accuracy of ground finding in the current version of ATL08 (ATL08 release 005) to not fully capture the terrain in steep hills and valleys due to over-smoothing. Subsequently, these ground finding algorithm blunders further contribute to the perceived relative canopy height measurement errors. In addition, LVIS facility datasets are spatially-dispersed and thereby prohibit our implementation of geolocation correction from ATL08 segments. Though ICESat-2 is featured with relatively high geolocation accuracy (~3.5 m), the sampling error can still result in high uncertainties on canopy height estimations in regions occupied by high-relief terrain. Note, however, that due to spatial configurations of ATL08 and LVIS footprints, our derived uncertainties of ATL08 measurements in high-relief terrain might be overestimated once the footprint grazes the edge of segments. Despite that, we caution the use of ATL08 rel005 relative canopy height data in regions occupied by high-relief terrain due to the challenges within the ground finding algorithm on ATL08. The ATL08 data product, however, is continuously updated since public release with new improvements offered in each version. A comprehensive validation for each version of the product is critical to inform the vegetation community about progress and performance (Malambo and Popescu, 2021). In this study, we carried out an in-depth evaluation of the latest ATL08 release 005 dataset; However we anticipate improvements with the ATL08 release 006 data products. Through the inter-version comparison, we did not observe significant improvements on terrain and canopy height estimates while the accuracies were consistent when compared with the version 2 & 4 of the product, carried out by Malambo and Popescu (2021) and Liu et al. (2021) respectively.
The extensive nature of LVIS-Facility dataset with total flightline of ~78,000 km provided us an unprecedented opportunity to examine the quality of multi-resolution ATL08 datasets across large latitudinal, species, and topographic gradients in the North America boreal. While LVIS provided a high-quality reference dataset, some limitations still exist when LVIS was adopted for the validation of ATL08. First, the geolocation uncertainties from spaceborne lidar cannot be well adjusted when LVIS-Facility was used as reference given the partial coverage of large footprints relative to ATL08 segments. Since ICESat-2/ATLAS have been addressed with beam-dependent geolocation error ranging from 2.5 m for beam 6 to 4.4 m for beam 2 (Luthcke et al., 2021), previous studies validating ICESat-2 ATL08 attempted to shift discrete-returned lidar (DRL) point cloud over a horizontal plane for the best match of terrain elevation (Liu et al., 2021;Malambo and Popescu, 2021;Neuenschwander et al., 2020) with assumption that the shifted distance is attributed to the geolocation error. In contrast, the spatial configuration of LVIS and ATL08 constrained our geolocation correction procedure while direct comparison between ATL08 and spatially matched LVIS-Facility footprints can thereby introduce additional uncertainties. In short, some spatial discrepancy between LVIS and ICESat-2 likely added uncertainty to this comparison.
Secondly, LVIS is a sampling laser ranging system which provides consistent but spatially-dispersed observations with each footprint as an observational unit. While LVIS-Facility has a nominal footprint diameter of ~10 m, the exact footprint size is altitude-dependent, and can vary between 5 m and 20 m. Due to the lack of altitude information during time of operation, applying a fixed footprint diameter on LVIS-Facility footprints can either overestimate or underestimate the true canopy height within ATL08 segments. Additionally, the inter-footprint is a function of multiple factors (e.g., flying speed, scanning repetition rate, etc.), thereby can be highly variable. The assessment of ATL08 terrain and canopy heights without sufficient LVIS footprints can make the reference datasets less representative. The utilization of DRL as reference data can confront similar issues as well when the limited vertical sampling density of ALS was not high enough to characterize the canopy structure. Overall, each reference dataset has a respective trade-off and a strategic selection of high-quality validation data can potentially avoid uncertainties on our accuracy assessments.
The three-year operational period of ICESat-2 primary mission provides billions of ground samples with equatorial spacing of <2 km (Markus et al., 2017). The development of high-resolution ATL08 dataset with a spatial scale equivalent to Landsat imagery can facilitate our extension of lidar plots over space and time for the advancements of wall-to-wall forest attributes modeling and mapping. While some initial efforts have been made on the integration of ICESat-2 and optical data for the spatially continuous mapping of forest structural attributes, most approaches are either based on spatial interpolation (Liu et al., 2022) or the adoption of coarse-resolution imagery such as MODIS as training datasets (Xi et al., 2022) due to the unmatched spatial extent between passive imagery and lidar plots. Through comprehensive and transparent assessment, the high-resolution ICESat-2 ATL08 has demonstrated strong promise for terrain and canopy height estimation and the results found in this study can serve as a baseline for highquality lidar plot extraction prior to the integration of multi-sourced remote sensing observations.

Conclusion
The recent development of high-resolution ATL08 dataset at 30 m × 11 m spatial scale has the potential to broadly extend the suite of ICESat-2 applications for forest mapping and monitoring. In this study, we examined the accuracies of multi-resolution ATL08 terrain and canopy heights estimates, at both 100 m × 11 m and 30 m × 11 m spatial scale and across the North America boreal ecosystem. The LVIS-Facility data acquired during the Summer of 2019 enabled the validation of ICESat-2 products over broad geographic and structural gradients. Overall, the downscaling of ATL08 to 30 m × 11 m spatial resolution resulted in a slightly lower accuracy on terrain elevation and canopy height estimates, but accuracies were still generally high. The highest quality canopy height measurements were from strong beams, acquired during the night. We found that snow cover had limited influence on the accuracy of terrain and canopy estimates, but this finding may be latitudinally dependent and deserves further research. We speculate that snow increases the reflectivity of the ground, and thus the accuracy of terrain retrievals. However, deep snow would likely bias the height results, decreasing the observed forest height, and therefore data acquired late in the winter season at high latitudes should be used with caution. The RMSE and Bias were relatively consistent across reference height and land-cover gradients while some high RMSE and Bias were observed, mainly attributed to the presence of short vegetation within given reference height and land-cover intervals. The ATL08 product has demonstrated strong potential on canopy height measurements across the North America boreal ecosystem and the high-quality vegetation structure data can serve as a critical input for multi-scale biomass and carbon cycle modeling.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Data availability
Data will be made available on request.