Tree canopy and snow depth relationships at fine scales with terrestrial laser scanning

. Understanding the impact of tree structure on snow depth and extent is important in order to make predictions of snow amounts, and how changes in forest cover may affect future water resources. In this work, we investigate snow depth under tree canopies and in open areas to quantify the role of tree structure in controlling snow depth, as well as the controls from wind and topography. We use fine scale terrestrial laser scanning (TLS) data collected across Grand Mesa, Colorado, USA, to measure the snow depth and extract horizontal and vertical tree descriptors (metrics) at six sites. We apply the Marker- 15 controlled watershed algorithm for individual tree segmentation and measure the snow depth using the Multi-scale Model to Model Cloud Comparison algorithm. Canopy, topography and snow interaction results indicate that vegetation structural metrics (specifically foliage height diversity) along with local scale processes such as wind are highly influential on snow depth variation. Our study specifies that windward slopes show greater impact on snow accumulation than vegetation metrics. In addition, the results emphasize the importance of tree species and distribution on snow depth patterns. Fine scale analysis 20 from TLS provides information on local scale controls, and provides an opportunity to be readily coupled with airborne or spaceborne lidar to investigate larger-scale controls on snow depth.


Terrestrial Laser Scanning
We collected TLS data in snow-off (fall 2016) and snow-on (winter 2017) conditions at Grand Mesa at several sites ( Fig. 1, Table 1) (Glenn et al., 2019;Hiemstra et al., 2019). The winter 2017 data collection occurred over 16 days but without significant snowfall between days. Each site was scanned once during the duration. A Riegl VZ-1000 (1550 nm) and Leica Scan Station C10 (532 nm) were used. Multiple scans were obtained at each site and coregistered to produce a single point 100 cloud for each site and date. Coregistered scans were then georegistered using surveyed locations within the plots. The georegistered scans (i.e. area of analysis) for each site ranges from approximately 10,000 to 38,000 m 2 ( Table 1).
The TLS data were then utilized to derive snow depths, vegetation metrics, and topographic indices. From these data, we developed a series of statistical analyses to investigate relationships between the canopy and snow depths under the canopy at each of our sites. We also describe snow depths in open areas with no trees (of a least 0.5 m height). Methods on identifying 105 individual trees, under the canopy and in the open, are described below.

TLS Data Processing
We processed the TLS data into returns from ground and vegetation (fall 2016) or snow and vegetation (winter 2017), and estimated snow depths at each of the sites, using several sub-routines in CloudCompare (v2.11 alpha; retrieved from http://www.cloudcompare.org/). We also used the TLS data to perform individual tree segmentation and extract vegetation 110 parameters using R 3.5.3 (R Core Team, 2019), lidR (v3.1.1; Roussel) and rLiDAR (v0.1.1; Silva) packages. These steps are outlined in Fig. 3.

Ground, snow, and vegetation classification
We used the CANUPO method in CloudCompare to separate vegetation from ground and snow returns. This method includes training and classification. In the training step, we used 10,000 snow and vegetation samples to construct the classifier. We 115 trained the algorithm at 15 different scales to assign features related to each class and selected the 9 best combinations of scales (0.1m, 0.2m, 0.25m, 0.5m, 0.75m, 1m, 2m, 3m, 5m) to properly separate different classes. The combination of information from these scales helped the algorithm detect the dimension of each feature and assign snow and vegetation labels to the unclassified point clouds (Lague et al., 2013). We found that CANUPO misclassified snow data points near tree stems as vegetation, and thus we reclassified these points manually using the software TerraScan (Helsinki, Finland). 120

Snow depth estimation
To estimate under-canopy and open-area snow depths, we used the M3C2 algorithm (Lague et al., 2013) in CloudCompare. In this algorithm, for every single point in the ground point cloud we defined a cylinder with a range of different radii (projection scales) varying from 10 cm to 3 m and a length (height) of 3 m (see Lague et al., 2013, for details on these parameters). The orientation of the cylinder was along the normal vector of planes fitted on the ground points within a 10cm radius. We projected all points within the cylinder onto the cylinder axis, took the vertical distance between projected snow, and ground points as the snow depth estimation. Through iteration, we found a balance between including enough TLS points for subsequent analysis and the accuracy of the snow depths by using a 1 m projection scale. Our resulting snow depth measurement has a relative accuracy of approximately 2.5 cm based on the maximum standard deviation from M3C2. Utilizing these measurements, we compared snow depths under the canopy and in the open at each site. We also defined a transition zone of 130 10 m from the canopy to the open to identify any relevant differences within this zone.

Individual tree segmentation and vegetation metrics
We developed a canopy height model at 0.5 m resolution and identified tree tops to segment individual trees in the R package lidR. We detected a local maxima to identify tree tops using window sizes ranging from 1-3 m and minimum tree heights from 2-6m, depending upon the site. For areas with lower tree heights (0.5 -2m), we tiled the data that contained these trees and 135 segmented them in a similar approach. This allowed us to more accurately segment distinctly shorter and taller tree populations within sites, by adjusting segmentation parameters that worked better for those areas. Based on our preliminary analyses, we found that the Marker-controlled watershed segmentation algorithm was most accurate (compared to li2012, dalponte2016, and watershed, all available in the lidR package). We found that in cases where tall and short trees are close to each other, the algorithm could not detect shorter trees with large crown radii, and if a small crown radii is used, the branches far from the top 140 of the tree may be considered as an individual tree. We resolved this problem by tiling the las files and processing each separately and then combining the results. An example of the segmentation results from site F is shown in Fig. 4. A similar process was performed for all sites.
To define under the canopy and in the open, we first performed segmentation to identify individual trees. Under the canopy was defined by all snow depth points within the tree polygons. To define the open area, we merged individual tree polygons 145 that were less than 3 m from each other (patches of trees) and used the remaining areas as open. Site A was the only site dominated by shrubs ( Fig. 1, Table 1) and we considered the shrub area as open (we removed shrubs in the processing and kept the ground points below) at this site since the focus of our study was on tree-snow relationships.
We computed 22 vegetation metrics for each individual tree identified in the segmentation process (Tables 2, A1). We then used these metrics to predict snow depths at each site using a simple linear regression. 150

Influence of canopy edge on snow depth
We used the individual trees to predict snow depth at distances of 1 to 10 m away from the canopy edge. This represents how snow depth changes when we move from the edge of individual trees to the open within 10 m distance from the edge. We subsampled our data to only include trees that had good snow coverage (from TLS) under the canopy. This was determined based on an individual tree polygon containing a minimum of 10 snow pixels (pixel resolution of 0.5 m).

Tree height
We investigated the correlation of maximum tree height (at the individual tree level) with snow depth as a function of distance from the canopy. We used 1 m distance bins and calculated the correlation between the maximum tree height and mean snow depth at each bin.

Topography 160
Slope and aspect were derived for each site using a nearest neighbor method at 1m grid resolution. We evaluated snow depth changes related to plot-scale (1 m) slope and aspect variations. We did this in both under the canopy and open areas at each site.

Gap distribution
We explored whether any of our sites were suitable for understanding the role of forest gaps (i.e. shading, interception) on 165 snow depth distributions. While our study was not designed to analyze a range of gap distributions, the inherent forest density and distribution gradient that spanned our sites across Grand Mesa provided this opportunity. In particular, we sought to identify if sites had a dispersed tree pattern, such that the gaps were large enough to prevent canopy interception of snow, and thus accumulated deeper snow. Seyednasrollah and Kumar, (2014) used a relationship of tree height and gap radius for evaluation of net radiation. We derived a similar but simplified gap distribution approach (Equation 1). We calculated the 170 average distance of 10 nearest trees to each individual tree. This gave us a rough estimate of a gap size around each tree (D).
In the next step, we divided that average distance (D) by the average height (H) of those 10 nearest trees (D/H). This ultimately provides a ratio by which we can investigate the impact of shading from trees on gaps combined with gap size.
Equation (1), illustrates the gap distribution for an individual tree (j) where, & is the mean distance of the k closest trees to 175 tree j; & is the average height of k closest trees to tree j; k is the number of neighbors and '& and ℎ '& are the distance and height of tree i to tree j, respectively.
Secondly, we performed an average nearest neighbor analysis of the distribution of trees at each of the sites. In this analysis, we tested for tightly clustered trees in which gaps were minimal (clustered), randomly distributed trees where gaps could potentially lead to deeper snow accumulation (random), or dispersed trees where no particular pattern exists and thus gaps are 180 likely not prevalent (dispersed).

Directional analysis
We also investigated relationships between tree heights and snow depth based on direction. We did this using the 10 m transition zone (buffer) for each individual tree. We classified snow depths within each buffer in the four cardinal directions https://doi.org/10.5194/tc-2020-277 Preprint. Discussion started: 23 October 2020 c Author(s) 2020. CC BY 4.0 License. and computed the correlation between tree heights and mean snow depth per each direction. We compiled our tree height data 185 with all tree heights analyzed together, as well as subsets of data with tree heights < 10 m and > 10 m. We also performed a directional analysis with a Wilcoxon signed-rank test for comparing snow depth on the north and south sides (and east and west) for individual trees at each site.

Snow depths 190
Using our individual tree analysis, we found higher snow depths in open regions and lower snow depths in areas dominated with trees (see Table A2, Fig. B1). Snow depths were 12-28 % higher in the open than under canopy. Mean snow depth percent change between the 10m transition zone and under the canopy ranges from less than 1 % for sites A and K to a maximum of 7 % at site M. We found the lowest mean snow depths in our most westerly site (A), which is dominated by dense clusters of relatively rigid shrubs (Dasiphorafruticosa) and has the lowest tree cover of all sites. The standard deviation (SD) of snow 195 depths was similar between the transition zone and under the canopy for four sites (A, F, K, and O). We found a significantly lower SD of open area snow depths at three of the sites (M, N, O) compared with under-canopy and transition zones (Table   A2).

Influence of vegetation metrics on snow depth
Based on a linear regression, most sites had mid-to-high correlation between a specific vegetation metric (from Table 2) from 200 individual trees and snow depths, with foliage height diversity (FHD), which was the most influential vegetation metric at four of the six sites (Tables 1, 3). Figure 5 shows the distribution of FHD at each of the sites, with higher FHD demonstrating more evenly spaced foliar arrangement along an individual tree. Most of the sites had two distributions of FHD. At sites F, K, N, O, the FHD and snow depth had negative correlations ranging from 0.35-0.75; as the foliar arrangement was more evenly spaced along an individual tree, snow depth decreased. At site A, a high negative correlation (0.68) between percentage of returns 205 above 1 m (pzabov1) and snow depth occurred. FHD also had a relatively high negative correlation at this site (0.47) with snow depth.

Influence of slope and aspect on snow depth
We found that at site A, slope was influential on snow depths in the open. For example, slope explained 44 % of the variance of snow depth in the open at this site (Table 3)

Influence of canopy edge on snow depth 215
We found that snow depths increase with distance from the canopy edge into the open for the majority of individual trees (Fig.   7a). These correlations were typically above 0.60 and in many sites over 0.80. At some sites we found negative correlations between distance and snow depth ( Fig. 7b and Figs. B2-7 for individual trees at each site). For example, at site O, strong negative correlations (0.80) were apparent on the northwest side of the tree patches in the southeast portion of the study area. This is the area of site O where the north facing slope is likely the largest influence on snow depths. The high positive 220 correlations between canopy edge and snow depth occur in the north where snow depths are low (less than 1 m). Site A also had negative correlations in the north/north-western sampled region, and we propose this is likely due to northeast winds and deeper snow depths in the northeast facing slopes in the southern portion of the site. We also found that single trees not belonging to a patch or trees in small patches (e.g. <10 trees surrounding), had negative correlations between canopy edge and snow depth. 225

Influence of tree height and distance on snow depth
We found that the correlations of maximum tree height with snow depth are highest within 6 m, though several sites had the highest correlation at 3 to 4 m distance away from individual trees (Fig. B8). Site K had high correlation (0.56) at 1 m whereas site O had high correlation at 6 m (0.74). Overall sites M and N had relatively lower correlations (0.36 and 0.28, respectively) with distance than the other sites. These sites both had the highest percent tree cover, though site M had tall trees (mean height 230 of 21.6 m) and Site N had short trees (mean height of 10.5 m) (Table 1). Site A had a correlation of 0.6 at 4 m and site F had a correlation of 0.39 at 3 m.

Gap distribution
Our results show that site N has the largest median D/H ratio (0.74) compared to all other sites of <0.5 (Table 1). Site N is the only site with a dispersed tree pattern (Fig. B9) and thus the most likely site to experience lower interception, possibly resulting 235 in deeper snow.

Directional analysis
We found negative correlations between tree heights and snow depths based on direction at all sites (Table A3) and regardless of tree height for sites K, N, and O (Tables A4-5). We found positive correlations at site F in the south direction and site M in the south and east directions for trees < 10m and in the west and south directions for site A for trees > 10 m (Tables A4-5). 240 We found snow depths were different between the north and south sides of trees at sites A, K, and O but not for any other sites or directions (Table A6).

Discussion
We observed several interesting relationships between vegetation canopies, topography, wind and snow depths across our sites. As expected, snow depths were deeper in the open compared to under canopy. However, describing the relationships 245 between vegetation and snow is complicated by the structure, distribution (pattern), and type of vegetation. The relationship is further convoluted by local topography and wind speed/direction. For example, we found that slope, aspect, and wind (rather than vegetation) might control snow depths at local scales at two of the sites, A and O. This is not surprising, as site A was dominated by 0.4-0.6 m tall shrubs and wind exposed, and site O had a relatively low tree canopy cover. While site A had the lowest tree canopy cover in our dataset, we only sampled the edge of a much larger patch of trees (based on field observations). 250 High snow depths were found in open northeast facing slopes at site A where wind deposits snow. Thus, local topographic interactions with wind have a major influence on snow accumulation, especially when we do not consider the much larger landscape controls. While sites A and O have slopes within the overall range of all of our sites, the combination of local slope and aspect for site O, appear to be driving factors in snow depths. In fact, site O has the highest mean snow depth (1.44 m), likely due to these local site conditions. 255 When our analyses were confined to under the canopy of individual trees, we generally found a significant relationship between the vertical spatial arrangement of the foliage (based on the FHD) and snow depth, but this relationship did not hold across all sites. For example, FHD explained less than 20 % of the variance of snow depth at site M. This site has the highest mean tree heights of the study. Taken together with a negative correlation between Zq90 and snow depth, the height of the tree likely had higher control on snow depths than on the particular foliar arrangement of the trees at this site. Overall, sites M and N had 260 the lowest correlations between vegetation metrics and snow depth. This may be due to the vegetation patterns at these sites (under the canopy slope and aspect have no effect at both sites M and N (Table 3, Fig. 6a)). Site M is a relatively open area with mature Engelmann spruce and subalpine fir trees in the SW and NE areas of our site. Subalpine fir trees are generally more slender than Engelmann spruce, and thus their shape may not be as influential on accumulation of under canopy snow depths. Site N has the highest percent cover and the smallest trees (mean tree height 10.5 m, SD of 2.62 m, Table 1). This 265 second growth canopy is the only site dominated by lodgepole pine, which are also slender. While the mean FHD is similar to the other sites, the spatial arrangement of the trees may have a larger control on snow depths. This is the only site with trees in a dispersed pattern in which the size of the gaps likely prevents snow interception, and thus provides an opportunity for snow accumulation. In fact, site N had the second highest mean snow depth under the canopy (1.38 m, compared with 1.44 m at site O, Table A2, Fig. 6a). Testing for a dispersed tree pattern could be beneficial to future studies, especially because 270 previous research (e.g. Ning Sun et al., 2018) found gap size to be a control on snowmelt timing; however, our study was during the accumulation phase so we cannot draw similar conclusions for site N.
Similar to previous work (e.g. Moeser et al., 2015a;Mazzotti et al., 2019), we found that measures of distance from canopy are important correlates with snow depth. However, most of the previous work was performed with airborne lidar across larger spatial extents and lower vertical resolutions in the canopy. Our data are best suited to fine scale interactions between individual trees, or clusters of trees, and under the canopy or surrounding snow depths. Ultimately, understanding controls at the scale in which TLS provides is complementary to airborne lidar, which can help test larger scale features, such as gap area across space.
Our canopy edge analyses generally found that as the distance increases from the canopy edge, snow depths also increase. We also found that small patches or individual trees do not influence snow depths significantly with distance. In addition, on the 280 north side of site A, wind and slope may influence the relationship between canopy edge and snow depths (where we found negative correlations, Fig. B2).
We did not find high snow depth accumulation or variability within a transition zone similar to the findings of Broxton et al. (2015). While their study included similar tree species, wind speeds and elevation, their spatial scale of analysis was larger with the use of airborne lidar. Interestingly, correlations between maximum tree height and snow depth in Grand Mesa occurs 285 within 6 m of an individual tree at the sites, reinforcing local scale controls. Further, the negative correlations between snow depth and all cardinal directions from individual trees indicate that we should expect shallower snow within the 10 m transition zone from taller trees (Tables A3-5). In other words, two adjacent trees with different heights affect snow depth differently in any one direction and shorter trees keep deeper snow in all directions. We expected the opposite to occur i.e. taller trees should create larger shadows and provide more shading/sheltering. As our snow-on datasets are from the accumulation season, we 290 may not see shading effects of taller trees in the transition zone; negligible melt had occurred at the time of these surveys. If our datasets extended throughout the season, we might expect these relationships to change. Note that due to sampling extent, our transition zone analysis was performed at 1 m increments instead of at multiples of mean tree height as in previous literature (e.g. Currier and Lundquist, 2018).
Following previous studies that showed a directional relationship with snow depths (e.g. Mazzotti et al., 2019;Currier and 295 Lundquist, 2018), we found significantly different snow depths between the north and south sides of trees at site A, K, and O, but not other sites. This may be due to the local topography and wind at sites A and O. Site K had a maximum correlation between tree height and snow depth at 1m, and thus there may be local controls with distance and direction at site K.
Additionally, previous lidar-based canopy snow interaction studies (Trujillo et al., 2007(Trujillo et al., , 2009Deems et al., 2006) relied on more simple canopy models using maximum height. Our results show that in nearly all situations, structural information 300 contained in denser lidar point clouds have much more predictive capability.

Conclusions
Our study indicates that even with fine scale, individual tree observations from TLS, vegetation structural metrics are not enough to describe snow depth during the accumulation season. Local scale topography and wind should also be considered.
While our sites were not designed solely for intercomparison, we found notable trends in our site comparisons. The vertical 305 arrangement of foliage (e.g. FHD) of individual trees and tree height influences under canopy snow depths, and in some cases, quite strongly. Further studies should be designed to test this within and between species. For example, our sites were primarily https://doi.org/10.5194/tc-2020-277 Preprint. Discussion started: 23 October 2020 c Author(s) 2020. CC BY 4.0 License.
Engleman spruce, subalpine fir, and lodgepole pine, all of which have different canopy structural shapes. Further studies targeting samples of each of these with different foliar arrangements and heights should be undertaken to fully understand the implications of FHD and tree heights on snow depths at local scales. 310 We also found that topography had greater control than vegetation at sites where slopes favored wind conditions for increasing snow depths, or where vegetation presence was minimal. While the latter may be obvious, increased observations with varying vegetation cover, wind, and topography should be considered with TLS.
This study highlights the complementary nature of TLS observations to airborne lidar, where TLS can provide fine scale observations within the canopy and relationships with under the canopy snow depth. Data from TLS can also be used to validate 315 airborne lidar (e.g. Currier et al., 2019), and further studies should investigate how vegetation metrics such as FHD compare between TLS and airborne lidar in these snow-dominated forest ecosystems. Further, along with airborne lidar, TLS provides a complementary dataset for upscaling to similar types of vegetation structure and topography observed from satellites such as ICESat-2, or missions such as GEDI. Importantly, results from this study and others with TLS and airborne lidar for forestsnow observations can also be the foundation for the 2017 Decadal Survey designated observable, Surface Topography 320 Vegetation study (National Academy of Sciences, Engineering, Medicine, 2018).