Comparison of LiDAR and Digital Aerial Photogrammetry for Characterizing Canopy Openings in the Boreal Forest of Northern Alberta

Forest canopy openings are a key element of forest structure, influencing a host of ecological dynamics. Light detection and ranging (LiDAR) is the de-facto standard for measuring three-dimensional forest structure, but digital aerial photogrammetry (DAP) has emerged as a viable and economical alternative. We compared the performance of LiDAR and DAP data for characterizing canopy openings and no-openings across a 1-km2 expanse of boreal forest in northern Alberta, Canada. Structural openings in canopy cover were delineated using three canopy height model (CHM) alternatives, from (i) LiDAR, (ii) DAP, and (iii) a LiDAR/DAP hybrid. From a point-based detectability perspective, the LiDAR CHM produced the best results (87% overall accuracy), followed by the hybrid and DAP models (47% and 46%, respectively). The hybrid and DAP CHMs experienced large errors of omission (9–53%), particularly with small openings up to 20m2, which are an important element of boreal forest structure. By missing these, DAP and hybrid datasets substantially under-reported the total area of openings across our site (152,470 m2 and 159,848 m2, respectively) compared to LiDAR (245,920 m2). Our results illustrate DAP’s sensitivity to occlusions, mismatched tie points, and other optical challenges inherent to using structure-from-motion workflows in complex forest scenes. These under-documented constraints currently limit the technology’s capacity to fully characterize canopy structure. For now, we recommend that operational use of DAP in forests be limited to mapping large canopy openings, and area-based attributes that are well-documented in the literature.


Introduction
Canopy openings occur regularly and at various spatial scales throughout a forest, due to numerous exogenous and endogenous factors. In the boreal forest of Alberta, Canada, large openings are caused by fires [1], insect manifestation [2], windthrow [3], anthropogenic disturbances [4,5], and persistent suboptimal growing conditions associated with excessive soil moisture [6], water, or rock outcroppings [7]. Smaller openings might be the outcome of individual limb fall or the death of single trees [8]. In all cases, open areas can experience markedly different biotic and abiotic conditions, including a locally increased supply of photosynthetically active radiation (e.g., sunflecks) [9,10]. If left undisturbed, these openings constitute ecologically important patches of pronounced tree recruitment and plant growth. Thus, canopy openings support the process of regenerating forest vegetation, and help maintain meaningful biodiversity within the forest ecosystem [11][12][13]. The altered physical Table 1. A summary of recent (circa~10 years) studies comparing forest-structure information derived from light detection and ranging (LiDAR) and digital airborne photogrammetry (DAP) reveals the dearth of research assessing canopy openings directly. This is not intended to be an exhaustive list.

Motivation and Objectives
We require additional studies comparing LiDAR and DAP for mapping discrete canopy openings. In particular, the literature could benefit from work in other forest types, using high-density point clouds, high along-and across-track image overlap, and ground-based reference data. In this study, we aim to contribute such additional insights through experiments in the boreal forest.
Our motivation was to assess the capacity modern 3D remote-sensing data to characterize forest canopy structure, not to contribute to the disturbance-ecology literature. For that reason, we are interested in canopy openings of any size, including both very small openings and those that are part of large open-ended systems: both of which would fall outside of most classical definitions of "canopy gap" [59]. We did this for two reasons. First, our primary goal was to learn about the respective capabilities of DAP and airborne LiDAR data for capturing discrete forest attributes. As we have indicated, the vast majority of previous studies have focused on trees; relatively few have looked at the spaces between them. Second, it is our contention that the size and ecological impact of canopy openings must be viewed through the lens of an individual species in order to measure their effects. As a result, spatial limitations seemed too arbitrary to be employed in our work.
In order to differentiate our motivations from the those in the canopy-gap literature, we employ the term canopy opening rather than canopy gap, which we perceive to refer to the functional response of openings rather than their physical structure. Our work is designed to provide insights into the nature of emerging DAP technologies, and to help future practitioners make informed decisions on the use of LiDAR and DAP in forested environments. While our findings may inform future studies on canopy gaps, that subject is outside of our current scope.
This study is the first (that we are aware of) in the peer-reviewed literature to compare airborne LiDAR and DAP for characterizing discrete canopy openings in the boreal forest, and one of just a handful to have done so in any biome. In addition, this work employs high-density DAP data extracted from high-overlap (80% along-track and 60% across-track) aerial photography, which has not been evaluated in previous works

Materials and Methods
We delineated canopy openings in the boreal forest of northern Alberta, Canada, using canopy height models (CHMs) derived from three different data sources: (i) LiDAR, (ii) DAP, and (iii) a LiDAR/DAP hybrid. In addition to individual datasets, we used the combination of LiDAR and DAP to generate a hybrid model, theorizing that a canopy height model that combined the relative strengths of DAP (higher-density point clouds) and LiDAR (more reliable ground detection) could lead to better performance than the DAP model alone. A graphical summary of the workflow used to conduct our experiments is shown in Figure 1. Details are provided in the subsequent subsections. LiDAR/DAP hybrid. In addition to individual datasets, we used the combination of LiDAR and DAP to generate a hybrid model, theorizing that a canopy height model that combined the relative strengths of DAP (higher-density point clouds) and LiDAR (more reliable ground detection) could lead to better performance than the DAP model alone. A graphical summary of the workflow used to conduct our experiments is shown in Figure 1. Details are provided in the subsequent subsections.

Figure 1.
Workflow chart depicting the steps performed in this study, with numbers referring to the corresponding subsections of this paper. Two point clouds-digital aerial photogrammetry (DAP) and light detection and ranging (LiDAR) served as the basis for the derivation of three canopy height model (CHM) alternatives: LiDAR (LiDAR_CHM), DAP (DAP_CHM), and a LiDAR-DAP hybrid (Hybrid_CHM). Evaluation data was generated during a field campaign. Both an area-based as well as a point-based analysis of the resulting opening were performed. Detailed flowcharts for the data processing using LAStools are provided in the supplementary material.

Study Area
Our study took place in the boreal forest biome of northern Alberta, Canada: an environment that has been heavily disturbed by petroleum exploration and other forms of resource extraction. This region is in the boreal snow climate, and is characterized by moderately warm, humid summers and cold, dry winters [60]. The average air temperature is below -10 °C in January and below 20 °C in July, with annual precipitation about 445 mm [60]. The elevation of the study area is approximately 700 m above mean sea level.

Study Area
Our study took place in the boreal forest biome of northern Alberta, Canada: an environment that has been heavily disturbed by petroleum exploration and other forms of resource extraction. This region is in the boreal snow climate, and is characterized by moderately warm, humid summers and cold, dry winters [60]. The average air temperature is below −10 • C in January and below 20 • C in July, with annual precipitation about 445 mm [60]. The elevation of the study area is approximately 700 m above mean sea level.
The specific study area for this research was a 1 km × 1 km expanse situated 40 km south of the town of Conklin ( Figure 2). It includes natural upland-forest stands, fens, and human disturbances related to bitumen extraction. Due to varying growing conditions, tree density was observed to vary markedly, as were the size and shape of openings; especially between the wet lowlands and the drier uplands. Human-disturbance features including seismic lines (petroleum exploration corridors), pipelines, well sites, and roads occur throughout the study area.

Remote Sensing Data Collection
ORTHOSHOP Geomatics Ltd. (OGL) Engineering of Calgary used a piloted aircraft (Cessna 210T) to acquire remote-sensing data for this project over two flights in 2017: one in May (leaf-off conditions) and the other in August (leaf-on conditions). The platform carried a four-band (blue, green, red, and near infrared) optical sensor (Leica RCD30) and a LiDAR instrument (Leica ALS70). The average flying height for both missions was 850 m above ground level. The projection system/datum used for this survey was NAD83 (CSRS) Universal Transverse Mercator Zone 12N. We installed 250 clearly visible ground control points (GCPs) prior to the first mission. These were left in the field until the completion of the second mission. We recorded the geolocations of these GCPs using a Trimble R8 real-time kinematic (RTK) global navigation satellite system (GNSS) with 8 mm horizontal and 15 mm vertical precision.
The purpose of the spring flight (leaf-off conditions) was to obtain data that would deliver the best-possible digital terrain models (DTMs), while the summer flights (leaf-on conditions) were to capture full-canopy digital surface models (DSMs). Unfortunately, the LiDAR sensor experienced technical difficulty during the spring flight, so DAP was the only leaf-off data available for our analysis.
For the optical data, flight lines were designed to deliver a minimum of 80% endlap and 60% sidelap, with a ground resolution of 5.5 cm or less. The raw digital photos were processed into photogrammetric point clouds and four-band digital orthophotos (RGBN) using a structure-from-motion workflow [62] in Pix4D (Version 4.3.27). Both image-based point clouds (leaf-on and leaf-off) had mean point densities greater than 200 points/m 2 for the leaf-on data. Of the 250 GCPs surveyed, we used 100 for georeferencing purposes (point cloud and orthophoto) and reserved the remaining 150 for accuracy assessment. The X, Y, and Z accuracy (root mean square error) of the final leaf-off optical products were 5, 6, and 13 cm, respectively. The X, Y, and Z accuracy of the leaf-on optical products were 6, 6, and 19 cm root mean square error (RMSE), respectively. A 1 km × 1 km sample area was clipped from the larger mission footprint (5.5 km × 6 km) to match the specific study area for this project.
The discrete-return LiDAR data (2 returns/pulse) were converted to a raw point cloud using inertial and GNSS measurements. Mean point density of the LiDAR data was~40 points/m 2 with a mean point spacing of 15 cm.   Human disturbances related to petroleum extraction-roads, seismic lines (petroleum exploration corridors), and exploration wells-are easily distinguishable from the underlying undisturbed forest. The in-situ reference points are displayed using blue circles and image analysis-based reference points are displayed using red circles.

LiDAR-Based CHM
Leaf-on LiDAR data was used to obtain the LiDAR-based CHM. We used Bentley Microstation with Terrasolid software to remove noise from the raw LiDAR point cloud, and then classified it into ground and non-ground points using LAStools (Version 170802). We next rasterized (20 cm pixels) the ground points into a digital terrain model (DTM LiDAR ), and the first-return points into a digital surface model (DSM LiDAR ). Finally, a canopy height model (CHM LiDAR ) was obtained by subtracting the DTM LiDAR from the DSM LiDAR .

DAP-Based CHM
We used LAStools to classify the leaf-off DAP point clouds into ground and non-ground points. The ground points were then rasterized to obtain DTM DAP. The leaf-off data was used for this purpose, under the assumption that ground can be better modeled under leaf-off conditions. We then rasterized the leaf-on DAP points (all points) to obtain DSM DAP . Finally, we subtracted the DTM DAP from the DSM DAP to create a CHM DAP . All DAP-derived products had a spatial resolution of 20 cm, matching those from the LiDAR.

Hybrid CHM
LiDAR sensors have better canopy penetration abilities than DAP, and are therefore better suited to modeling the terrain [42]. Hence, we produced a hybrid CHM (CHM Hybrid ) using the DSM DAP and the DTM LiDAR . We theorized that this LiDAR-fixed DAP CHM might lead to better detection of canopy openings, given the reliance on fixed-height thresholds for delineating canopy openings, as described in Section 2.3, and the high-density point cloud produced by DAP.

Delineating Canopy Openings
We used a fixed-height approach to extract canopy openings from each of the experimental CHMs. The purpose of this processing step was to separate low-lying and regrowth vegetation from established canopy forming vegetation. While some studies have used variable-height approaches for this task, where the threshold varies as a function of stand height [e.g., 8], we elected for the more parsimonious fixed-height approach. This decision was primarily designed to simplify the job of field personnel.
Given the presence of short, stunted trees within the large wetland portions of our study area, we set a relatively low canopy threshold of 1.3 m. Thus, vegetation below 1.3 m was considered non-canopy vegetation. Vegetation above this threshold was considered established vegetation, as it is mostly safe from ungulates [54]. We made no distinction between trees and other forms of vegetation, such as shrubs. The three CHMs (LiDAR, DAP, and hybrid) were converted to binary opening/no-opening maps using this threshold (i.e., CHM > 1.3 = no-opening and CHM ≤ 1.3 = opening).

Validating Canopy Openings-Point-Based Analysis
We used a multi-stage sampling approach to create a reference dataset suitable for validating the accuracy of the three canopy-opening maps. First, we divided the study area into two strata, disturbed and undisturbed, since the patterns of openings within these two classes across the study area were observed to be fundamentally different. Disturbed areas included roads, seismic lines, and well sites, which we delineated manually using the 5.5-cm orthophotos. Undisturbed areas were defined as the natural areas outside of these anthropogenic features. It is important to note that many of the disturbed portions of our study area supported vegetation taller than 1.3 m, which would be subsequently classified as no-opening.
The second sampling stage involved delineating a stratum describing canopy-opening size. We defined five opening-size classes to provide categories that could be distinguished unambiguously by personnel in the field: Sample size (n) was calculated following Kershaw et al. [63], who described an equation designed to incorporate factors of variability (CV) and allowable error (E): where t is the probability level (Student's t value) implied by the associated E value, and CV is the coefficient of variation of canopy openings across the study area. We used an E value of 20% (confidence level 80%), which resulted in a t value of 1.3. The variability of openings and no-openings across the study area (CV) was estimated using the canopy-opening map produced with the CHM LiDAR layer. We selected CHM LiDAR for this purpose given its status in the literature and the experience of White et al. [8], who found it to be a good source of canopy-opening information.
Equation (1) produced a recommended sample size of 1835, which we distributed proportionally across the study area using ArcMAP (Version 10.5). The distribution of these sample points is summarized in Table 2. Table 2. Distribution of validation samples across opening size class and disturbance regime. Equation (1) suggested higher values for Class 0 (both in disturbed and in undisturbed areas) and a smaller number for Class 2 in undisturbed areas. All three values were adjusted for feasibility and statistical robustness, respectively.

Opening Class
Disturbed Undisturbed   0  100  300  1  116  102  2  36  30  3  69  104  4  396  582  Total  717  1118 We used a combination of in-situ measurements (879 points) and visual image analysis (956 points) to collect opening-class reference attributes associated with the 1835 locations in our sample. We elected to perform as many in-situ observations as possible, and used cluster sampling to maximize the efficiency of personnel in the field. For safety and logistical reasons, field sampling was limited to within 30 m of 84 previously existing study sites, which were being visited by researchers working on a related project. These study sites were carefully established (using a stratified-random sampling technique) to account for the structural and functional variability present within the study area. The field campaign was executed in July 2018, and represents leaf-on conditions.
Field personnel used an RTK GNSS with a horizontal accuracy of 10 cm or better to navigate to the in-situ sample points. Upon reaching a sample location, the surveyor checked if the sky was visible at 1.3 m above ground to determine if the location was an opening or no-opening. If it was an opening, the surveyor also assigned a class label to the opening to the best of her judgement. A photographic record of each site was made with a vertically oriented camera pointed at the sky at 1.3 m; representative photographs from each class are shown in Figure 3. Field calls were made by the same surveyor throughout the study for consistency. We visited 879 field locations for in-situ measurements.
Visual image analysis was performed on the remaining 956 points using orthophotography and the LiDAR CHM. Once again, image calls were made by the same surveyor for consistency.
We acknowledge that using LiDAR data in our visual assessment had the potential to bias our reference data. However, a preliminary assessment of LiDAR accuracy using the 879 in-situ points produced an overall accuracy of 97.8%. This provided us with confidence in being able to use LiDAR data as a visual reference-as well as orthophotography-for supplementary image-based samples. Our in-situ sampling design (cluster sampling) led to concerns regarding spatial autocorrelation in the reference data ( Figure 2). To address this, we inspected the range of a semi-variogram created from the DSMLiDAR, and determined that 20 m was the distance at which sample locations within our study area were no longer spatially autocorrelated. Consequently, we subdivided our full set of validation points into 100 random subsets. Within each subset, individual samples were kept least 20 m apart. This was achieved with the following steps: 1. Identify and isolate points that do not have a neighboring point within a 20 m distance. Our in-situ sampling design (cluster sampling) led to concerns regarding spatial autocorrelation in the reference data ( Figure 2). To address this, we inspected the range of a semi-variogram created from the DSM LiDAR , and determined that 20 m was the distance at which sample locations within our study area were no longer spatially autocorrelated. Consequently, we subdivided our full set of validation points into 100 random subsets. Within each subset, individual samples were kept least 20 m apart. This was achieved with the following steps:

1.
Identify and isolate points that do not have a neighboring point within a 20 m distance.

2.
Randomly select one point from the remaining points. Identify all remaining points within a 20 m distance of this point and remove them from the list.

3.
From the remaining points in the list, select another random point. Identify all the points within a 20 m distance of this point and remove them from the list. 4.
Repeat step 3 until the list is empty.

5.
Add the randomly selected samples (steps 2-4) and the isolated samples (step 1) together to obtain a complete set of spatially uncorrelated samples.
These steps were repeated 100 times to generate 100 subsets of random samples. The sample distribution of these 100 subsets across classes is presented in Table 3. This exercise allowed us to obtain an average of 358 unbiased samples in each sampling round. These samples were observed to be evenly distributed across opening classes. The low standard deviation values in Table 3 indicate that the number of samples generated in each run did not deviate from the mean by a large margin. Consequently, the average results of these 100 sample sets were used in the remaining sections of this paper for analyses and accuracy assessments.

Examination of Opening Characteristics-Area-Based Analysis
While the point-based validation exercise (Section 2.4) provided important information on the ability of the three candidate products to detect canopy openings and no-openings observed at single points, it did not explicitly compare the resulting scaled-up aerial estimates summed across the whole study area. As a result, we performed a second analysis that compared the total area of openings measured by LiDAR with those detected by each of the two DAP-based products.
To further investigate this, we also compared the area of individual openings produced by DAP/hybrid with that of LiDAR. To do so, we converted each LiDAR-derived opening into a polygon feature and attributed it with a corresponding class. Then, the LiDAR polygons were overlaid onto the DAP and hybrid CHMs (Figure 4). In this manner, we could quantify the percent area of the LiDAR-based opening polygons that were identified by DAP and hybrid. For visual analysis, we also produced a raster of agreement, which communicates the number of datasets documenting any given opening.

Examination of the Value of DAP in Measuring Height
In addition to experiments designed to compare the capacity of the experimental datasets to characterize canopy openings, we also conducted a simple analysis to assess the value of the hybrid model in estimating canopy height, compared to the DAP model. To do so, we used the LiDAR-measued heights as reference and created difference maps for DAP (CHM LiDAR -CHM DAP ) and hybrid (CHM LiDAR -CHM Hybrid ) models. Subsequently, we classified the difference maps into five classes (−30 to −10 m, −10 to −2 m, −2 to 2 m, 2 to 10 m, and 10 to 28 m) for a visual assessment.
To further investigate this, we also compared the area of individual openings produced by DAP/hybrid with that of LiDAR. To do so, we converted each LiDAR-derived opening into a polygon feature and attributed it with a corresponding class. Then, the LiDAR polygons were overlaid onto the DAP and hybrid CHMs (Figure 4). In this manner, we could quantify the percent area of the LiDAR-based opening polygons that were identified by DAP and hybrid. For visual analysis, we also produced a raster of agreement, which communicates the number of datasets documenting any given opening.

Examination of the Value of DAP in Measuring Height
In addition to experiments designed to compare the capacity of the experimental datasets to characterize canopy openings, we also conducted a simple analysis to assess the value of the hybrid model in estimating canopy height, compared to the DAP model. To do so, we used the LiDAR-  Figure 5 presents the three derived maps delineating canopy openings and no-openings from DSM LiDAR , DSM DAP , and DSM Hybrid . Substantial differences between the LiDAR and DAP/hybridyHyH maps are apparent, especially in terms of number of openings and the total area classified as opening. However, there is very little visible difference between DAP and hybrid.

Point-Based Analysis
When the opening maps were evaluated against our point-based validation data (Table 4), LiDAR was found to produce the highest overall accuracy: 94%. DAP and hybrid data had overall accuracies of 77% and 78%, respectively. A Z-test revealed that all these accuracies were significantly different from each other (p < 0.05) at the 95% confidence level. Table 4 further reports on errors of omission and commission, which can be translated into underestimation and overestimation of openings. LiDAR tended to slightly underestimate the true number of canopy openings (13% errors of omission; 1% errors of commission), and correspondingly overestimate the number of no-openings (0% errors of omission; 9% errors of commission). The same error pattern-underestimation of openings; overestimation of no-openings-was observed in the DAP and the hybrid datasets, though more exaggerated. Table 4. Summary of opening and no-opening detection confusion matrices: mean overall accuracies, and errors of omission and commission. One-hundred sets of uncorrelated samples were used to produce these statistics.  Figure 6 presents relative-detection accuracies within different opening-size classes. Accuracies were observed to rise for all three data sets as opening-size class increased. No-openings (Class 0) were detected at close to 100% by all three approaches. LiDAR-detection results were notably higher than those produced by DAP and hybrid data sets for all four classes of openings. While DAP and hybrid CHMs were found to be reasonably able to detect openings larger than 200 m 2 (Class 4; >70%), small openings <20 m 2 had very low (<20%) detection rates. LiDAR, on the other hand, showed detection rates >70% for all size classes and >90% for openings larger than 20 m 2 (Classes 3 and 4). On an average, overall opening detection accuracy of LiDAR, DAP and hybrid models were 87%, 46%, and 47% respectively (excluding no-opening samples).  1 Figure 6. Results of canopy opening and no-opening detection presenting the proportion of correctly detected openings, relative to their size class. Class 0 is no-opening, Class 1 is openings >0-4m 2 , Class 2 is openings >4-20m 2 , Class 3 is openings >20-200m 2 , and Class 4 is openings >200m 2 . The evaluation of size class-dependent detection rates is based on the average of 100 sample runs.

Area-Based Analysis
The total area of openings detected by the LiDAR dataset was 245,920 m 2 , which is equivalent to about 33% of the study area. DAP and hybrid data estimated considerably smaller total areas of openings: 152,470 m 2 and 159,848 m 2 , respectively. This represents just 62% (DAP) and 65% (hybrid) of the amount of openings detected by LiDAR. Table 5 provides a more detailed assessment of our area-based analysis, with openings grouped by class and percentage of area detected. The vast majority (96% and 97%) of the very small (0-4 m 2 ) openings recorded by LiDAR were not detected at all in the DAP or hybrid datasets. Of the~4% that were detected by DAP and hybrid, about half (1.8%) covered less than 20% of the polygon areas. Less than 0.5% of the very small openings detected by LiDAR were represented completely (80-100% of the area covered) by DAP and hybrid. This pattern of DAP and hybrid either missing canopy openings detected by LiDAR or delineating just a portion of their area improved steadily as opening size increased, to the point where all of the very large (>200 m 2 ) openings detected by LiDAR were also delineated by DAD and hybrid, though never as thoroughly (only~4-8% of the areas fully depicted).  Figure 7 shows the raster of agreement displaying the spatial overlap among all three opening maps. Large openings on roads, clear cuts, conventional seismic lines, and water bodies were consistently detected by all three approaches (pixels displayed in green). Higher disagreement was found among smaller openings (yellow, orange, and red pixels). Abundant orange dots found throughout the image indicate that many small openings were captured only by the LiDAR dataset. These graphical observations comply with the numeric results presented previously.

Value of the Hybrid Model-CHM Hybrid vs. CHM DAP Relative to CHM LiDAR
There was no major difference in results obtained by the DAP and hybrid datasets in identifying and delineating openings and no-openings (Sections 3.1 and 3.2), suggesting that there is little value in integrating LiDAR and DAP for this purpose. In retrospect this is understandable, since both DAP and hybrid techniques used the same DSM, which is the controlling layer for binary opening/no-opening detection. However, our assumption was that the added value of the hybrid model would be in estimating canopy height, as it uses a better DTM (from LiDAR). To examine this assumption, we subtracted CHM DAP and CHM Hybrid from CHM LiDAR and analyzed the resulting differential images (Figure 8). Figure 8a shows large, solid patches of high divergence (illustrated in shades of turquoise and dark brown) in the differential image between CHM DAP and CHM LiDAR . These patterns are the result of the overestimation of ground height and underestimation of vegetation height by digital photogrammetry [64]. While 70% of all pixels of CHM DAP were within the +/− 2 m range of CHM LiDAR canopy heights, this proportion increased to 81% when comparing CHM Hybrid and CHM LiDAR . Root mean square errors were 3.6 m and 3.1 m for CHM DAP and CHM Hybrid , respectively. This simple analysis provides clear evidence that canopy height estimation by the hybrid dataset is more consistent than the DAP dataset. However, this did little to assist in the detection of canopy openings.

Discussion
We examined the capability of three different 3D datasets to map canopy openings in the boreal forest: LiDAR, DAP, and hybrid. As expected, the LiDAR dataset exhibited the highest accuracy measures when compared to the validation data, both in terms of opening detectability (73%, 81%, 91%, and 93% for classes 1 through 4, respectively) and overall accuracy (94%, which includes the noopening class). The DAP and hybrid approaches produced much more modest results, both for opening detectability (DAP 8%, 21%, 46%, and 70%; hybrid 7%, 18%, 48%, and 73% for classes 1 through 4, respectively) and overall accuracy (77% and 78% for DAP and hybrid, respectively). Upon closer inspection, DAP and hybrid data sets were found to perform worst in high-density tree stands, where openings are generally small. This can be seen in Figure 9, where the DAP data set omits many small openings. When we stratified the results by opening size, very large (> 200 m 2 ) openings and no-opening areas were mapped reliably by both DAP and hybrid. As a result, we recommend these technologies as a good, low-cost solution for mapping large openings, such as those produced by windthrow, fire, and many anthropogenic disturbances. However, LiDAR's superior accuracy with smaller openings makes it the preferred technology for mapping the complete range of openings encountered in the boreal forest, which includes very small openings.

Discussion
We examined the capability of three different 3D datasets to map canopy openings in the boreal forest: LiDAR, DAP, and hybrid. As expected, the LiDAR dataset exhibited the highest accuracy measures when compared to the validation data, both in terms of opening detectability (73%, 81%, 91%, and 93% for classes 1 through 4, respectively) and overall accuracy (94%, which includes the no-opening class). The DAP and hybrid approaches produced much more modest results, both for opening detectability (DAP 8%, 21%, 46%, and 70%; hybrid 7%, 18%, 48%, and 73% for classes 1 through 4, respectively) and overall accuracy (77% and 78% for DAP and hybrid, respectively). Upon closer inspection, DAP and hybrid data sets were found to perform worst in high-density tree stands, where openings are generally small. This can be seen in Figure 9, where the DAP data set omits many small openings. When we stratified the results by opening size, very large (>200 m 2 ) openings and no-opening areas were mapped reliably by both DAP and hybrid. As a result, we recommend these technologies as a good, low-cost solution for mapping large openings, such as those produced by windthrow, fire, and many anthropogenic disturbances. However, LiDAR's superior accuracy with smaller openings makes it the preferred technology for mapping the complete range of openings encountered in the boreal forest, which includes very small openings.
Remote Sens. 2019, 11, x FOR PEER REVIEW 19 of 25 Figure 9. A visual comparison between LiDAR and DAP opening map segments in the study area. In the map segments, the purple areas represent pixels classified as opening by the respective approach, whereas the trees visibility is enhanced by using a green overlay. Black rectangles in each map segments represent an area for which corresponding vertical profiles of original elevation data (point cloud) are displayed below. Figure 6A displays a sample area in natural landscape where small openings are picked up by LiDAR (left) but not by DAP (right). Figure 6B represents a sample area in disturbed landscape where roads and seismic lines are more or less perfectly detected by both techniques. However, the LiDAR data identify other smaller openings in the surrounding landscape that DAP does not identify. In the map segments, the purple areas represent pixels classified as opening by the respective approach, whereas the trees visibility is enhanced by using a green overlay. Black rectangles in each map segments represent an area for which corresponding vertical profiles of original elevation data (point cloud) are displayed below. Figure 6A displays a sample area in natural landscape where small openings are picked up by LiDAR (left) but not by DAP (right). Figure 6B represents a sample area in disturbed landscape where roads and seismic lines are more or less perfectly detected by both techniques. However, the LiDAR data identify other smaller openings in the surrounding landscape that DAP does not identify.
This propensity of DAP-derived CHMs to miss smaller openings seems ill-matched to the very high density of the DAP point cloud we employed (>200 pts/m 2 ), though it mirrors the observations of White et al. [8] and Zielewska-Büttner et al. [65]. The omission of small openings might be due to a number of factors, including occlusions, shadows, and tree sway between photographs, as reported previously [8]. However, our photography was acquired in light-wind and strong, diffuse-light conditions specifically designed to reduce these known limitations. Even though our DAP data set was acquired with much higher image acquisition overlap and point density (80% along-track, 60% across-track and 270 points/m 2 ) than previous studies, e.g., White et al. [8] (60% along-track and 20% across-track and 12.27 points/m 2 ), the DAP data showed a significant shortcoming in detecting smaller openings. Interestingly, we found canopy openings-both large and small-to be easily observed using the same photoset that underlies the DAP data, when viewed in a softcopy environment. Our suspicion is that the photogrammetric data was degraded by the structure-from-motion workflow. Previous authors have noted the difficulty of detecting individual trees automatically, and the propensity of structure-from-motion workflows to confuse tie points in such complex scenes [40]. In addition to observing best practices related to flight parameters, image overlap, and sensor resolution when planning photogrammetric flights [58], we encourage additional enhancements to point-matching algorithms and other elements of the structure-from-motion workflow designed specifically to deal with complex forest scenes. We believe that there is more information in stereo photographs than what is currently being delivered by contemporary software such as Pix4D.
Most of the factors that currently limit DAP datasets have little impact on the data-collection process of LiDAR. As a result, our LiDAR data was able to detect very small openings, even though its point density was much lower (~40 pts/m 2 ) than that of our DAP data (>200 pts/m 2 ). As mentioned, this observation confirms the findings of White et al. [8], who compared lower-density LiDAR (11.6 points/m 2 ) and DAP (12.27 m 2 ) in the temperate rain forest of British Columbia. They found significant differences in opening sizes and numbers, as well as large differences in the overall accuracies of the detection methods. These large errors further commensurate the results produced by Zielewska-Büttner et al. [65], who used fixed-height thresholds on DAP data to detect canopy openings in southern Germany. They found errors of omission to be 48% in tall-growth forests.
Our results showed that there was little value in integrating LiDAR and DAP for the purpose of mapping canopy openings in our study area. This is because the same DSM is used for both the DAP and hybrid techniques, which is the most important factor to binary opening/no-opening detection. The real value of the hybrid dataset comes in estimating canopy height, which is improved substantially by the canopy-penetrating abilities of LiDAR.
It must be noted that our recommendations regarding the use of DAP and hybrid datasets is directed at mapping structural-opening features only. For other applications (e.g., canopy height mapping [36,40], plot metrics [3,36], and microtopography characterization [66]), photogrammetric techniques may produce much better results. However, evaluating those attributes was beyond the scope of this study.

Conclusions
Our investigation reveals DAP-derived products, either with or without the enhancement of a LiDAR-based DTM, to be poor alternatives to LiDAR for mapping the full range of canopy openings in the boreal forest. While large openings are well-mapped with photogrammetric data, we found these products to be hindered by occlusions, mis-matched tie points, and other optical challenges, making small openings difficult to detect and delineate. While we did not quantify it specifically, we suspect that there is more information in the underlying stereo photography than what is currently being derived by contemporary structure-from-motion software such as Pix4D.
Given that DAP is an attractive, affordable, and easily accessible data source, further research into the derivation of detailed CHMs from either DAP or hybrid data sets is recommended, and we look forward to future methodological enhancements designed to address the specific challenges posed by structure-from-motion in complex forest scenes. At this moment, DAP's and the hybrid model's strength is clearly the affordable detection and delineation of larger-sized (>200 m 2 ) forest openings. Expanding the capabilities of DAP would facilitate additional monitoring, modeling, and ultimately management applications in the boreal forest. However, currently, our observation is that DAP cannot produce the accuracies needed for the complete characterization of complex canopy structures, which include very small (<10 m 2 ) openings. However, the results produced in this and other studies justify optimism that with further development, enhanced 3D modeling of forest structure with DAP is forthcoming.