A new method for error modeling in the kinematic calibration of redundantly actuated parallel kinematic machine

: This paper presents a new method for error modeling and studies the kinematic calibration of redundantly actuated parallel kinematic machines (RA-PKM). First, a n-DOF RA-PKM is split into several n-DOF non-redundantly actuated sub-mechanisms by removing actuators in limbs in an ergodic manner without changing the DOF. The error model of the sub-mechanisms is established by differentiating the forward kinematics. Then, the complete error model of the RA-PKM is obtained by a weighted summation of errors for all sub-mechanisms. Finally, a kinematic calibration experiments are performed on a 3-DOF RA-PKM to verify the method of error modeling. The positioning and orientation error of the moving platform is replaced by the positioning error of the tool center point, which was reduced considerably from 3.427 mm to 0.177 mm through kinematic calibration. The experimental results demonstrate the improvement of the kinematic accuracy after kinematic calibration using the proposed error modeling method.


Introduction 1
Recently, the redundantly actuated parallel kinematic machines (RA-PKMs) have received much more attention because actuation redundancy may eliminate singularities, enlarge the usable workspace, and improve the stiffness in different configurations [1][2][3][4]. The actuation redundancy in parallel mechanisms can be classified into two types: (1) in-branch redundancy where passive joints are actuated, and (2) additional branch redundancy where extra actuated limbs are added to the minimal number of limbs needed. The latter was preferred by many researchers because of better force distribution and higher stiffness. This paper focuses on the latter type of RA-PKMs.
Kinematic calibration of RA-PKMs has been proved to be an efficient method to improve accuracy [5], which can be divided into three categories, namely, constrained calibration, self-calibration, and external calibration. Constrained calibration methods keep some kinematic parameters of the mechanism constant in the calibration [6][7][8] and therefore has less need for additional measuring devices. Self-calibration uses internal sensors to measure the pose of the moving platform [9,10], which is cheap in cost and can be achieved online. However, neither the self-calibration method nor the constraint calibration method can identify all parameters. External calibration, which employs devices such as laser trackers [11], ball bars [12,13], and vision systems [14][15][16] to acquire the pose of a parallel mechanism, is widely applied in calibration. For kinematic calibration of RA-PKMs, external calibration is a more suitable method as it does not require mechanical constraints or extra sensors on passive joints. Normally, the external calibration procedure consist of four steps [17,18]: error modeling, error measurement, parameter identification, and parameter compensation. The error modeling is to establish the mapping between geometric errors and pose errors of a manipulator or machine, which is fundamental in kinematic calibration. There are two main ways to establish the error model of the mechanism. The first method integrates error models of each limbs considering geometry constraints. The error modeling of each limb can be generally established based on Denavit-Hartenberg (D-H) convention [19], product-of-exponentials (PoE) formula [20], and screw theory [21,22]. The method is simple and has clear physical meaning. The second method uses the differential-to-kinematic equations, which may yield error models expressed in analytical form for specific PKMs. Hollerbach and co-workers [23,24] successfully obtained the error model of a mechanism with the implicit loop method. However, the redundant limbs with actuatoin of a RA-PKM introduces more couplings and error sources, which cause great difficulties. Thus, the error modeling of a RA-PKM is still an open problem that has not been well addressed [25,26]..
The main contribution of this work is to propose a simple and effective error modeling method of RA-PKMs. The approximate errors of the RA-RAM are calculated by using the errors of all non-redundantly actuated sub-mechanism. A 3-degrees-of-freedom (3-DOF) RA-PKMs is taken as an example to verify the error modeling method proposed for kinematic calibrations. The rest of paper is organized as follows: In Section 2, a systematic process for kinematic error modeling of the non-redundant PKM and RA-PKMs is first introduced. Section 3 gives a case study for error modeling of a 3-DOF 2UPR-2RPU RA-PKM. Kinematic calibration ·3· experiments of the 2UPR-2RPU RA-PKM are described in Section 4. Finally, conclusions are drawn.

Kinematic error modeling
The errors of RA-PKM are derived from two parts: geometric errors and non-geometric errors. The Geometric errors include deviations in the length of the connecting rod and the joint axis, and zero errors in the actuated joint; while the non-geometric errors are generated from heat, backlash, et.al. Generally, the non-geometric errors are too complicated to model analytically, and the kinematic errors in geometric errors are the most important error source for the RA-PKM. Therefore, only the kinematic errors were taken into account in error modeling in this paper.

Kinematic error modeling of non-redundantly actuated PKMs
A PKM without redundantly actuated limbs ( Fig. 1) has n (n ≤ 6) DOF and n actuated limbs. Two hinge points for each limb, denoted Ai and Bi (i = 1, 2, …, n), are attached to the fixed base and the moving platform, respectively. The geometric centers of the platforms are the points labeled OA and OB. Normally, to derive the error model of a mechanism, the kinematic model of the mechanism should be established first. The fixed frame {OA-xyz} is located at the fixed base, and the moving frame {OB-uvw} is attached to the moving platform. For n actuated limbs, n constraint equations can be obtained as ( ) , , 0 1, 2,... , x y (1) where qi and xi represent the driving parameter and kinematic parameters of the i-th limb, respectively, and y denotes the pose of the moving platform. Combining the n constraint equations yields the kinematic model of the parallel mechanism   ( )  Because the driving errors are determined by the actuation system, usually they can be ignored to simplify the error model [30]. Expanding both sides of Eq.(2) to first-order in perturbation, the error transfer matrix can be obtained as represent the position error and orientation errors of the moving platform, and x represents the kinematic errors of the mechanism. The matrices Jr and JΩ are the Jacobian matrices of r and Ω relative to the kinematic errors x.
Normally, to facilitate the measurement of parameters during calibration, a tool is mounted on the moving platform. The position vector of the tool center point P (TCP) is denoted by p and written , c =+ p r Rl (4) where lc represents the position vector of the TCP w.r.t. the moving frame, and R denotes the rotation matrix of the moving frame w.r.t. the fixed frame, which is always expressed as Introducing the first-order terms in a perturbation on both sides of Eq. (4) yields the error model of the PKM . cc  =  +  +  p r Rl R l (6) Lei-Ying He et al.

·4·
Hence ΔR can be written as 0 - where S(ΔΩ) is an antisymmetric matrix about the vector ΔΩ. Thus Eq. (6) can be simplified as (8) where Δlc is the position vector error of the tool w.r.t moving frame, xl . The error model of the mechanism is finally obtained.
The constraint equations similar to the Eq. (1) can be built for each sub-mechanism. Assuming each sub-mechanism has only one pose in a specific configuration, N poses of the moving platform can be calculated by solving all normal constraint equations. Therefore, in a specific configuration, the N position vectors of the TCP can also be calculated from (10) where rj and Rj are respectively the position vector and the rotation matrix of the moving frame w.r.t. the fixed frame for the j-th sub-mechanism. The nominal driving parameter is calculated according to the pose of mechanism and kinematic model. The position of the TCP calculated through different direct ·5· kinematic model will be different due to the existence of kinematic error. In order to solve this problem, the driving parameters of redundant limbs are adjusted accordingly to make the mechanism meet the constraint. As a result, it is assumed that there is no conflict between kinematic limbs and no deformation of parts in the process of calibration. Since the sub-mechanism is a part of RA-PKM indeed, the position of the TCP calculated from different sub-mechanism is supposed to be consistent.
The N position vectors Pj is supposed to be And the actual position of TCP can be expressed as  (14) where Jj and Δεj are the error transfer matrix and kinematic error of the j-th sub-mechanism, and T TT j By substituting Eq. (14) into Eq. (13), the complete error model of the RA-PKM is written (15) wherein Δε includes kinematic errors Δx of the m limbs and the geometric error Δlc of the tool, Δlc represent the position vector error of the tool w.r.t moving frame. Therefore, before using Eq. (15), all Δεj should be augmented with Δε, and correspondingly, the Jacobian matrices Jj also need augment by inserting zero-row vectors. In this way, the error modeling method applicable to non-redundant PKM can also be used in the error modeling of RA-PKM.

Example
In this section, the kinematic calibration of a 3-DOF RA-PKM is studied to verify the proposed method of error modeling. A 2UPR-2RPU RA-PKM (Fig. 4) is composed of a fixed base, a moving platform, a tool, two UPR limbs, and two RPU limbs. The tool is connected to the center of the moving platform. Each limb has one universal (U) joint, one prismatic (P) joint, and one revolute (R) joint, the prismatic (P) joints with underline denote actuated joints. The difference in these two types of limbs is the way they are installed.

Descriptions of the 3-DOF RA-PKM
In Fig. 5, the geometrical centers of the universal joints are the points A1, A2, B3, B4, whereas the intersection points of the axes of the revolute joints and the axes of the prismatic joints are the points A3, A4, B1, B2. The TCP of the mechanism is denoted as P. Points Ai (i = 1, 2, 3, 4) and Bi (i = 1, 2, 3, 4) form two separate squares, the centers of which are points OA and OB, respectively. The fixed frame {OA-xyz} is established at the origin OA with its z-axis pointing downward and perpendicular to the fixed base; the x and y axes are parallel to OAA3 and OAA2, respectively. The moving frame{OB-uvw} is assigned to point OB with the w-axis vertical to the moving platform; the u and v axes are parallel to OBB3 and OBB2, respectively. The axes of the revolute joints and the universal joints in the fixed base (Fig. 5) are parallel to either x or y axis, whereas the axes of the joints in the moving platform are parallel to either u or v axis.
Lei-Ying He et al.
Generally, the length of link OAAi and OBBi (i = 1, 2, 3, 4) are supposed to be a and b, respectively. However, the lengths of all links may not be equal to their ideal values in practice because of manufacturing and assembly errors. In the following sections, the actual lengths of OAAi and OBBi (i = 1, 2, 3, 4), the zero position errors of the four P-pairs, and the position vector of the TCP w.r.t. the moving frame are considered to be identified through a kinematic calibration.

Kinematic Modeling
The 2UPR-2RPU RA-PKM is constructed by adding a RPU extra-actuated limb on the PKM 2UPR-RPU (Fig. 6). Like the original mechanism, the 2UPR-2RPU RA-PKM has 3 DOF: translation along the z-axis, rotation relative to the y-axis, and rotation relative to the u-axis [31]. The angle of rotation relative to the y-axis (u-axis) is denoted by θ (ψ). Because the constraints of R-pairs on the UPR limb, vector OAOB is always perpendicular to vector B3B4, the angle of rotation relative to the y-axis is equal to the angle of rotation relative to the v-axis, as shown in Fig. 6. So that the rotation matrix R between the two coordinate systems can be written as Because of the constraints provided by the RPU limbs, observing the machine in the direction of A2OA, the vector r which represents the translation of the moving platform (Fig. 6) is written   T tan 0 .
zz  = r (17)  Considering OAAiBiOB (i=1,2,3,4) as a closed-loop vector kinematic chain, the constraint equations of the machine can be expressed as , 1, 2,3, 4 , (18) where ai and qi denote the link OAAi and AiBi w.r.t. the fixed coordinate system, bi denotes the link OBBi w.r.t. the moving coordinate system. Taking into consideration the manufacturing and assembling errors, the actual length of link OAAi and OBBi is not always equal to designed values, the constraint equations are then written as where δqi 0 (i = 1,2,3,4) denote the zero position errors of the four P-pairs, and where δai and δbi denote the error of the lengths of OAAi and OBBi (i = 1, 2, 3, 4). For the 2UPR-2RPU RA-PKM, the pose of the moving platform can be described using three quantities, namely, z, θ, and ψ. As mentioned above, four sub-mechanisms can be constructed by selecting three of four limbs as actuated limbs (Fig. 7). From the figure, limb 4, 1, 3, and 2 are set as the passive limb corresponding to the four sub-mechanisms labeled as SMj, (j = 1,2,3,4). There are two types of sub-mechanisms, specifically 2UPR-RPU and 2RPU-UPR. Both SM1 and SM3 are of type ·7· 2UPR-RPU PKM, whereas SM2 and SM4 belong to type of 2RPU-UPR PKM.
The kinematic model of a sub-mechanism can be solved theoretically. Hereinafter, SM1 is used as an example to briefly describe the procedure in solving the kinematic model. For SM1, the first three constraint equations in Eq. (19)   Similar to SM1, the kinematic model of the remaining SMj, (j = 2, 3, 4) can be deduced by the same procedure, and their corresponding position of TCP, denoted p2, p3, and p4, can also can be determined.
As introduced in the section 2.2, the actual position of the TCP can be written Here all the coefficients ωj are set to 0.25 to simplify the calculation according to the constraint condition in Eq. (12).

Error modeling
Given the new proposed kinematic model, the error model  (27) As mentioned above, the error models of the four sub-mechanisms can be obtained firstly by taking a first-order perturbation expansion on both sides of Eq. (25); that is, Considering the first-order terms on both sides of Eq. (16) and (17) The first sub-mechanism SM1 is used as an example to outline briefly how the solution of the error transfer matrix is obtained. The Eq. (19) can be rewritten as   ; ; ( 1, 2,3, 4) ; ; ; The error model of SM1 can be obtained by solving the first three equations of Eq. (32), (34) where J1 0 denote the mapping relation between   T z     and the kinematic error of SM1.
, and Δqi 0 is the zero position errors of the four P-pairs; Δai and Δbi are the position errors of the U joint and R joint in the fixed base and moving platform, respectively. Note that the error parameter Δqi 0 stems from the assembly precision of the prismatic joint rather than the actuated system. Substituting Eq. (34) into Eq. (30), the error model of the j-th sub-mechanism is simplified giving 1 1 1 ,  =  pJ ε (35)

Kinematic calibration experiments
To validate the error model of the 2RPU-2UPR RA-PKM, a kinematic calibration experiment was performed (Fig. 8).
A laser tracker (Leica-AT901-LR, Leica Geosystems AG) was used to measure the 3D coordinate values of the TCP. The measurement uncertainty for the laser tracker is ±15 μm + 5 μm/m within a 2.5 m × 5 m × 10 m volume. The repeatability of the prototype of the 2RPU-2UPR RA-PKM was measured less than 0.15 mm.

Measurement of errors of the TCP
The nominal and measured positions of the TCP w.r.t. the fixed frame ( respectively, where Rt and rt are the rotation matrix and the translation vector of the measurement frame {Om-xmymzm} w.r.t. the fixed frame, and lm is the measured position vector of the TCP w.r.t. the measurement frame. Given the kinematic parameters and driving parameters, the nominal position of the TCP p N was calculated by the kinematic model as mentioned above. The measured value of the TCP p M was obtained using the laser tracker. The measurement error at one configuration was defined as the difference between the nominal and measure values of the TCP. For example, for the k-th calibration point, its measurement error was expressed as ( )

Parameters identification
During error identification, the object function and the fitting method need to be established. The error identification of the PKM can normally be converted to a minimization problem, The minimization is a nonlinear least square problem with an iterative relation. To solve this problem, various methods have been proposed, such as the Levenburg-Marquardt algorithm [32], the Kalman filtering approach [33], and the Ridge estimation method [34,35]. Here it given as where t represents the iteration times, M denote the identification matrix which is defined in Eq. (41). The termination condition of the iteration is that the error residual ||Δe|| or the change in Δε between two adjacent iterations is sufficiently small. For the selection of the measurement configurations, the calibration is well performed in a sensitive area of the entire workspace, for example, the boundary of the workspace. The workspace of the 2UPR-2RPU machine with a tool of [0 0 125] T is shown in Fig. 9. Three ellipses on the surface that fit the shape of the boundary of the workspace were constructed with the intention to take a certain number of calibration points distributed evenly on these ellipses.
To avoid parameter identification failure caused by the singularity of the identification matrix, the number of calibration points also has certain requirements. Since there are in total 15 kinematic errors considered in the error model, here 16 calibration points on each ellipse, totally 48 calibration points in the workspace (Fig. 9), were selected. That is, the number of calibration points is more than twice the number of error parameters to be identified, thus ensuring the validity and accuracy of the calibration [36].

·10·
Workspace Calibration point Fig. 9. Distribution of the calibration points within the workspace.

Experiments
A calibration experiment on the parallel manipulator using a laser tracker was performed (Fig 10) which is obtained through the fitting of plane, and the unit of rt is millimeter. Table 1 lists the values of the nominal kinematic parameters for this mechanism. The initial kinematic errors are set to zero.  The calibration comprises four steps: error modeling, error measurement, parameter identification, and error compensation, and its procedure is shown in Fig. 11.

Results and discussion
Since the machine has three DOFs: two rotational DOFs and one translational DOF, three independent motion parameters ψ, θ and z are sufficient to describe the configuration of the moving platform. However, it is impossible to compare these three parameters in a space uniformly, because their dimensions are not the same. Alternatively, the position vector of the TCP p = [px py pz] T is something should be paid more attention to when calibrating, whose error denotes the position error and orientation error of moving platform. In this paper, the position vector p, which can be expressed by three motion parameters ψ, θ and z using Eq. (4), replaces these three parameters to analyze the accuracy of the machine. The position error and orientation error of moving platform can be denoted by the position error of TCP. Fig. 12 shows the errors of 48 calibration points (position errors of TCP) before and after calibration. The mean, maximum, and standard derivation of the position residual decrease from 3.427 mm to 0.177 mm, from 4.740 mm to 0.384 mm, and from 0.547 mm to 0.084 mm, respectively. In other words, the calibration significantly reduces the residuals. For an actual mechanism, Δa and Δb are the deviations of distances from the U-pair and R-pair to the center of the fixed base and moving platform, respectively, Δq 0 is the error of the origin position of the P-pair, and Δlc the position deviation of the TCP w.r.t. the moving platform. The calibration result shows that Δa1 and Δa2 were smaller than Δa3 and Δa4 because the accuracy of the length of link A1A2 is easier to ensure than that of A3A4 in the assembly process. Also, Δb3 and Δb4 were larger than Δb1 and Δb2 because the former represent the position error of the U joint and the latter are those of the R joint. Recall that, in terms of manufacturing and assembly, the U-pair is more complex than the R-pair.
To verify the effect of compensation with the identified kinematic parameters in the whole workspace, some test points, either inside or outside the calibration area, are chosen. The distribution of test points (Fig. 13) show that their distribution characteristics are basically the same as those of the calibration points. The difference is the distance from the point to the boundary of the workspace. Of course, these test points are still in the workspace. The driving parameter is calculated with the compensated kinematic model; the result is presented in Fig. 14   In the test area outside the calibration area, the mean and maximum values of the deviation in position decreased from 3.570 mm and 4.767 mm to 0.440 mm and 0.845 mm, respectively. With the test area inside the calibration area, the mean and maximum values of the deviation in position decreased from 2.832 mm and 3.606 mm to 0.366 mm and 0.764 mm, respectively. The conclusion drawn from Fig. 14 is that, regardless of the location of the test point, the deviation in position is significantly reduced after compensation, thereby demonstrating the effectiveness of this calibration in the whole workspace and showing that in practice the proposed error model is suitable for the 2UPR-2RPU RA-PKM. Nevertheless, the compensated result from the test area is seen to be not as ideal as in the calibration area. Also, the compensation results of the test points inside the calibration area are slightly better than those outside the calibration area. The reason is that the closer the test point is to the boundary of the workspace, the more significant the effect of sources of error on the deviation in position.
Although, the error of the RA-PKM will decrease dramatically after kinematic calibration using the proposed error modeling, some challenges and limitations still exist. The accuracy of the error model is dependent heavily on the distribution of weight coefficients ω in Eq. (13). Actually, the weight coefficients ω are affected both by the stiffness and geometric errors of the links. However, their clear relationship is difficult to establish. For simplification, equal weight coefficients were determined in the paper, which is an approximate and simple solution, but still effective for kinematic calibration. To further improve the accuracy of the error model, a more approximate influence of weight coefficients should be studied in the future work.

Conclusions
A new error modeling method for kinematic calibration of RA-PKMs is presented considering the kinematic errors. The complete geometric error model of a RA-PKM can be calculated as the weighted summation of error models of all non-redundantly actuated sub-mechanisms which keep the same mobility as the PA -PKM. The proposed method is verified to be effective through a kinematic calibration experiment using a 3-DOF 2UPR-2RPU RA-PKM. It is shown that the mean positioning error of the TCP of the RA-PKM is decreased from 3.427 mm to 0.177 mm after kinematic calibration, which shows the improvement of the accuracy of the RA-PKM. Additionally, the result of the experiment also indicates the effectiveness of the parameters compensation after kinematic calibration in the whole workspace.