Development of Sensing Algorithms for Object Tracking and Predictive Safety Evaluation of Autonomous Excavators

: This article presents the sensing and safety algorithms for autonomous excavators operating on construction sites. Safety is a key concern for autonomous construction to reduce collisions and machinery damage. Taking this point into consideration, our study deals with LiDAR data processing that allows for object detection, motion tracking/prediction, and track management, as well as safety evaluation in terms of potential collision risk. In the safety algorithm developed in this study, potential collision risks can be evaluated based on information from excavator working areas, predicted states of detected objects, and calculated safety indices. Experiments were performed using a modiﬁed mini hydraulic excavator with Velodyne VLP-16 LiDAR. Experimental validations prove that the developed algorithms are capable of tracking objects, predicting their future states, and assessing the degree of collision risks with respect to distance and time. Hence, the proposed algorithms can be applied to diverse autonomous machines for safety enhancement.


Introduction
Excavation is a major task that is performed in the construction sector. It is the process of earthmoving, such as creating building foundations, land grading, pile driving, highway maintenance, etc. Therefore, excavators are considered essential equipment for any construction. With the increasing number of fatal accidents on construction sites, researchers have conducted numerous studies to raise the issue of safety concerns. The construction sites for excavation are dynamic as they include interactions between different objects such as machines and workers, which can increase the chance of potential collisions due to poor visibility, blind spots, and operator carelessness. Thus, there has been a need for more efficient and advanced safety management to tackle such high-risk accidents.
To fully achieve autonomous operations, advanced safety systems are required to protect the workers and construction infrastructure by gathering information from machine dynamics (proprioceptive) and the external environment (exteroceptive) [1]. Various safety management systems have been proposed for construction equipment to detect dangerous proximities, alert operators and workers about these, and identify blind spots by using sensing data obtained from RF (radio frequency) ID techniques, GPS (global positioning system), cameras, and laser scanners [2].
During the past decade, RFID (radio frequency identification) has been the most widely used sensor for safety management, followed by wireless, camera-based, and UWB (ultra-wide-band) sensors [3]. The RTLS (real-time location system) developed by Lee et al. [4], was intended to enhance worker safety using RFID tags. The implementation of this method requires tags to be linked with the target objects and the leaving-out of objects without tags. However, the study lacks in providing a strategy to prevent collisions. safe margins. Therefore, quantitative and predictive metrics are highly necessary to assess unsafe conditions and potential risks in fully automated excavators.
To tackle the above issue, this study proposes 3D Lidar-based safety algorithms for autonomous excavators that can effectively detect, track, and predict object motion around the excavator, and then evaluate the severity of collision risks using the defined risk metric.
Section 1 describes the introduction and literature review of state-of-the-art methodologies. Section 2 briefly explains a framework of the proposed methodologies that covers object tracking and safety evaluation. Section 3 presents the sensory data processing and safety index calculations. Section 4 provides experimental results and analysis, and the paper is concluded in Section 5.

Framework of Object Tracking and Safety Evaluation
The primary goal of this research is to develop a safety module with 3D Lidar to evaluate the severity of potential risks for an autonomous excavator with respect to time and distance by predicting the motions of the detected objects. Figure 1 presents the framework of this study to achieve our objectives. The framework comprises four components that include object detection and tracking (C1), calculation of excavation working area through kinematic analysis (C2), calculation of potential collision severity with respect to time and distance (C3), and safety evaluation (C4).
The majority of prior studies in the area of safety management focused on proximitydetection and utilization of current states of the detected objects that may not be comprehensive in capturing potential safety risks. For example, the proximity between construction equipment and workers may not be used to indicate potential unsafe conditions if a manipulator or vehicle remains stationary (i.e., does not move) or moves rapidly despite sufficient safe margins. Therefore, quantitative and predictive metrics are highly necessary to assess unsafe conditions and potential risks in fully automated excavators.
To tackle the above issue, this study proposes 3D Lidar-based safety algorithms for autonomous excavators that can effectively detect, track, and predict object motion around the excavator, and then evaluate the severity of collision risks using the defined risk metric.
Section 1 describes the introduction and literature review of state-of-the-art methodologies. Section 2 briefly explains a framework of the proposed methodologies that covers object tracking and safety evaluation. Section 3 presents the sensory data processing and safety index calculations. Section 4 provides experimental results and analysis, and the paper is concluded in Section 5.

Framework of Object Tracking and Safety Evaluation
The primary goal of this research is to develop a safety module with 3D Lidar to evaluate the severity of potential risks for an autonomous excavator with respect to time and distance by predicting the motions of the detected objects. Figure 1 presents the framework of this study to achieve our objectives. The framework comprises four components that include object detection and tracking (C1), calculation of excavation working area through kinematic analysis (C2), calculation of potential collision severity with respect to time and distance (C3), and safety evaluation (C4). The first component deals with the estimation and tracking of object positions and velocities through adopting the IMM (interacting multiple model) coupled with a UKF (unscented Kalman filter) and JPDA (joint probabilistic data association). This algorithm facilitates object tracking in a robust manner and the track management can automatically assign track IDs to each detected object [16]. The reason why the IMM-UKF-JPDA approach was selected was to handle the various motion models of different moving objects along with non-linear states. The second part of this work is to calculate the actual and predicted working areas of a modified mini hydraulic excavator by applying a kinematic analysis. Using the position states of the detected objects and working areas and applying the defined safety indices, the time and distance to a collision are estimated as the next step. The first component deals with the estimation and tracking of object positions and velocities through adopting the IMM (interacting multiple model) coupled with a UKF (unscented Kalman filter) and JPDA (joint probabilistic data association). This algorithm facilitates object tracking in a robust manner and the track management can automatically assign track IDs to each detected object [16]. The reason why the IMM-UKF-JPDA approach was selected was to handle the various motion models of different moving objects along with non-linear states. The second part of this work is to calculate the actual and predicted working areas of a modified mini hydraulic excavator by applying a kinematic analysis. Using the position states of the detected objects and working areas and applying the defined safety indices, the time and distance to a collision are estimated as the next step.
Lastly, the collision risks at the current moment were categorized into safe, warning, and emergency based on the computed safety index values. The aforementioned algorithms were developed in the MATLAB environment. Detailed explanations on each step are provided in the following section.

Data Processing and Kinematic Analysis
The extraction and processing of sensor data are described in this section, as shown in Figure 2. The VLP-16 LiDAR [18] (VLP-16, Velodyne Lidar, San Jose, CA, US, 2020) was selected as the sensor used for this research as it provides a horizontal detection range of 360 degrees and a vertical range of 30 degrees, respectively. It also generates approximately Appl. Sci. 2021, 11, 6366 4 of 21 300,000 points per second, which allows for fast and high-resolution detection. In addition, the point cloud of each Lidar includes the information of the x, y, and z coordinates that can be used to identify the position of the detected objects. In the study, MATLAB (R2020a, MathWorks, Natick, MA, US, 2020) software was used to store, process, and visualize the point cloud data. The data processing for object tracking and track management was taken through the steps discussed below.
are provided in the following section.

Data Processing and Kinematic Analysis
The extraction and processing of sensor data are described in this section, as shown in Figure 2. The VLP-16 LiDAR [18] (VLP-16, Velodyne Lidar, San Jose, CA, US, 2020) was selected as the sensor used for this research as it provides a horizontal detection range of 360 degrees and a vertical range of 30 degrees, respectively. It also generates approximately 300,000 points per second, which allows for fast and high-resolution detection. In addition, the point cloud of each Lidar includes the information of the x, y, and z coordinates that can be used to identify the position of the detected objects. In the study, MATLAB (R2020a, MathWorks, Natick, MA, US, 2020) software was used to store, process, and visualize the point cloud data. The data processing for object tracking and track management was taken through the steps discussed below.

Extraction of Ground Point Cloud
The raw point cloud data represents many points containing ground and non-ground points ( Figure 3a). Since our focus is on the detection of moving objects around the excavator, the point clouds representing a ground place (earth surface) firstly need to be eliminated from the point cloud. This is because the moving objects belong to the non-ground points and the ground plane can act as unwanted points during the object detection. Therefore, the extraction of these unwanted points is an essential step, and it can also help reduce the computation load by removing them from the large quantity of unprocessed raw data of LiDAR measurements.
In this research, the point cloud data representing the ground plane were extracted by adopting the MSAC (M-estimator sample and consensus) algorithm as a variant of RANSAC (random sample consensus) [19]. This is a repetitive method to estimate a model from a given dataset containing outliers. To extract the ground points, a plane model was fitted to the raw points cloud by the MSAC algorithm. The RANSAC algorithm can determine the number of outliers by selecting a subset of data and fitting a model to it.
The point clouds representing non-ground points are shown in Figure 3b. After detecting the non-ground points, the next step was to cluster the non-ground points and detect the obstacles (Figure 3c), which is discussed in the next section.

Extraction of Ground Point Cloud
The raw point cloud data represents many points containing ground and non-ground points ( Figure 3a). Since our focus is on the detection of moving objects around the excavator, the point clouds representing a ground place (earth surface) firstly need to be eliminated from the point cloud. This is because the moving objects belong to the non-ground points and the ground plane can act as unwanted points during the object detection. Therefore, the extraction of these unwanted points is an essential step, and it can also help reduce the computation load by removing them from the large quantity of unprocessed raw data of LiDAR measurements.
In this research, the point cloud data representing the ground plane were extracted by adopting the MSAC (M-estimator sample and consensus) algorithm as a variant of RANSAC (random sample consensus) [19]. This is a repetitive method to estimate a model from a given dataset containing outliers. To extract the ground points, a plane model was fitted to the raw points cloud by the MSAC algorithm. The RANSAC algorithm can determine the number of outliers by selecting a subset of data and fitting a model to it.
The point clouds representing non-ground points are shown in Figure 3b. After detecting the non-ground points, the next step was to cluster the non-ground points and detect the obstacles (Figure 3c), which is discussed in the next section.

Obstacle Points
After the removal of non-ground points, the remaining point clouds were formed into clusters. However, these forming clusters are still computationally high, as the LiDAR has a long detection range, thus a region of interest (ROI) was restricted by the defined radius to solve this problem. The point cloud points inside the ROI were considered obstacle points [20]. Clustering was carried out with the obstacle points.

Clustering and Object Detection
The Euclidean distance required for the clustering process was set to 0.4 m (minimum) 2 m (maximum) that correspond to the averages of a human (worker) chest width and vehicle size, respectively. Then, the only objects detected within this range were considered for clustering. In this way, the objects that are smaller than a human and bigger than a vehicle were ignored. As the next step, a bounding box was assigned to fit into a single cluster that indicates a detected object. The bounding box generates the center point of each detected object in terms of x, y, and z coordinates. This center point was used for a tracking filter to estimate and predict the object states.

Obstacle Points
After the removal of non-ground points, the remaining point clouds were formed into clusters. However, these forming clusters are still computationally high, as the LiDAR has a long detection range, thus a region of interest (ROI) was restricted by the defined radius to solve this problem. The point cloud points inside the ROI were considered obstacle points [20]. Clustering was carried out with the obstacle points.

Clustering and Object Detection
The Euclidean distance required for the clustering process was set to 0.4 m (minimum) ~2 m (maximum) that correspond to the averages of a human (worker) chest width and vehicle size, respectively. Then, the only objects detected within this range were considered for clustering. In this way, the objects that are smaller than a human and bigger than a vehicle were ignored. As the next step, a bounding box was assigned to fit into a single cluster that indicates a detected object. The bounding box generates the center point of each detected object in terms of x, y, and z coordinates. This center point was used for a tracking filter to estimate and predict the object states.

Object Tracking and Management
To track multiple objects around the excavator, we integrated several filters to design our tracking filter composed of IMM (interacting-multiple-model), UKF (unscented Kalman filter), and JPDA (joint probabilistic data association filter) techniques. The reason for choosing the integrated filters was to handle multiple models representing different mo-

Object Tracking and Management
To track multiple objects around the excavator, we integrated several filters to design our tracking filter composed of IMM (interacting-multiple-model), UKF (unscented Kalman filter), and JPDA (joint probabilistic data association filter) techniques. The reason for choosing the integrated filters was to handle multiple models representing different motion behaviors (IMM), estimate their non-linear states (UKF), and provide a tracking management system that allows for multitarget tracking in clutter at the stage of the data association (JPDA). The IMM-UKF-JPDA approach is computationally effective for recursive state estimation and the mode probabilities of targets. For the implementation of this filter, we consider model set M = M j r j=1 with r = 2 meaning that the M 1 , constant velocity (CV), and M 2 , constant turn (CT) models are considered. These two models share a common state vector x = p x , v x , p y , v y , p z , v z , θ, ω, l, w, h T consisting of bounding box center points in the x, y, and z directions p x , p y , and p z , their corresponding velocities (v x , v y , and v z ), and the yaw angle (θ), turn rate (ω), and the length (l), width (w), and height (h) of a bounding box. The discrete time-state equation for the CV model is given as: and the state for the CT model is: where T is the sampling time, u k is the control input of the excavator, and both models include a white Gaussian noise W 1,k , and W 2,k having a mean at zero and the covariance matrices Q 1,k and Q 2,k , respectively. Such models are common for maneuvering multiobject tracking and are extensions of the models found in [21]. Among the two models, the common measurement model is given as: with the measurement noise V k having the covariance matrix R k . The next step is to estimate the probability to identify the measurements according to the CV and CT models. The progress of the system among the models is operated on the Markovian model transition probability matrix.
where p i,j represents the probability of mode transition. Several studies have addressed the integration of the above filter, but this paper follows the five-step process introduced in [21,22]. These steps include (1) interaction, (2) state prediction and measurement validation, (3) data association and tracking management, (4) mode probability update, and (5) combination. Each step can be explained in mathematical representations as follows.

Interaction
In the interaction step, the initial state and covariance of the model filter are probabilistically mixed with the state and covariance of the previous time step.
The conditional mode probabilities denoted by µ (i|j),k−1 indicate that the system's transition from mode i (previous time step) to mode j (current time step), and can be calculated as: The conditional mode probabilities are dependent on the priori mode probabilities of the current frame denoted by As seen in Equation (6), the conditional mode probabilities of the current frame can be expressed using the prediction previous frame's probabilities, and elements of the matrix from Equation (4).

State Prediction and Measurement Validation
The next step is to predict the current states using the stochastic non-linear model. To deal with nonlinearities, UKF is implemented to predict the states and covariances. x Sigma points, X i with i = 0, . . . , 2n are selected with a square root decomposition of the mixed initial covariance of each filter. Each sigma point has scalar weights w m i and w c i that are associated, and can be calculated as: where α U , β U , λ, and γ are the scalar parameters. The sigma points are propagated with the system function, f. As seen in Equation (7), the sigma points are propagated using the transition function f, and the weighted sigma points are used to calculate the predicted state (x − k ) and covariance (P − k ). Similarly, the sigma points are projected by the observation function h in Equation (8) and then the weighted sigma points are involved in generating the predicted measurement (ẑ − k ) and predicted measurement covariance (S k ). A measurement is regarded valid if it falls within the following elliptical validation gate G(γ G ). The validation gate of an object q is the same for all the models in M.
where i q is the index of the model in M, γ G is the gate threshold, and y k is the set of validated measurements.

Tracking Management and Data Association
Data association is the step of labeling each detected object with a unique ID and associating the tracking information obtained from a LiDAR sensor that includes the estimation/prediction of the detected object's position/velocity, with each label. The data association in the tracking filter was dealt with through JPDA [23], which follows the procedure of initialization, confirmation, correction, prediction, and deletion of tracks. Among various data association algorithms, JPDA enables the calculation of all possible target association probability of each observation (see Equations (14)- (17)).
The 'tentative' label was given to a newly detected track through the data association process; it was labeled as 'confirmed' if a specific track was sufficiently detected to record; the identical ID was assigned if the same track as the current one was detected. The logic flow [24] of JPDA starts with the division of detections from multiple sensors into multiple groups (however, a single LiDAR was used in this study). For each sensor, the distances from detections to existing tracks and a validation matrix were calculated by turns. After the calculation of the validation matrix, tracks and detections were separated into clusters.
The following steps outline how to update each cluster: • S1: Generation of feasible joint events. To formulate JPDA, we assume to have x 1 t , . . . , x N t and z 1 t , . . . , z M t as the states of N objects and M measurements at time t, respectively. Here, the state vector x j t means the position and velocity for the j th detected object and the measurements are the cluttered and noisy detected positions observed by a sensor.
The probability of the data association, p t d j i which represents the measurement index i [M] 0 {0, 1, . . . , M}, is defined as follows [24].
where p D denotes the detection probability, β represents the density of the false detection, x j t is the prediction state of the object, ∑ S is the innovation covariance matrix of the Kalman filter, and N is the normal distribution. p t d j i is assumed to be a linear Gaussian model and its 0 value means missed or dummy detection. The joint data association space, Θ, defined in Equation (15) comprises all possible combinations of measurement-to-object assignments in order that each measurement is assigned to at most one object (A1) and each object is uniquely assigned to a measurement (A2).
where |Θ| =: n h is the total number of joint assignments and θ represents the binary vector that is one possible solution to the data association problem. θ ∈ Θ ⊆ B N×(M+1) . Then, the marginalized JPDA probability, q t d j i over Θ can be calculated as where As the last step, the normalized q t d j i is used to update the state of the j th detected object.

Mode Probability Update
The mode probabilities are updated in this step using the best-fitting measurements and the uniform mixture model likelihoods λ j,k .

Combination
In the last step, the updated states of each filter are used to generate the final states and covariances using the model probabilities of tracks as follows.
3.5. Calculation of the Excavator's Working Area Figure 4 illustrates the positions of the boom, arm, and bucket links. Using the measurements with the stroke sensors and rotational encoder, the position of each link and the rotation of the main body were identified. Stroke sensors and the encoder were mounted on the links of the excavator and its bottom as mentioned in Section 5. Figure 5 presents the measurement data from each sensor. Using the information of the measured stroke and Equations (21)-(23), the angle of each link can be computed. Then, the application of a kinematic analysis with the angles of links enables the calculation of the predicted and actual working areas (see Figure 6).
where µ = α + β. Figure 4 illustrates the positions of the boom, arm, and bucket links. Using the measurements with the stroke sensors and rotational encoder, the position of each link and the rotation of the main body were identified. Stroke sensors and the encoder were mounted on the links of the excavator and its bottom as mentioned in Section 5. Figure 5 presents the measurement data from each sensor. Using the information of the measured stroke and Equations (21)-(23), the angle of each link can be computed. Then, the application of a kinematic analysis with the angles of links enables the calculation of the predicted and actual working areas (see Figure 6).
where = +   In Figure 6, 1 denotes the arm joint, 2 represents the bucket joint, 3−2 is the midpoint of the bucket, and 3−1 shows the bucket tip, shows the distance between the body center (origin of the defined coordinate system) and the boom link. The maximum working area can be achieved by the full extension of the links that can be denoted by   In Figure 6, 1 denotes the arm joint, 2 represents the bucket joint, 3−2 is the midpoint of the bucket, and 3−1 shows the bucket tip, shows the distance between the body center (origin of the defined coordinate system) and the boom link. The maximum working area can be achieved by the full extension of the links that can be denoted by  In Figure 6, P 1 denotes the arm joint, P 2 represents the bucket joint, P 3−2 is the midpoint of the bucket, and P 3−1 shows the bucket tip, d shows the distance between the body center (origin of the defined coordinate system) and the boom link. The maximum working area can be achieved by the full extension of the links that can be denoted by x 3−1 (in this case, bucket tip position from the origin).
l boom cos(θ boom ) + d l boom cos(θ boom ) + l arm cos(θ boom + θ arm ) + d l boom cos(θ boom ) + l arm cos(θ boom + θ arm ) + l bucket cos(θ boom + θ arm + θ bucket ) + d l boom cos(θ boom ) + l arm cos(θ boom + θ arm ) + l bucket,2 cos(θ boom + θ arm + θ bucket + θ bucket,in ) + d where x pre,p,i is the maximum value among the predicted longitudinal displacements, x p is the current longitudinal displacement, p = 1, 2, 3-1, and 3-2 are the current main points of the excavator, and i denotes the prediction step (1,2, . . . , N); ∆x p,j = A Jacobian ∆θ i . The predicted angle of each part was calculated using Equation (24), which follows the Gaussian distribution, having a zero mean and a standard deviation of the noise in angular velocity. The maximum possible value of noise is 3σ. ∆t is the discrete-time increment. i is the prediction step. By using Equations (24) and (25), the x (horizontal position) of the point p (P 1~P3-2 ) is calculated as shown in Equation (26). The predicted horizontal displacements, x pre,p,i can be generated using Equation (27). The maximum radius for the predicted working area is chosen from the maximum value among the predicted displacements. Figure 6 illustrates the horizontal position of p (left subfigure) and the predicted and maximum working areas (right subfigure). One notes that the maximum working area is defined as the maximum working limit, and therefore it does not change since it can be measured by extending the three links fully. The actual working area is determined by calculating the maximum values among all x components.

Safety Indices and Predictive Safety Evaluation
As the required safety indices for the predictive safety evaluation, time to collision (TTC) is defined as the time required for an object to collide with the excavator. The warning index, x, is the normalized parameter to evaluate excavator safety in terms of distance information. To calculate the safety indices, the angle to collision, clearance, and distance between an object and the excavator should be calculated first. The safety indices are calculated based on the excavator's current states and the object's predicted states.

Time to Collision (TTC)
TTC is computed using the current swing acceleration and maximum swing velocity. The following equations were used to calculate TTC. ||r obs || = x obj 2 + y obj 2 . (28) If θ ATC < 0, TTC 2nd ≤ t m , where r obj in Figure 7 denotes the object's distance to the coordinate origin, the object's estimated coordinates are denoted by x ob and y obj , ϕ sw is the swing angle of the excavator, θ ATC shows the angle to a collision, B c is the bucket clearance, t m denotes the predicted time to reach maximum velocity, . θ sw,k is the swing velocity, . θ obj,k presents the angular velocity of a detected object, TTC 2nd is the second-order TTC that represents TTC while considering the current swing acceleration. between an object and the excavator should be calculated first. The safety indices are calculated based on the excavator's current states and the object's predicted states.

Time to Collision (TTC)
TTC is computed using the current swing acceleration and maximum swing velocity. The following equations were used to calculate TTC.
where in Figure 7 denotes the object's distance to the coordinate origin, the object's estimated coordinates are denoted by and , is the swing angle of the excavator, shows the angle to a collision, is the bucket clearance, denotes the predicted time to reach maximum velocity, , is the swing velocity, , presents the angular velocity of a detected object, 2 is the second-order that represents TTC while considering the current swing acceleration.

Warning Index (x)
The warning index, x to evaluate the excavator's working safety is defined in terms of distance [18] and computed as: where c br is the braking distance, α decel,max represents the maximum rotational deceleration of the main body, c sm is the safety clearance, r obj,current is the current object radius, u sensor is the sensor signal noise, u control is the control uncertainty, d min is the minimum clearance.
x is the warning index, c is the distance (clearance) between an object and the excavator. A smaller value of x indicates a higher collision risk and the excavator is required to stop its motion when x becomes one in the defined risk metric (i.e., safety level) as shown in Figure 8.

Warning Index (x)
The warning index, x to evaluate the excavator's working safety is defined in terms of distance [18] and computed as: where is the braking distance, , represents the maximum rotational deceleration of the main body, is the safety clearance, , is the current object radius, is the sensor signal noise, is the control uncertainty, is the minimum clearance. is the warning index, c is the distance (clearance) between an object and the excavator. A smaller value of x indicates a higher collision risk and the excavator is required to stop its motion when x becomes one in the defined risk metric (i.e., safety level) as shown in Figure 8.

Safety Levels
For this study, three safety levels, i.e., safe, warning, and emergency braking were considered, as shown in Figure 8. The safety level is a two-dimensional parameter that is based on the calculated values of and . Primarily, the safety level is the time and distance to a collision. The excavator can continue working without any risk under the safe region. Alarm signals are generated in the warning region for an operator to recognize a hazardous situation. Finally, immediate braking needs to be applied at maximum deceleration to avoid any collision in the emergency region.

Safety Levels
For this study, three safety levels, i.e., safe, warning, and emergency braking were considered, as shown in Figure 8. The safety level is a two-dimensional parameter that is based on the calculated values of TTC and x. Primarily, the safety level is the time and distance to a collision.
The excavator can continue working without any risk under the safe region. Alarm signals are generated in the warning region for an operator to recognize a hazardous situation. Finally, immediate braking needs to be applied at maximum deceleration to avoid any collision in the emergency region.

Experimental Results
A LiDAR sensor (Velodyne VLP-16) in Figure 9a was used to detect objects around the excavator during experimental tests. Data obtained from the LiDAR was post-processed with MATLAB whose libraries allow us real-time processing.

Experimental Results
A LiDAR sensor (Velodyne VLP-16) in Figure 9a was used to detect objects around the excavator during experimental tests. Data obtained from the LiDAR was post-processed with MATLAB whose libraries allow us real-time processing. The algorithms were tested on the autonomous mini excavator at AVEC lab (Ontario Tech University) shown in Figure 9b that was modified from an existing manual type excavator (RMD 500, Yantai Rima Machinery, Yantai, SD, China, 2020) by implementing a HW controller (dSPACE), valve blocks, and several types of sensors.
The developed test platform excavator was comprised of three subsystems, i.e., hydraulic, electronic, and mechanical systems. The hydraulic subsystem consisted of electrohydraulic proportional valves (EHPVs), directional control valves (DCVs), a hydraulic reservoir, hydraulic actuators, and a hydraulic pump. The electronic subsystem was composed of a power supply, relay box, pressure sensors, LVDT stroke sensors (Figure 9c), and a rotary encoder (Figure 9d). The mechanical subsystem was made up of three links (boom, arm, and bucket) and the main body.

Experimental Test Scenario
Experiments were carried out in the presence of multiple static objects, as well as moving ones that were walking around the excavator to represent a collision scenario. The algorithms were tested on the autonomous mini excavator at AVEC lab (Ontario Tech University) shown in Figure 9b that was modified from an existing manual type excavator (RMD 500, Yantai Rima Machinery, Yantai, SD, China, 2020) by implementing a HW controller (dSPACE), valve blocks, and several types of sensors.
The developed test platform excavator was comprised of three subsystems, i.e., hydraulic, electronic, and mechanical systems. The hydraulic subsystem consisted of electrohydraulic proportional valves (EHPVs), directional control valves (DCVs), a hydraulic reservoir, hydraulic actuators, and a hydraulic pump. The electronic subsystem was composed of a power supply, relay box, pressure sensors, LVDT stroke sensors (Figure 9c), and a rotary encoder (Figure 9d). The mechanical subsystem was made up of three links (boom, arm, and bucket) and the main body.

Experimental Test Scenario
Experiments were carried out in the presence of multiple static objects, as well as moving ones that were walking around the excavator to represent a collision scenario.
In Figure 10, one monitored static object T2, and two moving objects, T1 and T3, are indicated. The objects present in the ROI carry their unique labels. If an object leaves the ROI, its corresponding label is deleted, and a new label will be assigned upon returning to the ROI. Therefore, some labels in Figure 10 are missing because the objects that had left the ROI did not present in the ROI anymore. As the object enters the ROI and continues to move closer to the excavator, object detection, tracking, and safety evaluation are performed.
In Figure 10, one monitored static object T2, and two moving objects, T1 and T3, are indicated. The objects present in the ROI carry their unique labels. If an object leaves the ROI, its corresponding label is deleted, and a new label will be assigned upon returning to the ROI. Therefore, some labels in Figure 10 are missing because the objects that had left the ROI did not present in the ROI anymore. As the object enters the ROI and continues to move closer to the excavator, object detection, tracking, and safety evaluation are performed. Note that some points in Figure 10 are not clustered at the current time frame due to the issue of processing time. Therefore, these points may be in a cluster formed after a few frames. Figure 10 shows the raw point cloud data for detected objects. Note that every object has its unique label that helps track management to associate data to each detected object ( 1, 2, and 3). This could also help in counting the number of objects around the excavator. Figure 11 presents the operator's view to see the excavator and moving objects ( 1, 2, and 3) around it. The box (orange color) and fan shape (purple color) in the figure represent the excavator and the operator's visible sight, respectively. The bounding boxes denote the objects around the excavator. The track IDs in Figures 10 and 11 are identical to maintain unique labels for the objects. Note that some points in Figure 10 are not clustered at the current time frame due to the issue of processing time. Therefore, these points may be in a cluster formed after a few frames. Figure 10 shows the raw point cloud data for detected objects. Note that every object has its unique label that helps track management to associate data to each detected object (T1, T2, and T3). This could also help in counting the number of objects around the excavator. Figure 11 presents the operator's view to see the excavator and moving objects (T1, T2, and T3) around it. The box (orange color) and fan shape (purple color) in the figure represent the excavator and the operator's visible sight, respectively. The bounding boxes denote the objects around the excavator. The track IDs in Figures 10 and 11 are identical to maintain unique labels for the objects. The detailed tracking results are presented in Figure 12 that includes the information of the , , and coordinates of each object. These coordinates indicate the current and predicted (1 s ahead) positions of the detected objects. The tracking details were utilized for object tracking and safety evaluation. In Figure 13, the maximum and the actual working areas of the excavator are shown. It also shows the rotation angle of the excavator, the heading angle of tracked moving objects, and their current and predicted states. The detailed tracking results are presented in Figure 12 that includes the information of the x, y, and z coordinates of each object. These coordinates indicate the current and predicted (1 s ahead) positions of the detected objects. The tracking details were utilized for object tracking and safety evaluation. In Figure 13, the maximum and the actual working areas of the excavator are shown. It also shows the rotation angle of the excavator, the heading angle of tracked moving objects, and their current and predicted states. Figure 11. Object detection and tracking around the excavator in the operator's view.

Results
The detailed tracking results are presented in Figure 12 that includes the information of the , , and coordinates of each object. These coordinates indicate the current and predicted (1 s ahead) positions of the detected objects. The tracking details were utilized for object tracking and safety evaluation. In Figure 13, the maximum and the actual working areas of the excavator are shown. It also shows the rotation angle of the excavator, the heading angle of tracked moving objects, and their current and predicted states.  The current area of work represented by the inner circle was obtained by conducting the kinematic analysis, while the maximum working area was achieved through a full extension of three links. The safety indices were calculated based on the actual working area. Note that 2 is the static object, and thus its predicted and current positions did not change.  Figure 14 shows the results of the safety evaluation results computed with the safety indices presented in Figure 8. The current states of T1 and T3 with a circle mark are in the warning and safe regions, respectively. However, the predicted state (1 s ahead) of T1 with the × mark indicates that it is entering the emergency region, and thus the chance of The current area of work represented by the inner circle was obtained by conducting the kinematic analysis, while the maximum working area was achieved through a full extension of three links. The safety indices were calculated based on the actual working area. Note that T2 is the static object, and thus its predicted and current positions did not change. Figure 14 shows the results of the safety evaluation results computed with the safety indices presented in Figure 8. The current states of T1 and T3 with a circle mark are in the warning and safe regions, respectively. However, the predicted state (1 s ahead) of T1 with the × mark indicates that it is entering the emergency region, and thus the chance of collision has highly increased. Therefore, immediate braking should be taken in 1 s to prevent a collision accident for T1 since TTC is less than 1 s. All the detected static objects were far from the maximum working area in the test scenario (i.e., no potential hazard exists in this case), and thus only moving objects are included for the safety evaluation, as seen in Figure 14. Monitoring in Figure 13 and safety evaluation in Figure 14 can work synchronously (see the attached video clip). Finally, the performance of the developed safety algorithms was evaluated in terms of prediction accuracy using the same experimental data from the test scenario in Figures  13 and 14. For this evaluation, we compared the predicted positions of the moving objects T1 and T3 at each time frame with the actual positions 1 s after prediction. Since the minimum Euclidean distance for clustering (i.e., the minimum size of a bounding box to detect moving objects) was set as 0.4 m (average of a worker's chest size), this value was determined as a threshold to judge the success or failure of the prediction. Specifically, the prediction is recognized as 'success' if the error between the predicted and actual positions at each time frame (0.1 s) per object is smaller than the threshold value. Table 1 includes evaluation results to show a successful operation (prediction) rate for each object achieved by the developed algorithms. During the considered test scenario, the successful operation rate for T1 was 85.27% to indicate that the error is less than the threshold in 220 time frames out of the total, 258. The average error across all the time frames (258) was 0.29 m for T1. However, the successful operation rate for T3 was slightly lower (83.59%) than the case of T1, and its prediction accuracy was more or less increased with the decreased average error of 0.28 m. Therefore, it is noted that the proposed algorithms provide satisfactory performance on position prediction in most of the given time frames. Table 1. Evaluation of prediction accuracy of the developed safety algorithms.

T1
T3 No. of time frames with suc- Finally, the performance of the developed safety algorithms was evaluated in terms of prediction accuracy using the same experimental data from the test scenario in Figures 13  and 14. For this evaluation, we compared the predicted positions of the moving objects T1 and T3 at each time frame with the actual positions 1 s after prediction. Since the minimum Euclidean distance for clustering (i.e., the minimum size of a bounding box to detect moving objects) was set as 0.4 m (average of a worker's chest size), this value was determined as a threshold to judge the success or failure of the prediction. Specifically, the prediction is recognized as 'success' if the error between the predicted and actual positions at each time frame (0.1 s) per object is smaller than the threshold value. Table 1 includes evaluation results to show a successful operation (prediction) rate for each object achieved by the developed algorithms. During the considered test scenario, the successful operation rate for T1 was 85.27% to indicate that the error is less than the threshold in 220 time frames out of the total, 258. The average error across all the time frames (258) was 0.29 m for T1. However, the successful operation rate for T3 was slightly lower (83.59%) than the case of T1, and its prediction accuracy was more or less increased with the decreased average error of 0.28 m. Therefore, it is noted that the proposed algorithms provide satisfactory performance on position prediction in most of the given time frames.

Conclusions
This study proposed sensing algorithms for object tracking and predictive safety evaluation for an autonomous excavator that effectively enables us to avoid potential collision risks between the excavator and the workers around it. Firstly, the proposed method extracted non-ground points from raw point cloud data with 3D LiDAR. The non-ground points were then clustered, and each cluster was fit into a bounding box that was considered as a detected object. The positions and velocities of the detected objects were then estimated and predicted using the integrated tracking filter where the JPDA for multi-target tracking was coupled with UKF and IMM to improve the performance of the tracking management by dealing with the non-linear state estimation and different motion behaviors (i.e., constant velocity and constant turn models), respectively.
The developed algorithm can also monitor the current working areas of the excavator and predict them to conduct safety evaluations. The defined safety indices, TTC and x were computed through the prediction of the motion of the detected objects and calculation of the working areas of the excavator based on kinematic analysis. The severity of potential collision risks was defined using the safety indices.
Validations of the developed algorithms were performed with a modified mini excavator having a VLP-16 3D LiDAR, LVDT stroke sensors, and a rotational encoder. Test results indicate that the proposed predictive safety algorithms can successfully conduct tracking of objects, predict their motion, and provide safety evaluations to prevent a collision accident.
A novel feature of the developed safety strategy is the ability to predict the future states of both multiple objects and an excavator within a prediction horizon (i.e., 1 s in advance) during operations. Since previous studies of construction equipment have focused on monitoring its present state, the predictive attribute of this study, conceptualized through a risk metric with time and distance to potential collisions, is the most significant advance, and it will be used to control autonomous excavation equipment that requires more elaborate and precise safety strategies. As part of the evaluation of the developed safety algorithms, this work provides a unique methodology for defining a threshold for judging success or failure predictions and analyzing the success rates of operations with respect to time frame. Therefore, the study proposes an advanced safety management framework that can be applied to many industrial applications, such as construction equipment in transition to full automation and/or autonomy, in which predictive safety is a prerequisite component.
As the constant velocity and constant turn models were assumed for the prediction of the moving objects in the tracking filter, the application of the variable velocity and random motion models will assist in improving the performance of the safety algorithms in predicting the nonlinear motions of the objects in future work. System dynamics [25,26] can also be considered for an extension of the current work to enhance the predictive capability of the developed algorithms by investigating the interrelations among the factors causing safety incidents, such as worker-machine interaction.
Excavator safety can be improved by pairing the proposed safety algorithms with the development of sophisticated sensing algorithms or platforms that can avoid the restricted FOV and occlusions caused by obstacles in excavation spaces. Finally, a further advantage of state prediction is that the actuators of an autonomous excavator can be controlled in an adaptive manner depending on the degree of potential collision risk, thereby providing maximum productivity and safety over simply stopping the excavator.