A Novel Multi-LiDAR-Based Point Cloud Stitching Method Based on a Constrained Particle Filter

Abstract

In coal-fired power plants, coal piles serve as the fundamental management units. Acquiring point clouds of coal piles facilitates the convenient measurement of daily coal consumption and combustion efficiency. When using servo motors to drive Light Detection and Ranging (LiDAR) scanning of large-scale coal piles, the motors are subject to rotational errors due to gravitational effects. As a result, the acquired point clouds often contain significant noise. To address this issue, we proposes a Rapid Point Cloud Stitching–Constrained Particle Filter (RPCS-CPF) method. By introducing random noise to simulate servo motor rotational errors, both local and global point clouds are sequentially subjected to RPCS-CPF operations, resulting in smooth and continuous coal pile point clouds. Moreover, this paper presents a coal pile boundary detection method based on gradient region growing clustering. Experimental results demonstrate that our proposed RPCS-CPF method can generate smooth and continuous coal pile point clouds, even in the presence of servo motor rotational errors.

1. Introduction

With the rapid growth in demand for electricity supply, the need for coal is also steadily increasing. In the operation of coal-fired power plants, coal piles serve as the fundamental units of management. Obtaining point cloud data of coal piles allows for the easy measurement of various parameters such as volume [1,2,3], density, and boundaries, which are crucial for ensuring safe operations at coal yards and effectively managing coal combustion efficiency.

1.1. 3D Reconstruction of Coal Piles

Traditional methods for obtaining point clouds of coal piles involve the tedious process of manual scanning using handheld scanners [4]. This approach demands significant time and effort. Recently, researchers have explored more efficient 3D reconstruction techniques [5,6], such as employing GNSS-RTK technology, utilizing multiple fixed LiDARs for scanning coal heaps [7], and deploying drones. Each of these methods offers their own advantages and disadvantages depending on specific circumstances and applications.

Handheld devices, typically employing laser scanners or cameras, are utilized by operators to manually scan coal piles for data collection. However, the scanning efficiency of this method is relatively low, making it unsuitable for rapid scanning and continuous monitoring of large coal piles. Moreover, measuring larger coal piles presents significant challenges for workers.

Drone scanning involves the use of unmanned aerial vehicles (UAVs) equipped with cameras or laser scanners to perform aerial scans of coal piles [8,9,10,11]. This system can rapidly capture surface data of extensive coal piles within a short timeframe, thereby offering the advantages of improved resolution and comprehensiveness. However, drone scanning does have certain limitations. A study by Alsayed et al. [9] revealed that drones equipped with LiDAR sensors may have blind spots during data scanning, particularly in enclosed coal yard environments, which pose challenges for optimizing flight trajectories and may require further refinement. Additionally, research by Davis et al. [12] demonstrated that increasing the altitude of drone flights can result in greater errors and reduced reliability.

Mahlberg et al. [13] pioneered a portable LiDAR device mounted on a pole for scanning point clouds in expansive granaries. However, this method relies on complex point cloud registration techniques, leading to considerable evaluation complexity. Farhood et al. [14] proposed using smartphone cameras to extract material point clouds by capturing moving images. However, it is important to note that this approach is restricted to smaller materials and faces challenges in ensuring the accuracy of point cloud reconstruction.

While GNSS-RTK can provide precise location data, it often faces challenges like signal obstruction and multipath interference in indoor enclosed environments. These issues negatively affect its positioning accuracy and reliability, making it unable to provide continuous point cloud data [15]. Raevaa et al. [8] utilized GNSS-RTK for measurements in an open-pit quarry, but they found that the measurement speed was significantly slow, resulting in decreased work efficiency.

In contrast, using laser scanning [16,17,18] allows for the acquisition of complete point cloud data for coal piles, but, for large-scale coal piles, a single laser scanner operates at a slow speed and produces sparse point clouds [19,20]. Therefore, this study employs a system of multiple fixed LiDARs to construct a coal pile scanning system, which offers several advantages. Firstly, fixed LiDARs exhibit higher scanning precision and stability [21,22,23]. By adjusting the installation height and angle appropriately, they can adapt to coal piles of varying heights, thus better covering the entire surface of the coal pile and obtaining more comprehensive and accurate point cloud data [12,24]. Secondly, fixed LiDARs enable long-term continuous scanning, allowing for the continuous acquisition of point cloud data for coal piles. This makes them suitable for long-term monitoring of changes in coal pile morphology and real-time volume calculations. Therefore, this study adopts the approach of using multiple laser scanners to obtain coal pile point clouds.

However, employing multiple LiDARs for scanning coal piles does pose certain limitations. Applying pressure to actuators can induce instability, leading to deviations in the generated point clouds. These discrepancies may cause non-smooth surfaces when stitching point clouds using the existing coordinate system. To address this challenge, this study introduces the RPCS-CPF method as a solution for coal pile point cloud stitching.

1.2. Point Cloud Edge Detection

During coal pile operations, vehicles often need to climb to the top of coal piles for tasks such as loading, leveling, or measurement. However, the surface terrain of coal piles is complex and variable. When nearing the edges of the pile or encountering depressed areas, vehicles are susceptible to sliding hazards, resulting in potential casualties. Therefore, detecting boundaries of coal piles and providing early warnings are crucial for ensuring the safety of personnel.

By analyzing and meticulously processing dense point cloud data, it becomes possible to accurately detect the boundaries and depressions on the surface of the coal pile, thus revealing potential safety risks. Chen et al. [25] introduced a 3D boundary identification technique that utilizes DBSCAN clustering. This method demonstrated favorable outcomes when applied to point clouds exhibiting conventional local shapes. However, it was found to be unsuitable for point clouds representing complex coal piles. Furthermore, the controllability of the parameters in this clustering method is limited, posing a challenge in identifying the Pareto optimum. In their study, Mineo et al. [26] introduced a unique algorithm called BPD for boundary point identification, along with a spatial FFT-based filtering approach. This approach is effective in fitting surfaces with polynomials and is particularly suitable for smooth coal pile surfaces. However, its performance is suboptimal when dealing with complex coal heaps that have pits and tunnels. Additionally, Runge’s phenomenon arises when the order of the polynomial exceeds 10. In their study, Yang et al. [27] introduced an algorithm that utilizes multi-scale directional curvature to extract and quantify the borders of accumulations from 3D point cloud data. This algorithm is capable of detecting the boundaries of various wave peak materials, but it faces difficulties in detecting pits and trenches. Hu et al. [28] introduced a boundary identification technique that utilizes semantic segmentation. This method demonstrated favorable outcomes when used on the S3DIS and ScanNet datasets. However, it encounters challenges in accurately differentiating the intricate distribution of coal pile surface terrain.

Given these challenges, it is crucial to devise a robust and flexible strategy for surface identification and assessment of coal piles, ensuring the safety and efficiency of coal handling operations. The adoption of such technology holds promise for significantly reducing accident rates and optimizing the performance of coal handling facilities.

The innovations in this work are as follows:

(1)

A rapid point cloud stitching algorithm grounded in the Constrained Particle Filter (CPF) is presented, which addresses the stochastic rotational errors of servos through mathematical modeling and has undergone algorithmic validation on a large coal pile. Utilizing multiple LiDAR–servo units, we scanned the coal pile and initially processed the point cloud generated by a single LiDAR scan with the CPF. Following this, we applied the CPF to the point cloud resulting from the stitching of multiple LiDAR scans. Experimental results have confirmed that our stitching algorithm not only ensures a smooth transition at the junction points but also maintains the surface integrity of the coal pile’s point cloud.

(2)

We propose a complex coal pile surface edge detection algorithm based on gradient region growing clustering. Initially, we estimate the normal vectors and calculate the gradients of the stitched point cloud. Subsequently, clustering is performed using the slope and gradient magnitude of the coal pile. By setting specific slope and magnitude intervals, we extract the boundaries of the coal pile. Experimental results indicate that our method is capable of detecting the boundaries of hazardous terrains such as pits, aisles, and ridges within the coal pile, thereby enhancing the safety of coal pile operations. This approach holds broad application value.

2. Method

When acquiring point cloud data for a large coal pile in an enclosed coal yard, using a handheld scanner is impractical for scanning the entire coal pile, and employing drones requires planning complex scanning routes while avoiding obstacles within the coal yard. Therefore, we opted to use multiple fixed LiDARs to scan the coal pile. Each LiDAR was mounted on a servo, which rotated the LiDAR to scan a local area of the coal pile. Ultimately, the point clouds obtained from the multiple LiDAR scans were stitched together. With the known spatial relationships of the LiDARs’ coordinates, we could conveniently merge the coal pile point clouds, thereby circumventing the complex process of point cloud registration.

The point cloud processing workflow employed in this study is illustrated in Figure 1. Initially, a series of LiDAR scans were utilized to acquire point cloud data from various sections of the coal stockpile. To mitigate noise arising from random errors in servo motor rotation angles, the RPCS-CPF algorithm was applied to each cloud data point for refined filtering. Subsequently, coordinate transformations were applied to the filtered point clouds prior to fusion, followed by another round of RPCS-CPF algorithmic treatment to obtain a smooth point cloud representation of the entire coal pile. Building upon this framework, a gradient-based region growing clustering algorithm is proposed for identifying boundaries on the coal pile surface, with the effectiveness of the algorithm corroborated through the meticulous construction of point cloud models.

2.1. Point Cloud Coordinate Transformation

The original point cloud’s coordinates are aligned with the three-dimensional coordinate system of the LiDAR. To analyze the 3D point clouds collected by the LiDAR at various angles, it is essential to calculate the coordinate transformation between them. As illustrated in Figure 2, the actuator coordinate system rotates around point O, which is the center of rotation. Assuming the distance from point O to the center of the LiDAR is r, the actuator captures a 3D point cloud of the coal pile with each rotation of ${\Delta }_{\phi }$. At each distinct rotational position (n = 0, 1, …), the 3D points are expressed in the actuator’s coordinate system ($X_n^s$, $Y_n^s$, $Z_n^s$) with coordinates ($x_n^s$, $y_n^s$, $z_n^s$). The actuator’s rotation is around the $X_n^s$-axis.

Then, the coordinates of point cloud in servo coordinate system can be obtained as follows:

$$P_n^{\mathrm{S}}=\left[\begin{array}{c}{x}_n^{\mathrm{S}}\\ y_n^{\mathrm{S}}\\ z_n^{\mathrm{S}}\end{array}\right]=\left[\begin{array}{c}{x}_n^{\mathrm{L}}\\ y_n^{\mathrm{L}}\\ z_n^{\mathrm{L}}+r\end{array}\right]=P_n^{\mathrm{L}}+\left[\begin{array}{c}0\\ 0\\ r\end{array}\right]$$

(1)

where $P_n^L$ represents the coordinates of the point in the LiDAR coordinate system ($X_n^L$, $Y_n^L$, $Z_n^L$). During the rotation of the actuator, since the rotation axis is parallel to the X-axis, the x-coordinate of the point remains unchanged.The relationship of transformation between the coordinates of the point after the n-th rotation and the coordinates in the initial state (n = 0) can be expressed as follows:

$$P_n^{\mathrm{S}}=\left[\begin{array}{c}{x}_n^{\mathrm{S}}\\ y_n^{\mathrm{S}}\\ z_n^{\mathrm{S}}\end{array}\right]=R_n^{\mathrm{S}}\left[\begin{array}{c}{x}_0^{\mathrm{S}}\\ y_0^{\mathrm{S}}\\ z_0^{\mathrm{S}}\end{array}\right]$$

(2)

where the rotation matrix $R_n^S$ is:

$$R_n^{\mathrm{S}}=\left[\begin{array}{ccc}1& 0& 0\\ 0& cos(n·\Delta \phi )& -sin(n·\Delta \phi )\\ 0& sin(n·\Delta \phi )& cos(n·\Delta \phi )\end{array}\right]$$

(3)

According to Formula (1), Formula (2) can be expressed as follows:

$$P_n^{\mathrm{L}}+\left[\begin{array}{c}0\hfill \\ 0\hfill \\ r\hfill \end{array}\right]=P_n^{\mathrm{S}}=\left[\begin{array}{c}{x}_n^{\mathrm{S}}\\ y_n^{\mathrm{S}}\\ z_n^{\mathrm{S}}\end{array}\right]=R_n^{\mathrm{S}}\left[\begin{array}{c}{x}_0^{\mathrm{S}}\hfill \\ y_0^{\mathrm{S}}\hfill \\ z_0^{\mathrm{S}}\hfill \end{array}\right]=R_n^{\mathrm{S}}\left(P_0^{\mathrm{L}}+\left[\begin{array}{c}0\hfill \\ 0\hfill \\ r\hfill \end{array}\right]\right)$$

(4)

Therefore, we derive the subsequent transformation relation:

$$P_n^{\mathrm{L}}=R_n^{\mathrm{S}}{P}_0^{\mathrm{L}}+\left(R_n^{\mathrm{S}}-1\right)\left[\begin{array}{c}0\hfill \\ 0\hfill \\ r\hfill \end{array}\right]$$

(5)

$$P_0^{\mathrm{L}}={\left(R_n^{\mathrm{S}}\right)}^{-1}{P}_n^{\mathrm{L}}+\left({\left(R_n^{\mathrm{S}}\right)}^{-1}-1\right)\left[\begin{array}{c}0\hfill \\ 0\hfill \\ r\hfill \end{array}\right]$$

(6)

Based on Equation (6), it is feasible to transform the point cloud coordinates from various scanning angles into the initial state within the LiDAR coordinate system. This makes it easier to integrate point clouds from various scanning angles on a single LiDAR. A thorough point cloud of a sizable coal pile can be obtained by using many LiDARs and matching their coordinate systems into a single global coordinate system.

Given the large area of a coal pile, the comprehensive acquisition of its point cloud data generally requires the deployment multiple LiDARs positioned above the pile to facilitate scanning. Following this, it is necessary to integrate the point clouds acquired from several LiDAR scans. Nevertheless, the accuracy of rotation cannot be assured as a result of the gravitational impact exerted by the LiDAR on the actuator. Consequently, following the process of coordinate transformation, the point clouds exhibit an inability to be effectively integrated. Given the aforementioned concern, it is evident that Equation (3) is no longer capable of accurately depicting the rotation matrix of the point clouds. Consequently, it is restructured in the following manner:

$$R_n^{{\mathrm{S}}^{\prime }}=\left[\begin{array}{ccc}1& 0& 0\\ 0& cos\left(n·\Delta \phi +{\theta }_{\mathrm{n}}\right)& -sin\left(n·\Delta \phi +{\theta }_{\mathrm{n}}\right)\\ 0& sin\left(n·\Delta \phi +{\theta }_{\mathrm{n}}\right)& cos\left(n·\Delta \phi +{\theta }_{\mathrm{n}}\right)\end{array}\right]$$

(7)

where ${\theta }_n$ denotes the error in the nth rotation of the actuator. The rotation deviation is assumed to follow a uniform distribution with a range of [0, ${△}_{\theta }$], where ${△}_{\theta }$ represents the range of rotation error from 0 to 5 degrees. Consequently, the new point cloud coordinates can be expressed as follows:

$$P_0^{{\mathrm{L}}^{\prime }}={\left(R_n^{{\mathrm{S}}^{\prime }}\right)}^{-1}{P}_n^{\mathrm{L}}+\left({\left(R_n^{{\mathrm{S}}^{\prime }}\right)}^{-1}-1\right)\left[\begin{array}{c}0\hfill \\ 0\hfill \\ r\hfill \end{array}\right]$$

(8)

2.2. RPCS-CPF Method

This study introduces a rapid point cloud registration approach called RPCS-CPF to tackle the problem of inadequate robustness in point cloud stitching caused by LiDAR rotation deviation. The algorithm assumes that the rotation plane remains constant. The formula can be expressed as follows:

$$P_{\mathrm{m}}={\left\{MLS\left(G_s\left(G_v\left(P_i\right)\right)\right)\right\}}_{i=1}^N$$

(9)

where $P_i$ represents the original point cloud, $P_m$ denotes the point cloud after undergoing particle filtering, $G_V$ signifies voxel filtering, $G_S$ indicates statistical filtering, and MLS represents the moving least square.

2.2.1. Voxel Filtering

The expanded number of laser emission lines in a 16-line LiDAR system enhances the spatial resolution of the resultant point set, yielding a densely populated point cloud. However, this results in an escalated computational load and a reduction in processing efficiency. To address this, we employed voxel filtering to downsample the point cloud of the coal pile, while preserving its geometric integrity.

Voxel filtering is a method used for downsampling point clouds, wherein the density of the point cloud is reduced by partitioning it into individual cubic cells. This partitioning process ensures that just one representative point is retained within each cell. The process of voxel filtering involves dividing an original point cloud, represented as $P={(x_i,y_i,z_i)}_{i=1}^N$, into cubic cells of size 0.3 m³. The centroid of each cell is used to replace the points within that cell. The point cloud obtained following the application of voxel filtering can be expressed as follows:

$$P_v=\left\{\left(x_i^{\prime },y_i^{\prime },z_i^{\prime }\right)\right\}={\left\{\left(\frac{1}{n}\sum _{j=1}^n{x}_{ij},\frac{1}{n}\sum _{j=1}^n{y}_{ij},\frac{1}{n}\sum _{j=1}^n{z}_{ij}\right)\right\}}_{i=1}^{M_v}$$

(10)

where $P_v$ represents the filtered point cloud after the application of the statistical method, $M_v$ denotes the number of points remaining after the filtering process.

2.2.2. Statistical Filtering

The point cloud obtained from LiDAR scanning displays slight positional discrepancies attributed to imprecisions in the rotation angle of the actuator during each rotation. As a result, the stitched point cloud contains an increased number of outlier points, which adversely affects the reconstruction of the coal pile surface. To alleviate the presence of outlier points and improve the overall smoothness of the point cloud surface, this study employs a statistical filtering technique for point cloud processing.

The process of statistical filtering entails assessing the point cloud by calculating the statistical attributes of individual points and their neighboring points to identify any outliers. The method commences by computing the mean ${\mu }_i$ and standard deviation ${\sigma }_i$ for each point, assuming that the set of the nearest 100 points for each point is represented as $N_i$. Subsequently, it ascertains if a given data point is classified as an outlier. After applying statistical filtration, the resulting point cloud may be expressed as follows:

$${\mathrm{P}}_{\mathrm{s}}={\left\{\left(x_i,y_i,z_i\right)\left|\right|z_i-{\mu }_i\mid <k·{\sigma }_i\right\}}_{i=1}^{{\mathrm{M}}_{\mathrm{s}}}$$

(11)

where $P_s$ represents the filtered point cloud after the application of the statistical method, $M_s$ denotes the number of points remaining after the filtering process, and k is the threshold value that controls the criteria for identifying outlier points.

2.2.3. Moving Least Squares

Following voxel and statistical filtering, the point cloud tends to exhibit sparsity and local point cloud voids. Hence, a two-step moving least squares method was employed to upsample the point cloud and fit smooth surfaces, thereby enhancing the continuity and smoothness of the coal pile point cloud surface and achieving a desirable point cloud density. For each point i in the point cloud $P_s$, a neighbor set $N_i^s$ is selected within a 0.5 m radius, and a plane is fitted using the moving least squares method:

$$z_i=a{x}_i+b{y}_i+c$$

(12)

The projection coordinates on the plane for each point i can be determined as $({\widehat{x}}_i,{\widehat{y}}_i,{\widehat{z}}_i)$. To fill gaps in the point cloud, we incrementally move 5 cm along the x- and y-axes, using Formula (12) to generate new points for the surface. The resultant point cloud is commonly represented as follows:

$$P_m={\left\{{\widehat{\mathrm{x}}}_i,{\widehat{y}}_i,{\widehat{z}}_i\right\}}_{i=1}^{M_m}$$

(13)

where $P_m$ represents the point cloud after being fitted by the least squares method, and $M_m$ denotes the number of points in the fitted $P_m$.

2.3. Edge Detection Algorithm Based on Gradient Clustering

This study introduces a novel approach to address the complex terrain of coal piles by employing a gradient-based region-expanding clustering algorithm for edge detection. The method begins with the calculation of local normal vectors for the point cloud representing the coal pile surface. Subsequently, the local slope is determined based on these normal vectors.

The initial step involves constructing a KD-tree structure, followed by loading the point cloud data into the KD tree to facilitate nearest neighbor searches. Within a radius of $r=0.5\mathrm{m}$, the nearest neighbors for point p are determined.

The covariance matrix is then computed using Equation (14). Utilizing Equation (15), the covariance matrix is decomposed, leading to the determination of eigenvalues and eigenvectors. The eigenvector corresponding to the minimum eigenvalue provides an approximation of the normal vector at the given point.

$$\mathcal{C}=\frac{1}{k}\sum _{i=1}^k·\left({\mathbf{p}}_i-\overline{\mathbf{p}}\right)·{\left({\mathbf{p}}_i-\overline{\mathbf{p}}\right)}^T$$

(14)

$$\mathcal{C}·\overrightarrow{{\mathbf{v}}_j}={\lambda }_j·\overrightarrow{{\mathbf{v}}_j},j\in \{0,1,2\}$$

(15)

where C denotes the covariance matrix, k represents the size of the set of nearest neighbor points, $P_i$ denotes the i-th nearest neighbor point, $\overline{P}$ denotes the centroid of the set of nearest neighbor points, which is the average of all nearest neighbor points, ${\lambda }_j$ represents the j-th eigenvalue, and $\overrightarrow{v_j}$ corresponds to the eigenvector associated with ${\lambda }_j$.

In order to maintain consistency in the acquired normal directions, all normals are oriented towards the view vector $V_p$, as depicted in Figure 3.

$${\overrightarrow{\mathbf{n}}}_i·\left({\mathrm{v}}_p-{\mathbf{p}}_i\right)>0$$

(16)

where $\overrightarrow{n_i}$ represents the normal vector of the i-th point, $V_p$ denotes the view vector, which is aligned with the positive Z-axis direction in this paper, and $P_i$ indicates the coordinates of the i-th point.

After obtaining the normal vectors for each point, the approximate gradient for each point can be calculated using the finite difference method:

$$\nabla f\left(P_i\right)=(\frac{f\left(P_i+h\Delta x\right)-f\left(P_i\right)}{h},\frac{f\left(P_i+h\Delta y\right)-f\left(P_i\right)}{h})$$

(17)

where h is the step size parameter, $\Delta x$ and $\Delta y$ are the unit vectors along the x- and y-axes, respectively. $f(P_i+h\Delta x)$ represents the corresponding height at point $P_i$ after moving h along the x-axis direction.

Once the gradient for each point is obtained, we proceed with clustering using a method that relies on the growth of the gradient region. The first step involves establishing a neighborhood connection, where each point is designated as the center and a neighborhood radius r is defined. Points that are located within a distance of less than r from a certain point are regarded as constituents of the neighborhood of that point. The procedure proceeds to cycle through all seed points, denoted as $S_i\in S$, within the set S. For every seed point $S_i$, it verifies if the points $P_j$ in its vicinity meet the following criteria:

$$cos\left({\theta }_{ij}\right)=\frac{G_i·G_j}{∥G_i∥∥G_j∥}>cos\left({\theta }_{\mathrm{threshold}}\right)$$

(18)

$$\left|\frac{M_i-M_j}{M_i}\right|<\epsilon ,M_i=\left|G_i\right|$$

(19)

where $G_i$ and $G_j$ represent the gradient direction vectors of points $G_i$ and $G_j$, respectively, ${\theta }_{ij}$ is the angle between the two direction vectors, ${\theta }_{\mathrm{threshold}}$ is the maximum allowable angle of gradient direction variation which is set to 15 degrees, $M_i$ and $M_j$ are the magnitudes of the gradients for $P_i$ and $P_j$, respectively, and $\epsilon$ is the maximum allowable relative change in gradient magnitudes, which is set to 2.

Point $P_j$ is incorporated into the same region as $S_i$ if it meets the stated conditions, allowing for a smooth merging process. By means of consecutive iterations, this methodology utilizes gradient information to efficiently propagate the process of region growth.

3. Experiments

3.1. Hardware System

To validate the efficacy of the RPCS-CPF and the gradient-based region-growing clustering algorithm, experimental validation was conducted at Guoneng Ningxia Lingwu Power Generation Co., Ltd., Yinchuan, China. Four LiDARs were deployed across a significant coal pile measuring 200 m by 50 m, ensuring extensive coverage of the coal pile area, as depicted in Figure 4.

We employed 16-line 3D LiDARs as our primary sensors. LiDARs are instrumental in acquiring precise three-dimensional positional data, allowing for the determination of an object’s position, size, external morphology, and even material composition. The LiDARs chosen for this study possess a horizontal field of view spanning 360^∘ and a vertical field of view of 30°. To facilitate the scanning of large coal piles, we mounted the LiDARs on servo actuators, thereby expanding their vertical scanning range. These actuators rotate the LiDARs to scan the entire coal pile.

The installation setup of the LiDARs is depicted in Figure 5. Initially, the servo actuators were securely affixed beneath the I-beams of the coal storage facility’s canopy, while the LiDARs were mounted on the rotating structures of these actuators. With a maximum rotation angle of 180^∘, the LiDARs could capture a comprehensive point cloud of the coal pile.

3.2. Single LiDAR Particle Filter Results

The coal pile under investigation is illustrated in Figure 6. However, owing to the restricted field of view of the camera utilized, the representation is confined to a segment of the entire coal pile.

The point clouds of the coal pile acquired from the four lidar scans are shown in Figure 7.

To achieve seamless and continuous point cloud data, we initially employ the RPCS-CPF algorithm on the point cloud generated from a single LiDAR scan. Subsequently, the same RPCS-CPF algorithm was applied to the point cloud 3 of the coal pile, yielding the point cloud depicted in Figure 8.

3.3. Point Cloud Registration

We first conducted a coordinate transformation on the point cloud depicted in Figure 7 to align it with the coordinate system of the coal yard. Subsequently, a straightforward stitching operation was performed on the point cloud. Following the stitching process, we filtered the point cloud corresponding to the greenhouse section based on their coordinates. The resulting stitched point cloud is illustrated in Figure 9.

Figure 9 indicates that the point cloud of the coal pile obtained from the 16-line LiDAR scan is highly dense, witch leads to a significant computational load for subsequent point cloud processing. The objective of our study is to perform downsampling on the point cloud while maintaining the integrity and smoothness of the surface. To achieve this, we proceeded with voxel filtering. The point cloud after voxel filtering is presented in Figure 10a.

Statistical filtering. Figure 10a reveals that the merged point cloud contains a significant number of outlier points. This is attributed to the overlapping regions scanned by two adjacent LiDARs, where the presence of servo actuator rotation errors prevented the proper alignment of the point clouds within these areas. To mitigate this issue, statistical filtering was implemented to eliminate the outliers, with the post-filtering results displayed in Figure 10b.

Moving least squares. After using voxel and statistical filtering procedures, a significant decrease in point cloud density was observed, resulting in a notably smoother surface characterized by a minimal presence of outlier points. The utilization of a moving least squares technique was employed for surface fitting in order to improve the continuity and smoothness of the resultant point cloud. The results of this fitting procedure are depicted in Figure 10c,d. The point cloud of the coal pile surface obtained after the application of the moving least squares approach demonstrates enhanced smoothness, as there are no observable outliers in the merged sections. Nevertheless, following the process of voxel downsampling and statistical filtering, the point cloud underwent a reduction in density, resulting in the presence of several regions with localized point loss. In order to tackle this issue, an additional iteration of moving least squares upsampling was implemented, and the outcomes are illustrated in Figure 10e.

As depicted in Figure 10, the RPCS-CPF algorithm demonstrated noteworthy outcomes. The combined point cloud preserves the coherence and uniformity while precisely representing the topographic arrangement of the coal pile’s surface. In contrast to alternative techniques like point cloud registration, this particular approach exhibits reduced computing complexity, thereby indicating its wide-ranging potential for various applications.

Smoothness comparison. To compare the point cloud stitching effects before and after the application of the particle filtering algorithm, we extracted the point cloud at the stitching areas and voxelized them with a size of 0.5 × 0.5 along the x- and y-coordinates. The standard deviation of the height values and the maximum height difference were subsequently computed for each point cloud. The mathematical expression used to get the standard deviation is as follows:

$$\sigma =\sqrt{\frac{{\sum }_{\mathrm{i}=1}^{\mathrm{n}}\left(z_i-{\mu }^2\right)}{n}}$$

(20)

where $\sigma$ represents the standard deviation of the Z-values within the same voxel, $z_i$ is the height value of the i-th point within the voxel, and $\mu$ is the average height value of the points within the same voxel. The average standard deviation and average maximum difference of all voxels near the same stitching boundary are used to measure the smoothness at the stitching site, as shown in

Table 1. Comparison of smoothness at the stitching before and after the RPCS-CPF algorithm.

Standard Deviation—Before	Maximum Difference—Before	Standard Deviation—After	Maximum Difference—After
0.47	1.75	0.19	0.82
0.35	1.17	0.20	0.85
0.47	1.48	0.21	0.90

.

Figure 11 displays the point cloud that was extracted in close proximity to the stitching line. The point cloud’s range in the X direction has been extended to enhance clarity. Figure 11a,c depict the point cloud from various viewpoints prior to the implementation of the particle filtering technique. The presence of layering in the up and down regions of the stitching site is evident, indicating a lack of smooth flow. Figure 11b,d depict the point cloud pictures subsequent to the implementation of the particle filtering technique. Based on the presented data, it is apparent that the stitching site exhibits minimal layering and a seamless transition, suggesting that the RPCS-CPF algorithm, as proposed, has yielded satisfactory outcomes.

3.4. Edge Detection

We conducted algorithm validation within an enclosed coal yard, as depicted in Figure 12. This paper introduces a region-growing clustering algorithm that leverages preliminary edge cues derived from gradient information. The algorithm initiates with seed points exhibiting strong gradient responses and incrementally incorporates neighboring pixels into the same region, adhering to a predetermined consistency in gradient direction, until a specified termination condition is met.

The experimental results affirmed the algorithm’s proficiency at defining the boundaries of unique features within the coal pile. The algorithm was first utilized to extract the boundary of the coal pile as shown in Figure 13. Following this, visual inspection revealed that the algorithm could also identify the edges of critical terrain features, including depressions (Figure 14a), corridors (Figure 14b), and ridges (Figure 14c). This observation indicates that, although our algorithm is primarily intended for precise extraction of the coal pile boundary, it effectively captures the boundaries of specific terrain features as well.

4. Discussion

Coal piles, as fundamental units of management, allow for the convenient measurement of daily coal usage and the efficiency of coal combustion when their point clouds are obtained. The current popular methods for acquiring point clouds of large coal piles involve using drones equipped with LiDARs for scanning or employing movable 2D LiDAR devices. Both methods have certain drawbacks. For example, scanning with drones requires planning complex flight paths, and, if the internal environment of the coal yard changes, the original path may no longer be suitable. Although 2D LiDAR can capture the entire point cloud of a coal pile, the scanning speed is slow, and the resulting point cloud is relatively sparse. Therefore, we propose that using multiple LiDAR devices to scan large coal piles and then stitching together the obtained segments can both increase the scanning speed and ensure the completeness of the point cloud.

In the experimental portion of our study, we deployed four LiDAR devices, each affixed to a servo to allow for scanning across the breadth of the coal pile. However, when the servos operated the LiDAR devices to scan the coal pile, there was a rotational discrepancy due to the effects of gravity, which introduced noise into the scanned point cloud. This generated two primary issues: firstly, the point cloud derived from the LiDAR scan was non-continuous and contained numerous outliers; secondly, the post-stitching point cloud did not merge seamlessly, as illustrated in Figure 9, where a distinct demarcation line is visible at the juncture. This demarcation is a direct consequence of the servo’s rotational error. Consequently, the RPCS-CPF method introduced in this paper is of significant importance. In this section, we will be comparing the stitching method predicated on CPF with the original method, which relies on coordinate transformation, from the perspectives of point cloud smoothness and the number of parameters involved.

4.1. Smoothness Comparison

To verify the effectiveness of our CPF (Conditional Point Filtering) filter, we conducted comparative experiments by applying the CPF filter to both individual point clouds and the entire point cloud, as shown in

Table 2. Comparison of smoothness between CPF and the original method.

	Standard Deviation	Maximum Difference	SP	EP
Original method	0.43	1.47
CPF (ours)	0.31	1.22	√
CPF (ours)	0.26	1.05		√
CPF (ours)	0.2	0.85	√	√

“SP” denotes the application of CPF to an individual point cloud, whereas “EP” signifies its application to the concatenated point cloud.

. We have compiled and presented the mean values of the standard deviation and the maximum difference at three coal pile point cloud stitching locations. It is evident from the table that, prior to the CPF process, the point clouds were merely stitched together through coordinate rotation transformations. However, due to the presence of servo rotational errors, it was nearly impossible to ensure a successful stitch at every junction. After applying our CPF method to filter the individual point clouds followed by the entire point cloud, this approach substantially lessened the influence of servo rotational errors during the stitching process of point clouds. Therefore, it is essential to apply CPF filtering to individual point clouds as well as to the entire point cloud ensemble.

4.2. Comparison of Parameters

As can be observed from Figure 9, the point cloud immediately after scanning is extremely dense, containing a significant amount of noise and outliers. To reduce the density of the point cloud and enhance computational speed, downsampling was initially applied to the point cloud. This was followed by upsampling to achieve a better presentation. To facilitate a more comprehensive demonstration, statistics on the quantity and density of the point clouds at each step were compiled, as shown in

Table 3. Comparison of density between CPF and the original method.

	Number	Density	Voxel Filtering	Statictical Filtering	MLS-Upsample
Original method	1,057,584	105
CPF (ours)	100,510	10	√
CPF (ours)	99,689	9	√	√
CPF (ours)	797,512	80	√	√	√

. It can be seen from the table that, after downsampling the initial point cloud, the number of parameters was reduced by 90%, and the processing time was also reduced by the same percentage. Upon final upsampling, the density of the point cloud was similar to that before processing, yet the smoothness and integrity of the point cloud were significantly improved, which indicates that our method has yielded satisfactory results.

4.3. Limitations

This study still has several limitations in the acquisition and processing of coal pile point cloud data. Due to the mechanical structure, it is challenging for multiple LiDAR devices to achieve fully synchronized scanning, which may lead to temporal errors during point cloud stitching. These temporal discrepancies could potentially impact the consistency and accuracy of the point cloud data. Although the CPF proposed in this study has demonstrated satisfactory performance in practical applications, there is scope for improvement in terms of algorithmic precision and the extent of surface detail restoration of the coal pile point cloud. The current approach may not fully recover all the details of the coal pile surface, particularly in instances where significant discontinuities are present on the coal pile surface. We intend to explore more efficient and accurate algorithms in future research to replicate the surface characteristics of the coal pile point cloud as closely as possible.

5. Conclusions

This paper focused on addressing the efficient stitching of coal pile point cloud data and the precise detection of boundaries, with the aim of enhancing safety and efficiency in the coal mining industry. The primary contributions of this study are as follows:

(1)

The RPCS-CPF (Rapid Point Cloud Stitching–Constrained Particle Filter) algorithm is proposed, specifically optimized for integrating point cloud data in large-scale coal pile environments. Experimental validations conducted on real large-scale coal piles demonstrated the unique advantages of this algorithm. It not only facilitates smooth transitions in stitched areas, but also ensures the consistency and integrity of the overall point cloud data while preserving the detailed geometric features of the coal pile surface.

(2)

Proposal of an edge detection algorithm based on gradient region expanding clustering to address the complex surface characteristics of coal piles. Experimental results validated the capability of this method to accurately identify boundaries, thereby significantly contributing to safety assessments and guidance in coal mining operations.

In conclusion, this study has not only presented the novel RPCS-CPF algorithm for point cloud stitching and edge detection, but also provided a foundation for future research. The outputs of this research have immediate practical applications in the coal mining industry, improving both the safety and efficiency of operations. For future work, we suggest exploring the integration of these algorithms with real-time monitoring systems and expanding the study to include a wider range of coal pile environments and conditions. Additionally, we will continue to explore more efficient and precise algorithms to replicate the surface characteristics of coal pile point clouds as closely as possible.