A method for estimating tree crown porosity based on porous media theory and computer graphics

Abstract

The present invention discloses a tree crown porosity estimation method based on porous media theory and computer graphics, including collecting tree laser point cloud data; classifying leaf and branch point clouds, fitting branches by cylinders with adaptive radius allocation, extracting point cloud data of a single leaf; performing denoising on the point cloud data of a single leaf; determining the edge points of each leaf by a cubic polynomial curve fitting method; calculating the equivalent thickness of each leaf, and then enveloping the point cloud of each leaf by a hexagonal prism formed by adaptive parameter allocation to obtain the equivalent volume of each leaf; calculating the total equivalent volume of all leaves in the tree crown and the total volume of all branches, and estimating the porosity of the tree crown in combination with the tree crown volume. The GF global solution obtained by the present invention does not change with the change of the observation angle, and overcomes the limitation of generating local data of HP by projection transformation, and the limited visual coverage of omitting other information from a fixed position.

Description

Technical Field

[0001] The invention belongs to the technical field of tree crown porosity analysis, and in particular relates to a tree crown porosity estimation method based on porous medium theory and computer graphics. Background Art

[0002] The porosity (GF) of tree crowns, also known as canopy porosity, is very important for characterizing forest habitats. GF ​​is not only a key parameter for the accurate quantification of the leaf area index (LAI), but also an important factor for a better understanding of plant-environment exchange processes, such as photosynthesis, plant transpiration efficiency, and canopy precipitation interception. The lack of air movement in the plant canopy produces a more humid microclimate, but increasing GF favors air movement within the canopy, which ensures that leaf elements within the canopy are exposed to solar radiation, reduces shadowed areas with poor lighting conditions, and prevents the occurrence of fungal infections within the forest canopy. In addition, GF can be used as a bio-indicator to detect leaf fall and leaf area density. Therefore, accurate retrieval of canopy GF is of central importance for various physiological, climatological, and biogeochemical studies.

[0003] Existing methods were built with the help of modern optical imaging and laser scanning technologies. Existing methods for calculating GF can be divided into two categories. The first category is based on image processing methods; these methods consist of capturing fisheye photographs (HP) of the studied forest upward or projecting the canopy elements into a two-dimensional image and calculating GF as the number of sky pixels (excluding pixels covered by the canopy) to the total image pixels. The other category is based on methods based on light transport models; these methods define GF based on the probability that a emitted ray will pass through the canopy and travel a certain distance in the canopy without encountering a vegetation element.

[0004] For GF calculation methods based on image processing, there are large differences in sky and canopy pixels in the collected images under different sky conditions, forest stand appearance, and shaded vegetation elements with variable solar irradiation angles, which are likely to make the reference readings and canopy readings of imaging techniques erroneous. In addition, some image quality degradation is caused by camera equipment specifications and the imaging environment of the forest, such as optical aberrations caused by changes in camera exposure values, whether the intrinsic and extrinsic parameters of camera calibration are aligned, airborne organisms moving within the canopy, and dust or water vapor in the air can generate noise in the image. GF calculation methods based on light transport models are always affected by ghost points or multi-scan registration deviations from wind on the target tree. In addition, laser scanning sensors exhibit nonlinear scattering effects depending on atmospheric conditions and canopy geometry, and the reflective properties of different plant species materials add additional interference and lead to instability in the estimation of GF.

[0005] In addition, since the crown is a special geometric target and the location of measurement is not affected by the local crown pattern, the calculation results of the above two methods show different results with the change of observation angle and observation position. For example, when analyzing HP based on image processing methods, different GF magnitudes are obtained from different annular segments, corresponding to different zenith angles in a forest area, and there is no unified standard or mechanism to explain these differences. Similarly, the method based on optical transmission model is also sensitive to the light emission angle. At a fixed emission angle, the partial interception and gap identification obtained from ranging LiDAR data are related to the horizontal and vertical distribution of crown mass to stem base position along this emission angle, which always depends on the statistical sample of a tree, and it is difficult to quantify the crown with asymmetric phenotypic characteristics and allometric growth characteristics. In addition, some statistical methods are only applicable to large-scale forests with a large amount of data analysis, or standardize the distance and direction of sample data collection to compensate for the deviation of individual specimens, but these methods are not suitable for single trees without the guidance of prior analysis of larger forests with a large number of homogeneous populations.

[0006] A porous medium is a material that contains pores (voids). The skeletal part of the material is often referred to as the "matrix" or "framework". The pores are usually filled with a fluid (liquid or gas), while the skeletal material is usually solid. Many natural materials such as rocks, soils, and biological tissues such as bone, wood, cork, and cotton can be considered porous media due to their inherent properties similar to porous media. These materials all have their own porosity size, which is mentioned in the discussion section. As an important biological form on Earth, the canopy of a tree consists of stems, leaves, and reproductive structures that remain almost stationary above the ground as air currents and sunlight flow through the canopy. Therefore, if both the photosynthetic and non-photosynthetic parts are considered as skeletal materials, and the pores between the plant materials in the canopy are considered as pores that allow air currents and light to penetrate, then the canopy can be considered as a multiphase and porous medium.

[0007] Based on the porosity index in porous media theory, where porosity is defined as the portion of the channel volume occupied by the gas phase, the concept of porosity can be derived to define the volume-based GF of the crown by evaluating the volume occupied by the branches and leaves aggregation of the crown, and deriving the remaining volume ratio of the pores (gaps) between the vegetation elements in the crown. This volume-based GF is a normalized variable that can describe the overall spatial volume of the crown, transcending the limitation of obtaining information from a fixed perspective. However, the bottleneck of implementing this method is the need to finely describe the spatial distribution and volume of each vegetation element in the crown, and enable computer machines to graphically and geometrically represent it in a way that combines analytical and visual capabilities. Three-dimensional modeling in computer graphics is the process of making a three-dimensional representation of any surface or object based on mathematical coordinates through specialized algorithms. If the surface of each leaf and branch can be finely simulated into a geometric primitive based on computer graphics technology, the converted geometric primitive has a known volume size and constitutes a crown surrounded by pores, making it possible to estimate the volume-based GF in the porous matrix of the crown filled with vegetation materials. Summary of the invention

[0008] The technical problem to be solved by the present invention is that the existing GF estimation method based on fisheye photo (HP) processing and light transmission model ignores the intrinsic geometric structure of the crown, that is, the natural morphological characteristics of the spatial branch and leaf components of the crown, and is limited by the observation angle and can only obtain local information. A crown porosity estimation method based on porous media theory and computer graphics is provided. The GF global solution obtained by the crown porosity estimation method based on porous media theory and computer graphics does not change with the change of observation angle, and overcomes the limitations of the local data expansion of HP generated by projection transformation, and the limited visual coverage of omitting other information from a fixed position.

[0009] In order to achieve the above technical objectives, the technical solution adopted by the present invention is:

[0010] A method for estimating tree crown porosity based on porous media theory and computer graphics includes:

[0011] Step 1: Collect tree laser point cloud data;

[0012] Step 2: Use the branch-leaf separation algorithm to realize the point cloud classification of leaves and branches, fit the branches by cylinders with adaptive radius allocation to model the branches, and extract the point cloud data of a single leaf;

[0013] Step 3: Use the moving least squares method to denoise the point cloud data of a single leaf;

[0014] Step 4: Using a cubic polynomial curve fitting method to determine the edge points of each blade;

[0015] Step 5: Calculate the equivalent thickness of each leaf, and then use the hexagonal prism formed by adaptive parameter allocation to envelop the point cloud of each leaf to obtain the equivalent volume of each leaf;

[0016] Step 6: Calculate the total equivalent volume of all leaves in the crown and the total volume of all branches based on cylinder fitting, and combine them with the crown volume to estimate the porosity of the crown.

[0017] As a further improved technical solution of the present invention, the step 1 is:

[0018] Tree laser point cloud data was collected using a terrestrial laser scanner.

[0019] As a further improved technical solution of the present invention, the step 3 is specifically:

[0020] In the local subdomain of the fitting area, let point p i is the i-th point cloud of the leaf, For point p i coordinates, where i = 1, 2, 3..., n, in order to eliminate the noise caused by blade vibration interference, point p i of The value needs to be recalculated, and the new value obtained by recalculation is recorded as Retrieving New Values The fitting function is expressed as follows:

[0021]

[0022] in It is p i There are m neighboring points p in the circular neighborhood centered at j x,y coordinates, where j = 1,2,3...,m; basis function is with p i ξ×1 column vector related to the neighboring points; is the uncertainty coefficient. When formula (2) is the minimum, Take the best value:

[0023]

[0024] where p j Yes i Neighborhood points, j = 1, 2, ... m; ω (||p i -p j || 2 ) is based on two points p i and p j The weight function with compact support property of the Euclidean distance of , where s = ||p i -p j ||2 and where r is the radius of the collection of circular adjacency domains; then,

[0025] As a further improved technical solution of the present invention, the step 4 is specifically as follows:

[0026] (4.1) The two endpoints on the midrib of the leaf are denoted as p e (x e ,y e ,z e ),p s (x s ,y s ,z s ), where (x e ,y e ,z e ) indicates p e The coordinates of (x s ,y s ,z s ) represents point p s The coordinates of the two endpoints form a central axis L 1 , p e (x e ,y e ,z e )∈P MLS ,p s (x s ,y s ,z s )∈P MLS , where P MLS This is the fitted leaf point cloud after denoising;

[0027] (4.2) Project the leaf surface onto the XY plane, the central axis L 1 After projection, it becomes L′ 1 , L′ 1 The point p on is: p = p e (x e ,y e )+t×(p s (x s ,y s )-p e (x e ,y e )), where t is the slope of The independent variable of the straight line;

[0028] (4.3) In L′ 1 Sample n-1 points uniformly on the top and define multiple lines perpendicular to L′ 1 And the vertical line L passing through n-1 points respectively2 , vertical line L 2 The slope of Together with L 1 The two endpoints p e and p s Together, through the vertical line L 2 Divide the blade point cloud into multiple segments at equal intervals;

[0029] (4.4) Through the back-projection transformation from the two-dimensional plane to the three-dimensional space, we can determine the 2 The left endpoint of each segment of the direction The corresponding x, y values, the corresponding z value and the right endpoint The corresponding x, y values ​​and the corresponding z values, where the left endpoint The corresponding x, y values ​​and the corresponding z values ​​are recorded as Right endpoint The corresponding x, y values ​​and the corresponding z values ​​are recorded as Where k = 1, 2, 3...n-1;

[0030] (4.5) Using the cubic polynomial curve fitting method, a cubic polynomial is used to fit the independent variable to all left edge points. With dependent variable and The nonlinear relationship between them is as follows:

[0031]

[0032] in The left edge point Fitting The mean square error, The left edge point Fitting The mean square error of

[0033] Minimize the mean square error in formula (3), expressed as:

[0034] and Where b = 0, 1, 2, 3;

[0035] Then find the coefficient and Therefore, after fitting and The values ​​are expressed as:

[0036] and

[0037] Similarly, calculate the fitting values ​​corresponding to all right edge points according to the method in step (4.5) and

[0038] Finally, the accurate x, y and z values ​​of the left edge point are fitted and the exact x, y, and z values ​​of the right edge point Where k = 1, 2, 3...n-1.

[0039] As a further improved technical solution of the present invention, the step 5 is specifically as follows:

[0040] (5.1) For each leaf with n-1 left and right edge points, select round((n-1) / 4) left edge points and round(3(n-1) / 4) left edge points from the n-1 left edge points in top-down order, and select round((n-1) / 4) right edge points and round(3(n-1) / 4) right edge points from the n-1 right edge points in top-down order. Place the two selected left edge points, the two right edge points, and the endpoint p e and endpoint p s As the six vertices of a hexagonal base, round() represents the number in the brackets to the nearest integer;

[0041] (5.2) Calculate the normal vector of each leaf surface using the covariance matrix according to the fitted points, and select the highest point along the direction of the calculated leaf normal vector. and the lowest point highest point and the lowest point The distance between normal is taken as the equivalent thickness of the current blade;

[0042] (5.3) dist normal It is used as the height and combined with the above hexagonal base to form a hexagonal prism, which in turn wraps the fitting point cloud of a single leaf;

[0043] (5.4) Calculate the volume of the hexagonal prism, which is the equivalent volume of the current leaf

[0044] As a further improved technical solution of the present invention, the calculation formula of the porosity of the tree crown in step 6 is:

[0045]

[0046] Where N represents the total number of leaves in the crown, M represents the total number of branches in the crown, represents the equivalent volume of a leaf in the crown, Represents the volume of a branch in the crown, Vcanopy Represents the tree crown volume.

[0047] The beneficial effects of the present invention are:

[0048] The porosity (GF) of vegetation canopies is an important growth property that determines the transmission of light and airflow through the canopy. Existing GF estimation methods based on image processing and light transport models ignore the intrinsic geometric structure of the canopy, including branch and leaf elements, and are limited by the limited field of view of a fixed emission position. The method of the present invention innovatively applies the information of individual leaves and branches at the scanning point to perform GF inversion calculations. In addition, porous media theory is applicable to the evaluation of most natural and artificial materials, such as porous stones or sponges, allowing airflow or light to pass through the pores in these materials. The canopy can be represented as a porous entity with a mixture of pores. The porous structure of the canopy has a decisive influence on the ventilation and light transmission of the tree body, allowing light and airflow to propagate through the pores between leaves. Therefore, it is reasonable to regard the canopy as a porous medium for retrieving GF. On the other hand, the curling and twisting of leaves is a common phenomenon in most deciduous tree species, which is caused by many intrinsic or external factors, such as natural characteristics, insect damage, abiotic disorders and even herbicides. For most tree species, the front side of their leaves is not as flat as the back side due to the presence of rugged leaf veins. From the side view of any leaf or the view parallel to the midrib, the equivalent thickness caused by the twisting or drooping behavior of the leaf blocks the transmission of light and the passage of air from the side. Therefore, the present invention proposes the concept of equivalent thickness of the leaf, which is mainly based on the distance between the highest and lowest fitting points along the normal vector of the leaf surface, plus the refined edge points to form the bottom surface of the hexagon to define the shape parameters of the hexagonal prism, enveloping the scanning points of each leaf. Therefore, the equivalent leaf volume is approximated to the corresponding value of the hexagonal prism, while the wood volume and crown volume can be evaluated by fitting cylinders and Alphashape algorithms, respectively, to reversely calculate the plant wood volume in the current crown. Therefore, the method of the present invention takes into account the natural morphological characteristics of the spatial branch and leaf components of the crown, and the resulting GF global solution does not change with the change of the observation angle, and overcomes the limitations of generating local data of HP through projection transformation, and the limited visual coverage of omitting other information from fixed positions.

[0049] The present invention uses LiDAR data derived from porous media theory and creates a volume-based GF calculation method using computer graphics techniques. Individual leaves are extracted from the point cloud based on our previous segmentation algorithm, and then the noise is removed and the complete edge of each leaf is extracted by using moving least squares and cubic polynomial fitting. The newly defined equivalent thickness depicts the degree of curling and drooping of the leaf flakes, and the fitted leaf and branch points are enclosed by hexagonal prisms and cylinders, respectively, allowing the total volume ratio of nutrients in the crown to be the measured canopy volume. All effects allow for volume-based GF calculations, which we apply to three different canopies and compare with existing methods. BRIEF DESCRIPTION OF THE DRAWINGS

[0050] Figure 1 (a) is a schematic diagram of the scanning point cloud of Crape Myrtle collected by ground laser scanning.

[0051] Figure 1 (b) is a schematic diagram of the scanning point cloud of cherry blossoms collected by ground laser scanning.

[0052] Figure 1 (c) is a schematic diagram of the scanning point cloud of Magnolia collected by terrestrial laser scanning.

[0053] Figure 2 (a) is a schematic diagram of the branch and leaf separation results of the ground laser scanning data of Crape Myrtle.

[0054] Figure 2 (b) is a schematic diagram of the branch segment reconstruction results of the ground laser scanning data of Crape Myrtle.

[0055] Figure 2 (c) is a schematic diagram of the branch and leaf separation results of the ground laser scanning data of cherry blossoms.

[0056] Figure 2 (d) is a schematic diagram of the branch segment reconstruction results of the ground laser scanning data of cherry blossoms.

[0057] Figure 2 (e) is a schematic diagram of the branch and leaf separation results of the ground laser scanning data of Magnolia.

[0058] Figure 2 (f) is a schematic diagram of the branch segment reconstruction results of the ground laser scanning data of Magnolia.

[0059] Figure 2 (g) is a schematic diagram of the data extraction results of a single leaf from the ground laser scanning data of Crape Myrtle.

[0060] Figure 2 (h) is a schematic diagram of the data extraction results of a single leaf from the ground laser scanning data of cherry blossoms.

[0061] Figure 2 (i) is a schematic diagram of the data extraction results of a single leaf from the ground laser scanning data of Magnolia.

[0062] Figure 3 (a1) is a schematic diagram of relative thickness shown by the arrows as viewed from the side of the crape myrtle leaf.

[0063] Figure 3 (a2) is a schematic diagram of relative thickness indicated by arrows parallel to the midrib of the Lagerstroemia leaf.

[0064] Figure 3 (b1) is a schematic diagram of relative thickness shown by the arrow when viewed from the side of a cherry blossom leaf.

[0065] Figure 3 (b2) is a schematic diagram of relative thickness indicated by arrows parallel to the midrib of a cherry leaf.

[0066] Figure 3 (c1) is a schematic diagram of relative thickness shown by the arrow from the side of the magnolia leaf.

[0067] Figure 3 (c2) is a schematic diagram of relative thickness indicated by arrows parallel to the midrib of a magnolia leaf.

[0068] Figure 4 (a) is a schematic diagram of fitting the crape myrtle leaves using the moving least squares method to remove noise and smooth the surface.

[0069] Figure 4 (b) is a schematic diagram of fitting the cherry blossom leaves using the moving least squares method to remove noise and smooth the surface.

[0070] Figure 4 (c) is a schematic diagram of fitting magnolia leaves using moving least squares fitting to remove noise and smooth the surface.

[0071] Figure 4 (d) is the vertical line L perpendicular to the midrib of the crape myrtle leaf 2 Schematic diagram of dividing the leaf surface into many segments.

[0072] Figure 4 (e) is a vertical line L perpendicular to the midrib of the cherry leaf 2 Schematic diagram of dividing the leaf surface into many segments.

[0073] Figure 4 (f) is the vertical line L perpendicular to the midrib of the magnolia leaf 2 Schematic diagram of dividing the leaf surface into many segments.

[0074] Figure 4(g) is a schematic diagram of the hexagonal prism constructed on the crape myrtle leaf.

[0075] Figure 4 (h) is a schematic diagram of the hexagonal prism constructed on the cherry blossom leaf.

[0076] Figure 4 (i) is a schematic diagram of the hexagonal prism constructed on the magnolia leaf.

[0077] Figure 5 (a) is a schematic diagram showing that each crape myrtle leaf is wrapped by a hexagonal prism.

[0078] Figure 5 (b) is a schematic diagram showing that each cherry blossom leaf is wrapped by a hexagonal prism.

[0079] Figure 5 (c) is a schematic diagram showing that each magnolia leaf is wrapped by a hexagonal prism.

[0080] Figure 5 (d) is a schematic diagram of the crape myrtle crown enclosure based on the Alphashape algorithm.

[0081] Figure 5 (e) is a schematic diagram of the cherry tree crown enclosure based on the Alphashape algorithm.

[0082] Figure 5 (f) is a schematic diagram of the magnolia crown enclosure based on the Alphashape algorithm.

[0083] Figure 6 (a) is a schematic diagram of the HP generated by transforming the Ziwei scanning data using Lambert azimuthal conformal projection.

[0084] Figure 6 (b) is a schematic diagram of the HP generated by converting the cherry blossom scan data using Lambert azimuthal conformal projection.

[0085] Figure 6 (c) is a schematic diagram of the HP generated by converting the Yulan scan data using Lambert azimuthal conformal projection.

[0086] Figure 6 (d) is a schematic diagram of the three-dimensional visualization of the crape myrtle tree based on cubic voxels.

[0087] Figure 6 (e) is a schematic diagram of the 3D visualization of the cherry blossom tree based on cubic voxels.

[0088] Figure 6 (f) is a schematic diagram of the three-dimensional visualization of the magnolia tree based on cubic voxels.

[0089] Figure 6 (g) is a schematic diagram comparing the calculation result matrix of the porosity of each crape myrtle crown based on the optical transmission method and the calculation result of the crape myrtle crown by the method of the present invention.

[0090] Figure 6 (h) is a schematic diagram comparing the calculation result matrix of the porosity of each cherry tree crown based on the optical transmission method with the method of the present invention.

[0091] Figure 6 (i) is a schematic diagram comparing the calculation result matrix of the porosity of each magnolia tree crown based on the optical transmission method with the method of the present invention.

[0092] Figure 7 (a) is a schematic diagram showing the result of applying the method of the present invention to the leaf segmentation of rubber tree.

[0093] Figure 7 (b) in Figure 7 Schematic diagram of three local magnification results in (a) and schematic diagram of a single leaf enclosed by a hexagonal prism in three local magnification areas.

[0094] Figure 7 (c) is a schematic diagram showing the result of applying the method of the present invention to the segmentation of Ginkgo biloba leaves.

[0095] Figure 7 (d) in Figure 7 Schematic diagram of three local magnification results in (c) and schematic diagram of a single leaf enclosed by a hexagonal prism in three local magnification areas.

[0096] Figure 8 (a) shows an image of a dense forest plot of rubber tree clone PR 107 and a schematic diagram of scan points.

[0097] Figure 8 (b) is a schematic diagram of converting the scanning points into a synthetic fisheye photo (HP) with the center of the plot as the projection center of the spherical projection.

[0098] Figure 8 (c) is a schematic diagram of using the algorithm of the present invention to magnify part of the tree crown, highlight and identify individual leaves with different colors, and enclose them by a semi-transparent hexagonal prism. DETAILED DESCRIPTION

[0099] The specific embodiments of the present invention are further described below according to the accompanying drawings:

[0100] The porosity (GF) of vegetation canopies is an important growth property that determines the transmission of light and airflow through the canopy. Existing GF estimation methods based on image processing and light transport models ignore the intrinsic geometry of the canopy, including branch and leaf elements, and are limited by a limited field of view with a fixed emission position. Our work builds on porous media theory and computer graphics techniques to treat vegetation elements in the canopy as a solid matrix and the gaps between them as pores for volume-based GF calculation. By extracting individual leaves from terrestrial laser scanning data and using a fitting algorithm to remove noise and depict leaf edges, we innovatively define the concept of equivalent thickness that describes the degree of leaf curling and drooping, representing the interception of light and airflow from the side, and calculate it from the fitting points. Next, we implement wrapping of the scan points of each leaf with a hexagonal prism with adaptive size allocation, and build a cylindrical model that fits each branch to calculate the equivalent leaf and branch volumes of the canopy. Finally, the volume-based GF of the crown is obtained by subtracting the ratio of the total volume of plant leaves and trunk to the crown volume from 1, which is consistent with the porosity in porous media theory. This method was tested on three trees with different crown structures, and the estimated total equivalent leaf volume and volume-based GF were 0.018 m3 and 0.989 for crape myrtle, 0.051 m3 and 0.978 for magnolia, and 0.186 m3 and 0.968 for cherry tree. This work proposes a unified standard framework for volume-based crown GF calculation, which is the same as the porosity retrieval of most materials in the world.

[0101] The specific research is as follows:

[0102] 1. Materials and methods

[0103] 1.1. Research location and data collection:

[0104] The selected experimental trees are located on the campus of Nanjing Forestry University (32°08'N, 118°81'E, WGS-84), including individual cherry (Cerasus sp.), crape myrtle (Lagerstroemia indica L), and magnolia (Michelia L) trees. The point clouds of the three trees, crowns and skeletons, were acquired using a Leica C10 terrestrial laser scanning (TLS) system.

[0105] In October 2019, a complete point cloud of the three trees was obtained using a Leica C10 terrestrial laser scanning (TLS) system. The Leica C10 instrument is a 532 nm phase-based scanner with a 360° × 270° upward field of view and a laser rate of 5000 points per second. The instrument was set at a height of 1.5 m and a target measurement accuracy of ±1.5 mm (at a distance of approximately 3.0 m). The circular laser beam diameter and beam divergence at the scanner exit were 3 mm and 0.22 mrad, respectively, and the minimum distance between the resulting continuous beams was approximately 0.4 mm at a distance of 3 m from the instrument, with a scanning speed of up to 50,000 points per second. Based on the different heights of the three experimental trees, the distance between the experimental tree and the TLS was set to 3 m, and two lateral scans were collected from opposite sides of each crown and pointed to the center of the crown to achieve a complete phenotypic feature scan of the target tree.

[0106] The point cloud structures of the three experimental trees scanned by TLS are shown in Figure 2. Figure 1 In addition, the tree height, diameter at breast height (DBH), crown width, and leaf size of each tree were measured on site. The tree height was measured using a Vertex IV moisture meter, and the DBH was measured on the outside of the bark at 1.37 m above the ground using a diameter ruler. Table 1 lists the measured height, crown width, number of scanning points, and leaf size of the three experimental trees.

[0107] Table 1 shows the growth properties of the three studied trees and the number of scan points collected:

[0108]

[0109]

[0110] 1.2 Data preprocessing:

[0111] The raw point cloud data were processed using the branch-leaf separation algorithm to achieve the classification of scanned points into leaves and wood, and branches were modeled by a sequence of fitted cylinders with adaptive radius assignment (Huang et al., 2020) for wood volume estimation. Since the extracted leaf point cloud contained numerous leaves, our single leaf extraction algorithm (Xu et al., 2019) was applied. The results of wood-leaf separation, branch segment reconstruction, and single leaf data extraction are shown in Figure 2 middle. Figure 2 Results of branch-leaf separation, branch reconstruction, and single leaf extraction from terrestrial laser scanning data of three experimental trees. Figure 2 (a) is a schematic diagram of the branch and leaf separation results of the ground laser scanning data of Crape Myrtle. Figure 2 (b) is a schematic diagram of the branch segment reconstruction results of the ground laser scanning data of Crape Myrtle. Figure 2 (c) is a schematic diagram of the branch and leaf separation results of cherry blossom ground laser scanning data. Figure 2 (d) is a schematic diagram of the branch segment reconstruction results of the ground laser scanning data of cherry blossoms. Figure 2 (e) is a schematic diagram of the branch and leaf separation results of Magnolia terrestrial laser scanning data. Figure 2 (f) is a schematic diagram of the branch segment reconstruction results of Magnolia terrestrial laser scanning data. Figure 2 (g) is a schematic diagram of the data extraction results of a single leaf of Crape Myrtle from ground laser scanning data. Figure 2 (h) is a schematic diagram of the data extraction results of a single leaf from the ground laser scanning data of cherry blossoms. Figure 2 (i) in the figure is a schematic diagram of the data extraction results of a single leaf of Magnolia terrestrial laser scanning data. Figure 2 As shown in (a), (c) and (e) in Figure 1, dark grey and light grey represent the classified leaves and branches, respectively. The reconstruction results of the branches fitted using the cylinder model are shown in Figure 2 (b), (d) and (f) in Fig. 3. The extracted individual leaves are coded with the same color, e.g. Figure 2 As shown in (g), (h) and (i) in Figure 3, enlarged images of some extraction results are shown on the right side of each figure. Figure 2 (a), (b) and (g) are crape myrtle trees. Figure 2 (c), (d) and (h) are cherry trees. Figure 2 (e), (f) and (i) are magnolia trees.

[0112] 1.3. Physical and equivalent blade thickness:

[0113] The physical thickness of leaves is usually several millimeters, which is related to the number of anatomical features (such as cell walls and chloroplasts), biomass allocation, net productivity and crop response to water shortage. In addition, leaves can reflect, absorb or transmit light, perform photosynthesis to utilize light energy, or release water vapor through the movement of airflow around them for transpiration. However, sunlight or airflow always changes with different light emission angles, incoming airflow directions and time of day. In addition, the zenith angle and azimuth angle of leaves vary with different growth stages, environmental conditions and physiological characteristics. For a single leaf, the efficiency of photosynthetic energy conversion or water evaporation is related to the included angle between the incident solar beam or airflow and the normal vector of the leaf surface. In principle, when the incident light is perpendicular to the leaf surface, there is the largest illuminated leaf area or most stomata join the interaction with the airflow, and the photosynthesis or transpiration efficiency is the highest. Conversely, when the incident sunlight or airflow is parallel to the leaf surface, the photosynthesis or transpiration efficiency is the highest. So a question raised here is, when the sunlight direction or airflow is parallel to the leaf surface, only a few millimeters of leaf thickness is covered by sunlight or against the airflow. At that time, does the leaf provide enough photosynthesis area or too few stomata to absorb light and release water vapor?

[0114] In fact, there are different degrees of distortion on the surface of leaves. The surface of leaves is often twisted, curled, drooped, or even spiraled. There can be many reasons for these phenomena, such as water loss due to evaporation, continuous rainless weather, sunburn in summer or hot weather, pests eating leaves, physiological diseases, excessive spraying of pesticides on leaves, and uneven veins. Here, we infer that these distortions may be due to natural selection that causes the leaf surface to absorb more ambient light and optimize respiration and transpiration under the action of airflow from different directions. Therefore, when viewed from the side parallel to the plane of the leaf surface, that is, the side ( Figure 3 (a1), (b1) and (c1) in the figure), or parallel to the midrib ( Figure 3 When observing the blade in (a2), (b2) and (c2) in Fig. 2, the observed thickness (here called equivalent thickness) is greater than the physical thickness of the blade, such as Figure 3 shown. Figure 3 Different leaves of three tree species, namely (a) Crape Myrtle, (b) Sakura, and (c) Magnolia, are shown, where the equivalent thickness caused by leaf curling and drooping can be seen. We randomly sampled 100 leaves uniformly from the crown of each tree, and the range of manually measured equivalent thickness of leaves belonging to the three target trees from side views and parallel views to the midrib is shown in the last column of Table 3.

[0115] 1.4. Calculate the equivalent thickness and volume of each leaf:

[0116] TLS can obtain dense point clouds that can well describe the phenotypic characteristics of spatial objects. By registering two opposite scans, a comprehensive point cloud of each leaf in the canopy is obtained; here, a computer graphics-based method is used to estimate the equivalent thickness and volume of each leaf.

[0117] 1.4.1. Leaf surface point cloud fitting based on MLS:

[0118] During the scanning process, trees are blown by the wind, causing the vibration of leaves and branches, which will produce ghost points of position deviation of scanned leaves, inevitably leading to error propagation in the final leaf area estimation. Therefore, it is necessary to constrain each leaf to the best fitting surface to eliminate ghost points. Here, the moving least squares (MLS) method is used to process each leaf point after the single leaf point extraction process in Section 1.2.

[0119] MLS is a method for reconstructing a continuous function from an unorganized set of point samples by computing a weighted least squares solution around the point that is biased towards the desired reconstruction value. MLS is used to perform an efficient fit by introducing a subdomain of the local circular neighborhood to preserve the detailed geometric features of the leaf. On the local subdomain of the fitting area, we assume is the i-th scanning point of the blade. In order to eliminate the noise caused by the vibration interference of the blade, each point p i of Needs to be recalculated. Retrieve the new value The fitting function is expressed as follows:

[0120]

[0121] in It is p i There are m neighboring points p in the circular neighborhood centered at j x,y coordinates, where j = 1,2,3...,m. Basis function is with p i The ξ×1 column vector related to the neighboring points of , usually has two forms, for example, a linear basis with ξ = 3 elements and a quadratic basis with ξ = 6 elements Each point The undetermined coefficients of the ξ×1 column vector The best choice of coefficients is always to minimize the function (2) defined below:

[0122]

[0123] where p j ,j=1,2,...m is p i Neighborhood point, ω(||p i -p j || 2 ) is based on two points p i and p j The weight function with compact support property of the Euclidean distance of , where s = ||p i -p j || 2 and where r is the radius of the collection of circular adjacency domains; then,

[0124] In order to minimize formula (2) and obtain The best value of J relative to The derivative of is processed as follows:

[0125]

[0126] in is a and The ξ×ξ matrix, is a and A ξ×m matrix. is a matrix of size m×1. The size of is m×ξ, The size is m×m.

[0127] Finally, the coefficient is Finally, each point p i The fitted value of It can be calculated as Regression results After noise reduction and elimination of wind direction fluctuations, the true morphological structure of each leaf is restored.

[0128] Figure 4 The 3D modeling process of some leaves of different tree species is shown, using a hexagonal prism envelope to calculate the corresponding equivalent volume of each leaf. Figure 4 (a), (d), and (g) show the process of modeling the crape myrtle leaf. Figure 4 The original scan points represented in green in (a) are fitted using the moving least squares method to remove noise and smooth the surface, and the fitting results are represented by the red points. Figure 4 The black line in (d) indicates the midrib L 1 , green line L 2 Perpendicular to the midrib 1 , green line L 2 The leaf surface is divided into many segments. The edge points of each segment on both sides refined using polynomial curve fitting are shown in blue. Figure 4 The hexagonal prism envelope constructed in (g) is modeled based on the refined edge points, and its height is equal to the distance between the highest and lowest fitted leaf points along the leaf surface normal vector. It properly encloses the fitted leaf points and represents the equivalent volume of the leaf. Figure 4 (b), (e), and (h) show the same modeling process for cherry blossom leaves. Figure 4 (c), (f), and (i) show the same modeling process of magnolia leaves.

[0129] 1.5. Determination of precise blade edge points based on polynomial curve fitting:

[0130] After using MLS to remove noise and obtain points on the smooth leaf surface, the depiction of the true left and right edges of each leaf surface is necessary for leaf model construction.

[0131] The two endpoints p of the leaf midrib e (x e ,y e ,z e ),p s (x s ,y s ,z s )∈PMLS Forming the central axis L 1 , which has the longest distance through the leaf surface, can be easily identified, where p e is the point at the tip of the petiole, p s It is the top vertex of the petiole (i.e. the top of the midrib). Since infinite orthogonal lines can be drawn for points on a straight line in 3D space, the dimensionality reduction strategy of projecting points is adopted here to place each leaf on the XY plane, and only the x and y values ​​of each point are used to define the central axis L perpendicular to the XY plane. 1 The purpose is to find the corresponding left and right edge points of each leaf. Therefore, the leaf surface is projected onto the central axis L of the XY plane. 1 Described as follows: L′ 1 :p=p e (x e ,y e )+t×(p s (x s ,y s )-p e (x e ,y e )), where t is the slope of The independent variable of the straight line. 1 Uniformly sampling n-1 points on the top can define a series of vertical lines L 2 To L′ on the XY plane 1 , the slope is The intercept values ​​are different. Together with L 1 The two endpoints p e and p s Together, through the vertical line L 2 Divide the leaf point cloud into multiple segments at equal distances.

[0132] Then, by back-projection transformation from the two-dimensional plane to the three-dimensional space, we can determine the 2 The x, y and z values ​​corresponding to the left and right endpoints (scan points) of each segment of the direction are represented by the left endpoint and the right endpoint Where k = 1, 2, 3...n-1, such as Figure 4 Similarly, due to the influence of the environment and the occlusion in the canopy, the quality of the scanned data is aggravated. The connection of the endpoints (scanning points) on both sides of the leaf usually forms an irregular curve, resulting in a rough edge contour. Therefore, the cubic polynomial curve fitting method is used to fit the independent variable y to all edge points on the left. k With the dependent variable x k and z k The nonlinear relationship between them. The specific formula is as follows:

[0133]

[0134] in The left edge point Fitting The mean square error, The left edge point Fitting In order to minimize the mean square error, the partial derivative of equation (4) is obtained, which is recorded as: and Where b = 0, 1, 2, 3, the normal equation is obtained according to the method of finding the extreme value as follows:

[0135]

[0136]

[0137] Then find the coefficient and Therefore, after fitting and The values ​​are expressed as: and Similarly, the same method is used to fit the edge points on the right. Finally, the fitting value is obtained. and Accurately express left and right edge points exist Figure 4 In (d), (e) and (f), the optimized edge points are highlighted in blue, which gives a better representation of the leaf edge morphology.

[0138] 1.6. Use a hexagonal prism of equal thickness to surround each leaf:

[0139] For each leaf with n-1 refined left and right endpoints, we select two points on each edge with counts of round((n-1) / 4) and round(3(n-1) / 4) in top-down order, denoted as and The two endpoints of the midrib e and p s The six vertices that together form a hexagonal base, such as Figure 4 As shown in (g), (h) and (i) in the figure, round() represents the number in the brackets and is expressed to the nearest integer.

[0140] Therefore, we use the covariance matrix to calculate the normal vector of each leaf surface based on the fitted points. Then, along the direction of the calculated leaf normal vector, the highest point and the lowest point The distance betweennormal is used as the equivalent thickness of the current leaf to describe the curling, twisting, and drooping of the leaf surface. normal It is used as the height (thickness) and combined with the above hexagonal base to form a hexagonal prism, which can properly wrap the fitting point cloud P of a single leaf. MLS The hexagonal base of the hexagonal prism can be divided into three parts, namely the upper triangle, the middle trapezoid and the lower triangle, plus the equivalent thickness dist normal Therefore, the equivalent volume of different degrees of leaf deformation based on the wrapped hexagonal prism can be calculated Figure 4 (g), (h) and (i) are experimental images of point clouds of different individual leaves surrounded by hexagonal prisms, where the black arrow, blue star and red star represent the normal vector, the highest point and the lowest point along the normal vector direction of the leaf, respectively.

[0141] 1.7. Calculation of porosity of tree crown:

[0142] The equivalent volume of each leaf in the crown after the hexagonal prism with automatic adaptive parameter assignment is appropriately enclosed for each leaf in the crown. can be quantitatively calculated. In addition, the volume of wood It can be obtained from a set of fitted radii of cylinders along the branch direction, such as Figure 2 As shown. Crown volume V canopy It can be obtained using the 3D Alphashape algorithm based on scanned points of vegetation elements above the lowest branch height. According to the definition of void fraction (porosity) of porous media, it describes the ratio of pore space in the material, where pores may contain air or water. The volume-based GF of the canopy is defined as:

[0143]

[0144] The meaning of formula (7) shows that the volume-based GF can be calculated as 1 minus the ratio of the volume occupied by leaves and branches to the total volume of the crown, where and represent the total number of leaves and branches in the study crown, respectively.

[0145] 2. Results:

[0146] 2.1. Porosity retrieval using the volume-based method of the present invention:

[0147] Figure 5(a), (b), and (c) in Figure 3 show the results of each processing step for three trees, namely Crape Myrtle, Cherry Blossom, and Magnolia. This includes individual leaf point segmentation, MLS fitting to remove noise, and generating hexagonal prisms to envelop the leaf points. Hexagonal prisms are fitted to the point cloud of each leaf, so that the equivalent volume of each leaf can be calculated. Brown cylinders along the growth direction of the branches are assembled to represent the skeleton of the tree, approximating the wood volume of the crown. The overall volume of the crown is calculated based on the point set of the crown using the Alphashape convex hull algorithm. Figure 5 Illustrative results for the three tree species are shown in (d), (e), and (f) in Figure 2. Equation (7) uses these volumes to calculate the volume-based porosity at the scale of individual trees based on porous media theory. Table 2 lists the detailed results of our calculations for the three studied trees.

[0148] Figure 5 Each leaf in (a), (b) and (c) is wrapped by a hexagonal prism. Figure 5 The crown enclosures in (d), (e) and (f) are based on the Alphashape algorithm and our method is used to assess the porosity of each crown. Figure 5 (a) and (d) in the figure represent Crape Myrtle trees. Figure 5 (b) and (e) in the figure represent cherry trees. Figure 5 (c) and (f) in Fig. 3 represent magnolia trees, where the results in the dashed boxes are enlarged to illustrate the effectiveness of our algorithm.

[0149] Table 2 shows the inversion results of crown porosity and related branch and leaf parameters of the three experimental tree species:

[0150]

[0151] P. crape myrtle is a small shrub with small leaves (5.16 ± 2.05 cm long and 3.64 ± 1.32 cm wide). Each side of the lamellae always bends inward from the midrib, resulting in a small equivalent thickness (1.02 cm on average). Due to the sparse distribution of leaves and the spread of the crown, the total leaf space volume occupied (0.018 m3) and the estimated wood volume (0.014 m3) were the smallest among the three trees, resulting in a small LAI (0.81) and a higher volume-based GF (0.989). The studied cherry trees were in the mid-deciduous stage, preparing for dormancy; therefore, a small number of medium-sized leaves (2870 leaves, 9.98 ± 3.45 cm in length, and 5.36 ± 1.83 cm in average width) and some leafless branches were found in the crown. Many curled and senescent leaves were yellow and drooping, exhibiting a medium equivalent thickness of 1.10 cm and a low volume-based GF of 0.968. In the early stages of Magnolia growth, the crown is cylindrical and the leaves are large (average length 14.28 ± 5.06 cm, average width 6.32 ± 2.47 cm). The larger Magnolia leaves, which are reinforced or supported by fiber bundles, are composed of overlapping thick-walled cells, forming hard and rough leaves that are relatively twisted, resulting in a greater equivalent thickness (average 1.89 cm). Due to the small number of leaves present in the crown, the total equivalent leaf volume is 0.051 m3, so the GF value based on the medium volume is 0.978.

[0152] 2.2 Comparison with existing methods:

[0153] Two typical GF retrieval methods based on HP processing and light transport model are performed on the same three experimental trees to provide a quantitative evaluation of our results with existing methods. Figure 6 .according to Figure 6 The simulated HPs generated from the scan data of the Lambert azimuthal conformal projection shown in (a), (b) and (c) above calculate the GF of each annular sector of the HP corresponding to different zenith angles, i.e., the ratio of the number of sky pixels to the total number of pixels. Figure 6(a), (b), and (c) in Figure 3 show the GF calculated at 20° intervals. Generally, the GF at a zenith angle equal to 57.3° is considered to be a representative indicator for describing the GF of target trees. Here, the GF values ​​obtained from the HP with a zenith angle of 57.3° (0.693 for Crape Myrtle, 0.590 for Sakura, and 0.762 for Magnolia) were smaller than those of our volume-based method. A plausible explanation for this is that the GF retrieved from the HP is derived from the overlapping leaf and branch materials in the observation from the root to the treetop after the Lambert azimuthal equal-area projection transformation, and the original spatial pore distribution is reduced; in addition, the intrinsic spatial geometric information of the crown is compressed into a two-dimensional domain, with considerable loss of three-dimensional traits, and limited representation of plant phenotypic characteristics. Although qualitatively similar results were obtained from the HP method, such as the ranking of tree GF remained unchanged, the difference in methods produced different results at different scales.

[0154] Another approach to GF estimation is based on voxelization and a model of light transport through the canopy. The point cloud of the canopy is subdivided into numerous cubic voxels. A ray tracing algorithm that calculates the distance each ray travels through the canopy before being intercepted by plant elements is used to calculate the total leaf area. Figure 6 As shown in (d), (e) and (f) in the figure, the GF at the voxel scale, i.e., the canopy volume minus the total volume of all voxels containing the tree scan points, is 0.932 for the crape myrtle tree, 0.903 for the cherry tree, and 0.965 for the magnolia tree, which is the canopy volume represented by the semi-transparent red convex hull. In addition, the color of each voxel is assigned according to the ratio of the number of points contained in the current voxel to the number of the largest points in the current canopy voxel. Most voxels are blue or light blue, indicating that the voxel with the largest number of points only represents 10-30% of the points contained. This shows that the leaf elements are evenly and sparsely distributed in the canopy space, and a large proportion of light is not intercepted by the leaf element points in each voxel, resulting in a higher GF value based on the light transport model. Figure 6 (g), (h) and (i) in the figure show the GF sizes retrieved based on the light transport model for the three trees, and their values ​​are roughly the same as those obtained by the volume-based method of the present invention. The reason is that the light transport model takes into account two key factors, namely the vegetation material contained in each spatial voxel and the average distance of the vegetation elements hit or missed when the emitted laser beam passes through the canopy, which is related to the geometric exploration of the canopy and reflects the leaf element indicators in each voxel. Since the heterogeneous canopy produces different degrees of occlusion from different observation angles, the GF calculation method based on the light transport model produces different evaluation results as the azimuth and zenith angle of the laser scanning emission position relative to the target tree changes. Therefore, the results show a two-dimensional matrix form for each tree, such as Figure 6As shown in (g), (h) and (i) in the figure, both the zenith angle and the azimuth angle vary. The GF based on the light transmission model obtained from the side of each tree, i.e., the vertical section with a larger zenith angle and more spatial information, still shows that the GF value of the cherry tree is the lowest (0.933), while the GF value of the crape myrtle tree is higher (0.953), which is qualitatively consistent with the results of the present invention.

[0155] Figure 6 The results of the method of the present invention are compared with those of existing porosity retrieval methods for two studied trees, namely, crape myrtle, cherry blossom, and magnolia trees. Figure 6 The porosity in (a), (b) and (c) is calculated from the generated HP by Lambert azimuthal conformal projection of the scan tree points. Figure 6 (d), (e) and (f) are 3D visualizations of the three studied trees based on cubic voxels, where the different colors of each voxel represent the ratio of the number of scan points in each voxel to the maximum number of scan points in a voxel within the canopy. The figure shows the disparity score at the voxel scale for each canopy, i.e., the disparity between the canopy volume minus the total volume of all voxels containing the tree scan points and the canopy volume. Figure 6 (g), (h) and (i) show the calculation result matrix of the porosity of each crown based on the light transmission method, which varies with different azimuth and zenith angles of the emitted light. When the zenith angle of the emitted light is large, its value is roughly equal to that of the volume-based method of the present invention.

[0156] 3. Discussion:

[0157] 3.1. Different types of GF estimation methods:

[0158] The pioneering research on GF retrieval of tree crowns mainly focuses on tree crown image processing technology or light transmission algorithm. The image processing-based methods derive the below-canopy images, for example, captured digital hemispherical images or synthetic HP generated from scanned points, combined with image segmentation algorithms to calculate the proportion of sky pixels in the image, which provides a convenient way for GF calculation applied to single trees or forest maps, but the dimensionality reduction from spatial plant material identification to 2D projection image segmentation cannot preserve the original geometric characteristics in the reduced space and cannot provide an ideal recovery to the intrinsic 3D dimension of the plant.

[0159] The GF method based on the light transport model is built on the Bouguer-Beer-Lambert law, which relates the attenuation of light propagating through the properties of a homogeneous material. The original Beer's law is often applied to chemical analysis measurements or to understand attenuation in physical optics, but the performance of Beer's law is always compromised when the medium that attenuates the light is inhomogeneous and turbid, resulting in multiple scattered photons and requiring a correction term to the Lambert-Beer law. Similarly, when Beer's law (8) is used to calculate the GF P(θ) of the canopy for turbid media, there are uncertainties and other obstacles as described below.

[0160] P(θ)=exp[-G(θ)·Ω·LAI / cosθ] (8).

[0161] In formula (8), the extinction coefficient G(θ) depends on the leaf inclination and azimuth range of the incident laser beam, with a magnitude of 0.5 and a standard deviation within 0.1. The clumping index Ω describes the deviation of the leaf arrangement from random distribution, and its magnitude ranges from 0.3 to 0.8 for coniferous forests and 0.3 to 0.8 for broadleaved forests. In general, the phenotypic traits of trees vary greatly, with irregular and discontinuous crowns, different leaf orientations, and multi-factor driven leaf clumping states, which exacerbates the accurate magnitude estimation of the two key parameters G(θ) and Ω in formula (8), suggesting that biased results may appear in the final estimated GF P(θ). In addition, the GF results obtained by the light transport model are affected by the light emission position and incident angle of the laser scanner with a limited field of view. Any fluctuations in the laser beam divergence, scanning resolution, and reflectivity of different object surfaces will lead to perturbations and error propagation in the final crown GF estimation.

[0162] From the perspective of porous media theory, GF (porosity or void fraction) is a concept that describes the ratio of the volume of pores or void space to the total volume of the solid. It is widely used in fields such as geology, hydrogeology and soil science. The volume-based GF is a fraction between 0 and 1, usually ranging from less than 0.005 for solid granite to more than 0.5 for peat and clay (Dutta and Tarafdar, 2003). Other earth materials also have their measurable volume GF or porosity. For example, sand has spaces between its particles, and its average volume GF is 37.7%, 42.3% and 46.3% in packed, natural (in situ) and loose packed conditions, respectively. Expanded polystyrene foam used for insulation usually contains 3% polystyrene and 97% fixed air. Tree crowns grow on the earth. The crown volume is relatively large, containing thin leaves and not much woody material. According to the unified framework of volumetric GF of porous materials in the world, tree crowns should have a high porosity magnitude compared with other materials on Earth because of their relatively large volume and many voids that allow light, air, and liquid to pass through. At the same time, unlike underground rocks or sediments that require sonic and neutron logging to measure the porosity of clay-rich coal rocks, and unlike industrial and daily materials that use scanning electron microscopy for porosity quantification, a framework based on visual trait heuristics and computer graphics driven vegetation element modeling should be explored to explore insights about volume-based canopy GF determination. Quantitative assessment of spatial occupancy of plant materials in canopies and refined geometric polyhedral representations will promote the exploration of canopy volume-based GF and be consistent with existing evaluation standards for porosity of most materials on Earth.

[0163] 3.2. Application of the method of the present invention in large trees and forest land:

[0164] To test the applicability of the method of the present invention to large trees and forest land, we applied the method of the present invention to two larger trees: a rubber tree located in a rubber tree plantation in Danzhou City, Hainan Island, China, with a tree height of 16.15m, a breast diameter of 0.21m, a north-south width of 6.86m and an east-west width of 6.24m, and a ginkgo tree on campus with a tree height of 7.34m, a breast diameter of 0.11m, and a crown width of 4.43m in the north-south direction and 3.31m in the east-west direction. Compared with small tree crowns, large tree crowns contain more leaves, which will produce more data loss and shadow effects when collecting scan data using ground-based LiDAR, resulting in lower accuracy in the separation of individual leaves and longer program execution time (6-10 minutes).

[0165] Rubber trees are tall perennial tropical trees in the upper layers of tropical rainforests with compound leaves, each of which is usually composed of three leaflets connected to the middle vein but with their own petiole. We randomly sampled 100 leaflets from the crown of rubber trees of clone CATAS 7-20-59 and manually measured the equivalent thickness of the leaflets, which was 0.3-1.4 cm in the side view and 0.8-3.0 cm in the view parallel to the midrib, as well as the length of the leaflets of 11.1-24.6 cm and the width of 3.9-10.2 cm. The results of the segmentation of a single leaf of rubber tree are shown in Figure 7 In (a), the local magnification result is the same as Figure 7 The constructed hexagonal prisms shown in (b) of FIG. 1 are paired to surround each leaf point, where the relatively highly accurate segmentation results and the geometrically clear leaflet layers presented by the hexagonal prisms are illustrated. Although there are slight over-segmentation and under-segmentation problems with individual leaflets, the overall visual effect of the rubber tree leaflets is satisfactory. The results of the method of the present invention show that in the larger crown, a total of 11,633 leaflets with an equivalent thickness of 0.42-2.97 cm occupy an equivalent volume of 1.71 cubic meters (volume of 121.17 cubic meters), and the wood volume is 0.56 cubic meters, which results in a volume-based GF of 0.981 for the rubber tree crown.

[0166] Figure 7 The experimental results of applying the method of the present invention to two large trees are shown. Figure 7 (a) and (c) show the leaf separation results of rubber tree and ginkgo tree, respectively, and the segmented individual leaves are highlighted in different colors. Figure 7 Close-ups of four local regions of each tree crown in (b) and (d), showing extracted individual leaf points and the corresponding enclosed hexagonal prisms, shown separately and together.

[0167] The fan-shaped leaves of Ginkgo are drooping with irregular notches. They usually have deep grooves in the middle of the leaves, producing two distinct lobes. Ginkgo trees have many forked and parallel fine veins, which cannot support the leaves well, causing the leaves to fall from the middle to both ends, making the Ginkgo leaves appear soft and drooping. Similarly, 100 leaves in the crown of the Ginkgo tree were randomly sampled, and the manually measured equivalent thickness of the leaves was 0.2-1.6 cm in the side view and 0.9-3.3 cm in the view parallel to the midrib of the leaf. The leaf length was 2.8-6.6 cm and the leaf width was 4-10.5 cm. The crown of the Ginkgo tree is conical, and the leaf clusters grow near the first-order or second-order branches, which leads to the dominant occlusion effect in the upper part of the crown, and many branch data cannot be correctly collected from the ground view. Therefore, the quality of the scanned data of the crown deteriorates as the tree height increases. Due to the unique fan-shaped structure of the Ginkgo leaf, we used a hexagonal prism base consisting of a larger upper triangle, a smaller lower triangle, and a middle trapezoid connecting them to align the leaf shape, e.g. Figure 7 As shown in (c) and (d) in . Figure 7 As can be seen from (c) and (d) in the figure, the point density in the upper part of the tree canopy is reduced to a certain extent, which has some adverse effects on the leaf segmentation and visualization results using the method of the present invention. Many leaves that are identified as different colors due to incomplete point clouds caused by occlusion are also appropriately surrounded by hexagonal prisms ( Figure 7 The results obtained by the method of the present invention show that in a medium-sized crown (volume of 18.86 cubic meters), there are 12,654 leaves with a thickness of 0.23-3.2 cm occupying an equivalent volume of 0.52 cubic meters, and the wood volume is 0.101 cubic meters, which leads to a volume-based GF of 0.967 for the Ginkgo canopy.

[0168] A dense forest consisting of 24 rubber trees of clone PR 107 was selected ( Figure 8 (a) in the figure, with a row spacing of approximately 3 m and an inter-row spacing of 7 m, is used to calculate the volume-based GF. The calculated average height and crown volume of each tree are 19.51 m and 110.47 m3, respectively. The observation position is set at the ground center of the plot. After the branch and leaf separation and individual leaf extraction processes, the synthetic HP of the forest area is generated using the spherical projection of the scan points, from which the HP-based GF retrieved is 0.284 ( Figure 8(b) in the figure). The average number of detected leaves, the average calculated total equivalent leaf volume, and the wood volume of each tree are 13120, 2.02 cubic meters, and 0.53 cubic meters, respectively. Then, the volume-based GF retrieval of the entire plot using the method of the present invention is 0.976. The point cloud of the entire plot is obtained by registering many ground scans at different scanning positions, which leads to an uneven point density distribution due to the influence of mutual occlusion or multiple scans of local areas, the anisotropic crown shape resulting in a weakened ability to record details, and the laser beam is not fully covered. Some ghost points originate from the shaking of leaves in the breeze. Therefore, from Figure 8 As can be seen in (b) and (c), some leaves in the upper part of the forest canopy have low density and incomplete scan points, resulting in a fragmented geometric approximation, leading to a slight overestimation of our volume-based GF, while in some leaf clusters at lower tree heights, leaves are highly concentrated and ghost points exist, producing spurious larger-sized hexagonal prisms, leading to a subtle underestimation of the volume-based GF. Therefore, due to the limitations of the scanner's technical specifications and the actual complex forest environment, we expect that the volume-based GF of dense forest plots will have slight deviations.

[0169] Figure 8 The calculation of volume-based porosity for the rubber tree plot is shown. Figure 8 (a) Image and scan points of a dense forest plot of Hevea brasiliensis clone PR 107 with interlaced branches and clumps of leaves. Figure 8 (b) The plot center is used as the projection center of the spherical projection, and the scan points are converted into a synthetic fisheye photo (HP) for the estimation of crown openness, i.e., the gap score based on the HP. Figure 8 (c) in the figure uses the algorithm of the present invention to magnify part of the tree crown, highlight the identified individual leaves with different colors, and enclose them by a translucent hexagonal prism.

[0170] 3.3. Advantages and limitations of the method of the present invention:

[0171] The purpose of the method of the present invention is to move from two-dimensional image analysis or coarse geometric description to fine three-dimensional spatial element rendering for crown GF retrieval. In recent years, researchers have made efforts to describe the phenotypic characteristics of vegetation and detailed crown structure. Some progress has been made, such as using regular geometric models such as ellipsoids, umbrellas, cylinders or cones to comprehensively simulate the crown, using the crown as a pseudo-turbulent medium to describe the internal content and evaluate solar radiation, using voxel-based representation for spatial vegetation tissue distribution, and capturing images from different observation angles to obtain multi-angle profiles to analyze material agglomeration within the crown. It has been proposed to design robust algorithms using computer graphics technology to achieve segmentation of individual leaves from point clouds, enhance the ability to outline the detailed features of the real crown, and enable computers to analyze individual leaves and branches separately. In addition, the crown of a tree growing on the earth has a porous structure, which has a decisive impact on the ventilation and light transmission of the tree body. According to the unified framework of porous media theory to evaluate the natural and artificial materials in the world, it is reasonable to treat the crown as a porous medium to calculate the GF size of its volume. At the same time, the concept of leaf equivalent thickness is proposed here to describe the degree of curling and drooping of leaves, and a tightly closed hexagonal prism is used to quantify the leaf volume. Wood volume and crown volume can also be evaluated using the fitted cylinder and convex hull algorithms, respectively. Therefore, our method based on porous media theory and computer graphics technology can reversely calculate the crown volume GF, which overcomes the limitations of existing crown volume estimation methods, such as the GF method based on image processing, which compresses the original large-area canopy space (2651.3 cubic meters) onto a two-dimensional HP image ( Figure 8 (a) and (b) in the figure, the undesirable overlapping effect of plant materials is greatly exaggerated, or the performance of GF methods based on light transport models is always limited by the limited field of view of local canopy data at fixed locations, resulting in narrow information acquisition. With further development, our method will show the possibility of revealing the biophysical and ecological relationships between volume-based GF and other spatially explicit growth characteristics, such as the topological skeleton of trees, crown shape, spatial leaf distribution, and leaf area density.

[0172] The method of the present invention has certain defects. The accuracy of the algorithm for separating individual leaves will produce certain deviations in the final equivalent leaf volume calculation. A higher scanning resolution is required to ensure that each leaf is covered with at least 40 scanning points, that is, a point cloud sampling spacing of about 0.5 cm is preferred. If the density of the collected crown scan points is very sparse, making it difficult for humans to distinguish each leaf by visual inspection, the method of the present invention cannot cope with the segmentation of each leaf in this case. Fortunately, many well-known brands of fixed terrestrial laser scanners, such as Leica, RIEGL or FARO, currently have mature technologies to meet our density requirements. In addition, users can adjust the scanning mode by reducing the scanning distance or aligning multiple scans in the parent coordinate system to improve the resolution of the tree point cloud. In addition, as described in Section 3.2, when the leaf area density of large trees or forests is high, the occlusion effect is emphasized in the dense crown. Ground observations always miss the information in the middle and upper parts of the crown of large trees, resulting in the impairment of the ability to identify leaves and accurate hexagonal prisms using the method of the present invention, and resulting in unfavorable deviations in the final volume-based GF calculation results. Meanwhile, the method of the present invention is very time-consuming compared to the GF methods based on image processing, which only carry a hemispherical camera to capture some HPs in the forest, and combine any existing software to segment the sky pixels from the HPs for HP-based GF calculation. In contrast, the extraction of individual leaves and the calculation process of the appropriate size of the corresponding hexagonal prism envelope will take more time, and finally the display of numerous individual leaves and the corresponding envelope rendering need to utilize the graphics processing unit (GPU) to accelerate the creation and operation of the generated fine forest scene model, such as Figure 8 (c) In addition, for leaves with irregular and complex shapes, such as golden disks and conifer needles, it is not convenient to choose the right geometric primitives to envelop the leaves. If the crown is under different climatic conditions, such as hurricanes or precipitation, the compression of the tree structure caused by wind or rainstorm will cause our volume-based GF estimation to decline, but it is not easy to quantify this because it is difficult to collect scanning data that record the instantaneous compression state of the crown in a storm environment or during the vibration behavior of the crown in response to wind excitation. Recently, we are building a large fan simulation wind tunnel to provide continuous and stable wind flow mechanics under different speed loads for rubber trees, and combined with the laser scanning data collected at that time, to explore the geometric knowledge of the crown deformation and track the changes in the pore volume between the crowns. In general, with the upgrading of computer hardware, achieving excellent graphics processing capabilities, and the advancement of artificial intelligence algorithms empowering the application of forest surveys, we believe that in the future, the tasks of fine tree organ identification and plant phenotyping will become easier and more accurate, which will promote the development of precision forestry and explain the volume GF of trees in a coherent perspective.

[0173] 4. Conclusion:

[0174] This work describes a new measurement of volume-based porosity (GF) of tree crowns from terrestrial laser scanning (TLS) data based on computer graphics and porous media theory. In addition, the method produces several simultaneous measurements, such as the equivalent thickness and equivalent volume of leaf lamellae, which stem from the common physiological phenomenon that drooping and curling leaves produce considerable lateral area. By individual leaf segmentation and moving least squares (MLS) fitting, the distance between the highest and lowest fitting points of the normal vector along each leaf surface is taken as the equivalent thickness, coupled with a hexagonal prism formed by adaptive parameter assignment to properly enclose each leaf point and obtain the equivalent volume of each leaf. Then, the total volume of leaves in the crown and the trunk volume based on cylinder fitting are calculated, considering the crown as a multiphase porous medium, and then combined with the crown volume evaluated by the convex hull algorithm, the porosity of each crown, i.e., the volume-based GF, is estimated according to porous media theory. This calculation provides a uniform evaluation standard for crown GF, which, like the porosity assessment of most materials in the world, is considered an important property of materials in materials science, chemical engineering, soil mechanics and biomechanics. Five trees (three small trees and two large trees) and one forest plot were selected for the experiment, and the calculated volume-based GF ranged from 0.96 to 0.99. This range is approximately the porosity of expanded polystyrene foam, which has a similar porous structure that allows for large amounts of air to be stored and contains many small polystyrene beads, like tiny leaves within the crown. Overall, this study provides more evidence that the combined multidisciplinary approach is widely applicable to volume-based crown GF estimation from fine-scale phenotypic features inherent in the visual content of plant materials.

[0175] The protection scope of the present invention includes but is not limited to the above embodiments. The protection scope of the present invention shall be based on the claims. Any replacement, deformation, and improvement of the technology that can be easily thought of by technicians in this field shall fall within the protection scope of the present invention.

Claims

1. A method for estimating tree crown porosity based on porous media theory and computer graphics, characterized by: include: Step 1: Collect tree laser point cloud data; Step 2: Use the branch-leaf separation algorithm to realize the point cloud classification of leaves and branches, fit the branches with cylinders with adaptive radius allocation to model the branches, and extract the point cloud data of a single leaf; Step 3: Use the moving least squares method to denoise the point cloud data of a single leaf; Step 4: Using a cubic polynomial curve fitting method to determine the edge points of each blade; Step 5: Calculate the equivalent thickness of each leaf, and then use the hexagonal prism formed by adaptive parameter allocation to envelop the point cloud of each leaf to obtain the equivalent volume of each leaf; Step 6: Calculate the total equivalent volume of all leaves in the crown and the total volume of all branches based on cylinder fitting, and combine them with the crown volume to estimate the porosity of the crown.

2. The method for estimating tree crown porosity based on porous media theory and computer graphics according to claim 1, characterized in that: The step 1 is: Tree laser point cloud data was collected using a terrestrial laser scanner.

3. The method for estimating tree crown porosity based on porous media theory and computer graphics according to claim 2, characterized in that: The step 3 is specifically as follows: In the local subdomain of the fitting area, let point p i is the i-th point cloud of the leaf, For point p i coordinates, where i = 1, 2, 3..., n, in order to eliminate the noise caused by blade vibration interference, point p i of The value needs to be recalculated, and the new value obtained by recalculation is recorded as Retrieving New Values The fitting function is expressed as follows: in It is p i There are m neighboring points p in the circular neighborhood centered at j x,y coordinates, where j = 1,2,3...,m; basis function is with p i ξ×1 column vector related to the neighboring points; is the uncertainty coefficient. When formula (2) is the minimum, Take the best value: where p j Yes i Neighborhood points, j = 1, 2, ... m; ω (||p i -p j ||2) is based on two points p i and p j The weight function with compact support property of the Euclidean distance of , where s = ||p i -p j ||2 and where r is the radius of the collection of circular adjacency domains; then, 4. The method for estimating tree crown porosity based on porous media theory and computer graphics according to claim 2, characterized in that: The step 4 is specifically as follows: (4.1) The two endpoints on the midrib of the leaf are denoted as p e (x e ,y e ,z e ),p s (x s ,y s ,z s ), where (x e ,y e ,z e ) indicates p e The coordinates of (x s ,y s ,z s ) represents point p s The coordinates of the two endpoints form a central axis L1, p e (x e ,y e ,z e )∈P MLS ,p s (x s ,y s ,z s )∈P MLS , where P MLS It is the fitted leaf point cloud after denoising; (4.2) Project the leaf surface onto the XY plane. The central axis L1 becomes L′1 after projection. The point p on L′1 is: p = p e (x e ,y e )+t×(p s (x s ,y s )-p e (x e ,y e )), where t is the slope of The independent variable of the straight line; (4.3) Sample n-1 points uniformly on L′1, define multiple vertical lines L2 that are perpendicular to L′1 and pass through n-1 points respectively, and the slope of the vertical line L2 is Together with the two endpoints p of L1 e and p s Together, the blade point cloud is equidistantly divided into multiple segments by the vertical line L2; (4.4) Determine the left endpoint of each segment along the L2 direction by back-projection transformation from the two-dimensional plane to the three-dimensional space The corresponding x, y values, the corresponding z value and the right endpoint The corresponding x, y values ​​and the corresponding z values, where the left endpoint The corresponding x, y values ​​and the corresponding z values ​​are recorded as Right endpoint The corresponding x, y values ​​and the corresponding z values ​​are recorded as Where k = 1, 2, 3...n-1; (4.5) Using the cubic polynomial curve fitting method, a cubic polynomial is used to fit the independent variable to all left edge points. With dependent variable and The nonlinear relationship between them is as follows: in The left edge point Fit x′ k L The mean square error, The left edge point Fitting z′ k L The mean square error of Minimize the mean square error in formula (3), expressed as: and Where b = 0, 1, 2, 3; Then find the coefficient and Therefore, the fitted x′ k L and z′ k L The values ​​are expressed as: and Similarly, calculate the fitting value x′ corresponding to all right edge points according to the method in step (4.5) k R and z′ k R ; Finally, the accurate x, y and z values ​​of the left edge point are fitted and the exact x, y, and z values ​​of the right edge point Where k = 1, 2, 3...n-1.

5. The method for estimating tree crown porosity based on porous media theory and computer graphics according to claim 2, characterized in that: The step 5 is specifically as follows: (5.1) For each leaf with n-1 left and right edge points, select round((n-1) / 4) left edge points and round(3(n-1) / 4) left edge points from the n-1 left edge points in top-down order, and select round((n-1) / 4) right edge points and round(3(n-1) / 4) right edge points from the n-1 right edge points in top-down order. Place the two selected left edge points, the two right edge points, and the endpoint p e and endpoint p s As the six vertices of a hexagonal base, round() represents the number in the brackets to the nearest integer; (5.2) Calculate the normal vector of each leaf surface using the covariance matrix according to the fitted points, and select the highest point along the direction of the calculated leaf normal vector. and the lowest point highest point and the lowest point The distance between normal is taken as the equivalent thickness of the current blade; (5.3) dist normal It is used as the height and combined with the above hexagonal base to form a hexagonal prism, which in turn wraps the fitting point cloud of a single leaf; (5.4) Calculate the volume of the hexagonal prism, which is the equivalent volume of the current leaf 6. The method for estimating tree crown porosity based on porous media theory and computer graphics according to claim 2, characterized in that: The calculation formula of the porosity of the tree crown in step 6 is: Where N represents the total number of leaves in the crown, M represents the total number of branches in the crown, represents the equivalent volume of a leaf in the crown, Represents the volume of a branch in the crown, V canopy Represents the tree crown volume.

Images