SEGMENTATION AND CROWN PARAMETER EXTRACTION OF INDIVIDUAL TREES IN AN AIRBORNE TOMOSAR POINT CLOUD

The analysis of individual trees is an important field of research in the forest remote sensing community. While the current state-of-theart mostly focuses on the exploitation of optical imagery and airborne LiDAR data, modern SAR sensors have not yet met the interest of the research community in that regard. This paper describes how several critical parameters of individual deciduous trees can be extraced from airborne multi-aspect TomoSAR point clouds: First, the point cloud is segmented by unsupervised mean shift clustering. Then ellipsoid models are fitted to the points of each cluster. Finally, from these 3D ellipsoids the geometrical tree parameters location, height and crown radius are extracted. Evaluation with respect to a manually derived reference dataset prove that almost 86% of all trees are localized, thus providing a promising perspective for further research towards individual tree recognition from SAR data.


INTRODUCTION
The analysis of individual trees in remote sensing data until now has mainly focused on the exploitation of aerial imagery or Li-DAR point clouds.In this framework, many studies have been published about the detection and localization of individual trees (Pollock, 1996, Wulder et al., 2000, Leckie et al., 2005, Chen et al., 2006, Chang et al., 2013) as well as the delineation of their tree crowns (Culvenor, 2002, Pouliot et al., 2002, Erikson, 2003, Koch et al., 2006, Jing et al., 2012).In contrast to that, the analysis of forested areas on the single-tree level by means of synthetic aperture radar (SAR) remote sensing has not yet met the interest of the community, although modern sensors have reached submeter resolutions down to the decimeter-range in recent years.However, recently it has been shown that it is possible to generate almost fully layover-and shadow-free point clouds by means of airborne single-pass SAR tomography using millimeterwave sensors (Schmitt and Stilla, 2014).Therefore, in analogy to the approaches based on 3D LiDAR point clouds, in this paper an unsupervised segmentation of the TomoSAR point cloud aiming at the reconstruction of individual trees is proposed.While, for example, the studies of (Morsdorf et al., 2004) or (Gupta et al., 2010) are suggesting to employ k-means clustering for tree segmentation in LiDAR point clouds, in the presented work the unsupervised mean-shift clustering algorithm is used.This way the need to know the number of expected clusters and an initialization of their centers a priori is avoided and a fully automatic procedure is enabled (Comaniciu and Meer, 2002), which has already been proven for the reconstruction of buildings in To-moSAR point clouds (Shahzad and Zhu, 2015).After clustering, rotational ellipsoids are used to model the individual segments in order to approximate the tree crown shapes.From these ellipsoids the tree positions, heights and crown diameters can be extracted.This tree reconstruction strategy is evaluated using a 3D TomoSAR point cloud, which was generated from airborne millimeterwave InSAR data acquired from multiple aspect angles.

POINT CLOUD SEGMENTATION BY MEAN SHIFT CLUSTERING
The basis of the method proposed in this paper is the clustering of the 3D TomoSAR point cloud by the mean shift algorithm as described by (Comaniciu and Meer, 2002).Since the tree crowns generally show a comparably high point density, the points are clustered in the spatial domain, i.e. the feature space is comprised of spatial coordinates in Euclidian space.The kernel density estimate at any point pi of the n 3D points is given by the expression where b is the bandwidth parameter and g (x) is a non-negative, non-increasing, piecewise continuous function with definite integral, i.e. ∞ 0 g (x) dx < ∞.Based on the concept of kernels discussed by (Cheng, 1995) and (Comaniciu and Meer, 2002), the function g (x) is defined as the profile of the radially symmetric kernel G (x) satisfying where c is a normalization constant ensuring that G (x) integrates to 1. Different kernels, such as the Epanechnikov kernel and the Gaussian kernel can be used to define the density Dp i .Mean shift clustering essentially seeks modes of the kernel density estimates and works iteratively by shifting every data point toward the weighted mean of points within its neighborhood (defined to be cylindrical in the presented case).The shift vector m (pi) always points towards the direction of the maximum increase in the density Dp i and is computed as The iteration process continues until there is no or only little shift in m (pi) anymore, i.e. the length of the shift vector m (pi) is close to 0. Due to the gradient ascent nature, the mean shift algorithm returns clusters using the concept attraction of basin, i.e. those points whose trajectories lead to the same mode form the basin of attraction for that mode and are clustered into one group.The clustering procedure is repeated until all points are assigned to their respective modes.
Clustering via mean shift can be considered an unsupervised procedure since it does not require the number of clusters a priori, The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XL-3/W2, 2015 PIA15+HRIGI15 -Joint ISPRS conference 2015, 25-27 March 2015, Munich, Germany This contribution has been peer-reviewed.doi:10.5194/isprsarchives-XL-3-W2-205-2015nor does it need any pre-defined model for the shape of the resulting clusters.Nevertheless, it still does require a bandwidth parameter (corresponding to the radius of the kernel), which affects the number of clusters, i.e. the number of modes, that are returned by the algorithm.However, unlike other clustering algorithms such as k-means, fuzzy c-means, expectation maximization etc., the bandwidth parameter exhibits a physical meaning for variables in spatial coordinates and can be set based on prior knowledge as, e.g., the expected average radius of the tree crowns in the scene.

ROBUST MODELING OF ELLIPSOIDIC TREE CROWNS
Since three-dimensional rotational ellipsoids can be seen as a good approximation of deciduous tree crowns, the individual tree clusters are modeled using generalized tri-axial ellipsoids.For this purpose, parameters of an arbitrarily oriented minimum volume enclosing ellipse (MVEE) are estimated by first projecting points belonging to individual tree clusters onto the xy-plane followed by extruding the 2D xy-ellipse in z-direction to form a 3D ellipsoid.The motivation for expanding the ellipsoid along the z-axis is based on geometrical considerations: It is assumed that correct tree models may have an arbitrary orientation in the xyplane, but remain upright or vertical with respect to the ground (cf.Fig. 1).This is based on the light prior that tree trunks are modeled to be vertical to the ground surface.

Computation of the MVEE
If K = {ki|i = 1, . . ., m} denotes m clusters returned by the mean shift algorithm, and Q = {qu|u = 1, . . ., r} denotes the set of r points qu belonging to a particular cluster k f (f ∈ i), then any arbitrarily oriented ellipse ε can be a candidate for MVEE(Q), if and only if all points in Q lie on or inside its boundary, i.e. if the following condition is satisfied (Kumar and Yildirim, 2005): In this equation, A is a d × d positive definite matrix, where d refers to the dimension 2 in the presented case, and c k f is the center of the ellipse surrounding the clustered points Q.The semiaxes si of such an ellipse are given as where vi denote the eigenvectors of A, which correspond to the directions of the semi-axes.λi denotes the eigenvalues of A, which are related to the length of these axes: The length of each axis is equal to 1 √ λ i . The area of an ellipse or volume of an ellipsoid, respectively, is thus directly proportional to det 1 √ A .Therefore, in order to obtain an MVEE(Q), det 1 √ A has to be minimized such that (4) is satisfied in conjunction with A being positive definite.In order to solve this minimization, Khachyan's first order algorithm is used, which formulates the problem as optimization using Lagrangian duality (Khachiyan, 1996).
The computed MVEE(Q) is extended to the third dimension by extruding it in z-axis in order to form a 3D ellipsoid.The zcoordinate of the ellipsoid center and its semi-axis length s3 in where h min,k f ,i and h max,k f ,i (i = 1, . . ., N ) are the N lowest heights and the N largest heights of all points in the cluster k f , respectively.
Once this modeling is complete, the tree parameters tree height, crown diameter, and trunk location can directly be extracted from the ellipsoid model: The tree height is the maximum height of the ellipsoid in z-direction, the tree crown radii are given by the xand y-semi-axes of the ellipsoid, and the xy-coordinates of the ellipsoid center point provide the location of the tree trunk.Of course, this is a simplifying model only valid for deciduous trees of approximately ellipsoidal shape, but an extension toward a more general tree model as, e.g.described by (Sheng et al., 2001) basically seems possible.

Test and Reference Data
The input data for the experiments presented in this paper is a multi-aspect TomoSAR point cloud generated from an airborne millimeterwave SAR dataset acquired during a flight campaign over the city of Munich, Germany, in 2013 (Schmitt and Stilla, 2014).Using two opposing viewing directions, four images from four simultaneously receiving antennas per pass, and a discretization of 50 cm × 50 cm × 50 cm during processing, the resulting 3D point cloud consists of about 1.66 mio.points, corresponding to an average point density of about 22 points/m 2 .As comparison to helicopter-borne LiDAR point cloud showed (point density ca. 3 points/m 2 ), the three-dimensional localization accuracy of the TomoSAR points lies between 0.7 m and 1.4 m -negatively biased by systematic errors in the comparison due to different point densities.
The test scene consists of the "Alter Nordfriedhof", an abandoned cemetery, which is used as a public park today.As can be seen in Fig. 2, it is mainly characterized by a light planting of deciduous trees, resembling a grove or little wood.The TomoSAR point cloud is displayed in Fig. 3 (a).
Figure 2. Orthophoto of the test scene "Nordfriedhof" in Munich, Germany.
As a reference dataset, a helicopter-borne LiDAR point cloud containing approximately 0.16 mio.points (i.e. 3 points/m 2 ) in conjunction with a co-registered orthophoto was analyzed by a human operator, who extracted tree positions, diameters and heights manually.The result can be seen in Fig. 4. In total, the reference consists of 570 trees with average height of 14.56 m (median) or 12.33 m (mode), respectively, and an average radius of 3.70 m (median) or 3.32 m.

Clustering Results
Fig. 5 compiles the clustering results for varying kernel bandwidth parameter.Obviously, the optimal bandwidth parameter is 3.2 m, giving an optimal detection of 70.88% of the trees, plus oversegmented detections (i.e. more than one cluster center for a reference tree) at the rate of 14.91%.Thus, in total 85.79% of all reference trees are discovered, only 14.21% are missed.In this context, it is interesting to note that while the median tree radius of the reference trees is 3.70 m, the mode of the tree radii is 3.32 m.That means that only light prior knowledge about the expected tree radii of the scene of interest is sufficient to tune the clustering process, while keeping it otherwise fully unsupervised.
The result of the mean shift clustering of the point cloud with the thus-determined bandwidth parameter of 3.2 m is displayed in

Ellipsoid Modelling Results
The final result of the ellipsoid modeling process can be assessed in Fig. 3 (c), including tree crowns of different shape and hypothetical stem positions.A projection of the ellipsoids onto the 2D reference data is shown in Fig. 6.A summary of the tree parameter reconstruction errors is given in Tab. 1.In addition, the error distributions for tree heights and crown radii are shown in Fig. 7.

DISCUSSION
Although the results of this study are already very promising, they also show there is still room for further improvement: First of all, it is obvious that the clustering is dependent on choosing the optimal bandwidth parameter.Although this can be handled by using some light prior knowledge, an adaptive setting of the bandwidth parameter could possibly enhance the segmentation accuracy, in particular concerning the case of oversegmented trees.Secondly, the ellipsoid model of course is only a coarse approximation of real-life tree crowns, and only useful for deciduous trees at that.
Here, e.g. a generalized ellipsoid model also accounting for varying crown curvature could help to create a more universal approach and more detailed results.In addition, a more robust estimation of the tree heights and the crown radii is expected to reduce the over-estimation bias in these parameters significantly.
Concerning the number of missed trees, there is unfortunately always the sensor-inherent limitation: If a small tree is surrounded by large trees on all sides, not even multi-aspect SAR data will  help to avoid missing that tree due to the side-looking nature of the SAR imaging process.In such a case, only approaches based on volume tomography might provide a viable solution.
Last but not least, it has to be mentioned that the reference data also provides some potential for erroneous modeling, since no analysis of any kind of data can replace in-situ observations.For example, the smallest tree in the reference data is only 0.25 m high, i.e. in a real ground truth dataset, it would possibly not have been included at all.

CONCLUSION
In this article, an unsupervised approach for localization and reconstruction of individual trees from multi-aspect TomoSAR point clouds has been described.The point cloud is first segmented by unsupervised mean shift clustering.Then for every cluster a three-dimensional ellipsoid is modeled to the contained points.Since these ellipsoids are supposed to serve as satisfying approximations of deciduous tree crowns, three important tree parameters are extracted from each ellipsoid: tree location, tree height and tree crown diameter.Experiments based on a TomoSAR point cloud derived from an airborne millimeterwave dataset of two opposing aspects acquired over a cemetery in the city of Munich, Germany, have shown that about 86% of all trees can be localized and reconstructed by the presented technology.Although the side-looking SAR imaging geometry serves as a system-inherent limitation and leads to the fact that particularly small trees fully surrounded by large trees will always be missed, the results presented in this paper are expected to further stimulate the research interest in exploiting SAR imagery for forest remote sensing on the individual tree level.

Figure 1 .
Figure 1.Illustration of the ellipsoid modeling: (a) MVEE computed using 3D points projected onto the xy-plane; s1 and s2 are the computed semi-axes of the MVEE.(b) The MVEE of (a) is extruded in z-direction both upwards and downwards forming a 3D ellipsoid with a third semi-axis denoted as s3.x , y and z in (b) represent axes of the local coordinate system aligned to the ellipsoid semi-axes.The red points in both (a) and (b) represent the ellipsoid centers c k f .

Figure 3 .
Figure 3.The scene shown in different processing stages: (a) The 3D point cloud as derived by multi-aspect TomoSAR data fusion; (b) the clustered point cloud; (c) the reconstructed tree models.

Figure 4 .
Figure 4. Reference data of the test scene, created from a LiDAR point cloud and a co-registered orthophoto.Every circle indicates one manually extracted reference tree.

Figure 6 .
Figure 6.Ellipsoid models projected onto the 2D reference dataset for one-to-one comparison.

Figure 7 .
Figure 7. Distributions of (a) the tree height errors and (b) the crown radii errors.It can be seen that both tree heights and crown radii tend to be slightly overestimated.

Table 1 .
Mean absolute errors (MAEs) of reconstructed tree parameters.