Modelling and Comparing Shading Effects of 3D Tree Structures with Virtual Leaves

: Reduced solar radiation brought about by trees on agricultural land can both positively and negatively affect crop growth. For a better understanding of this issue, we aim for an improved simulation of the shade cast by trees in agroforestry systems and a precise estimation of insolation reduction. We present a leaf creation algorithm to generate realistic leaves to be placed upon quantitative structure models (QSMs) of real trees. Further, we couple it with an enhanced approach of a 3D model capable of quantifying shading effects of a tree, at a high temporal and spatial resolution. Hence, 3D data derived from wild cherry trees ( Prunus avium L.) generated by terrestrial laser scanner technology formed a basis for the tree reconstruction, and served as leaf-off mode. Two leaf-on modes were simulated: realistic leaves, fed with leaf data from wild cherry trees; and ellipsoidal leaves, having ellipsoids as leaf-replacement. For comparison, we assessed the shading effects using hemispherical photography as an alternative method. Results showed that insolation reduction was higher using realistic leaves, and that the shaded area was greater in size than with the ellipsoidal leaves or leaf-off conditions. All shading effects were similarly distributed on the ground, with the exception of those derived through hemispherical photography, which were greater in size, but with less insolation reduction than realistic leaves. The main achievements of this study are: the enhancement of the leaf-on mode for QSMs with realistic leaves, the updates of the shadow model, and the comparison of shading effects. We provide evidence that the inclusion of realistic leaves with precise 3D data might be fundamental to accurately model the shading effects of trees. solar exposure. ground plane. The six leaf-geometry points were presented in XYZ coordinates relative to the input QSM. The leaf edge coordinates output dataset was used as input within the shadow model, jointly with the tree-cylinder model.


Introduction
Established plants act as physical ecosystem engineers that directly or indirectly control the availability of resources to other organisms [1], for example, trees casting shade over food crops. Although shade is frequently treated simply in its literal sense as low light conditions [2], in a more general and ecological sense, it involves a wider range of environmental factors (altered atmospheric and substrate conditions) with various effects on plants (e.g., biotic interactions). As plant growth is driven by these factors and effects, i.e., resources (light, water, and nutrients) and conditions (temperature and shading), aside from specific plant traits, a change in quantity or quality of any of these factors triggers a change in plant growth. Thus, increasing the knowledge on resources availability, control, and modulation effects of the tree component, such as with a precise assessment of distribution of incoming solar radiation around a tree canopy, can aid a better understanding of resource dynamics, crop performance, and organisms' interrelations.
Ecological habitats and production systems can be described in terms of their structural components. Structural diversity is unarguably poor in monocropping agricultural systems. In any productive system, an increase in structural diversity generally correlates with an increase in the complexity of the interdependencies of the individual elements. The latter Table 1. Useful concepts, terminologies and definitions used within the framework of this study.

Solar Radiation Energy
The radiant energy from the sun that drives photosynthesis, the energy source of tree growth, absorbed by leaves [10]. It is the short-wave electromagnetic radiation [11].
Insolation Amount of solar radiation received on a given surface in a given time period [10], or simply the power per unit area received from the sun. Synonyms: solar irradiation/irradiance, solar exposure.

Reduced Insolation
Amount of solar radiation received on a given surface in a given time period after discounting the shading effect of obstacles (e.g., trees and leaves). This is assessed by comparing the insolation of one-unit area in comparison to another, under full radiation conditions.

Light
Sometimes is referred as a synonym for solar radiation energy. Other times, referred as to visible light, which is the radiation within the portion of the electromagnetic spectrum that can be perceived by the human eye. Moreover, a synonym for Photosynthetically Active Radiation (PAR).

Shadow
Absence of light; the dark area created by blocking the visible light from a light source (i.e., the sun)

Shading effect
Area subjected to reduced insolation. In this paper, often described by the "shaded area" (the surface not under full radiation conditions) and "insolation reduction" (the energy reduction in a given time period).
The measurement of light has a long tradition and is of high interest in environmental physics as well as in plant ecology [11]. To date, several techniques have been used to assess the share of incoming solar radiation under plant canopies. For the direct measurement of the light conditions different types of pyranometers were applied, while for indirect measurements devices such as the LAI-2000 Plant Canopy Analyser can be used to gather information about the leaf area index (LAI) of trees and forests [12].
The use of hemispherical photography is probably one of the most widely used indirect methods for gaining information on light conditions within wooded landscapes utilising a fisheye lens for whole sky photographs, a technique proposed nearly 100 years ago [13]. Evans and Coombe [14] were some of the first researchers applying it to a forest stand for measuring the direct and indirect radiation. Several follow-up studies appeared during the last decades and its applicability within other fields has also been tested. Comparative studies between hemispherical photography and other light sensors have been conducted by, e.g., [15][16][17][18]. This fact together with the ease of handling, the low price and the advantage of capturing all hemispherical directions simultaneously, are the Remote Sens. 2021, 13, 532 3 of 22 main reasons for the broad application of this method. Light transmission measurements have already been undertaken using hemispherical photography in AFS [19].
A more advanced, and arguably more precise approach is the modelling of the light or radiation regime around trees based on three-dimensional (3D) structures. Dupraz and Liagre [20] replaced the crown with an ellipsoid for modelling the light distribution in AFS. This approach was further refined by Talbot and Dupraz [21] including additional parameters such as a leaf clumping coefficient and the shading effect of stem and branches. Alternative 3D models have also been developed by other groups, but they are all constrained by the scarce richness of details on the assessment [22][23][24][25][26]. Schmidt et al. [26] used a block with dimensions of 30 × 30 m and height of 20 m as a replacement of a forest edge to simulate shading in the field-to-forest transition zones of agricultural fields and estimate crop yields. Though an effort to represent the shading of a cluster of trees, the approach is lacking details of actual tree structures.
During the last decade the use of terrestrial laser scanners (TLS) has opened up new possibilities for the 3D modelling of the shade cast by trees [27][28][29][30][31][32]. While the other abovementioned 3D models rely on simplified 3D structures, TLS is supplying users with realistic and detailed 3D data of trees. Likewise, TLS approaches are superior to several conventional canopy analyses with optical methods, since they enable detailed coverage of vegetation surface with properly designed multi-scan sampling [33].
On the modelling of shade cast by trees, a recurrent strategy has been scanning trees with a TLS during the vegetation period and processing voxel data with ray tracing algorithms. At the attempt to increase the accuracy of voxel-based models arose difficulties in handling occlusion within the 3D data, due to insufficient information of the inner crown structure and vegetation clumping [28,31]. On the use of the voxel-based approach, an auspicious technique was presented by van der Zande et al. [27], where leaf area density could be directly derived by transforming TLS data into small leaf-sized voxels.
An alternative for modelling the shade cast by trees is by scanning deciduous trees under leafless conditions (e.g., during winter dormancy), reconstructing tree architecture with cylinder-based algorithms and adding leaf substitutes in a computer environment. Such approaches allow for the simulation of the shading effect of trees in a leaf-on mode. This technique was chosen by Rosskopf et al. [32] using a tree-cylinder model created with the software "SimpleTree" [34,35] coupled with ellipsoids as replacement for leaves. Other models reconstructing tree architecture with cylinders as geometric primitives exist [36][37][38][39]. They are often named quantitative structure models (QSMs) and contain accurate topological, geometrical and volumetric tree properties. "TreeQSM" [37,40,41] has been a popular method for modelling tree structures. Meanwhile, algorithms to insert leaves or needles as virtual structures into QSMs (produced with TreeQSM), and according to an arbitrary distribution exist [42]. In addition, the use of QSMs could offer the option to include TLS information concerning tree growth [43][44][45] and integrate it within the modelling approach.
We present a leaf creation algorithm (LCA) aiming at an improved simulation of the shade cast by individual trees for the estimation of its shading effect within a defined time-period. In this study, to virtually create leaves, the LCA is fed with species-specific leaf parameters derived for wild cherry trees (Prunus avium L.) growing in an AFS within southern Germany. The novelty of this approach is to use QSMs augmented with realistic leaves to compute shadow projections. We apply the new LCA to six trees and use an enhanced 3D modelling approach (named "shadow model") capable of quantifying tree shading effect at high temporal and spatial resolution. Additionally, we simulate the shading effects for other leaf modes and compare the results with those achieved with hemispherical photography, as an alternative approach.
Our hypothesis is that a more realistic leaf-on mode (created by applying the LCA) does affect the simulated insolation on the surroundings of single-trees. Furthermore, our method is expected to be more accurate than concurrent methods for estimating spatial and temporal distribution of the shade cast by trees, and the correspondent insolation on the ground. The future validation of our shadow model will require the conduction of Remote Sens. 2021, 13, 532 4 of 22 field assessment and tuning of model parameters. We showcase tree shading effects and gained insights on the accuracy of our method. The evidence gathered in this paper can be strategically used for planning and designing integrated tree-crop agricultural systems, in any landscape where light distribution is of concern, and for future development of specific and/or generalist methods on this study field.

Site Description and Trees
We selected six wild cherry trees (Prunus avium L.) growing in AFS in the vicinity of Eichstetten (48 • 06 26.3 N 7 • 43 53.6 E, 193 m a.s.l.), in the Kaiserstuhl region of southwest Germany. The trees have deciduous character, are growing free of competition of other trees and are scattered in agricultural fields cultivated with different crops; it is likely that the individual trees have received various silvicultural treatments during their lifespan, thus altering the growth, development, and actual structure. At the time of measurements (April 2019), the trees had heights ranging between 6 and 9 m and diameters at breast height (DBH) of 18 to 27 cm. The slope of the terrain at the scanned tree locations on the tree locations was considered flat (0-2%).
Within the Kaiserstuhl region, average annual direct normal radiation is 4158 MJ m −2 (11.39 MJ m −2 day −1 ) and annual global horizontal irradiation is 4402 MJ m −2 , or 12.06 MJ m −2 day −1 [46]. The local climate is temperate and mild with a mean annual air temperature of 11.2 • C and a mean annual precipitation sum of 710 mm ensuring suitable growing conditions for trees [47]. The local soils consisting of periglacial loess which show a high-water permeability and low water storage capacity. As a consequence, the upper soil layers can be dry in summers experiencing low precipitation.

Overview of the Modelling Steps
The development of the shadow model can be divided in four parts (Parts A-D, see Figure 1):

•
Part A: Collection of 3D TLS data, processing and tree segmentation; • Part B: Creation of quantitative structure models (QSMs), representative of tree structures; • Part C: Introduction of the leaf creation algorithm (LCA) and its application; • Part D: Application of the shadow model using solar radiation data.
Our hypothesis is that a more realistic leaf-on mode (created by applying the LCA) does affect the simulated insolation on the surroundings of single-trees. Furthermore, our method is expected to be more accurate than concurrent methods for estimating spatial and temporal distribution of the shade cast by trees, and the correspondent insolation on the ground. The future validation of our shadow model will require the conduction of field assessment and tuning of model parameters. We showcase tree shading effects and gained insights on the accuracy of our method. The evidence gathered in this paper can be strategically used for planning and designing integrated tree-crop agricultural systems, in any landscape where light distribution is of concern, and for future development of specific and/or generalist methods on this study field.

Site Description and Trees
We selected six wild cherry trees (Prunus avium L.) growing in AFS in the vicinity of Eichstetten (48°06'26.3"N 7°43'53.6"E, 193 m a.s.l.), in the Kaiserstuhl region of southwest Germany. The trees have deciduous character, are growing free of competition of other trees and are scattered in agricultural fields cultivated with different crops; it is likely that the individual trees have received various silvicultural treatments during their lifespan, thus altering the growth, development, and actual structure. At the time of measurements (April 2019), the trees had heights ranging between 6 and 9 m and diameters at breast height (DBH) of 18 to 27 cm. The slope of the terrain at the scanned tree locations on the tree locations was considered flat (0-2%).
Within the Kaiserstuhl region, average annual direct normal radiation is 4,158 MJ m -2 (11.39 MJ m −2 day −1 ) and annual global horizontal irradiation is 4402 MJ m −2 , or 12.06 MJ m -2 day −1 [46]. The local climate is temperate and mild with a mean annual air temperature of 11.2 °C and a mean annual precipitation sum of 710 mm ensuring suitable growing conditions for trees [47]. The local soils consisting of periglacial loess which show a highwater permeability and low water storage capacity. As a consequence, the upper soil layers can be dry in summers experiencing low precipitation.

Overview of the Modelling Steps
The development of the shadow model can be divided in four parts (Parts A-D, see Figure 1): • Part A: Collection of 3D TLS data, processing and tree segmentation; • Part B: Creation of quantitative structure models (QSMs), representative of tree structures; • Part C: Introduction of the leaf creation algorithm (LCA) and its application; • Part D: Application of the shadow model using solar radiation data.  To update our shadow model, the six wild cherry trees were scanned with a TLS Faro Focus 3D S120 (FARO Technologies, Inc.; Lake Mary, FL, USA) in early spring (April 2019), under windless and leaf-off conditions, for a better visibility of the woody compartments of the trees [48]. The scanned tree surface is representative of the last growing season. We used a multiple scan approach with four scan positions around each tree, 90 • degrees apart (azimuth angle) and at a distance of 10 m from the base of the tree. TLS device sampling parameters were set to 1 4 for "resolution" and 4x for "quality". A minimum of five reflective targets were set out around each tree to merge the multiple scans. The scans were registered using the software Faro Scene 6.2.4.30 (FARO Technologies, Inc.; Lake Mary, FL, USA) with a reported mean point error below 4 mm. Duplicated points were eliminated with a "low" search radius. In the following step, tree point clouds were manually extracted; this step was straightforward as trees had no direct neighbours or hindering obstacles. Lastly, we filtered the point clouds to remove noise points and outliers in CloudCompare v2.10.2 [49].
Part B: Creation of Quantitative Structure Models (QSMs) The processed and filtered single-tree point clouds served as basis for reconstructing the architecture of trees using the MATLAB implementation of TreeQSM version 2.3 [37]. In the 3D reconstruction of QSMs, trees are modelled as a hierarchical collection of cylinders (e.g., geometric primitives) fitted to local details of the single-tree point cloud.
To optimise QSMs, we tested 32 combinations of key model input parameters (among them, cover patch diameter d, and relative cylinder length lcyl) and produced 15 models for each possible combination of inputs to select the best model, define the optimised parameters, and to assess the robustness of the reconstruction method and uncertainty of results [41,50]. The mean point-cylinder-model-distance was used as suitable metric for the optimisation [51]. All QSMs were reconstructed with the same optimised input parameters: in the first cover set d was 15 cm; for the second cover, the minimum and maximum d were 1 cm and 5 cm, respectively; lcyl was 3.5. The uncertainty of QSM-parameters stayed below 10% (CV; coefficient of variation) for each of the six trees. The regular, cylindrical-like, stem base of the trees did not require additional triangulation.
The QSM-derived tree parameters for the scanned trees are presented in Table 2. Trees were similar in terms of DBH and height, while total tree and branch volume, cylinder and branch counts, revealed the more complex structures hold by trees Pa_2, Pa_4 and Pa_5. We identified Pa_6 as the least elaborate tree structure. In addition, the count of branches and accumulated branch length are important as realistic leaves are placed along the branch-cylinders by the LCA (details in Section 2.2 Part C). The QSMs play a central role in our shadow model, as they provide the basic topological structure to which virtually created leaves are attached. The processing steps are illustrated in Figure 2 (Figure 2f-h). The full description of the new leaf-on mode and the details of the enhanced shadow algorithm are presented in Section 2.2 Part C and Part D, respectively; the use of these 3D structures for simulating the shading effects is described in Section 2.3. A video-visualisation of the 3D structure of tree Pa_5, the QSM and leaf modes is found in Supplementary File S1.
Remote Sens. 2020, 12, x FOR PEER REVIEW 6 of 22 four models as a worm's eye view ( Figure 2e); bird's eye views of the QSM, QSM with ellipsoidal leaf-replacement and with realistic leaves, respectively (Figure 2f-h). The full description of the new leaf-on mode and the details of the enhanced shadow algorithm are presented in Section 2.2 Part C and Part D, respectively; the use of these 3D structures for simulating the shading effects is described in Section 2.3. A video-visualisation of the 3D structure of tree Pa_5, the QSM and leaf modes is found in Supplementary File S1. In July 2018, leaves were collected from felled wild cherry trees for generating data to serve as basis for the virtual creation of leaves on overlaying the QSMs in a computer environment. The realistic leaves were an attempt to represent carefully the leaf topology of wild cherry trees, and were used as the RL leaf-on mode for the shadow simulations.
To gain insights in leaf patterns (i.e., shape distribution, size distribution, and spatial distribution within-trees), we sampled 58 branches from eight cherry trees growing in analogous conditions as the six laser-scanned trees. Branches were randomly selected from different tree heights, starting from the crown base to the tree top, and in multiple cardinal directions. We cut these branches at the branch collar, manually removed their leaves, and stored them in separate bins, which were weighed. A subset of 20 branches has been photographed with a scale to determine branch parameters, such as branch collar diameter, branch length and branch specific mass. From the leaf bins, 630 fresh leaves were randomly drawn and scanned with a ScanMaker 9800XL plus (Microtek International, Inc.; Hsinchu City, Taiwan) in a resolution of 600 dpi, to gather information of the leaf dimensions with the software ImageProPlus Version 7.0.1 (Media Cybernetics, Inc.; Rockville, MD, USA). Leaf length, leaf width and leaf area (one-sided area of the leaf blade) were measured. We divided the leaf area observations into five sets of equal length, every 15.1 cm² apart from minimum to maximum values ( Figure S1). In addition, we determined the leaf area distribution, defined five leaf size classes by taking the mean leaf area of each range (extra small, small, medium, large, extra-large) and obtained their proportions, as shown in Figure 3. Lastly, we computed hexagonal leaf-like shapes with the geometric parameters derived from leaf measurements. In July 2018, leaves were collected from felled wild cherry trees for generating data to serve as basis for the virtual creation of leaves on overlaying the QSMs in a computer environment. The realistic leaves were an attempt to represent carefully the leaf topology of wild cherry trees, and were used as the RL leaf-on mode for the shadow simulations.
To gain insights in leaf patterns (i.e., shape distribution, size distribution, and spatial distribution within-trees), we sampled 58 branches from eight cherry trees growing in analogous conditions as the six laser-scanned trees. Branches were randomly selected from different tree heights, starting from the crown base to the tree top, and in multiple cardinal directions. We cut these branches at the branch collar, manually removed their leaves, and stored them in separate bins, which were weighed. A subset of 20 branches has been photographed with a scale to determine branch parameters, such as branch collar diameter, branch length and branch specific mass. From the leaf bins, 630 fresh leaves were randomly drawn and scanned with a ScanMaker 9800XL plus (Microtek International, Inc.; Hsinchu City, Taiwan) in a resolution of 600 dpi, to gather information of the leaf dimensions with the software ImageProPlus Version 7.0.1 (Media Cybernetics, Inc.; Rockville, MD, USA). Leaf length, leaf width and leaf area (one-sided area of the leaf blade) were measured. We divided the leaf area observations into five sets of equal length, every 15.1 cm 2 apart from minimum to maximum values ( Figure S1). In addition, we determined the leaf area distribution, defined five leaf size classes by taking the mean leaf area of each range (extra small, small, medium, large, extra-large) and obtained their proportions, as shown in Figure 3. Lastly, we computed hexagonal leaf-like shapes with the geometric parameters derived from leaf measurements. edges were placed equidistant within the minimum diameter. The area of the convex hull of the six points was calculated as leaf area. Finally, the geometric edges were adjusted by a ratio, so that the mean leaf, representative of a size class, keeps the corresponding average leaf area (Table S1). The resulting irregular hexagons can be seen in Figure 3. The final realistic leaves representatives of the five classes were defined by six geometric points (Table S2). Leaf geometry, area, length and width, as well as insights on the distribution and proportion of leaves, were used to feed the leaf creation algorithm (LCA). The LCA was developed to overlay realistic leaves on the top of a QSM generated with TreeQSM. The distribution of leaves follows the pre-set information on leaf geometry and leaf size classes proportion, the 3D tree-cylinder model properties and the user-defined input parameter "leaf spacing". On a tree basis, for a range of virtually created leaves with increasing leaf spacing (0.5, 1.0, 1.5, 2.0, 2.5, and 4.0 cm), preliminary results [52] revealed the leaf spacing of 2 cm being the most suitable to match the estimated total leaf area in first order branches (mean percentual overestimation of 7%; Figures S2 and S3). The algorithm outputs consist of two data frames: the leaf edge coordinates dataset (six points defining leaf geometry and orientation), in a coordinate system respective to the input QSM, and; an informational table containing leaf attributes (leaf number, size class, area, origin).
The steps imbedded in the LCA are described in Algorithm 1. First order branches assumed a central role in the leaf generation process, while stem-cylinders were prohibited from receiving leaves. According to the insights obtained on leaf distribution and from the starting position of a first order branch, branch-cylinders were foliated if distanced in magnitude equal or greater to 8.47 % of the length of the respective first order branch ( Figure S4). Selected branch-cylinders were given leaf positions according to leaf spacing and an alternate distichous phyllotaxy (right and left arrangement). A virtual petiole (the stalk that attaches the leaf blade to the branches) of 2 cm was set from the initial leaf position. Furthermore, a random leaf area class was selected and the leaf-bottom point was matched with the extended leaf position. The other five leaf-geometry points were calculated by keeping leaf azimuth angle perpendicular to the cylinder direction (left or right direction) and leaf blade parallel to the ground plane. The six leaf-geometry points were presented in XYZ coordinates relative to the input QSM. The leaf edge coordinates output dataset was used as input within the shadow model, jointly with the tree-cylinder model. The realistic leaves were initially designed as regular six-sided polygons. We calculated leaf length and leaf width for the leaf classes by averaging minimum and maximum values for each leaf in the scanning procedure. Hence, the average lower leaf edges were at 38% of the leaf length, while the average upper leaf edges at 68%. Both lower and upper edges were placed equidistant within the minimum diameter. The area of the convex hull of the six points was calculated as leaf area. Finally, the geometric edges were adjusted by a ratio, so that the mean leaf, representative of a size class, keeps the corresponding average leaf area (Table S1). The resulting irregular hexagons can be seen in Figure 3. The final realistic leaves representatives of the five classes were defined by six geometric points (Table S2). Leaf geometry, area, length and width, as well as insights on the distribution and proportion of leaves, were used to feed the leaf creation algorithm (LCA).
The LCA was developed to overlay realistic leaves on the top of a QSM generated with TreeQSM. The distribution of leaves follows the pre-set information on leaf geometry and leaf size classes proportion, the 3D tree-cylinder model properties and the user-defined input parameter "leaf spacing". On a tree basis, for a range of virtually created leaves with increasing leaf spacing (0.5, 1.0, 1.5, 2.0, 2.5, and 4.0 cm), preliminary results [52] revealed the leaf spacing of 2 cm being the most suitable to match the estimated total leaf area in first order branches (mean percentual overestimation of 7%; Figures S2 and S3). The algorithm outputs consist of two data frames: the leaf edge coordinates dataset (six points defining leaf geometry and orientation), in a coordinate system respective to the input QSM, and; an informational table containing leaf attributes (leaf number, size class, area, origin).
The steps imbedded in the LCA are described in Algorithm 1. First order branches assumed a central role in the leaf generation process, while stem-cylinders were prohibited from receiving leaves. According to the insights obtained on leaf distribution and from the starting position of a first order branch, branch-cylinders were foliated if distanced in magnitude equal or greater to 8.47 % of the length of the respective first order branch ( Figure S4). Selected branch-cylinders were given leaf positions according to leaf spacing and an alternate distichous phyllotaxy (right and left arrangement). A virtual petiole (the stalk that attaches the leaf blade to the branches) of 2 cm was set from the initial leaf position. Furthermore, a random leaf area class was selected and the leaf-bottom point was matched with the extended leaf position. The other five leaf-geometry points were calculated by keeping leaf azimuth angle perpendicular to the cylinder direction (left or right direction) and leaf blade parallel to the ground plane. The six leaf-geometry points were presented in XYZ coordinates relative to the input QSM. The leaf edge coordinates output dataset was used as input within the shadow model, jointly with the tree-cylinder model. define branch section to be foliated (e.g., the first 8.47% of the length of a first order branch have no leaves). 5: for each branch cylinder to be foliated do 6: Establish positions for leaves along the directional axis between cylinder start and end, according to the parameter leaf spacing, and evenly distribute them, alternating between left and right side.

7:
Expand the preliminary leaf position with 2 cm (adding a virtual petiole), perpendicular to cylinder direction and horizontal to the ground, respecting the right or left orientation (+90 • or −90 • from cylinder direction) 8: for each leaf position do

9:
Randomly select a leaf-size class and match the lower leaf-geometry point with the established leaf position.

10:
Propagate the other five leaf-geometry points by keeping leaf oriented perpendicular to cylinder direction and horizontal to the ground. 11: end for (step 8) 12: end for (step 5) 13: end for (step 3) 14: return leaf edge coordinates dataframe and leaf attributes table A descriptive summary of the leaves created for the shadow simulations (presented in Section 2.3) is found within Table 3. The total leaf area correlates well with branch volume and other tree properties (i.e., cumulative branch length), and varied considerably between trees. Likewise, total leaf count ranged from 12,539 (Pa_6) to 35,194 (Pa_5). The proportions of leaves found for each leaf-size class confirmed our LCA is not biased. Table 3. Total leaf area, leaf count and proportion of realistic leaves in each size class. Leaves created with the Leaf Creation Algorithm (LCA), with leaf spacing of 2 cm, for each wild cherry tree (Pa_1 to Pa_6).

Tree
Total Leaf Area (m 2 ) Leaf Count We refined the shadow model proposed by [32] in two major aspects: (1) we integrated the output of the LCA to simulate the shading effects of trees (QSMs with realistic leaves); and (2) we expanded the functionalities of the algorithm to be compatible with the dataframe structure and building-logic of TreeQSM (previously based on SimpleTree), including the module related to the creation of ellipsoidal leaves (EL). The updated shadow model includes the realistic leaves following the same principles of the initial method, where cylinders and ellipsoids (as leaf-replacements) are taken as structures to project shadows on the ground under a given sun position.

Proportion of Leaves per Size Classes
For simulations of the EL leaf-on mode, we modified the cylinder-radius threshold for ellipsoids, proposed by Rosskopf et al. [32], from 0.5 to 1.0 cm (algorithm parameter defining where leaf-replacements are to be built upon) and adjusted the length of the longer/major ellipsoids axis, while the minor axis was set to 5.0 cm. Apart of the size limitation, ellipsoids were created for each set of cylinders in a branch order, where the Remote Sens. 2021, 13, 532 9 of 22 threshold applied. At the given level of implementation, a change in size to simulate leaf growth was not needed as the shade cast was modelled for the leaf attributes given in the month of July 2019.
For modelling the insolation on the ground plane, we used measured data of solar irradiance (global radiation and diffuse radiation) provided by the German Meteorological Service (DWD) [53] and chose the nearest meteorological station located in Freiburg (48 • 01 12.0 N 7 • 49 48.0 E, 236.5 m a.s.l.), about 15 km from the location of the scanned trees. The solar radiation data for July 2019, with a temporal resolution of 10 min, was integrated within the shadow model.
In summary, the updated shadow model was implemented in the open source language R, version 3.5.3 [54], and is based on key functions of the R packages "sp" [55,56] and "insol" [57]. In addition, we utilised functions of the package "rgl" [58] for visualisation purposes.

Shadow Simulations
The shading effects of trees were modelled for the period of 1 July to 31 July 2019, with three leaf modes:

•
QSMs with leaves created with the LCA (realistic leaves, RL); • QSMs with ellipsoids as leaf-replacements (ellipsoidal leaves; EL); • QSMs without leaves (no leaves; NoL): we were also interested in the effect of leaves on the shading in relation to a tree outside the vegetation period, under leaf-off conditions.
To quantify the loss of insolation due to the physical obstacles (the tree structure and leaves), we simulated the insolation for a raster grid size of one third of a hectare (70 m × 50 m), with cell size of 100 cm 2 (10 cm × 10 cm), and centralised on the tree's position. Insolation was determined for each raster cell at a 10-min interval. Shaded grid cells received the actual diffuse radiation value, while global radiation values were attributed to cells under full light conditions.
In order to compare shadow model results, we simulated insolation with an alternative established method named hemispherical photography (HP). We took photos of the same wild cherry trees in leaf-on conditions in July 2018 using a Canon EOS 700D camera (Canon Inc.; Tokyo, Japan) with a Sigma 4.

Comparisons and Analysis of Shading Effects
We explored the heterogeneity of the modelled shading effects in terms of monthly insolation, and compared it to full radiation conditions to estimate the shading effect: the total shaded area, and the absolute and relative insolation reduction.
The shaded area was defined as the area sum of grid cell receiving less than 98% of the maximum possible insolation (cells with insolation reduction), whereas the remaining cells were defined under full light conditions (no decrease in insolation). On each shaded area, total insolation reduction was calculated as the sum of differences between the maximum radiant energy potentially available on a grid cell and the actual radiant energy found on it. Furthermore, we investigated the shaded area proportion under different shade intensities (in percent), as defined: where x1 = sum of 10-min incoming radiation energy per grid cell over simulation period (31 days), under possible shaded conditions (insolation <98% of the maximum possible insolation); x2 = sum of 10-min incoming radiation energy per grid cell over simulation period (31 days), under full light conditions (insolation ≥98% of the maximum possible insolation).
Shade intensity values were split in classes of relative insolation reduction, as the intensity of shading have varied effects on the establishment and productivity of agricultural crops [9]. The five shade intensity classes were 2-5%; 5-10%; 10-20%; 20-30%; and >30% insolation reduction. We calculated the results for each tree under each of the four treatments (RL, EL, NoL, and HP). Lastly, we investigated the spatial heterogeneity for pairs of observations taking RL as references (RL against EL, NoL and HP).
In order to parameterise the bivariate spatial dependence, and to test the similarity of the spatial patterns of the shading effects, we calculated the bivariate association measure L [60,61]. For these spatial analyses, we used a subset data corresponding to a grid area of 128 m 2 (rectangle of 16 m × 8 m), encompassing the majority of the insolation reduction on ground, from the tree trunk 2 m to the south, 6 m to the north, and 8 m towards west and east directions.
The L measure, the univariate spatial association measures (SSSx and SSSy) and correlation coefficients associated to it, were estimated with the "spdep" library [55] in the R environment v 3.5.3 [54]. For the insolation reduction grid data, neighbours were created with the "cell2nb" function and the type of sharing boundary connectivity was set to "queen"; weights were given with the "nb2listw" function, and the globally standardised coding scheme style ("C") was chosen [62]. A permutation test for the L measure (400 random permutations, "lee.mc" function) of x and y (RL and EL, NoL or HP shading effects, respectively) for the given spatial weighting scheme established the rank of the observed statistic and calculated "pseudo p-values".

Results
We modelled the shade cast based on QSMs of six wild cherry trees for July 2019, at 10 min intervals, retrieved and compared the shading effects in the ground for the leaf modes (RL, EL, and NoL) and HP.

Modelling Insolation and Shading Effect with RL
Using our approach to simulate tree canopy with realistic leaves (RL), we modelled the insolation reduction for six wild cherry trees. In Figure 4, the shade cast by tree Pa_6 with RL is displayed for four timestamps (at 08:00, 12:00, 16:00, and 20:00 real local time; CEST) for 1 July 2019. In addition to the differing orientation of the shadow projections, the size and extent for the different daytimes can be observed. According to the sun movement for the northern hemisphere, the shadow projection is moving clockwise starting in a south-westerly direction early in the early morning and ending closer to the south-east direction in the nightfall. From the base of the tree, the projected shadow length was: 20 m towards the west, at the beginning of the day; 6 m towards the north, at midday and early afternoon; and 30 m to the east, prior to sunset. The number of active sun positions for the simulation period of July was 2865 (a total of 477.5 h of sunlight), meaning that the same number of shadow projections were modelled, for the six trees under different leaf modes, to estimate the shading effects: monthly insolation, insolation reduction, and shaded area. m towards the west, at the beginning of the day; 6 m towards the north, at midday and early afternoon; and 30 m to the east, prior to sunset. The number of active sun positions for the simulation period of July was 2865 (a total of 477.5 h of sunlight), meaning that the same number of shadow projections were modelled, for the six trees under different leaf modes, to estimate the shading effects: monthly insolation, insolation reduction, and shaded area. For the simulation period of July 2019, the maximum insolation was 691.98 MJ m −2 , or 22.32 MJ m −2 day −1 , on a daily basis. In Figure 5, the insolation around tree Pa_6 with RL is shown: the darker the colours, the greater the shading effect. Likewise, it is possible to observe areas of more intense energy reduction. For Pa_6, the mean insolation was 637.95 MJ m −2 , while the mean insolation reduction was 54.0 MJ m 2 .
In all trees with RL, the spatial distribution of the shading effect had a hyperbolic shape around the trees, with the strongest insolation reduction towards the north direction at a distance of ±7 m to the tree trunk. It spread roughly 16 m to the east and west direction, 8 m to the north and 4 m towards south-west and south-east. We found no shading effects towards the south. In addition, the mean insolation reduction was 60.7 MJ m −2 (±5.4) for RL, while the average shaded area was 175.91 m² (±42.56) per tree. On a daily basis, this resulted in an insolation reduction of 1.96 MJ m −2 day −1 (±0.17).   For the simulation period of July 2019, the maximum insolation was 691.98 MJ m −2 , or 22.32 MJ m −2 day −1 , on a daily basis. In Figure 5, the insolation around tree Pa_6 with RL is shown: the darker the colours, the greater the shading effect. Likewise, it is possible to observe areas of more intense energy reduction. For Pa_6, the mean insolation was 637.95 MJ m −2 , while the mean insolation reduction was 54.0 MJ m 2 .

Model Comparison
In all trees with RL, the spatial distribution of the shading effect had a hyperbolic shape around the trees, with the strongest insolation reduction towards the north direction at a distance of ±7 m to the tree trunk. It spread roughly 16 m to the east and west direction, 8 m to the north and 4 m towards south-west and south-east. We found no shading effects towards the south. In addition, the mean insolation reduction was 60.7 MJ m −2 (±5.4) for RL, while the average shaded area was 175.91 m² (±42.56) per tree. On a daily basis, this resulted in an insolation reduction of 1.96 MJ m −2 day −1 (±0.17).  In all trees with RL, the spatial distribution of the shading effect had a hyperbolic shape around the trees, with the strongest insolation reduction towards the north direction at a distance of ±7 m to the tree trunk. It spread roughly 16 m to the east and west direction, 8 m to the north and 4 m towards south-west and south-east. We found no shading effects towards the south. In addition, the mean insolation reduction was 60.7 MJ m −2 (±5.4) for RL, while the average shaded area was 175.91 m 2 (±42.56) per tree. On a daily basis, this resulted in an insolation reduction of 1.96 MJ m −2 day −1 (±0.17).

Model Comparison
We investigated the shading effects of all trees and leaf modes. Monthly insolation correlated strongly with the shaded area (~0.99), meaning the larger the shadows, the greater the shading effects. For all leaf modes, total monthly insolation ranged from 37,777.8 MJ (Pa_6 NoL) to 144,360.7 MJ (Pa_4 RL), and the total insolation reduction varied between 1782.7 MJ (Pa_6 NoL) and 14,753.2 MJ (Pa_4 RL). Furthermore, the shaded areas were greater for RL, and varied greatly within trees and leaf modes: from 57.17 m 2 (Pa_6 NoL) to 229.94 m 2 (Pa_5 RL).
The mean shading effects are presented in Table 4. Shaded area for HP were greater in size compared to the QSM-based leaf modes. The respective minimum and maximum shaded areas were: 181.15 m 2 (Pa_6 HP), and; 316.93 m 2 (Pa_5 HP). The shading effects of each tree and treatments is exemplified in Figure 6: the shaded area in comparison with the mean relative insolation reduction (colours specify the simulation set and symbols represent the sample trees). An effect of the differing leaf modes on the relative insolation reduction can be observed. RL showed more intense shading effects compared to concurrent modes, with a mean relative insolation reduction of 8.77% (±0.78). The mean relative insolation reduction in ascending order was: NoL < EL < HP < RL. The shading effect on a tree basis can be seen in Table S3.
We investigated the shading effects of all trees and leaf modes. Monthly insolation correlated strongly with the shaded area (~0.99), meaning the larger the shadows, the greater the shading effects. For all leaf modes, total monthly insolation ranged from 37,777.8 MJ (Pa_6 NoL) to 144,360.7 MJ (Pa_4 RL), and the total insolation reduction varied between 1782.7 MJ (Pa_6 NoL) and 14,753.2 MJ (Pa_4 RL). Furthermore, the shaded areas were greater for RL, and varied greatly within trees and leaf modes: from 57.17 m² (Pa_6 NoL) to 229.94 m² (Pa_5 RL).
The mean shading effects are presented in Table 4. Shaded area for HP were greater in size compared to the QSM-based leaf modes. The respective minimum and maximum shaded areas were: 181.15 m² (Pa_6 HP), and; 316.93 m² (Pa_5 HP). The shading effects of each tree and treatments is exemplified in Figure 6: the shaded area in comparison with the mean relative insolation reduction (colours specify the simulation set and symbols represent the sample trees). An effect of the differing leaf modes on the relative insolation reduction can be observed. RL showed more intense shading effects compared to concurrent modes, with a mean relative insolation reduction of 8.77% (±0.78). The mean relative insolation reduction in ascending order was: NoL < EL < HP < RL. The shading effect on a tree basis can be seen in Table S3.  The shaded area and mean proportion of each shade intensity classes can be observed in Figure 7. Over all simulations, we found a predominant shading effect in the two lowest intensity classes: 51.1% of shaded area are covered by the first class (2-5%), and 25.8% by the second class (5-10%). For RL, it was found an increased proportion of shaded area under the shade intensity classes 20-30% (9.9%) and ≥30% (1.3%), compared to EL and NoL. The portion of shaded area under 10-20% insolation reduction was similar between RL and EL, 20.3% and 20.1%, respectively. The shading effects of HP had a particular pattern: they varied in content for all intensity classes, and their proportions were more evenly distributed between classes with a negative exponential trend. intensity classes: 51.1% of shaded area are covered by the first class (2-5%), and 25.8% by the second class (5-10%). For RL, it was found an increased proportion of shaded area under the shade intensity classes 20-30% (9.9%) and ≥30% (1.3%), compared to EL and NoL. The portion of shaded area under 10-20% insolation reduction was similar between RL and EL, 20.3% and 20.1%, respectively. The shading effects of HP had a particular pattern: they varied in content for all intensity classes, and their proportions were more evenly distributed between classes with a negative exponential trend. The spatial distributions of EL and NoL shading effects were similar to RL: a hyperbolic shape around trees towards the north, where the strongest shading effect is manifested; with shade-free zone towards the south. However, they had reduced dimensions on the east-west and north-south axis compared to RL.
Hemispherical photography shadows spread differently. Shaded area was more regularly distributed around the tree trunk position (at 0,0,0), in a circular fashion, with a varying radius of 8 m to 14 m. For all trees, a more intense insolation reduction effect was found within the first 4 m around tree stem and directed to the north (similarly to the QSM-based shadows). Compared to RL, the shaded area had a shorter range in the east and west directions and, surprisingly, a strong insolation reduction around 4 m towards south.
The modelled insolation reduction is shown in Figure 8 for trees Pa_4 and Pa_6: the shading effects are zoomed to details limited at 8 m around tree trunk and presented in blocks comparing simulation sets. We observed areas of greater shading intensity on the shaded area (see Pa_4 with RL), which may correspond to branches of outstanding size, higher branchiness and/or higher leaf density in certain parts of the tree canopy.
Within the subset area focused on the insolation reduction of RL mode, we investigated the aggregation and similarity of the shadows in terms of spatial and aspatial association measures. All measures are reported in Table 5, and key measures are displayed in Figure 8. The spatial distributions of EL and NoL shading effects were similar to RL: a hyperbolic shape around trees towards the north, where the strongest shading effect is manifested; with shade-free zone towards the south. However, they had reduced dimensions on the east-west and north-south axis compared to RL.
Hemispherical photography shadows spread differently. Shaded area was more regularly distributed around the tree trunk position (at 0,0,0), in a circular fashion, with a varying radius of 8 m to 14 m. For all trees, a more intense insolation reduction effect was found within the first 4 m around tree stem and directed to the north (similarly to the QSM-based shadows). Compared to RL, the shaded area had a shorter range in the east and west directions and, surprisingly, a strong insolation reduction around 4 m towards south.
The modelled insolation reduction is shown in Figure 8 for trees Pa_4 and Pa_6: the shading effects are zoomed to details limited at 8 m around tree trunk and presented in blocks comparing simulation sets. We observed areas of greater shading intensity on the shaded area (see Pa_4 with RL), which may correspond to branches of outstanding size, higher branchiness and/or higher leaf density in certain parts of the tree canopy.
Within the subset area focused on the insolation reduction of RL mode, we investigated the aggregation and similarity of the shadows in terms of spatial and aspatial association measures. All measures are reported in Table 5, and key measures are displayed in Figure 8.
Values for the spatial smoothing scalars SSSx and SSSy confirmed that shading effects are spatially clustered. Often the clustering effect was stronger in EL than in RL, as ellipsoids project longer and more uniform fractions of shadow to the ground than RL, producing more similar neighbouring values on the ground. These measures of univariate spatial association offer a reference point (roof value) for the interpretation of the bivariate spatial association measure L.
For all trees, RL shading effects correlated strongly with EL (0.98-0.99), less strongly with NoL (0.89-0.95), and moderately with HP (0.55-0.80). The global L measure presented similar trend for the paired comparison of spatial distribution of shadows (RL shadows versus other modalities), which revealed a strong spatial similarity of the shading patterns RL-EL, and RL-NoL, and dissimilarity of RL-HP. Moreover, L statistic worked as a superior alternative to Pearson's correlation coefficient, as it captured the magnitude of the spatially clustering of the bivariate spatial associations [63].
The permutation test showed that all paired shading effects were spatially autocorrelated (pseudo p-value of 0.0025 unanimous; rank 400 of 400 simulations), and we therefore, accepted the alternative hypothesis that there is a greater variability in the shadows, which are not explained by the similarities of the spatial pattern of the insolation reduction area.
In the case of leaf modes, this is an evidence that the leaf pattern has an effect on the shading effects, since the woody structure defined by the QSM did not vary.

Discussion
The presented LCA and updated shadow model constitute a refined approach for the use of QSMs of real trees augmented with realistic leaves (as mean of canopy reconstruction) for modelling tree shading effects on the tree surroundings in high temporal and spatial resolution.

Leaf Creation Algorithm in Comparison to Others
To date several researchers have worked towards estimating tree and forest structural attributes and foliage properties from TLS point clouds. Côté et al. [64] worked with needle modelling for four conifer species to define the total amount of foliage in the crown and build the tree branching structure. Zheng and Moskal [65], presented a new method that indirectly and non-destructively retrieves foliage elements from TLS point clouds, assuming that tree leaf orientation, is an important attribute of forest canopy architecture and is critical in determining the within and below canopy solar radiation regimes. The separation of foliage and wood in TLS points clouds has been studied by many [28,50,66]. In all these studies the research goals have been diffuse, in cases towards reconstructing perfectly the woody architecture of trees or simply for deriving leaf attributes and insights on the light regime. Åkerblom et al. [42] presented a non-intersecting leaf insertion algorithm for tree structure models, its open MATLAB implementation, and conducted initial testing for English Oak (Querus robur L.). Similarly, to our LCA, the algorithm implementation is a tool aiming for realistic canopy reconstruction of trees and allows users to set the leaf shape (triangle-based geometry like) and distribution of leaf location, size and angles.
Many different parameters can be used to define leaf distribution on a tree to fully reconstruct canopy architecture. In this context, the LCA was designed for a simplified approach having few assumptions embedded, and parameters set to constant values: leaf inclination angle is 0 • , parallel to a theoretical horizontal ground plane; phyllotaxy is defined as alternate distichous, which delivers leaves arranged alternatively to the left and right on a same plane; leaf orientation is +90 • or −90 • (right or left) apart of the specific branch-cylinder directional axis; leaf disposition on the branch is ruled by the allometric relationship, establishing a leafless zone for first order branches. A main concern with the leaf inclination angle is that steeper angles decrease light captured when the sun is at high viewpoints in the sky (midday and summer), whilst increasing light captured from lower sun angles in the morning and afternoon, as well as in winter [67]. Moreover, we do not account for leaves that are clumped together, which is observed in many wild cherry individuals. All these aspects could be source of uncertainty when modelling shade cast by trees. The only user-defined parameter is the leaf spacing, which is used for controlling the total leaf area of a tree. The allometric equation for specifying the expected leaf area on a tree was based on the measurements of thirteen branches. As a limited dataset, the total leaf area equation does not account for complex branch structures (i.e., wind-damaged branches or from pruning treatments), and a more detailed and increased sampling size would be beneficial. Although all wild cherry trees include in this study belonged to AFS, our laser-sampled trees could have received different silvicultural treatments and were exposed to diverse growing conditions (i.e., in more open landscapes, isolated trees), than the trees utilised for the manual leaf sampling. Our sample trees had total leaf area values similar to the findings of Urban et al. [68] for P. serotina with comparable DBH and height ranges. Likewise, Miljković et al. [69] presented corresponding data on leaf area, leaf width, leaf length, and petiole length of P. avium individuals from five different natural populations in Bosnia and Herzegovina.
Another source of uncertainty is the tree reconstruction with QSMs, which was controlled by the specific optimisation process inside TreeQSM, testing many model parameters over repeated model runs [41]. In addition, the parameter "cylinder length" could play a major role in the EL leaf mode, as the construction of the ellipsoidal leaves is determined by the branch-cylinder diameter (i.e., 1 cm), so that an increased number of cylinders per branch is preferred.
Even if we used specific plants traits for model parameterisation and for defining assumptions, this could still be further adapted to any tree species by changing specific algorithm parts to user-defined options and values. An alternative to parameterise the leaf inclination angle would be by using the dataset provided by Chianucci et al. [70], which offers species-specific parameters and cover a great range of broadleaved tree and shrub species in temperate and boreal regions. Dataset observations are derived with the established levelled digital photography (LDP) approach, which is known of sensing only leaves oriented approximately parallel to the viewing direction of the camera, and are restricted to the source locations (Sweden, Estonia, USA, central Italy). Furthermore, LiDAR based approaches for assessing leaf angle distribution already overcame some of the shortcomings of LDP as they are not limited by leaf curvature [71], and might be a new data source for parameterisation of the LCA.

Shadow Simulations-Performance and Comparison
A full investigation of interrelations of QSM-derived tree parameters and shading effects was beyond the scope of this study. Nevertheless, we noticed greater shading effects with an increased tree structure (e.g., branch volume), the case for Pa_4 and Pa_5. Conventionally, tree height and crown size and density are noticeably affecting the shaded area and intensity of the shading effect. There was a trend associated with the insolation reduction which was clearly an effect of the leaf mode on a per tree basis, as the tree structure (QSM) of each individual was the same. We demonstrated the reduction of insolation on different leaf modes (RL, EL, and NoL) in the comparative blocks. For all shading effects modelled by means of the shadow algorithm, the visualised heterogeneity was the true effect of leaves (RL or EL).
Hemispherical photography is a 2D approach, for this reason the comparison with highly detailed 3D point cloud data is biased. In practice, using the hemispherical photography methodology a short tree with a dense crown in a photograph could be comparable to a tall tree with a lighter crown. Photographs taken in different distances to the tree, as carried out in this study, improve the situation, but it is still less precise than 3D data. As expected, we demonstrated that hemispherical photography reports different shading effect (intensity and size) and spatial distribution than simulations with our shadow model.
Within this study, we used factual weather data, corresponding to the time-period for the simulations, for RL (EL and NoL), while HP is based on average data which may suppress extreme values from July 2019 as applied to the simulations, possibly leading to differences in the shading effect of HP.
Dupraz and Liagre [20] presented examples of temporal shade cast by trees in France, their results displayed similar spatial distribution than our modelled shadows, but lacked the important fine spatial resolution. This approach of modelling the radiation regime around trees based on 3D structures used a replacement of the complete tree crown by an ellipsoid for modelling the light availability in AFS. The method was further developed by Talbot and Dupraz [21], which proved the efficiency of the model to predict the light pattern around an average tree, but not suitable for simulating the variability of individual trees.
It was observed that the shadows modelled with RL had a greater shading effect compared to the adjusted ellipsoidal leaf-replacement proposed by Rosskopf et al. [32]; unsurprisingly, RL and EL distributed similarly on space, and the shading effect of EL was often more clustered than RL. EL is a simplified leaf mode and is species independent. An additional calibration of the cylinder-radius threshold for the suggested ellipsoids (i.e., greater than 1.0 cm) could improve the results of EL in this study. It was noted that spherical-like leaf structures were formed with the EL approach, as well as few ellipsoids of distinguished length. On a specific branch order, single branch-cylinders with radius ≤ 1.0 cm, or a short set of cylinders, produced and ellipsoidal leaf akin to a sphere, or a coin shape. The extremely large ellipsoids are associated to the branch order of a branch ramification, and/or the uncertainty of the QSM. In the first case, where branch-cylinder radius were ≤ 1.0 cm, the lowest branch order could have had a large total length, and consequently a noticeably large ellipsoidal leaf. Secondly, due to noise or occlusion in the tree point clouds, unrealistic branch ramification could have been produced within the QSMs, which led to the creation of few ellipsoids of outstanding length. Though EL could be applied openly for any tree species, it is case-sensitive, since plant traits (i.e., leaf traits), tree architecture and reconstruction with QSMs are likely to influence the choice of the optimal parameters' values, what also justifies the optimisation of QSMs (conducted in this study).
The RL scenario assumption of horizontal leaf angle could have slightly overestimated the shading effects during the summer period (shaded area, and thus, a reduction of insolation). In reality, such leaf orientation is different to the simplified rules applied herein. Leaves could have a focal point on the sun's sky path, being horizontal when shaded or with varied angles to the incoming sunlight [72]. Moreover, leaves could be sparser, but bigger on the northern side to ensure adequate light interception [73]. Further efforts must be devoted to improve simulations and models based on the orientation of leaf planer surfaces affixed to the QSM. Represented by NoL, the leaf-off condition of trees remains relevant for AFS, since crops such as winter wheat could benefit from broadleaved trees without leaves in late autumn and winter and may be disproportionally impacted by shade cast due to reduced insolation at this time of year. NoL had reliable results as the spatial distribution of the shadows is very similar to the RL (in the focused shaded area), however, the shading effect is reduced to approximately 40% (area and intensity) in relative terms. This information could be used to simplify our approach, as adding leaves is time consuming, shadows produced with QSMs (alone) could already give precise location and a rough estimation of the shaded area and insolation reduction over a time period, possibly correcting it with an expansion factor.
Shadow simulations and comparisons of shading effect of leaf modes were carried out only for July, and these results agreed with the general spatial pattern presented by Rosskopf et al. [32]. Therefore, we would expect similar results independent of the month of the year.
Finally, we acknowledge a few limitations of our approach. Our shadow algorithm, implemented in R environment, requires high computational power for the leaf-on modes (≥32 Gb RAM, in many cases) and is still time consuming (~52 h of simulations for one tree, at the highest spatial and temporal resolution for one month). Simulations for the leaf-off mode were faster (~16 h). While a full assessment of algorithm performance was not targeted within this study, we recognise the need for further algorithm improvement and efficiency optimisation. Up to now, the shadow model does not have an openly available implementation. Once accessible, a fully standardised implementation of the shadow model should consider a standard way of delivering outputs, including analytical simulation reports. First-hand laser data obtained with TLS remains costly and inaccessible for many research groups around the word. In contrast, point clouds of forest stands and trees, as well as QSMs, are becoming increasingly available, enabling users to search for 3D data to be applied with the scope of modelling shading effects. The use of specific weather data makes the simulations of insolation case-specific, but for an open adoption of our method, it could offer a standardised protocol for including any region on earth. A full validation of the shadow model will require more efforts at field level and would have to include a thorough assessment for calibrating model parameters; such a solution is sought in the visible future.

Conclusions
We presented an updated shadow model that uses factual 3D data in combination with a leaf creation algorithm (LCA) to estimate individually the shading effect of six wild cherry trees growing under agroforestry conditions in southern Germany. The coupling of the shadow model with factual climate data added to the accuracy of the model. This approach was able to refine simulations of the shading effects of scattered trees represented by QSMs with realistic leaves (RL). We also evaluated the shading effects of other leaf moves (EL and NoL). As a simplified leaf-on mode, EL had similar distribution of the shading effect with slight underestimation in total insolation reduction and shaded area. Though RL might be closer to reality, the EL mode could be an alternative for the (at the time) exhaustive time consumption of simulations with RL. Likewise, the shading effects of NoL leaf-off mode were also similarly distributed to RL, though reduced by almost half of the intensity. The shading effects of NoL gave us insights on the potential of leafless woody structures to act as obstacles mediating the radiation regime on ground. The alternative method had dissimilar shading effects and distribution.
The shadow model is a suitable tool for the detailed quantification of the shading effects of single trees. We estimated insolation reduction and shaded area values focused on the particularities of the tree structures. Our approached is limited to work with QSMs, what implies in sample trees to be scanned outside the vegetation period. Derived results can be utilised to support decision-making by land managers, in agricultural and agroforestry systems, or even to facilitate the inclusion of trees in urban areas, as a strategic public health measure.

Outlook
We fulfilled our aim of developing an improved simulation of the shading effects of single-trees, but further improvements of the shadow model are foreseen. First, increasing algorithm capabilities for handling data from multiple trees [74]. Secondly, enabling the leaf inclination angle to be set as an input parameter and/or to follow a distribution function. Moreover, increasing capabilities to address phyllotaxy and leaf geometry of other tree species. Lastly, the inclusion of growth simulations to account for the growth of the trees during the vegetation period.
Supplementary Materials: The following are available online at https://www.mdpi.com/2072-429 2/13/3/532/s1, Figure S1: Leaf area count of the five leaf size classes (extra small, small, medium, large, extra-large), establishing the proportions of leaves to be created by the Leaf Creation Algorithm (LCA). Figure S2: Linear trend between branch length and first order branch diameter (N = 12). Figure S3: Branch leaf area as a function of first order branch diameter (N = 13). Figure S4: Relative leaf area according the relative branch position. Equation root at 8.47% relative branch position defined the initial position for insertion of leaves. Table S1: Measured and adjusted leaf parameters per leaf class. Table S2: Leaf attributes that define the leaf geometry points, from which six geometric points are derived: length, width, lower and upper edges. Table S3: Shaded area and insolation reduction of wild cherry trees (Pa_1 to Pa_6, in symbols) for different leaf modes (RL, EL, and NoL) and HP. Video S1: Wild cherry tree Pa_5 reconstructed as a QSM and augmented with the leaf modes; profile and top views.  Data Availability Statement: The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.