Point cloud models are a type of geometric model used to represent discrete points in 3D space. Each point in the model contains 3D coordinates and RGB values (Nirala and Agrawal, 2020). These models can be classified into scattered point clouds, array point clouds, and scanned point clouds. The points of the scanned point cloud all lie on the scan line, which is not fully present. The tracing line is not fully present means that the scanning process, the scanning line (or scanning path) does not completely cover or pass through the entire surface of the workpiece, or there are some missed parts. This may be due to the limitations of the scanning equipment, obstruction during the scanning process, scanning angle or positional problems, and so on. Array point cloud data is uniformly distributed with points at each intersection. Scattered point clouds, on the other hand, exhibit a completely random distribution. Regardless of the type of point cloud, during actual instrument measurements, the data are susceptible to noise and outliers due to the inherent characteristics of the measurement system and its precision. Therefore, in order to avoid data errors on the subsequent path generation, the acquisition and noise reduction processing of the point cloud data is particularly important. First, the study uses a high-performance 3D scanning equipment to scan the data. The process ensures a sufficiently high sampling density to capture the minute details of the workpiece surface. In addition, the sampling spacing was set to a fixed value, which ensured that the point cloud data were adequately captured. By using this scanning equipment, the study obtained accurate and detailed point cloud data models, which provided a reliable basis for subsequent algorithmic studies. Common denoising methods include radius filtering, Gaussian filtering, interactive denoising and the four-point denoising method (He et al., 2021). While removing noise, the four-point noise reduction method is more capable of retaining the detail information in the point cloud data compared to radius filtering and Gaussian filtering, and at the same time realizes a better smoothing effect, which makes the surface of the point cloud smoother without losing the key shape features (Zheng et al., 2022). Most importantly, the four-point noise reduction method has a certain degree of adaptivity, which can better adapt to the different densities and distributions in the point cloud (Leal et al., 2020). The four-point denoising method has been introduced in this study. The four-point denoising method uses a single scan line as a reference and calculates the spatial distances between consecutive points P i , j , P i , j + 1 , P i , j − 1 , P i , j + 2 . The distance calculation between any two points is given by Formula (Chen et al., 2021). d i s t = x 2 − x 1 2 + y 2 − y 1 2 + z 2 − z 1 2 (1)

In Formula 1 (Chen et al., 2021), x 1 , y 1 , z 1 represents the 3D spatial position of the first point, and x 2 , y 2 , z 2 represents the 3D spatial position of the second point. If points P i , j and P i , j + 1 are situated on the same side of the line connecting P i , j − 1 and P i , j + 2 , the four-point judgment process is described in Formula (Li and Tang, 2021). d 1 − d 2 > ε , P i , j d 1 − d 2 < ε , P i , j ∪ P i , j + 1 (2)

In Formula 2 (Li and Tang, 2021), ε is a given threshold, i represents the i -th scan line, and j represents a point on the scan line. The four-point denoising process is illustrated in Figure 1.

From Figure 1, it can be observed that the process involves setting values for i and j , selecting four points on the i -th scan line, and calculating the distances between pairs of points. If the difference in distances is greater than the threshold, point P i , j is identified as noise and removed. If the difference is less than the threshold, the noise point could be P i , j or P i , j + 1 . In this case, a smoothing denoising operation is performed using the arithmetic mean of PL, P i , j , and PI. If points P i , j and P i , j + 1 are on opposite sides of the line connecting P i , j − 1 and P i , j + 2 , then P i , j is identified as noise and is removed, and the arithmetic mean of P i , j − 1 and P i , j + 1 is used for interpolation. This process is repeated until the number of points on the scan line is less than n − 2 , and the number of scan lines is less than m , at which point the results are outputted. After denoising, down-sampling of the point cloud data is necessary to reduce the computational load of the toolpath algorithm. The study employs an octree structure to perform 3D grid processing on the point cloud model. The octree structure is illustrated in Figure 2.

In Figure 2, each face of the cube is subdivided in a fixed format. There is one point for each node, representing a cubic model. This process is carried out for each face of the node cube, and eight initial child nodes are yielded. This process is repeated until the volume of the cubic network or the point cloud data is insufficient to support further subdivision, at which point the operation ceases. In the actual scanned point cloud, a bounding box of length and width is created and divided into KW grids. Grids with data points are considered real grids, while those without data points are virtual grids. The center points of real grids are used to replace other point cloud data in the scan line, and this process is repeated until all scan lines are processed.

The toolpath refers to the trajectory of the tool’s movement relative to the workpiece in CNC machining. This process involves tool approach, machining, milling, and the path traveled back to the tool approach point. The toolpath comprises concepts such as tool position points, tool contact points, step size, pitch, and processing bandwidth (Choudhuri et al., 2023). General trajectory generation methods typically go through steps such as point cloud data modeling, surface data reconstruction, and surface path planning. Due to the time-consuming and labor-intensive nature of the reconstruction step, and considering that rough machining is primarily used to remove excess metal material without strict requirements on toolpaths, the study adopts a layered and iterative cutting approach. The schematic diagram of the process for slicing scattered point clouds is illustrated in Figure 3.

As depicted in Figure 3, the complete point cloud data set is divided into multiple sections through multi-layer slicing operations. By intersecting the point cloud data of each section, an overview of the entire point cloud model can be obtained, and the computational expression of this process is shown in Formula (Huang et al., 2020). E 1 = E ⁡ max E 2 + Δ d i s t = E ⁡ max E 3 + 2 Δ d i s t = E ⁡ max E i + i + 1 Δ d i s t = E ⁡ max E m + m − 1 Δ d i s t = E m a x (3)

In Formula 3 (Huang et al., 2020), Δ d i s t represents the distance between adjacent layers of slices, and E 1 , E 2 , E 3 , ⋯ E i , E m represents the number of layers of slices. As scattered point clouds are discretely distributed, to further optimize the contour of the point cloud model, the study introduces the concept of slicing thickness. The relationship between the slicing plane and the slicing thickness is expressed in Formula (Praveen et al., 2020). δ = k ` ρ (4)

In Formula 4 (Praveen et al., 2020), δ represents the slicing thickness, ρ represents point cloud density, and k ` represents the adjustment coefficient. The slice thickness is a critical parameter whose numerical determination is inextricably linked to computational efficiency, point cloud density and machine performance. First, small slice thickness increases the number of tool trajectories, which may lead to a decrease in computational efficiency. Therefore, a trade-off between computational efficiency and machining accuracy is required. Secondly, a high density point cloud can use a relatively large slice thickness, while a low density point cloud may require a smaller slice thickness. Finally, choosing the right machine tool for machining can improve the efficiency and accuracy of slicing thicknesses. At this point, selecting a point p within a section, the distance calculation formula between p and a nearby point m is expressed in Formula (Li et al., 2021). d = x − x i 2 + y − y i 2 + z − z i 2 (5)

In Formula 5 (Li et al., 2021), x i , y i , z i represents the 3D coordinates of point p , and x , y , z represents the coordinates of point m . The calculation of point cloud density at this point is expressed in Formula (He et al., 2022). ρ = ∑ i N d m , p / N (6)

In Formula 6 (He et al., 2022), N represents randomly selected N points, and d m , p represents the distance value between point m and point p . Additionally, considering that point clouds lack a topological structure, it is necessary to first perform a two-dimensional grid division of the data for each section, followed by equidistant processing of the point cloud network. The two-dimensional (2D) meshing is further divided into mesh generation, coordinate transformation and data interpolation. Firstly, the point cloud data is projected onto a plane by a specific mesh generation algorithm to form a 2D mesh. Secondly, the 3D coordinates of the point cloud data are mapped to the 2D mesh coordinate system by means of coordinate transformation. Finally, for the uncovered areas of the grid, interpolation can be used to fill them, ensuring that there are corresponding data points on the entire 2D grid. After obtaining the extreme values of the section data, a 2D bounding grid perpendicular to the Z-axis is constructed using the coordinate extreme values. The number of grids on the X and Y-axes is expressed in Formula (Liu et al., 2020). x _ n u m = f l o o r x max − x min / k ⋅ ρ + 1 y _ n u m = f l o o r y max − y min / k ⋅ ρ + 1 (7)

In Formula 7 (Liu et al., 2020), f l o o r ( x max − x min represents the number of layers on the X-axis, and f l o o r ( y max − y min represents the number of layers on the Y-axis. The grid position number of any point p is denoted as a , b , and the expression is given in Formula (Liang et al., 2020). a = x i − x min / l b = y i − y min / l (8)

In Formula 8 (Liang et al., 2020), l represents the side length of a square grid. The grid with the extreme values calculated is referred to as the boundary grid. After equidistant processing and sorting of the boundary grid, the tool position grid is obtained. The sorting is based on the center point of each cell, and assuming the grid coordinates for the a -th row and b -th column are S x s , y s , the coordinate calculation formula is expressed in Formula (Gawali et al., 2021). x s = x min + a × l − l 2 y s = y min + b × l − l 2 (9)

The algebraic meaning in Formula 9 (Gawali et al., 2021) is consistent with the previous explanation. Once the sorting process is completed, an equidistant processing can be performed. If S x s , y s belongs to one of the maximum grids in column b , the sequentially corresponding grid index after equidistant processing is marked as a − 1 , b , forming the set Q r . If S x s , y s belongs to one of the minimum grids in column b , the corresponding grid index is marked as a + 1 , b and denoted as Q r . The initial cutting grid then becomes the union of Q r and Q r . In summary, after evaluating and updating all cutting positions and combining them with row cutting, it forms the entire movement path of the rough machining tool positions. The schematic diagram of the tool path for rough machining at this point is depicted in Figure 4.

As depicted in Figure 4, a specific layer is treated as an example involves intersecting the cutting grid. This is achieved by drawing lines parallel to the X-axis and passing through the center of each grid. The outcome of this process is the generation of segmented cutting tool lines after truncation. The tool trajectory connects these lines in a zigzag pattern from left to right and bottom to top, starting from the bottom-left grid. Connecting the right-side tool trajectory using the same method and joining them together provides the complete trajectory for the rough machining tool. In summary, the integration of layer cutting method into the trajectory planning of roughing point cloud model in practical application is roughly divided into four steps. The first step is data preparation, i.e., selecting appropriate workpieces, ensuring that their point cloud data can contain complex geometric information, and collecting the point cloud data by high-quality 3D scanning equipment (e.g., laser scanner), setting appropriate resolution and sampling density. The next is data preprocessing, where the acquired information is denoised by four-point noise reduction to smooth the data. Then the layer cutting method is implemented to describe the layer cutting path generation of the point cloud data, set the thickness parameter, etc. Finally, the simulation test is conducted with specific engineering examples to verify the effect of the layer cutting roughing method.

Continuing the exploration into generating CNC milling machine finish machining tool paths, the research focuses on the field of five-axis machining using point cloud models. The specific path generation involves employing the equal residual height method. A typical schematic diagram of the equal residual height method is illustrated in Figure 5.

From Figure 5, it’s evident that initially, the longest scanning line within the point cloud boundary is selected as the first tool path. Subsequently, under the constraint of the residual height preset value, the tool contact points for adjacent rows are calculated and connected to form a complete tool path. In this case, the residual height represents the distance from the point to the surface of the workpiece and can be obtained by simple geometric calculations. The residual height threshold, on the other hand, is an important parameter that determines which residual heights are considered valid surfaces for machining. While this method can simulate tool path generation, the repetitive steps involved in practical operations increase the computational load of the algorithm. Hence, the study attempts to perform spatial grid processing on point cloud data, creating bounding boxes. The formula for calculating the point cloud model’s bounding box in polar coordinates within the spatial grid is shown in Formula (Liu et al., 2022). x ` _ n u m = f l o o r x max − x min / k ⋅ ρ + 1 y ` _ n u m = f l o o r y max − y min / k ⋅ ρ + 1 z ` _ n u m = f l o o r z max − z min / k ⋅ ρ + 1 (10)

In Formula 10 (Liu et al., 2022), x ` _ n u m , y ` _ n u m , and z ` _ n u m represent the number of grids along the X-axis, Y-axis, and Z-axis in the spatial grid, respectively. At this point, in the bounding box, a column is arbitrarily selected, and the data in that column are denoted as p i L . After selecting all the boundaries of this column and calculating the arch height error, the tool contact points for the first tool path can be derived. The computational formula for this process is shown in Formula (Sheng et al., 2022). x min + i − 1 ⋅ k ⋅ ρ ≤ p ⋅ x ≤ x min + i ⋅ k ⋅ ρ (11)

In Formula 11 (Sheng et al., 2022), the range of values for i is 1 ≤ i ≤ x ` _ n u m . The point with the minimum Y-coordinate in p i L is selected as the boundary point p i B D . Only when i = x ` _ n u m , all the boundary points can form a complete set. The calculation formula for this process is given in Formula (Sales et al., 2021). ε ` = p j B D ⋅ z − p j + n B D ⋅ z x k + p j B D ⋅ x − p j + n B D ⋅ x z k − p j B D ⋅ x ⋅ p j + n B D ⋅ z p j + n B D ⋅ z − p j B D ⋅ z 2 + p j B D ⋅ x − p j + n B D ⋅ x 2 (12)

In Formula 12 (Sales et al., 2021), p j B D represents the j th point in the boundary set, and p j + n B D represents the j + n th point in the boundary set. ε ` represents the general length of the line connecting p j B D and p j + n B D , and the midpoint of this line is p k B D . The definition of tool contact point judgment at this point is given by Formula (Nirala and Agrawal, 2020). p 1 = p j + n − 1 B D , ε ` > e p 2 = p j + n + 1 B D , ε ` < e ⋮ p 1 = p i , j C C ∪ p 2 = p i , j + 1 C C , ε ` = e (13)

In Formula 13 (Nirala and Agrawal, 2020), e represents the bow height error, p i , j C C and p i , j + 1 C C respectively represent the tool contact points of the j th and j + 1 th points in the i th column. The calculation of subsequent tool contact points continues based on the inference in the Formula, and finally, they are connected to form the first tool contact trajectory. This tool contact trajectory can be used to obtain the tool position trajectory for the first tool path. The tool contact point continues to update to p i , j + 1 C C until all tool contact points have undergone coordinate calculations. The tool position points and tool axis vectors during tool processing can be obtained from these coordinates. To summarize, the study preprocesses the point cloud data. Filtering, denoising and other operations are carried out before transforming the coordinate value mapping. Then the point cloud model is processed by the residual height method, the residual height value is calculated, and the finishing toolpath planning is obtained after filtering and threshold setting, as shown in Figure 6.

From Figure 6, it can be observed that the first step is to calculate the global coordinate extremes of the bounding box and select the Jth tool contact point of the Ith tool path. A local coordinate system is established based on this, and the contour point set and residual height value of this tool contact point are calculated. The adjacent row tool contact points are calculated using the residual height method. This process is repeated until all tool contact point calculations are completed. If, at this point, the tool path has covered the entire point cloud model, each path is connected to form the tool contact trajectory. Finally, the conversion of tool position points and tool axis vectors yields specific trajectory data. If the tool path has not yet covered the point cloud model, the next path is selected, and the process is continued. In summary, the integration of the residual height method into the trajectory planning of the finishing point cloud model in practice is also roughly divided into four steps. The first step is data preparation, i.e., selecting appropriate workpieces, ensuring that the point cloud data contains surfaces that need to be finely machined, and using high-precision 3D scanning equipment to obtain detailed surface information. Next is data preprocessing, where the collected point cloud data is denoised to ensure that the calculation of the residual height is not interfered by noise. Then comes the implementation of the residual height method, which calculates the residual height values and computes the path generation. Finally, an engineering example is cited for performance simulation test to evaluate the actual effect of residual height method finishing.

To verify the effectiveness of the proposed method for generating rough machining tool paths using a point cloud model, simulations were conducted on a portion of the saddle surface using MATLAB. The point cloud data dimensions were set to 108 mm × 45 mm × 25 mm, and rough machining was performed using a flat-bottom cutter D6 with a radius of 5 mm. The distance between adjacent layers of slices was set to 0.5 mm. Different slice thicknesses δ 1 = 3 m m , δ 2 = 4 m m , δ 3 = 5 m m , and δ 4 = 6 m m were considered as variables. The actual trajectory generated by the proposed method is shown in Figure 7.

Figure 7 illustrated the toolpaths under different thicknesses. Figure 7A represents the toolpath for a thickness of 3 mm, Figure 7B for 4 mm, Figure 7C for 5 mm, and Figure 7D for 6 mm. As observed from Figure 7, comparing the toolpaths for four different slicing thicknesses revealed that smaller slicing thickness led to an increase in the number of cross-sectional planes during point cloud data slicing, subsequently resulting in longer algorithm computation times. Moreover, if the slicing thickness was too large, it increased the number of tool contact points within each slice after conversion. Therefore, after comprehensive evaluation based on tool contact points and path length, the slicing thickness of 5 mm in Figure 7 yields the optimal trajectory generation, and subsequent tests were conducted using this thickness. To validate the proposed point cloud model layer slicing method against existing machining methods, such as equidistant cutting and unit cutting, the study conducted comparative tests using trajectory error rate and computation time as reference indicators, as depicted in Figure 8.

Figure 8A displays the error change rate curves for three different trajectory generation methods, while Figure 8B shows the computation time curves. It was evident from Figure 8 that the trajectory generation error for all three methods decreased with an increase in the number of tool contact points. Throughout this process, the proposed layer slicing method performed optimally, with the lowest error approaching 10%. Additionally, the layer slicing method exhibited the least computation time, with an average of approximately 35 s. These data indicate that 2D gridization of point cloud data using the layer slicing method can reduce the complexity of computation. In order to compare the performance effect of the three methods more accurately, the study tries to show the detail maps generated by the trajectories of these methods in X-O-Y planar projection as shown in Figure 9.

Figure 9A shows the trajectory projection of the unit cutting method, Figure 9B shows the trajectory projection of the isometric cutting method, and Figure 9C shows the trajectory projection of the layer cutting method. As can be seen from Figure 9, the number of trajectory paths of both unit cutting method and isometric cutting method was 6, and the number of tool contacts was obviously more than layer cutting method. The roughing of the workpiece by the layer cutting method produced only 5 trajectory paths, and the number of tool contacts was even less. The reason for this difference may be attributed to unit cutting dividing the workpiece into multiple units for sequential cutting, thereby generating more trajectory lines. Isometric cutting, on the other hand, mechanically made planning more difficult by calculating the distance between tool contacts. Furthermore, the layer slicing method had significantly fewer tool contact points than the unit cutting method, maintaining the overall machining structure unchanged. This can greatly reduce the data computation time for the point cloud model, ultimately optimizing the cutting process.

To validate the effectiveness of the proposed finish machining toolpath generation method, tests were conducted using MATLAB 2019a as the simulation software. The simulation involved 168,900 scattered points, and the effective cutting radius of the circular milling cutter was set to 4 mm. Residual height was set to 0.6 mm, bow height error to 0.05 mm, forward tilt angle to 10°, and lateral offset angle to 0°. Using a surface point cloud model as an example, the proposed method was tested against existing finish machining cutting methods such as polyhedron, projection, and offset surface methods, and the results are presented in Figure 10.

Figure 10A depicts the toolpath of the point cloud model under the polyhedron method, Figure 10B under the projection method, Figure 10C under the residual height method, and Figure 10D under the offset surface method. From Figure 10, it was observed that, in the segmentation of path data for the surface point cloud model, the polyhedron method exhibited complexity and multiple paths. Although the projection method showed improvement in trajectory generation, there were still instances of path crossings. The offset surface method enhanced the trajectory planning, but unevenness and bias persisted. Overall, the proposed residual height method performed optimally, demonstrating a certain superiority in model path generation. To further explore the machining effects of the residual height method, sampling along the residual height line was conducted, and the obtained residual heights are shown in Figure 11.

From Figure 11, it was evident that, after obtaining the sampling space coordinates with unchanged X and Y-axes, the Z-axis coordinates transformed by the residual height method mainly concentrated in the intervals of 0.2–0.25 and 0.25–0.3, while the remaining residual heights were scattered below 0.2. This indicates that the residual height method achieved an efficiency of over 60% in toolpath generation for finish machining workpieces, simultaneously avoiding the complex processes of surface fitting and triangulation, streamlining the steps of path generation. Finally, the study, in order to verify the practical application effect of the new roughing and finishing trajectory generation method, compared and verified the traditional finishing and improved finishing, and the traditional roughing and improved roughing, respectively, with H7 helical gears, HK1712 needle roller bearings, and M6 double-ended zinc-plated screws as the test parts, and the physical drawings are shown in Figure 12.

Figure 12A shows the physical drawing of H7 helical gear, Figure 12B shows the physical drawing of HK1712 needle roller bearing, and Figure 12C shows the physical drawing of M6 double-ended galvanized screw. In addition, the study set the roughing tool diameter to φ12 and finishing to φ6. The lathe speed settings were set to 3,500 m/s and 4,500 m/s, and the feeds were set to 1,000 and 800 mm/r. The cutting pitch was set to 1.2 mm and 4 mm, respectively. Furthermore, more machining methods of the same type were introduced, such as the grid of the same type as the layer cut roughing method cutting, contour cutting, and lithographic cutting, and adaptive cutting, topology optimization, and dynamic cutting, which are the same type as the residual height method for finishing. Trajectory generation error rate, computing time, trajectory smoothness and cutting efficiency were tested as reference indexes, and the test results are shown in Table 1.

As can be seen from Table 1, all 16 processing methods showed superior processing results in the multi-indicator test. Firstly, for the rough machining method comparison, the layer cutting method proposed in the study showed the best overall performance. This was followed by the photolithography cutting method, the contour cutting method and the grid cutting method. Comprehensively, the layer cutting roughing method had the lowest trajectory generation error rate of 10.17% in M6 double-ended galvanized screw, the lowest running time of 11.82s in H7 helical gear machining, the highest value of trajectory smoothness of 83.49% in M6 double-ended galvanized screw machining, and the highest cutting efficiency of 90.17% in M6 double-ended galvanized screw machining. On the contrary, the residual high finish machining method showed the lowest trajectory generation error rate of 7.68% in H7 helical gear machining, the lowest running time of 13.58 s in H7 helical gear machining, the highest value of trajectory smoothing of 84.17% in M6 double-head galvanized screw machining, and the highest cutting efficiency of 91.57% in M6 double-head galvanized screw machining. In summary, the proposed method of generating finishing tool trajectories using the residual height method to compute the point cloud model showed better performance and certain feasibility compared to the same type of arithmetic methods.