# Modeling and Prediction Method for CNC Machine Tools’ Errors Based on Spatial Feature Points

Authors: Guohua Chen; Lin Zhang; Hua Xiang; Yong Chen
Journal: Advances in Materials Science and Engineering
 (2020)

---

## Abstract

In order to improve the precision of CNC machine tools effectively, a method for modeling and predicting their spatial errors based on spatial feature points was proposed. Taking three-axis vertical CNC machine tools as the research object, we think that the whole space formed by machine tools’ working can be seen as the combination of a number of cubes, whose vertices are considered to be feature points, and others in the cubes are called nonfeature points. So, each nonfeature point’s errors can be predicted by the cube’s eight vertices’ errors. Based on the above ideas, an approach including the installing instrument for measuring any spatial feature point’s errors was put forward. In this way, all data of the feature points’ errors could be obtained. Moreover, according to these error data, the prediction model of nonfeature points’ errors was established by using the internal division ratio method. The method has the advantages of small interpolation operation, easy integration in the numerical control system, and high compensation precision. Finally, an example was used to prove its effectiveness and feasibility.


---

## Body

## 1. Introduction

With the rapid development of modern industry, CNC manufacturing technology is developing rapidly towards the direction of high precision and ultraprecision. As the main processing tool of precision manufacturing, CNC machine tools have a direct impact on the development of the whole manufacturing industry with their precision index [1]. Geometric errors and thermal errors are the two major errors of CNC machine tools [2, 3]. The geometrical errors of CNC machine tools are a large part of the overall errors of machine tools. At present, the research on the measurement and modeling of spatial geometric errors of machine tools has been going on for a long time. In particular, the measurement and modeling of spatial geometric errors of multiaxis machine tools has become a hot and difficult issue in recent years. In 2004, He et al. [4] used a laser interferometer together with a rotary encoder to measure roundness error and locate precision. This method can effectively solve the problem that the measuring range of the ballbar is limited by the length of the bar. In 2005, Schwenke et al. [5] put forward the error measurement method for machine tools and CMM by using the single-light laser tracking interferometer. In 2009, Schwenke et al. [6] studied the measurement method of continuous tracking of the laser tracker in the air and identified six error elements of the rotation axis and then analyzed the uncertainty of measurement. Li et al. [7] proposed a 13-line method to identify the geometric errors of machine tools. In 2012, Zhu et al. [8] established a spatial error model of TTTRR five-axis machine tools based on multibody system theory. In 2016, He et al. [9] used a laser interferometer to measure the errors of X, Y, and Z axes, as well as the diagonal of four individuals. Based on the nonuniform rational B-spline, an adaptive mathematical expression model was established to describe the spatial errors of machine tools. El Bechir et al. [10] suggested a simulation methodology for errors caused by the interpolations B-spline and C-spline in high-speed machining of warped shapes and developed analytical models expressing the basic paths of the interpolations. In 2017, Chen et al. [11] attempted to establish an error model of 12 terms of polynomial in total for the triaxial linkage of machine tools based on the error measurement data of the whole working space of machine tools and preliminarily discussed the error prediction method for the points at different spatial locations of machine tools. In 2018, Cao and Sun [12] offered the principle of the photodynamic laser Doppler measurement system by analyzing the geometric errors of triaxial machine tools, introduced the method of the space error test path, and compensated with the space error compensation function of i5 system. Based on the research in 2016, El Bechir et al. [13] developed a method of compensating for these errors based on the insertion of the nodes, while respecting the predefined tolerance. To do this, they modeled and simulated machining errors before and after compensation for each type of interpolation. In 2019, Zuo and Li [14] proposed a method to change the trajectory of the ballbar in its working space and carried out spatial error modeling based on this method, which could accurately solve the rotation error elements of the cutter coordinate system and the workbench coordinate system.Nowadays, taking average three-axis machine tools for an example in the research of CNC machine tool geometric error compensation [15–17], most of the researchers use the laser interferometer or other measuring instruments to measure geometric errors including perpendicularity errors at several points on three axes separately [18, 19], without considering the coupling effects between axes [20–22]. Then, according to the measured geometric error data of these points, they used all kinds of mathematical methods to calculate the machine tool error of the whole working space points. Obviously, the spatial error data obtained by this method are not accurate enough to represent the spatial errors of machine tools truly. In this paper, starting from the consideration of the three-axis combined superposition effect of machine tools, the error measurement model of a machine tools’ coordinate position and spatial feature points is established. This measurement method can realize the measurement of laser beam’s arbitrary movement in the whole working space of machine tools, so as to measure the position errors of each point on the predetermined spatial lattice. These data are used to set up a prediction model with the minimum interpolation algorithm and predict the spatial nonfeature points’ errors.


## 2. Modeling and Prediction Methods

### 2.1. Modeling Idea Based on Spatial Feature Points

The whole working space of three-axis vertical CNC machine tools can be seen as the combination of a number of cubes. All the vertices are considered to be feature points whose errors will be measured, and others are called nonfeature points. In fact, these measurement positions in the whole working space of machine tools are precisely the feature points. The physical model of the spatial structure and grids formed by these feature points are shown in Figure1.Figure 1
Hypothetical physical model.All measurement points can be viewed from a vertical view, and the measurement lattices formed in the space can be seen from Figure2.Figure 2
A lattice formed by feature points in the machine tools’ space.These lattices form the feature points of machine tools in the whole working space, which represent the basic attribute of machine tools’ spatial errors. These feature points have fixed spatial position coordinates in the machine tools’ spatial coordinate system, which can be described with a three-dimensional matrix (Figure3). The coordinate information of these feature points can be conveniently expressed by (x, y, z) and stored in the numerical control system. Among these feature points, eight adjacent points form a cube. Each vertex, namely, feature point, obtains the corresponding error vectors through a certain measuring method. These vectors contain six errors, but only three position errors’ data are actually needed for the research.Figure 3
Coordinate schematic of spatial feature points.In addition to the measurement feature points, the rest of machine tools’ spatial points are the nonfeature points. Any nonfeature point is included in the corresponding cube. Every spatial error of any nonfeature point can be seen as a function of the eight vertices in its cube. According to the data of feature points’ spatial errors, a prediction model of any nonfeature point can be established, which makes it possible to predict any point’s error in the whole working space. The prediction model established by this method has the advantages of higher reliability and precision than those based on the data of uniaxial errors.


### 2.2. Measurement Method

#### 2.2.1. Installation Method

According to the principle of laser measurement, through analyzing the laser interferometer installation method, machine tools’ spatial errors can be measured directly under the help of two optical steering mirrors with 90 degree angle and other measurement accessories.Meanwhile, according to the motion features of machine tools’ guideway and workbench, a method of measuring any point in all the working space was put forward, which can measure directly the deviation of any point’s errors in the 3 directions ofX-, Y-, and Z-axis in the machine tools’ moving space. It can be shown in Figure 4.Figure 4
Installation method of the laser interferometer on machine tools.


#### 2.2.2. Measuring Process

Step 1.
(planning measurement paths). According to theX-, Y-, Z-axis’ motion routes of machine tools, the requirements of the number of spatial points’ errors to be measured and the compensation accuracy and the distance between measurement points can be set. In addition, the number of measurement points can also be set according to the measurement interval. Moreover, the measurement paths can be planned according to the measurement sequence from points, lines, and planes to the space. As shown in Figure 5, the whole working space of machine tools was divided into 10 planes. In fact, it can be divided into more planes according to the needs of compensation precision.Figure 5
Planning method for measuring paths.Step 2.
(measuring all points’ errors on the first surface). Firstly, to measure three errors at each point ofX-axis, Z-axis is required to move to the vertical direction by one measuring distance of ΔL. Secondly, to get the three errors of each point at Z = ΔL, the head of the laser interferometer is required to move from the X-axis’ zero point towards its positive direction. Thirdly, to measure the 3 errors at each point at Z = 2ΔL by using the laser interferometer head, move Z-axis to another measuring distance of ΔL towards the positive direction. At last, what should be done is to repeat the above processes to complete measuring all points’ errors in the plane X0O0Z0.Step 3.
(measuring all points’ errors on other surfaces). In the same way, movingY-axis one measuring interval ΔL can get a new plane X1O1Z1. Meanwhile, we make the Z-axis move to the zero coordinate point. Next, the same measurement process can be used to complete measuring all positions’ errors in the X1O1Z1 plane. Secondly, moving Y-axis one measurement distance of ΔL gets a new plane X2O2Z2, and the spindle runs again to the zero coordinate point again. Next, we use the same steps to complete measuring all positions’ errors in the plane X2O2Z2. The optical path of the measurement process is shown in Figure 6.
Repeating the above processes can complete measuring all points’ error in all the planes in the whole working space of machine tools. In this way, the deviation data of feature points’ error at the three directions ofX, Y, and Z in the whole working space of machine tools are obtained.Figure 6
Optical path of measuring process.


### 2.3. Modeling and Prediction Method for Nonfeature Points’ Errors

All the measuring lines in the lattice divide the whole working space into many small cubes. All errors in the cubes are local, and the farther away from the center the point is, the less affected by the center it is. Therefore, it is necessary to divide the whole cube into small parts and to perform interpolation in the local part. Interpolation inside the small cubes maintains the continuity of the boundary. Therefore, it becomes possible to implement interpolation operations and establish a compensation prediction model in the whole working space.For a small cube, as shown in Figure7, there are eight vertices, represented by P1 (P1x, P1y, P1z), P2 (P2x, P2y, P2z), …, P8 (P8x, P8y, P8z), respectively, in which (Pnx, Pny, Pnz) (n = 1, 2, … …, 8) represents the three coordinates of the nth vertex. The coordinate of any internal point is expressed as P (Px, Py, Pz), and its spatial errors are expressed as (Cx, Cy, Cz). In Figure 7, (Cnx, Cny, Cnz) represents the three errors of the cube’s nth vertex.Figure 7
Schematic of calculating the internal ratios.The most direct interpolation algorithm is linear interpolation. Based on linear interpolation, the values ofCx, Cy, and Cz show the 3 errors of each vertex, which can be determined by all distances from point P to the 8 vertices. After calculating all these distances, their weights can be determined. The internal ratios of the X-axis, Y-axis, and Z-axis are defined as x, y, and z, respectively. As can be seen in Figure 7, they can be expressed as(1)x=px−p1xp2x−p1x,y=py−p1yp4y−p1y,z=pz−p1zp5z−p1z..Next, the pointP’s compensation prediction values Cx, Cy, and Cz can be obtained, which are expressed as(2)Cx=C1x,C2x,C3x,C4x,C5x,C6x,C7x,C8x⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z,Cy=C1y,C2y,C3y,C4y,C5y,C6y,C7y,C8y⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z,Cz=C1z,C2z,C3z,C4z,C5z,C6z,C7z,C8z⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z.After getting these values,Cx, Cy, and Cz, all points’ errors inside the small cube, including those points on the cube’s plane or line, have certain compensation predicted values. On each plane, the errors’ interpolation depends only on the four vertices on the plane, which are in coplanar with the neighboring cube. So, the interpolation is continuous on the neighboring cube’s plane. On each edge line, its points’ error interpolation only depends on the two endpoints, which are collinear with the four surrounding cubes. Therefore, for adjacent cubes, they are continuous on the edge line. This is especially true for vertices. It can be concluded that this interpolation method is continuous in the whole working space of machine tools.


## 2.1. Modeling Idea Based on Spatial Feature Points

The whole working space of three-axis vertical CNC machine tools can be seen as the combination of a number of cubes. All the vertices are considered to be feature points whose errors will be measured, and others are called nonfeature points. In fact, these measurement positions in the whole working space of machine tools are precisely the feature points. The physical model of the spatial structure and grids formed by these feature points are shown in Figure1.Figure 1
Hypothetical physical model.All measurement points can be viewed from a vertical view, and the measurement lattices formed in the space can be seen from Figure2.Figure 2
A lattice formed by feature points in the machine tools’ space.These lattices form the feature points of machine tools in the whole working space, which represent the basic attribute of machine tools’ spatial errors. These feature points have fixed spatial position coordinates in the machine tools’ spatial coordinate system, which can be described with a three-dimensional matrix (Figure3). The coordinate information of these feature points can be conveniently expressed by (x, y, z) and stored in the numerical control system. Among these feature points, eight adjacent points form a cube. Each vertex, namely, feature point, obtains the corresponding error vectors through a certain measuring method. These vectors contain six errors, but only three position errors’ data are actually needed for the research.Figure 3
Coordinate schematic of spatial feature points.In addition to the measurement feature points, the rest of machine tools’ spatial points are the nonfeature points. Any nonfeature point is included in the corresponding cube. Every spatial error of any nonfeature point can be seen as a function of the eight vertices in its cube. According to the data of feature points’ spatial errors, a prediction model of any nonfeature point can be established, which makes it possible to predict any point’s error in the whole working space. The prediction model established by this method has the advantages of higher reliability and precision than those based on the data of uniaxial errors.


## 2.2. Measurement Method

### 2.2.1. Installation Method

According to the principle of laser measurement, through analyzing the laser interferometer installation method, machine tools’ spatial errors can be measured directly under the help of two optical steering mirrors with 90 degree angle and other measurement accessories.Meanwhile, according to the motion features of machine tools’ guideway and workbench, a method of measuring any point in all the working space was put forward, which can measure directly the deviation of any point’s errors in the 3 directions ofX-, Y-, and Z-axis in the machine tools’ moving space. It can be shown in Figure 4.Figure 4
Installation method of the laser interferometer on machine tools.


### 2.2.2. Measuring Process

Step 1.
(planning measurement paths). According to theX-, Y-, Z-axis’ motion routes of machine tools, the requirements of the number of spatial points’ errors to be measured and the compensation accuracy and the distance between measurement points can be set. In addition, the number of measurement points can also be set according to the measurement interval. Moreover, the measurement paths can be planned according to the measurement sequence from points, lines, and planes to the space. As shown in Figure 5, the whole working space of machine tools was divided into 10 planes. In fact, it can be divided into more planes according to the needs of compensation precision.Figure 5
Planning method for measuring paths.Step 2.
(measuring all points’ errors on the first surface). Firstly, to measure three errors at each point ofX-axis, Z-axis is required to move to the vertical direction by one measuring distance of ΔL. Secondly, to get the three errors of each point at Z = ΔL, the head of the laser interferometer is required to move from the X-axis’ zero point towards its positive direction. Thirdly, to measure the 3 errors at each point at Z = 2ΔL by using the laser interferometer head, move Z-axis to another measuring distance of ΔL towards the positive direction. At last, what should be done is to repeat the above processes to complete measuring all points’ errors in the plane X0O0Z0.Step 3.
(measuring all points’ errors on other surfaces). In the same way, movingY-axis one measuring interval ΔL can get a new plane X1O1Z1. Meanwhile, we make the Z-axis move to the zero coordinate point. Next, the same measurement process can be used to complete measuring all positions’ errors in the X1O1Z1 plane. Secondly, moving Y-axis one measurement distance of ΔL gets a new plane X2O2Z2, and the spindle runs again to the zero coordinate point again. Next, we use the same steps to complete measuring all positions’ errors in the plane X2O2Z2. The optical path of the measurement process is shown in Figure 6.
Repeating the above processes can complete measuring all points’ error in all the planes in the whole working space of machine tools. In this way, the deviation data of feature points’ error at the three directions ofX, Y, and Z in the whole working space of machine tools are obtained.Figure 6
Optical path of measuring process.


## 2.2.1. Installation Method

According to the principle of laser measurement, through analyzing the laser interferometer installation method, machine tools’ spatial errors can be measured directly under the help of two optical steering mirrors with 90 degree angle and other measurement accessories.Meanwhile, according to the motion features of machine tools’ guideway and workbench, a method of measuring any point in all the working space was put forward, which can measure directly the deviation of any point’s errors in the 3 directions ofX-, Y-, and Z-axis in the machine tools’ moving space. It can be shown in Figure 4.Figure 4
Installation method of the laser interferometer on machine tools.


## 2.2.2. Measuring Process

Step 1.
(planning measurement paths). According to theX-, Y-, Z-axis’ motion routes of machine tools, the requirements of the number of spatial points’ errors to be measured and the compensation accuracy and the distance between measurement points can be set. In addition, the number of measurement points can also be set according to the measurement interval. Moreover, the measurement paths can be planned according to the measurement sequence from points, lines, and planes to the space. As shown in Figure 5, the whole working space of machine tools was divided into 10 planes. In fact, it can be divided into more planes according to the needs of compensation precision.Figure 5
Planning method for measuring paths.Step 2.
(measuring all points’ errors on the first surface). Firstly, to measure three errors at each point ofX-axis, Z-axis is required to move to the vertical direction by one measuring distance of ΔL. Secondly, to get the three errors of each point at Z = ΔL, the head of the laser interferometer is required to move from the X-axis’ zero point towards its positive direction. Thirdly, to measure the 3 errors at each point at Z = 2ΔL by using the laser interferometer head, move Z-axis to another measuring distance of ΔL towards the positive direction. At last, what should be done is to repeat the above processes to complete measuring all points’ errors in the plane X0O0Z0.Step 3.
(measuring all points’ errors on other surfaces). In the same way, movingY-axis one measuring interval ΔL can get a new plane X1O1Z1. Meanwhile, we make the Z-axis move to the zero coordinate point. Next, the same measurement process can be used to complete measuring all positions’ errors in the X1O1Z1 plane. Secondly, moving Y-axis one measurement distance of ΔL gets a new plane X2O2Z2, and the spindle runs again to the zero coordinate point again. Next, we use the same steps to complete measuring all positions’ errors in the plane X2O2Z2. The optical path of the measurement process is shown in Figure 6.
Repeating the above processes can complete measuring all points’ error in all the planes in the whole working space of machine tools. In this way, the deviation data of feature points’ error at the three directions ofX, Y, and Z in the whole working space of machine tools are obtained.Figure 6
Optical path of measuring process.


## 2.3. Modeling and Prediction Method for Nonfeature Points’ Errors

All the measuring lines in the lattice divide the whole working space into many small cubes. All errors in the cubes are local, and the farther away from the center the point is, the less affected by the center it is. Therefore, it is necessary to divide the whole cube into small parts and to perform interpolation in the local part. Interpolation inside the small cubes maintains the continuity of the boundary. Therefore, it becomes possible to implement interpolation operations and establish a compensation prediction model in the whole working space.For a small cube, as shown in Figure7, there are eight vertices, represented by P1 (P1x, P1y, P1z), P2 (P2x, P2y, P2z), …, P8 (P8x, P8y, P8z), respectively, in which (Pnx, Pny, Pnz) (n = 1, 2, … …, 8) represents the three coordinates of the nth vertex. The coordinate of any internal point is expressed as P (Px, Py, Pz), and its spatial errors are expressed as (Cx, Cy, Cz). In Figure 7, (Cnx, Cny, Cnz) represents the three errors of the cube’s nth vertex.Figure 7
Schematic of calculating the internal ratios.The most direct interpolation algorithm is linear interpolation. Based on linear interpolation, the values ofCx, Cy, and Cz show the 3 errors of each vertex, which can be determined by all distances from point P to the 8 vertices. After calculating all these distances, their weights can be determined. The internal ratios of the X-axis, Y-axis, and Z-axis are defined as x, y, and z, respectively. As can be seen in Figure 7, they can be expressed as(1)x=px−p1xp2x−p1x,y=py−p1yp4y−p1y,z=pz−p1zp5z−p1z..Next, the pointP’s compensation prediction values Cx, Cy, and Cz can be obtained, which are expressed as(2)Cx=C1x,C2x,C3x,C4x,C5x,C6x,C7x,C8x⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z,Cy=C1y,C2y,C3y,C4y,C5y,C6y,C7y,C8y⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z,Cz=C1z,C2z,C3z,C4z,C5z,C6z,C7z,C8z⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z.After getting these values,Cx, Cy, and Cz, all points’ errors inside the small cube, including those points on the cube’s plane or line, have certain compensation predicted values. On each plane, the errors’ interpolation depends only on the four vertices on the plane, which are in coplanar with the neighboring cube. So, the interpolation is continuous on the neighboring cube’s plane. On each edge line, its points’ error interpolation only depends on the two endpoints, which are collinear with the four surrounding cubes. Therefore, for adjacent cubes, they are continuous on the edge line. This is especially true for vertices. It can be concluded that this interpolation method is continuous in the whole working space of machine tools.


## 3. Application Verification

### 3.1. Basic Information

The research object of this paper is the XHK715 triaxial vertical machining center produced by Hubei Jiangshan Huake Digital Equipment Technology Co., Ltd. The CNC machining center adopts the HNC-818b/M bus central NC system, and the temperature sensors are installed on the bearing seat, nut seat, and spindle bearing seat of each feed shaft to detect the temperature change of the machine tools. Vibration sensors are also installed in the machine bed and spindle box to detect vibration frequency; the feed shaft of machine tools is equipped with the grating ruler, which can realize full closed-loop control. The feature dimensions of the overall structure of the machining center are as follows: the stokes ofX-axis, Y-axis, and Z-axis are 800 mm, 500 mm, and 550 mm, respectively, and the size of the worktable (length × width) is 500 mm × 1,050 mm.


### 3.2. Data Acquisition of Feature Points’ Errors

Considering the rationality of the experiment, in case of noncollision, the triaxial strokes are selected as follows:X-axis (−600, −250), Y-axis (−290, 60), and Z-axis (−120, −470), and their strokes are all 350 mm. Under the condition of meeting the relevant national requirements, the measurement interval is selected as 70 mm, and the X-axis, Y-axis, and Z-axis are all divided into 5 segments. Thus, there are 4 feature points inserted in the middle strokes. A large number of previous experiments proved that reading tended to be more stable after stopping the movement for 9 seconds. So, the measurement interval is selected as 9 seconds. According to the measurement method of lateral translation and the single-axis mode, each line with a point is measured in the order of x-y-z, y-z-x, and z-x-y. The location of the 18 measurement points selected in the space is shown in Figure 8.Figure 8
The location of the 18 measurement points selected in the space.The 6D sensor of the laser interferometer communicates with the laser head, which transmits data to the computer through Internet. Supporting software provided by API is used to collect data. After these measurement parameters are set, the interferometer will automatically sample, and the data will be automatically saved. The measurement partial data in our experiments are shown in Table1.Table 1
The partial data of feature points’ errors.
Measuring feature points
Error components (μm)
C
nx
C
ny
C
nz
P1 (−565, −255, −155)
−1.000
8.390
4.590
P2 (−425, −255, −155)
−6.040
−10.660
10.315
P3 (−285, −255, −155)
−35.900
−28.760
4.840
P4 (−565, −115, −155)
2.675
5.925
6.020
P5 (−425, −115, −155)
−2.900
−1.609
16.965
P6 (−285, −115, −155)
−26.530
−32.420
13.445
P7 (−565, 25, −155)
0.505
16.190
−3.750
P8 (−425, 25, −155)
−6.530
−2.440
−16.520
P9 (−285, 25, −155)
−32.350
−23.240
−41.650
P10 (−565, −255, −295)
−3.330
1.135
9.475
P11 (−425, −255, −295)
−18.710
−20.420
9.630
P12 (−285, −255, −295)
−42.530
−36.780
−5.190
P13 (−565, −115, −295)
0.235
−1.920
15.150
P14 (−425, −115, −295)
−11.630
−23.910
15.300
P15 (−285, −115, −295)
−36.200
−42.420
3.000
P16 (−565, 25, −295)
−3.750
7.910
−5.700
P17 (−425, 25, −295)
−16.520
−13.730
−4.380
P18 (−285, 25, −295)
−41.650
−34.610
−11.310The data of these 18 measurement points’ errors can also be shown in Figure9.Figure 9
The 18 measurement points’ errors.


### 3.3. Establishment of the Prediction Model for Nonfeature Points’ Errors

The nonfeature points’ errors in cubes are predicted by the method of internal ratio. Suppose that the errors’ data of all cubes’ eight vertices are all collected. In fact, every vertex has 6 errors, which can be expressed as (Cnx, Cny, Cnz, Cna, Cnb, Cnc); here, n = 1, 2, …, 8. This paper only needs to study the first three errors.In this paper, the length of the smallest cube in the measured lattice is 70 mm. In this section, for the convenience of calculation, it is assumed that the length is 1. Set up a rectangular coordinate system, and suppose P1 is located at the origin, and the lines from P1 to P2, P1 to P4, and P1 to P5 are in the positive direction of theX-axis, Y-axis, and Z-axis, respectively. In this way, their internal ratios can be completely calculated according to equation (1).Thus, theP point’s internal fraction is equal to the P point’s coordinate. The starting point and the ending point of interpolation are x0 = 0, x1 = 1, y0 = 0, y1 = 0, z0 = 0, z1 = 1, respectively. When the P point is at the center of the cube, the P point’s error values Cx, Cy, and Cz of the X-axis, Y-axis, and Z-axis can be obtained and expressed as(3)Cx=C1x,C2x,C3x,C4x,C5x,C6x,C7x,C8x⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z=C1x+C2x+C3x+C4x+C5x+C6x+C7x+C8x8,Cy=C1y,C2y,C3y,C4y,C5y,C6y,C7y,C8y⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z=C1y+C2y+C3y+C4y+C5y+C6y+C7y+C8y8,Cz=C1z,C2z,C3z,C4z,C5z,C6z,C7z,C8z⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z=C1z+C2z+C3z+C4z+C5z+C6z+C7z+C8z8.For the convenience of research, we chose a 3 × 3 × 3 lattice with a total of 27 points, which are, respectively, located in the center of the 27 smallest cubes. These 27 points are both nonfeature points and verification points, which are used to verify the effectiveness of the prediction method. The starting and ending coordinates of the lattice are (−600, −290, −120) and (−250, 60, −470), respectively, and the measurement step is 70 mm. In fact, the starting and ending coordinates of the measurement lattice are (−565, −255, −155) and (−285, 25, −435), respectively.Based on the above model and algorithm, we calculated the prediction values of these 27 points’ errors. The calculated results are shown in Table2.Table 2
The verification points’ prediction errors.
Predictive value
Y
1 (μm)
Y
2 (μm)
Y
3 (μm)
C
x
C
y
C
z
C
x
C
y
C
z
C
x
C
y
C
z
Z
1
X
1
2.140
0.210
0.290
−5.240
17.440
−23.230
−10.870
24.350
−28.510
X
2
−0.520
7.690
−9.920
−6.810
18.290
−26.170
−13.250
27.920
−31.860
X
3
−1.970
10.610
−20.340
−9.040
19.170
−28.420
−14.920
40.030
−35.750
Z
2
X
1
0.820
10.070
4.130
−3.120
19.130
−23.090
−11.990
28.660
−31.170
X
2
−1.160
11.540
−4.890
−5.130
19.870
−24.960
−12.810
40.190
−32.490
X
3
−1.980
16.120
−13.480
−8.730
24.380
−27.140
−14.950
42.810
−34.180
Z
3
X
1
3.170
15.120
−5.210
−3.620
20.240
−24.150
−9.770
38.040
−35.230
X
2
1.320
19.460
−13.780
−5.140
22.610
−27.210
−13.270
41.040
−36.890
X
3
−2.210
19.380
−22.940
−6.870
30.190
−31.360
−15.280
49.140
−37.540


### 3.4. Verification of the Prediction Model

To compare the difference between the prediction value and the actual value of the 27 points’ errors, we readjusted and installed the laser interferometer based on the method of measuring the whole working space points’ errors proposed in this paper. We measured directly the 27 verification points’ actual errors. These symbols ofDx, Dy, and Dz represent the 27 verification points’ errors of the X-axis, Y-axis, and Z-axis, respectively. The result can be seen from Table 3.Table 3
The verification points’ actual errors.
Measured value
Y
1 (μm)
Y
2 (μm)
Y
3 (μm)
D
x
D
y
D
z
D
x
D
y
D
z
D
x
D
y
D
z
Z
1
X
1
1.370
0.530
0.870
−4.610
16.980
−23.910
−11.420
23.930
−28.320
X
2
−1.270
8.150
−9.630
−6.030
18.540
−26.580
−12.490
28.260
−31.410
X
3
−2.850
10.280
−20.670
−9.480
18.950
−29.080
−14.530
40.410
−36.080
Z
2
X
1
1.570
10.210
3.690
−3.510
19.440
−22.580
−12.330
29.090
−30.820
X
2
−0.330
11.120
−5.150
−5.470
19.640
−25.210
−13.590
40.530
−32.730
X
3
−2.120
15.590
−13.170
−8.280
24.720
−27.580
−15.260
42.490
−34.710
Z
3
X
1
2.420
15.660
−4.890
−2.940
20.690
−24.730
−10.460
37.880
−34.880
X
2
1.140
19.870
−13.210
−4.530
22.850
−27.830
−12.910
41.370
−36.220
X
3
−1.860
19.150
−22.580
−7.450
29.770
−31.570
−14.420
49.660
−37.330From Tables2 and 3, we can get the comparison between the predicted value and the measured value of the 27 points’ errors. Its prediction accuracy reaches within ±1 μm, as it can be seen in Figure 10.Figure 10
The comparison between the predicted value and the measured value. (a) The comparison ofCx and Dx of X-axis error. (b) The comparison of Cy and Dy of Y-axis error. (c) The comparison of Cz and Dz of Z-axis error.
(a)
(b)
(c)


### 3.5. Evaluation of Verification Effect

In Tables2 and 3, we obtained the measurement value and the prediction value of the 27 verification points’ errors. The symbols of Δx, Δy, and Δz are defined as the difference between the prediction value and the measurement actual value. The formula is expressed as(4)Δx=Dx−Cx,Δy=Dy−Cy,Δz=Dz−Cz.The values of Δx, Δy, and Δz reflect the accuracy of the prediction model. The residual value is shown in Table 4. The comparison diagram between the prediction value and the measurement actual value of the 27 verification points’ errors is shown in Figure 11. The empirical analysis shows that the prediction accuracy of this method is only within 0.001 mm.Table 4
The residual value of the 27 verification points’ errors.
Residual value
Y
1 (μm)
Y
2 (μm)
Y
3 (μm)
∆x
∆y
∆z
∆x
∆y
∆z
∆x
∆y
∆z
Z
1
X
1
0.770
0.320
0.580
−0.630
−0.460
−0.680
0.550
−0.420
0.190
X
2
0.750
0.460
0.290
−0.780
0.250
−0.410
−0.760
0.340
0.450
X
3
0.880
−0.330
−0.330
0.440
−0.220
−0.660
−0.390
0.380
−0.330
Z
2
X
1
−0.750
0.140
−0.440
0.670
0.310
0.510
0.340
0.430
0.350
X
2
−0.830
−0.320
−0.260
0.340
−0.230
−0.250
0.780
0.340
−0.240
X
3
0.140
−0.530
0.310
−0.450
0.340
0.440
0.310
−0.320
−0.530
Z
3
X
1
0.750
0.540
0.320
−0.680
0.450
−0.580
0.690
−0.160
0.350
X
2
0.180
0.410
0.570
−0.610
0.240
−0.620
−0.360
0.230
0.670
X
3
−0.350
−0.230
0.360
0.580
−0.420
−0.210
−0.860
0.520
0.210Figure 11
The residual value of the verification points’ errors.


## 3.1. Basic Information

The research object of this paper is the XHK715 triaxial vertical machining center produced by Hubei Jiangshan Huake Digital Equipment Technology Co., Ltd. The CNC machining center adopts the HNC-818b/M bus central NC system, and the temperature sensors are installed on the bearing seat, nut seat, and spindle bearing seat of each feed shaft to detect the temperature change of the machine tools. Vibration sensors are also installed in the machine bed and spindle box to detect vibration frequency; the feed shaft of machine tools is equipped with the grating ruler, which can realize full closed-loop control. The feature dimensions of the overall structure of the machining center are as follows: the stokes ofX-axis, Y-axis, and Z-axis are 800 mm, 500 mm, and 550 mm, respectively, and the size of the worktable (length × width) is 500 mm × 1,050 mm.


## 3.2. Data Acquisition of Feature Points’ Errors

Considering the rationality of the experiment, in case of noncollision, the triaxial strokes are selected as follows:X-axis (−600, −250), Y-axis (−290, 60), and Z-axis (−120, −470), and their strokes are all 350 mm. Under the condition of meeting the relevant national requirements, the measurement interval is selected as 70 mm, and the X-axis, Y-axis, and Z-axis are all divided into 5 segments. Thus, there are 4 feature points inserted in the middle strokes. A large number of previous experiments proved that reading tended to be more stable after stopping the movement for 9 seconds. So, the measurement interval is selected as 9 seconds. According to the measurement method of lateral translation and the single-axis mode, each line with a point is measured in the order of x-y-z, y-z-x, and z-x-y. The location of the 18 measurement points selected in the space is shown in Figure 8.Figure 8
The location of the 18 measurement points selected in the space.The 6D sensor of the laser interferometer communicates with the laser head, which transmits data to the computer through Internet. Supporting software provided by API is used to collect data. After these measurement parameters are set, the interferometer will automatically sample, and the data will be automatically saved. The measurement partial data in our experiments are shown in Table1.Table 1
The partial data of feature points’ errors.
Measuring feature points
Error components (μm)
C
nx
C
ny
C
nz
P1 (−565, −255, −155)
−1.000
8.390
4.590
P2 (−425, −255, −155)
−6.040
−10.660
10.315
P3 (−285, −255, −155)
−35.900
−28.760
4.840
P4 (−565, −115, −155)
2.675
5.925
6.020
P5 (−425, −115, −155)
−2.900
−1.609
16.965
P6 (−285, −115, −155)
−26.530
−32.420
13.445
P7 (−565, 25, −155)
0.505
16.190
−3.750
P8 (−425, 25, −155)
−6.530
−2.440
−16.520
P9 (−285, 25, −155)
−32.350
−23.240
−41.650
P10 (−565, −255, −295)
−3.330
1.135
9.475
P11 (−425, −255, −295)
−18.710
−20.420
9.630
P12 (−285, −255, −295)
−42.530
−36.780
−5.190
P13 (−565, −115, −295)
0.235
−1.920
15.150
P14 (−425, −115, −295)
−11.630
−23.910
15.300
P15 (−285, −115, −295)
−36.200
−42.420
3.000
P16 (−565, 25, −295)
−3.750
7.910
−5.700
P17 (−425, 25, −295)
−16.520
−13.730
−4.380
P18 (−285, 25, −295)
−41.650
−34.610
−11.310The data of these 18 measurement points’ errors can also be shown in Figure9.Figure 9
The 18 measurement points’ errors.


## 3.3. Establishment of the Prediction Model for Nonfeature Points’ Errors

The nonfeature points’ errors in cubes are predicted by the method of internal ratio. Suppose that the errors’ data of all cubes’ eight vertices are all collected. In fact, every vertex has 6 errors, which can be expressed as (Cnx, Cny, Cnz, Cna, Cnb, Cnc); here, n = 1, 2, …, 8. This paper only needs to study the first three errors.In this paper, the length of the smallest cube in the measured lattice is 70 mm. In this section, for the convenience of calculation, it is assumed that the length is 1. Set up a rectangular coordinate system, and suppose P1 is located at the origin, and the lines from P1 to P2, P1 to P4, and P1 to P5 are in the positive direction of theX-axis, Y-axis, and Z-axis, respectively. In this way, their internal ratios can be completely calculated according to equation (1).Thus, theP point’s internal fraction is equal to the P point’s coordinate. The starting point and the ending point of interpolation are x0 = 0, x1 = 1, y0 = 0, y1 = 0, z0 = 0, z1 = 1, respectively. When the P point is at the center of the cube, the P point’s error values Cx, Cy, and Cz of the X-axis, Y-axis, and Z-axis can be obtained and expressed as(3)Cx=C1x,C2x,C3x,C4x,C5x,C6x,C7x,C8x⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z=C1x+C2x+C3x+C4x+C5x+C6x+C7x+C8x8,Cy=C1y,C2y,C3y,C4y,C5y,C6y,C7y,C8y⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z=C1y+C2y+C3y+C4y+C5y+C6y+C7y+C8y8,Cz=C1z,C2z,C3z,C4z,C5z,C6z,C7z,C8z⋅1−x×1−y×1−zx×1−y×1−zx×y×1−z1−x×y×1−z1−x×1−y×zx×1−y×zx×y×z1−x×y×z=C1z+C2z+C3z+C4z+C5z+C6z+C7z+C8z8.For the convenience of research, we chose a 3 × 3 × 3 lattice with a total of 27 points, which are, respectively, located in the center of the 27 smallest cubes. These 27 points are both nonfeature points and verification points, which are used to verify the effectiveness of the prediction method. The starting and ending coordinates of the lattice are (−600, −290, −120) and (−250, 60, −470), respectively, and the measurement step is 70 mm. In fact, the starting and ending coordinates of the measurement lattice are (−565, −255, −155) and (−285, 25, −435), respectively.Based on the above model and algorithm, we calculated the prediction values of these 27 points’ errors. The calculated results are shown in Table2.Table 2
The verification points’ prediction errors.
Predictive value
Y
1 (μm)
Y
2 (μm)
Y
3 (μm)
C
x
C
y
C
z
C
x
C
y
C
z
C
x
C
y
C
z
Z
1
X
1
2.140
0.210
0.290
−5.240
17.440
−23.230
−10.870
24.350
−28.510
X
2
−0.520
7.690
−9.920
−6.810
18.290
−26.170
−13.250
27.920
−31.860
X
3
−1.970
10.610
−20.340
−9.040
19.170
−28.420
−14.920
40.030
−35.750
Z
2
X
1
0.820
10.070
4.130
−3.120
19.130
−23.090
−11.990
28.660
−31.170
X
2
−1.160
11.540
−4.890
−5.130
19.870
−24.960
−12.810
40.190
−32.490
X
3
−1.980
16.120
−13.480
−8.730
24.380
−27.140
−14.950
42.810
−34.180
Z
3
X
1
3.170
15.120
−5.210
−3.620
20.240
−24.150
−9.770
38.040
−35.230
X
2
1.320
19.460
−13.780
−5.140
22.610
−27.210
−13.270
41.040
−36.890
X
3
−2.210
19.380
−22.940
−6.870
30.190
−31.360
−15.280
49.140
−37.540


## 3.4. Verification of the Prediction Model

To compare the difference between the prediction value and the actual value of the 27 points’ errors, we readjusted and installed the laser interferometer based on the method of measuring the whole working space points’ errors proposed in this paper. We measured directly the 27 verification points’ actual errors. These symbols ofDx, Dy, and Dz represent the 27 verification points’ errors of the X-axis, Y-axis, and Z-axis, respectively. The result can be seen from Table 3.Table 3
The verification points’ actual errors.
Measured value
Y
1 (μm)
Y
2 (μm)
Y
3 (μm)
D
x
D
y
D
z
D
x
D
y
D
z
D
x
D
y
D
z
Z
1
X
1
1.370
0.530
0.870
−4.610
16.980
−23.910
−11.420
23.930
−28.320
X
2
−1.270
8.150
−9.630
−6.030
18.540
−26.580
−12.490
28.260
−31.410
X
3
−2.850
10.280
−20.670
−9.480
18.950
−29.080
−14.530
40.410
−36.080
Z
2
X
1
1.570
10.210
3.690
−3.510
19.440
−22.580
−12.330
29.090
−30.820
X
2
−0.330
11.120
−5.150
−5.470
19.640
−25.210
−13.590
40.530
−32.730
X
3
−2.120
15.590
−13.170
−8.280
24.720
−27.580
−15.260
42.490
−34.710
Z
3
X
1
2.420
15.660
−4.890
−2.940
20.690
−24.730
−10.460
37.880
−34.880
X
2
1.140
19.870
−13.210
−4.530
22.850
−27.830
−12.910
41.370
−36.220
X
3
−1.860
19.150
−22.580
−7.450
29.770
−31.570
−14.420
49.660
−37.330From Tables2 and 3, we can get the comparison between the predicted value and the measured value of the 27 points’ errors. Its prediction accuracy reaches within ±1 μm, as it can be seen in Figure 10.Figure 10
The comparison between the predicted value and the measured value. (a) The comparison ofCx and Dx of X-axis error. (b) The comparison of Cy and Dy of Y-axis error. (c) The comparison of Cz and Dz of Z-axis error.
(a)
(b)
(c)


## 3.5. Evaluation of Verification Effect

In Tables2 and 3, we obtained the measurement value and the prediction value of the 27 verification points’ errors. The symbols of Δx, Δy, and Δz are defined as the difference between the prediction value and the measurement actual value. The formula is expressed as(4)Δx=Dx−Cx,Δy=Dy−Cy,Δz=Dz−Cz.The values of Δx, Δy, and Δz reflect the accuracy of the prediction model. The residual value is shown in Table 4. The comparison diagram between the prediction value and the measurement actual value of the 27 verification points’ errors is shown in Figure 11. The empirical analysis shows that the prediction accuracy of this method is only within 0.001 mm.Table 4
The residual value of the 27 verification points’ errors.
Residual value
Y
1 (μm)
Y
2 (μm)
Y
3 (μm)
∆x
∆y
∆z
∆x
∆y
∆z
∆x
∆y
∆z
Z
1
X
1
0.770
0.320
0.580
−0.630
−0.460
−0.680
0.550
−0.420
0.190
X
2
0.750
0.460
0.290
−0.780
0.250
−0.410
−0.760
0.340
0.450
X
3
0.880
−0.330
−0.330
0.440
−0.220
−0.660
−0.390
0.380
−0.330
Z
2
X
1
−0.750
0.140
−0.440
0.670
0.310
0.510
0.340
0.430
0.350
X
2
−0.830
−0.320
−0.260
0.340
−0.230
−0.250
0.780
0.340
−0.240
X
3
0.140
−0.530
0.310
−0.450
0.340
0.440
0.310
−0.320
−0.530
Z
3
X
1
0.750
0.540
0.320
−0.680
0.450
−0.580
0.690
−0.160
0.350
X
2
0.180
0.410
0.570
−0.610
0.240
−0.620
−0.360
0.230
0.670
X
3
−0.350
−0.230
0.360
0.580
−0.420
−0.210
−0.860
0.520
0.210Figure 11
The residual value of the verification points’ errors.


## 4. Conclusions

The current research on error modeling of machine tools ignores the problem that multiaxis linkage, coupling, and superposition will cause the spatial errors. A modeling and prediction method for CNC machine tools’ spatial errors was proposed based on feature points, which effectively overcomes the above shortcomings.In this method, through measuring the feature points’ errors to obtain the modeling data, the prediction model of nonfeature points’ errors was established by using the internal division ratio method. Thus, we can use this model to predict the compensation value of any point’s error in the whole space.The method has the advantages of small interpolation operation, easy integration in the numerical control system, and high compensation precision. Research results show that the smaller the cube is defined, the higher the prediction accuracy is. The example also proved this well, and its prediction accuracy reached within ±1μm, which shows this method is effective and feasible.The three key technologies to improve the spatial accuracy of machine tools are measurement, modeling, and compensation, respectively. In this paper, we solved the first. In fact, the study of compensation technology will be valuable and more complex. How to deeply integrate the error compensation model established in this paper with the CNC system will be challenging. We leave it as our further study.



---

*Source: 7254596-2020-05-12.xml*