Directionally Constrained Fully Convolutional Neural Network For Airborne Lidar Point Cloud Classification

TL;DR

This paper introduces a novel directionally constrained fully convolutional neural network (D-FCN) for semantic classification of airborne LiDAR point clouds, effectively handling unstructured data and improving accuracy without extra geometric features.

Contribution

The paper proposes a new D-Conv module that incorporates orientation information for local feature extraction and a multiscale FCN architecture for end-to-end semantic labeling of 3D point clouds.

Findings

Achieved state-of-the-art F1 score of 70.7% on ISPRS 3D labeling benchmark.

Significant performance improvements on categories with fewer points, like powerlines and facades.

Demonstrated effectiveness of orientation-aware convolution in LiDAR point cloud classification.

Abstract

Point cloud classification plays an important role in a wide range of airborne light detection and ranging (LiDAR) applications, such as topographic mapping, forest monitoring, power line detection, and road detection. However, due to the sensor noise, high redundancy, incompleteness, and complexity of airborne LiDAR systems, point cloud classification is challenging. In this paper, we proposed a directionally constrained fully convolutional neural network (D-FCN) that can take the original 3D coordinates and LiDAR intensity as input; thus, it can directly apply to unstructured 3D point clouds for semantic labeling. Specifically, we first introduce a novel directionally constrained point convolution (D-Conv) module to extract locally representative features of 3D point sets from the projected 2D receptive fields. To make full use of the orientation information of neighborhood points, the proposed D-Conv module performs convolution in an orientation-aware manner by using a directionally constrained nearest neighborhood search. Then, we designed a multiscale fully convolutional neural network with downsampling and upsampling blocks to enable multiscale point feature learning. The proposed D-FCN model can therefore process input point cloud with arbitrary sizes and directly predict the semantic labels for all the input points in an end-to-end manner. Without involving additional geometry features as input, the proposed method has demonstrated superior performance on the International Society for Photogrammetry and Remote Sensing (ISPRS) 3D labeling benchmark dataset. The results show that our model has achieved a new state-of-the-art level of performance with an average F1 score of 70.7%, and it has improved the performance by a large margin on categories with a small number of points (such as powerline, car, and facade).

Tables

Table 1. Table 1: Point category distribution of each scene.

Categories	Scene I	Scene II	Scene III	Scene II&III Total
Powerline	546	218	382	600
Low vegetation	180,850	38,104	60,586	98,690
Impervious surfaces	193,723	53,255	48,731	101,986
Car	4,614	2,184	1,524	3,708
Fence/Hedge	12,070	1,479	5,943	7,422
Roof	152,045	53,445	55,603	109,048
Facade	27,250	5,983	5,241	11,224
Shrub	47,605	1,414	23,404	24,818
Tree	135,173	23,915	30,311	54,226
Total	753,876	179,997	231,725	411,722

Table 2. Table 2: The overall accuracy (OA) and average F1 scores of the classification results of the models under different neighborhood partition strategies. From top to bottom are the model with no directionally constrained point convolution module and the models employing eight-directional search in 3D space, four-directional search in 2D space and eight-directional search in 2D space. The boldface entries indicate the model with the best performance.

Neighborhood partition strategy	OA	Average F1
no partition	0.817	0.663
3D space(8 directions)	0.818	0.668
2D space(4 directions)	0.820	0.669
2D space(8 directions)	0.833	0.699

Table 3. Table 3: The overall accuracy (OA) and average F1 score of model classification results under different numbers of sampling points ( N 𝑁 N ) and neighbor points ( K 𝐾 K ). The boldface text indicates the model with the best performance.

$N$	$K$	OA	Average F1
2,048	1	0.833	0.699
2,048	2	0.819	0.687
2,048	4	0.819	0.679
4,096	1	0.827	0.698
4,096	2	0.828	0.698
4,096	4	0.831	0.700
8,192	1	0.828	0.704
8,192	2	0.822	0.707
8,192	4	0.827	0.701

Table 4. Table 4: Classification results with different values of the class balance coefficient α 𝛼 \alpha . The overall accuracy (OA) and average F1 score are reported. The boldface text indicates the model with the best performance.

$\alpha$	OA	Average F1
NA	0.846	0.694
1.1	0.820	0.703
1.2	0.822	0.707
1.3	0.829	0.697
1.4	0.828	0.701
1.5	0.834	0.695

Table 5. Table 5: The classification confusion matrix of our directionally constrained fully convolutional neural network (D-FCN) model. Precision/correctness, recall/completeness, and F1 score are also reported. Our model obtained an overall accuracy of 82.2% and an average F1 score of 70.7%.

Categories	power	low_veg	imp_surf	car	fence_hedge	roof	facade	shrub	tree
power	0.690	0.003	0.000	0.000	0.000	0.143	0.053	0.003	0.107
low_veg	0.000	0.771	0.069	0.001	0.004	0.011	0.004	0.116	0.023
imp_surf	0.000	0.086	0.902	0.001	0.000	0.005	0.001	0.006	0.001
car	0.000	0.053	0.024	0.714	0.038	0.025	0.005	0.125	0.016
fence_hedge	0.000	0.094	0.015	0.008	0.264	0.013	0.005	0.483	0.119
roof	0.001	0.014	0.000	0.000	0.000	0.907	0.020	0.018	0.040
facade	0.002	0.061	0.004	0.001	0.004	0.116	0.575	0.138	0.099
shrub	0.000	0.104	0.006	0.005	0.016	0.037	0.017	0.626	0.189
tree	0.000	0.010	0.000	0.001	0.004	0.013	0.011	0.139	0.822
Precision/Correctness	0.718	0.836	0.927	0.862	0.621	0.954	0.637	0.364	0.767
Recall/Completeness	0.690	0.771	0.902	0.714	0.264	0.907	0.575	0.626	0.822
F1 score	0.704	0.802	0.914	0.781	0.370	0.930	0.605	0.460	0.794

Table 6. Table 6: Quantitative comparisons between our method and other state-of-art models on the ISPRS benchmark dataset. The numbers in the first nine columns of the table show the F1 scores for each category, and the last two columns show the overall accuracy (OA) and average F1 score (Average F1). The boldface text indicates the model with the best performance.

Categories	power	low_veg	imp_surf	car	fence_hedge	roof	facade	shrub	tree	OA	Average F1
UM	0.461	0.790	0.891	0.477	0.052	0.920	0.527	0.409	0.779	0.808	0.590
WhuY2	0.319	0.800	0.889	0.408	0.245	0.931	0.494	0.411	0.773	0.810	0.586
WhuY3	0.371	0.814	0.901	0.634	0.239	0.934	0.475	0.399	0.780	0.823	0.616
LUH	0.596	0.775	0.911	0.731	0.340	0.942	0.563	0.466	0.831	0.816	0.684
BIJ_W	0.138	0.785	0.905	0.564	0.363	0.922	0.532	0.433	0.784	0.815	0.603
RIT_1	0.375	0.779	0.915	0.734	0.180	0.940	0.493	0.459	0.825	0.816	0.633
NANJ2	0.620	0.888	0.912	0.667	0.407	0.936	0.426	0.559	0.826	0.852	0.693
WhuY4	0.425	0.827	0.914	0.747	0.537	0.943	0.531	0.479	0.828	0.849	0.692
Ours	0.704	0.802	0.914	0.781	0.370	0.930	0.605	0.460	0.794	0.822	0.707

Equations

$${\tilde{f}}_{ijc}=\sum _{k=1}^K{w}_{ck}\ast f_{ijk},c\in \{1,2,...,d_{out}\},j\in \{1,2,...,N_d\},$$

$$O_{ih}=\sum _{c=1}^{d_{out}}\sum _{j=1}^{N_d}{w}_{hj}\ast {\tilde{f}}_{ijc},h\in \{1,2,...,d_{out}\},$$

$$precision=\frac{TP}{TP+FP}$$

$$recall=\frac{TP}{TP+FN}$$

$$F1=2\ast \frac{precision\ast recall}{precision+recall},$$

$$W_c=\frac{1}{\ln (\alpha +\frac{N_c}{{\sum }_{c=1}^C{N}_c})},$$

$$\mathcal{L}=\sum _{i=1}^N{w}_i\ast \sum _{c=1}^C[y_{ic}log{p}_{ic}+(1-y_{ic})log(1-p_{ic})],$$

Full text

1 Introduction

The airborne light detection and ranging (LiDAR) system provides a new technical approach for acquiring 3D spatial data. Unlike commonly used satellite images, which mainly provide spectral information in a regular spatial grid, a LiDAR system captures data in 3D point cloud format, which provides reliable depth information and can be used to accurately localize objects and characterize their structure. Among various point cloud recognition tasks, the point cloud classification problem has been an active research area in remote sensing communities for decades (Rabbani et al.,, 2006; Vosselman et al.,, 2017) and serves as an essential building block for many applications, including topographic mapping, forest monitoring (Axelsson,, 2000; Mongus and Žalik,, 2013), power line detection (Andersen et al.,, 2005; Solberg et al.,, 2009; Zhao and Popescu,, 2009; Ene et al.,, 2017), road detection and planning, 3D building reconstruction (Kada and McKinley,, 2009; Yang et al., 2017a, ), and others. However, this process is still challenging due to the unstructured nature of 3D point clouds, their limited data resolution and high sensor noise, uneven density distribution, incompleteness and scene complexity (Zhu et al.,, 2017).

In general, the point cloud classification task can be divided into two basic steps: extracting representative point features both locally and globally and classifying each point into pre-defined semantic categories using the learned feature representations. Early works developed various hand-crafted point descriptors based on the geometric structure of local neighborhood points, e.g., density, curvature, roughness. In these methods, machine learning algorithms such as the Gaussian Mixture Model (Lalonde et al.,, 2005, 2006), Support Vector Machine (Zhang et al.,, 2013), AdaBoost (Lodha et al.,, 2007), and Random Forest (Babahajiani et al.,, 2017; Chehata et al.,, 2009) are commonly used to achieve point cloud classification based on hand-crafted features. However, these methods estimate the local features of each point independently and generate the label predictions without considering label consistency among neighborhood points. Therefore, the classification results tend to suffer from noise and label inconsistency.

Some studies have attempted to address this problem by incorporating contextual information into the classification models; the representative works include associative Markov networks (Munoz et al.,, 2009), non-associative Markov networks (Shapovalov et al.,, 2010) and conditional random fields (Niemeyer et al.,, 2014; Weinmann et al., 2015b, ; Niemeyer et al.,, 2012). Nevertheless, these machine learning-based point cloud classification methods comprehensively employ hand-crafted features to characterize each point in the input point cloud; thus, they have limited generalizability when applied to large-scale wild scenes.

In recent years, one of the most prevalent types of artificial intelligence algorithms (i.e., deep learning methods) has achieved great success in various real-world applications, including image classification, speech recognition, time series prediction, natural language processing and so on (LeCun et al.,, 2015; Chan et al.,, 2015; Hinton et al.,, 2012; Wen et al.,, 2019; Collobert and Weston,, 2008; Li et al.,, 2016, 2018). Following this trend, researchers have shifted their focus towards learning-based approaches and, more specifically, to deep neural networks (Qi et al.,, 2016; Su et al.,, 2015; Qi et al., 2017a, ; Qi et al., 2017b, ) for learning 3D local descriptors. For example, the authors of (Yang et al., 2017b, ; Zhao et al.,, 2018) proposed first transforming 3D point sets into regular 2D image representations and then employing a conventional 2D CNN for point classification. Such methods usually need to calculate additional local geometry features, such as planarity, sphericity and roughness, to enrich the 2D representations and further enhance the classification performance. However, the 3D to 2D conversion process results in information loss. More recent works have tried to apply convolution directly to unstructured point sets for classification (Yousefhussien et al.,, 2018; Wang et al., 2018b, ); however, these methods also require additional hand-crafted features (i.e., optical image features, local geometric features) as input to ensure the classification accuracy.

In this paper, we propose a novel directionally constrained fully convolutional network (D-FCN) that can be directly applied to unstructured 3D point clouds for airborne LiDAR point cloud classification. The raw 3D coordinates (i.e. X, Y, and Z) and intensity values obtained from the airborne LiDAR data are used as the input of our proposed model. Then, we introduce a novel directionally constrained point convolution (D-Conv) module to extract the locally representative features of 3D point sets. Based on this convolution module, we further design a multiscale fully convolutional neural network with downsampling and upsampling blocks to enable multiscale point feature learning. The resulting network can be trained in an end-to-end manner and can directly predict the semantic labels for all the input points. The main contributions of the proposed method are as follows:

This paper proposes a directionally constrained fully convolution neural network (D-FCN) that can be directly applied to raw point clouds for semantic labeling. The proposed D-FCN model can process input point cloud with arbitrary sizes and directly predict the semantic labels for all the input points in an end-to-end manner. 2. 2.

This paper introduces a novel directionally constrained point convolution (D-Conv) module to extract locally representative features of 3D point sets from the projected 2D receptive fields. Specifically, the proposed D-Conv module performs convolution in an orientation-aware manner by using a directionally constrained nearest neighborhood search. 3. 3.

This paper proposes a multiscale fully convolutional neural network with downsampling and upsampling blocks to enable multiscale point feature learning. 4. 4.

Without involving additional geometry features as input, the proposed method demonstrates superior performance on the ISPRS 3D labeling benchmark dataset.

The remainder of this paper is organized as follows. In Section 2, we provide a brief review of the airborne LiDAR point cloud classification methods. The proposed directional fully convolutional network (D-FCN) is described in detail in Section 3. In Section 4, we conduct experiments to verify the classification performance of the proposed method. Finally, the paper is concluded in Section 5.

2 Related Work

2.1 Point-based feature classification method

Traditional point cloud classification methods commonly regard each point in a point set as an independent entity, and classifying these points only depends on the local geometry features. Zhang et al., first calculated thirteen features of geometry, radiometry, topology and echo characteristics and then employed the Support Vector Machine (SVM) algorithm to classify airborne LiDAR point clouds of urban scene (Zhang et al.,, 2013). Lodha et al., adopted the AdaBoost algorithm to classify airborne LiDAR point clouds into four categories (i.e. road, grass, buildings, and trees) based on five features including height, height variation, normal variation, intensity, and image intensity (Lodha et al.,, 2007). Chehata et al., used the random forest algorithm to achieve the point cloud classification by selecting the most relevant features among the 21 features, which is composed of seventeen multi-echo and four full-waveform LiDAR features (Chehata et al.,, 2009).

However, the estimation of single-point local features is unstable due to the uneven distribution of point clouds, which can easily lead to classification noises and label inconsistency, especially for complex scenes (Weinmann et al., 2015a, ).

2.2 Context-based features classification method

To obtain smoother and more accurate classification results, many studies have introduced contextual information for point cloud classification based on single point local features. The associative Markov network (AMN) was presented by Munoz et al., to extract various features comprising the degree of scatter, linearity, planarity, etc.. These features contain multi-level contextual information. Then, contextual interactions were discriminatively considered to conduct 3D LiDAR point cloud classification (Munoz et al.,, 2009). Although the AMN model involving high-order interactions achieved more robust classification results than previous methods, it failed to detect both large and small objects due to oversmoothing, and it may lead to incorrect classifications of large region due to error propagation (Shapovalov et al.,, 2010). To address this issue, Shapovalov et al., further proposed the non-associative Markov network to overcome the weakness of the AMN model by using dynamic instead of constant pairwise potentials for pairs of different class labels based on 68-dimensional input features, including spectral and directional features, spin images, angular spin images and height distributions. In addition, Niemeyer et al., designed ten features, including amplitude, echo width, normalized echo number, distance to ground and so on, as the model input, and adopted the conditional random field (CRF) approach to classify airborne LiDAR point clouds. This proposed method enables contextual information to be incorporated (Niemeyer et al.,, 2012).

However, these classification methods require manually extracting context features in advance, which complicates the classification task. Moreover, their classification performances degrade when dealing with point clouds scanned from complex scenes (Zhao et al.,, 2018).

2.3 Deep learning-based classification method

As deep learning became increasingly prevalent, research on 3D point signature learning has followed the general trend in the machine vision community and shifted towards learning-based approaches for 3D point cloud classification.

In general, deep learning-based point cloud classification methods can be divided into two categories: feature image-based methods and point cloud-based methods.

2.3.1 Feature image based classification method

Convolutional neural networks (CNNs) have achieved great success on various 2D image recognition tasks, including scene classification, object detection, semantic segmentation, and many others. Following the success of 2D CNNs, many studies have tried to generate feature images from point clouds and then employ convolutional neural networks to achieve airborne LiDAR point cloud classification (Yang et al., 2017b, ; Yang et al., 2017b, ; Yang et al.,, 2018; Zhao et al.,, 2018). For example, Yang et al., 2017b proposed an image-generation method that transformed the local geometric features, global geometric features and full-waveform features of each point in the airborne LiDAR point clouds into 2D feature images and then applied conventional 2D convolutional neural networks to classify the feature images (Yang et al., 2017b, ). Yang et al., generated feature images of height, intensity, planarity, sphericity and variance of deviation angles across multiple scales and developed a multiscale convolutional neural network that learns the deep features of each point from the above feature images to enable airborne LiDAR point cloud classification. (Yang et al.,, 2018). Zhao et al., obtained a set of contextual images by calculating LiDAR point set features, including height, intensity and roughness, and then proposed a multiscale convolutional neural network to classify the airborne LiDAR point clouds into several object categories (Zhao et al.,, 2018).

However, these studies all involve generating 2D feature images from 3D point clouds, which can cause spatial information loss, and they still require handcrafted feature engineering before network training. In the feature engineering process, the selected features do not necessarily represent the most critical information of the original point cloud, and the manually defined combinations of different features can have a substantial influence on the final classification results.

2.3.2 Point cloud-based classification method

Due to the unordered and unstructured nature of point clouds, early efforts in the 3D point cloud recognition domain mainly attempted to convert 3D point clouds into more familiar representations and then apply conventional 2D/3D CNN to learn point features. Recently, Qi et al., 2017a proposed a deep learning framework called PointNet that can be applied directly to unstructured point clouds and learn representative point signatures from massive input point sets in a data-driven fashion. Experiments with various point cloud recognition tasks, including 3D shape classification, 3D ShapeNet part segmentation, 3D scene semantic segmentation and 3D object detection tasks, demonstrated the power of the PointNet model for point feature learning. To enable multiscale point feature learning, Qi et al., 2017b further designed a set abstract module and a feature propagation module for point feature downsampling and upsampling, respectively, and developed a hierarchical neural network that applies unit PointNet recursively on a nested partitioning of the input point set.

Following the great success of PointNet (Qi et al., 2017a, ) and PointNet++ (Qi et al., 2017b, ), studies have shifted their focus towards PointNet-like architectures for airborne LiDAR point cloud classification. For example, Yousefhussien et al., designed a 1D-fully convolutional network to achieve point cloud classification by inputting not only the raw coordinates but also three additional corresponding spectral features extracted from 2D georeferenced images for each point. Similarly, Wang et al., 2018b first extracted a 54-dimensional feature vector for each point in the input point set. Those features comprise a 6-dimensional eigenvalue feature vector and a 12-dimensional spin image feature vector using 3 different settings for the parameter K of a nearest-neighbor search. Then, they developed a deep neural network with spatial pooling to aggregate the above point-based features into cluster-based features, which were further fed into a multi-layer perceptron network for point cloud classification (Wang et al., 2018b, ). Although these works embody methods that learn 3D point clouds directly, they still require integrated handcrafted features to improve the classification accuracy on airborne LiDAR data, which makes the task tedious and complicated.

In this paper, we proposed a directionally constrained fully convolutional neural network (D-FCN) that can be directly applied to unstructured 3D point clouds and predict the end-to-end semantic labels for arbitrary input point sizes. Furthermore, compared to other existing point cloud-based deep learning methods, the proposed model can extract representative orientation-aware features. Further, it achieves a new state-of-the-art classification performance level on airborne LiDAR point clouds while taking as input only the raw 3D coordinates and intensity; it does not require any additional handcrafted features.

3 Methods

Our D-FCN model follows an encode-decode framework similar to the general semantic segmentation network (Badrinarayanan et al.,, 2017) for 2D image segmentation (illustrated in Figure 4).

In the downsampling stage, we recursively apply our proposed directionally constrained point convolution module (D-Conv) prior to each downsampling operation for hierarchical feature embedding. In the upsampling stage, our model achieves dense feature prediction by interleaving feature propagation with our proposed D-Conv module. The core component of our segmentation network is the D-Conv module, which enables orientation aware point feature extraction. We discuss the proposed directionally constrained convolution module in Section 3.1 and then introduce the overall architecture of our model in Section 3.2.

3.1 Directionally constrained point convolution

Given an input point cloud $P\in\mathbf{R}^{N\times 3}$ and a feature matrix $F\in\mathcal{R}^{N\times d_{in}}$ of one hidden layer, where $N$ represents the number of points in the input point set and $d_{in}$ denotes the dimension of the input feature map, our D-Conv module generates the output $O\in\mathcal{R}^{N\times d_{out}}$, and assigns a new $d_{out}$-dimension feature vector to each point. Here, similar to the conventional convolution operator, the proposed D-Conv module covers a small receptive field and is responsible for local feature extraction. We design our directionally constrained point convolution module following the construction patterns of the conventional convolution operator. Two steps are involved in the convolution module: we construct a local receptive field foreach point and then apply the convolution operation to the feature values within the receptive field using the kernel weights.

3.1.1 Directionally Constrained Nearest Neighbor

The conventional convolution operator starts by constructing a local receptive field to enable local feature extraction. Figure 1 (left) gives an example of conventional 3x3 convolution; given a central point (marked in red), the receptive field covers an area of 3x3 gridded pixels, and each pixel (marked in blue) can be regarded as the nearest neighbor in each of the 8 evenly divided directions, except for the central point itself.

Similarly, the proposed D-Conv module formulates its receptive field by finding the nearest neighbor points in $N_{d}$ (by default, we set $N_{d}$ to 8) evenly divided directions. More importantly, instead of dividing the original 3D point space into cone subspaces, we first project all the input 3D points onto an xy plane and then formulate the receptive field in 2D point space.

The reason why we do not divide the 3D space based on the original xyz coordinates is that the 3D point cloud captured by airborne LiDAR has a larger variance in the horizontal direction (xy coordinates) than in the vertical direction (z coordinates). Therefore, projecting points to the xy plane during the neighborhood search enables more representative neighbor points in the xy direction while paying less attention to points projected to nearby locations on the ground plane that share the same semantic labels. This is the reason why previous works such as (Zhao et al.,, 2018) and (Yang et al.,, 2018) directly project the original 3D point cloud into 2D images and use conventional 2D CNN to predict the semantic label for each point. Figure 2 shows an illustration of this statement. Building rooftop points such as $A,B$ have no upper neighbor points, while the ground points such as $C,D$ have no lower neighbor points. Therefore, searching directionally constrained nearest neighbors in the original 3D space causes many vacant subspaces and introduces fewer informative neighborhood points. In contrast, the proposed 2D search strategy guarantees fewer vacant subspaces and enables richer neighborhood information. Another advantageous property of this 2D search strategy is that dividing a 2D space into triangular areas is much easier and more computationally efficient than dividing a 3D space into cones. We verify the advantages of our 2D search strategy over the 3D search strategy in Section 4.4.

After projecting the neighborhood space onto the xy plane, we divide the 2D plane originating at each central point $p_{i}$ into $N_{d}$ triangle areas. Without special declarations, we set $N_{d}$ to 8 in this paper. Then, we search $K$-nearest neighbor points in each area within a radius of $R$ to formulate the receptive field of the central point (see Figure 1 (right) for illustration). In general, a larger $K$ enables a more reliable and stable representation of the captured features; however, a larger $K$ can also lead to overfitting. The effects of different values of $N_{d}$ and $K$ are discussed in Section 4.5. Moreover, the distant points tend to provide little information for describing the local patterns at the central point $p_{i}$. When no point exists within the search radius $R$ in a triangle area (e.g., area ‘1’ in the right-hand figure), we duplicate the central point $p_{i}$ as the nearest neighbor of itself. After searching for each central point, the proposed directionally constrained nearest neighbor generates $N_{d}$ groups of $K$ points in a fixed order.

Although we project all the 3D points onto the xy plane, the z coordinate feature is maintained for point feature encoding, as discussed in the next section.

3.1.2 Orientation-Aware Point Feature Encoding

In a conventional convolutional layer, the convolution operation for each kernel preserves the same order when sliding through the whole input feature maps, thus ensuring the translation invariance property. An important prerequisite for this operation is that in regular-gridded feature inputs, a well-defined order exists for organizing the neighbor points within the receptive field. Figure 1 (left) gives an example of the ordering of neighbor points starting from 1 to $N_{d}$.

However, 3D point clouds are presented in an unordered and unstructured format. No regular order exists for organizing the neighbor points. Early works partitioned the input 3D point space into gridded voxels to enable an order formulation. Recent PointNet-like works commonly apply unordered operations (e.g., max pooling (Qi et al., 2017b, ; Wang et al., 2018a, ; Shen et al.,, 2018; Jiang et al.,, 2018)) to aggregate the features of all the neighbor points when formulating the local feature for each central point. However, these unordered operations ignore order information that can capture the spatial structure of the input points; consequently, they also cause information loss (e.g., the max pooling operation discards all inputs except for the maximum).

To enable orientation-aware point feature learning, our D-Conv module is therefore expected to maintain a well-defined order when applied to each point, which can be easily achieved by adopting the order of each triangle area when formulating the receptive field (see Figure 1 (right) for illustration). In each of the $N_{d}$ triangle areas, the features of the $K$ neighbor points are first aggregated to generate the directionally aware features, resulting in a $1\times N_{d}$ feature vector for each central point. Given central point $p_{i}$ and $N_{d}$ groups of $K$ neighbor points, the directionally aware feature vector $\tilde{f}\in\mathcal{R}^{d_{out}}$ of the $j$th triangle area is calculated as follows:

[TABLE]

where $k\in\{1,2,...,K\}$ denotes the neighbor index in each triangle area, and $w_{ck}$ denotes the corresponding weight parameter of the $c$th output channel. This aggregation function can easily be implemented using a $1\times K$ convolution with a stride of $K$ and $d_{out}$ output channels on the grouped input features of dimension $B\times N\times(N_{d}\times K)\times d_{in}$, where $B$ denotes the batch size. The output tensor dimension is $B\times N\times N_{d}\times d_{out}$.

Then, our orientation-aware point convolution is formulated as a $1\times N_{d}$ convolution with a stride of $N_{d}$ and $d_{out}$ output channels. Given central point $p_{i}$ and the directionally aware feature vectors $\tilde{f}$, the output feature $O_{i}\in\mathbf{R}^{d_{out}}$ at central point $p_{i}$ is calculated as

[TABLE]

where $w_{hj}$ indicates the weight parameter for the $j$th triangle area of the $h$th output channel. The output tensor dimension is $B\times N\times 1\times d_{out}$.

In sum, our D-Conv module starts by finding the $K$-nearest neighbors in each of the $N_{d}$ triangle areas partitioned on the xy plane. Then, we aggregate the features in each triangle area by using a $1\times K$ convolution with a stride of $K$. Next, orientation-aware point convolution is applied to extract the point descriptors for each central point using a $1\times N_{d}$ convolution with a stride of $N_{d}$. Finally, we add the output features to the input using an element-wise summation operation. The overall architecture of our D-Conv module is illustrated in Figure 3. In Figure 3, we repeat the convolution block twice to enable more representative feature learning.

Note that our D-FCN module does not change the size of the input point set. In the next section, we introduce downsampling and upsampling blocks to enable multiscale point feature aggregation and propagation.

3.2 Network architecture

3.2.1 Downsampling and Upsampling blocks

To enable multiscale high-level feature learning, in this study, we adopt downsampling modules to extract multiscale features and upsampling modules to obtain point features for all the original points. These modules were introduced in PointNet++ (Qi et al., 2017b, ) to enable hierarchical point feature learning. Here, we provide a brief introduction to these two modules.

A downsampling module takes a $N_{l}\times d_{l}$ feature matrix as input, where $N_{l}$ represents the number of input points and $d_{l}$ indicates the feature dimension, and outputs a $N_{l+1}\times d_{l+1}$ feature matrix with a reduced number of points $N_{l+1}$ ($N_{l+1}\leq N_{l}$) and a feature dimension of $d_{l+1}$. The downsampling process is implemented by farthest point sampling to select a subset of points such that each selected point is the most distant point from the set of all the other points. After point downsampling, the embeddings for each selected point are calculated by passing the features of its neighbor points through a convolution module. A ball query or K-nearest neighbor search strategy can be used to define these neighbors, and the max-pooling function is used to aggregate the features in the convolution module. By default, we set the neighborhood size to 32 in this study.

An upsampling module propagates the point features from $N_{l}$ points to $N_{l-1}$ points (with $N_{l-1}\geq N_{l}$) using distance-based interpolation. Specifically, the features of the points dropped during the downsampling process are calculated using the inverse distance weighted average based on the k-nearest neighbors. As the default setting, for each point to be interpolated, we search for 32 neighbor points in the input points to produce the point feature for these dropped points. The interpolated features from the $N_{l-1}$ points are then concatenated with the features obtained during the downsampling process with the same number of points.

3.2.2 Overall architecture

Our point segmentation network is organized in a U-Net-like structure, as shown in Figure 4. It takes the 3D coordinates and the intensity feature as input and outputs a semantic label for each point. During downsampling, three successive downsampling operations are applied to shrink the size of the point set to 1,024, 256, 64. Then, during the upsampling stage, three successive upsampling operations are utilized to generate dense feature prediction. The point set size is increased to 256, 1,024, 8,192 points in each upsampling stage. We insert our D-Conv module is inserted before each downsampling operation and each upsampling operation. As a default setting, the search radius $R$ for each sampling level was set to 2, 5, 10 from bottom to top. Finally, the point features of the last upsampling layer are input into a fully connected layer to produce a semantic label for all the input points. Moreover, to incorporate the low-level information in the downsampling stage, we add skip connections between the downsampling and upsampling stages. The low level-level features from the downsampling stage are concatenated with the feature matrix of the same point set size in the upsampling stage. Previous works in 2D segmentation tasks have frequently reported that employing skip connections boost the convergence speed and improve performance (Ronneberger et al.,, 2015; Badrinarayanan et al.,, 2017).

3.3 Evaluation metrics

Following the standard convention of the ISPRS contest, precision, recall, F1 score, and overall accuracy (OA) are commonly used to evaluate the performance of 3D point cloud labeling. In general, overall accuracy is used to assess the overall classification accuracy for all categories, which is defined as the percentage of correctly classified points in the total test points. Moreover, the F1 score takes the precision and recall of the classification model into account and is generally more suitable for cases where the categories are unevenly distributed. The precision, recall and F1 score for each category are defined as follows:

[TABLE]

[TABLE]

[TABLE]

where TP (true positive), FP (false positive), and FN (false negative) respectively indicate the positive tuples that were correctly labeled by the classifier, the negative tuples that were incorrectly labeled as positive, and the positive tuples that were mislabeled as negative.

4 Experiments and Discussion

In this section, we present the extensive experiments conducted to evaluate the performance of the proposed model for Airborne LiDAR point cloud classification. In Section 4.1, we introduce the experimental dataset. In Section 4.2, we present two processing strategies, patch sampling and class balance. The experimental setup is given in 4.3. In Section 4.4, we evaluate the effect of the proposed directionally constrained point convolution module. In Section 4.5, we discuss the effect of different hyper-parameter configurations. The classification results of our model with optimal hyperparameters are given in Section 4.6. Finally, in Section 4.7, we compare our model with other state-of-the-art models for airborne LiDAR point cloud classification.

4.1 Experimental dataset

To evaluate the performance of the proposed method, we conducted experiments on the International Society for Photogrammetry and Remote Sensing (ISPRS) 3D labeling dataset (Niemeyer et al.,, 2014). This dataset includes three test areas that include reference data for various object classes. The airborne LiDAR point clouds were captured by a Leica ALS50 system, at an average flying height of 500 meters above the ground and a field of view of 45 degrees (Cramer,, 2010). The LiDAR points were labeled with 9 semantic categories, including powerline, low vegetation (low_veg), impervious surface (imp_surf), car, fence/hedge (fen/hed), roof, facade, shrub, and tree. The spatial 3D coordinates (XYZ), intensity values, number of returns and point labels were stored in LiDAR point files.

Following the standard setting of the ISPRS 3D labeling contest, the dataset was divided into two parts according to the scene areas. The first scene (see Figure 5 left) with 753,876 points was used as the training dataset, and the other two scenes (see Figure 5 right) with 411,722 points as the test dataset. The detailed category distribution of each scene is shown in Table 1. Note that the test dataset labels were used only for model evaluation, which is unknown during the model training stage.

4.2 Class balance and data augmentation strategies

As shown in Table 1, the category distribution of each scene was extremely uneven, as indicated by the number of points in each category. For example, in scene I, merely 546 points are labeled as powerline, while 193,723 points are labeled as impervious surfaces. Conducting training directly on this unbalanced dataset can potentially have a severe negative impact on the overall performance of deep networks, as discussed in a recent technique report (Masko and Hensman,, 2015). To address this issue, we added a category-specific weight coefficient for each category to the loss function of our D-FCN model. This strategy was designed to encourage our model to pay more attention to those points in the minority categories. The balance weight for each category was calculated by the logarithm function of the percentage of each category and then performing the reciprocal operation, formulated as Equation 6.

[TABLE]

where $W_{c}$ refers to the weight of the $c$th category, $N_{c}$ represents the number of points of the $c$th category, $C$ denotes the total number of categories, and $\alpha$ denotes the coefficient for class balance. The value of $\alpha$ is discussed in Section 4.5. After integrating the class balance weights, our final loss function is defined as follows:

[TABLE]

where $N$ denotes the number of sampling points in each training block, $y_{ic}$ and $p_{ic}$ denote the ground truth label and predicted probability of the $i$th point for the $c$th category, and $w_{i}$ denotes the balance weight for the $i$th sampling point, calculated as $w_{i}:=W_{c=y_{i}}$.

Moreover, to improve the generalizability, we introduce a data augmentation technique to further enhance the robustness of our proposed model. Considering that our training scene has an irregular boundary, we referenced the data augmentation algorithm introduced in PointNet++ (Qi et al., 2017b, ). During training, we randomly select a 30 m30 m40 m cuboid region from the whole scene and then randomly choose 8,192 points from the cuboid as the model input. To further reduce the risk of overfitting and make the model more robust, the selected 8,192 points are randomly dropped during the training stage. By default, the dropout ratio is set to 12.5% in the following sections. With respect to the testing dataset, the scenes were segmented into blocks of 30 m*30 m grids in the horizontal direction (see Figure 8). Note that small blocks generated at the edge of the scene were merged into the surrounding large blocks to ensure the integrity of each block. Also note that the test blocks can have substantially different numbers of points. Fortunately, due to the fully convolution nature of the proposed D-FCN model, all the points in each test block can be input directly into the model for point classification.

4.3 Experimental setup

We implemented our D-FCN method using the TensorFlow framework. During network training, we set the batch size to 6 and used the Adam optimizer with an initial learning rate of 0.01. We divided the learning rate by 2 every 3,000 steps. Training our model for 1,000 epochs until convergence required approximately 10 hours running on a TESLA K80 GPU.

During the test stage, we fed each test block directly into our model and generated the classification results in one network forward pass. We merged the prediction results of each block to obtain the final point cloud segmentation results, see Figure 9. Our code will be released at https://github.com/lixiang-ucas/D-FCN.

4.4 Effect of directionally constrained point convolutions

To explore the effect of the directionally constrained point convolutions module and further determine the optimal directional searching strategy of the neighbor points, we designed four models ( i.e., a model with no directionally constrained point convolution module, a model with an eight-directional search in 3D space, a model with a four-directional search in 2D space and a model with an eight-directional search in 2D space) and conduct point cloud classification using identical network configurations.

Here, as a default configuration, the number of sampling points $N$ was set to 2,048, and we sampled 1 point in each directionally constrained area (i.e., $K=1$). The $\alpha$ class balance coefficient was set to 1.2. The classification results of all four compared models are listed in Table 2. Here, we report both the overall accuracy and the average F1 score. The corresponding qualitative results are shown in Figure 7.

In Figure 7, the model without a directionally constrained point convolution module misclassified some roof points into the tree category. More importantly, the proposed models with directionally constrained searching in 2D space outperformed the model with directionally constrained search in 3D space regardless of the number of directions used, as indicated by the overall accuracy and average F1 score. This is probably because the projected 2D searching strategy enables richer informative neighbor points than does the 3D searching strategy, as discussed in Section 3.1. Moreover, comparing the results of the last two rows in Table 2, it is clear that instead of dividing the 2D space into 4 directionally constrained areas, using 8 directionally constrained areas leads to more representative neighbor points and achieves better classification results. Hereafter, our directionally constrained point convolution module adopts an 8-directional search in 2D space.

4.5 Model hyperparameters

In this section, we conduct extensive experiments to determine the optimal model hyperparameters, including the number of sampling points $N$ in each training block, the number of neighbor points $K$ in each triangle area and the class balance coefficient $\alpha$. Average F1 score was adopted as the evaluation metric to select the best hyperparameters.

4.5.1 Number of sampling points and neighbor points

As mentioned in Section 3.1, the number of neighbor points $K$ plays an important role in model performance. Selecting an appropriate number of these points enables a more reliable and stable representation of the captured features. Meanwhile, the number of sampling points $N$ in each training block controls the density of the input point set. In a sparse point set (i.e., where $N$ is small), each central point tends to have fewer candidate points that can be selected within the searching radius $R$. Therefore, we discuss the number of sampling points $N$ and neighbor points $K$ simultaneously to determine their optimal configurations.

In our experiments, we randomly sampled 2,048, 4,096, and 8,192 points from each training block as the input to our D-FCN model and searched 1, 2, and 4 neighbor points from each of the 8 directionally constrained triangle areas in our D-Conv module. In these experiments, the class balance coefficient $\alpha$ was 1.2. The classification performances were evaluated regarding both overall accuracy (OA) and average F1 score, as listed in Table 3.

As shown in Table 3, the model achieved a better classification performance by taking 1 nearest neighbor point in each direction (i.e., $K=1$) when fewer points were sampled in each training block (e.g., $N=2,048$). This result probably occurs because a larger neighborhood size involves more redundant information when the point set is sparse. As the number of sampling points increases (e.g., $N=8,192$), the increase in neighbor points ($K=2$) helps to obtain more stable features that improve the classification results. Note that using too many neighbor points ($K=4$) can have a negative effect on the classification performance when $N$ is sufficiently large. Our model achieved its best performance when $N$ and $K$ were set to 8,192 and 2, respectively, with an overall accuracy of 82.2% and an average F1 score of 70.7%. Hereafter, we use this configuration for $N$ and $K$ in our experiments.

4.5.2 Weight coefficient of data balance

In this section, we conducted comparative experiments with different class balance coefficient ($\alpha$) configurations; namely, a model with no data balancing, and models with $\alpha$ set to 1.1, 1.2, 1.3, 1.4 and 1.5. Table 4 shows the classification results under different $\alpha$ configurations.

As shown in Table 4, our model achieves the best overall accuracy (84.6%) with no class balancing strategy. However, its average F1 score is lower than the models with class balance. This result is probably caused by the uneven distribution of the input class categories (see Table 1). When using no class balancing strategy, our model paid more attention to the categories with more points (e.g., low vegetation, impervious surfaces) and thereby achieved higher accuracy for these majority categories. However, the categories with fewer points (e.g., powerline, car) may be ignored during the model training, resulting in lower accuracy for these categories.

When the class balancing strategy was used, the models improved the average F1 score. The highest average F1 score of 70.7% was achieved by the model where the class balance coefficient $\alpha$ was set to 1.2.

4.6 Classification results

After determining the best hyperparameters for the model, we employed Scene I of the ISPRS 3D Semantic Labeling Dataset to train the proposed D-FCN model until convergence. Then, each test block from Scene II and Scene III described in Section 4.2 was fed directly into the model without sampling or dropout. The final point cloud classification results and error maps are shown in Figure 8 and Figure 9. As shown in Figure 8 and Figure 9, the proposed D-FCN model successfully generates the correct label predictions for most of the points in the test scenes.

To quantitatively evaluate the classification performance, we calculated the classification confusion matrix, precision, recall and F1 score of each category and list the results in Table 5, where we can see that our proposed model obtained F1 scores above 70% for six of the categories, including powerline, low vegetation, impervious surfaces, car, roof, and tree. In addition, our model achieved a reasonable classification performance on the facade category, but its classification performance on the fence/hedge and shrub categories was relatively lower. As Table 5 shows, most of the fence/hedge points were incorrectly classified as shrubs. A primary reason for this result is that the fence/hedge category contains fewer points and presents spatial distribution and topological characteristics similar to the shrub category. These factors cause the model to be insufficiently trained, which prevents the model from being able to differentiate these two categories. Although the powerline and car categories have fewer points in the dataset, they present completely different characteristics from other categories; thus, they result in a better classification performance.

4.7 Comparisons with other methods

To demonstrate the advantages of the proposed model, we compared it with other models that participated in the ISPRS contest. The top eight models (i.e., UM (Horvat et al.,, 2016), WhuY2, WhuY3 (Yang et al., 2017b, ), LUH (Niemeyer et al.,, 2016), BIJ_W (Wang et al., 2018b, ), RIT_1 (Yousefhussien et al.,, 2018), NANJ2 (Zhao et al.,, 2018) and WhuY4 (Yang et al.,, 2018)) that achieved the best performance on this benchmark were selected for performance comparisons. The overall accuracy and F1 scores of our model and of all the compared models are listed in Table 6. As Table 6 shows, the proposed model achieved superior classification performances than did all the compared methods with regard to average F1 score. In particular, our model achieved substantially higher performances on the powerline, car and facade categories, for which it outperformed the state-of-the-art models by 8.4%, 3.4%, 7.3%, respectively.

However, analyzing the overall classification accuracy, our model failed to achieve the best performance compared to the other state-of-art models. This result occurs mainly because we adopted the average F1 score as the main indicator when determining the initial model parameters rather than overall accuracy. As mentioned in Section 4.5, many hyper-parameter settings exist that could lead to an overall accuracy above 82.2% (obtained with the optimal parameters setting aimed at obtaining the best average F1 score); see Table 3 and Table4. However, because adopting overall accuracy as the main evaluation metric will cause minority categories to be ignored, we focused primarily on the F1 score in this study.

To further compare the classification performance of different models, we selected the top three models from the ISPRS contest (i.e., RIT_1 (Yousefhussien et al.,, 2018), NANJ2 (Zhao et al.,, 2018) and WhuY4 (Yang et al.,, 2018)) as comparison models and show their classification results on a selected complicated scene area in Figure 10. Taking the points of the powerline category (marked in light yellow) as an example, our model clearly demonstrates higher classification accuracy compared to the other three models, which misclassify most of the powerline points into other categories. The main reason for this result is that only 600 points exist in the powerline category, accounting for merely 0.1% of the total 411,722 points in the test set. Therefore, this category is largely ignored by the three compared models.

5 Conclusions

In this paper, we proposed a directionally constrained fully convolution neural network (D-FCN) model that can be directly applied to unstructured 3D point sets to classify airborne LiDAR point clouds. Specifically, we proposed a novel directionally constrained point convolution (D-Conv) module to extract the local representative features of 3D point sets from projected 2D receptive fields and a multiscale fully convolutional neural network to enable multiscale point feature learning. Our model can process input point clouds of arbitrary sizes and directly predict the semantic labels for all the input points in an end-to-end manner. To demonstrate the advantages of the proposed model, we conducted experiments on the ISPRS 3D Labeling Benchmark Dataset and compared our model’s performance with those of other state-of-art methods. The results show that our model achieved a new state-of-art classification performance for average F1 score despite using only the original 3D coordinates and intensity as inputs. Our D-FCN model outperformed the state-of-the-art methods by 8.4%, 3.4%, 7.3% on the powerline, car, and facade categories, respectively.

ACKNOWLEDGEMENTS (Optional)

We thank the NYU MMVC Lab for providing GPU-equipped servers for these experiments. Additionally, we would like to gratefully acknowledge the ISPRS for providing airborne LiDAR data.