A Kinematic Calibration Method of a 3T1R 4-Degree-of-Freedom Symmetrical Parallel Manipulator

This paper proposes a method for kinematic calibration of a 3T1R, 4-degree-of-freedom symmetrical parallel manipulator driven by two pairs of linear actuators. The kinematic model of the individual branched chain is established by using the local product of exponentials formula. Based on this model, the model of the end effector’s pose error is established from a pair of symmetrical branched chains, and a recursive least square method is applied for the parameter identification. By installing built-in sensors at the passive joints, a calibration method for a serial manipulator is eventually extended to this parallel manipulator. Specifically, the sensor installed at the second revolute joint of each branched chain is saved, replaced by numerical calculation according to kinematic constraints. The simulation results validate the effectiveness of the proposed kinematic error modeling and identification methods. The procedure for pre-processing compensation on this 3T1R parallel manipulator is eventually given to improve its absolute positioning accuracy, using the inverse of the calibrated kinematic model.


Introduction
Compared with serial manipulators, parallel manipulators have the advantages of high rigidity and high loading capacity [1], but also have the disadvantages of complex structure and relatively small workspace. Recently, a 3T1R symmetrical parallel manipulator with a simple structure and a large working space has received widespread attention in academia and industrial applications [2]. This symmetrical parallel manipulator uses four open branched chains to connect the fixed platform with the moving platform. It can achieve three-degree-of-freedom translation along the X, Y, and Z axes and one-degree-of-freedom rotation around the Z axis (3T1R), for a total of four degrees of freedom. Its advantages such as large workspace and high speed make it be widely applicable to industrial automation scenarios such as pick-and-place and sorting, if its absolute positioning accuracy meets the particular requirements of these applications. Therefore, it is necessary to develop an effective calibration method to improve the absolute positioning accuracy of this symmetrical parallel manipulator [3,4].
The kinematic error of the manipulator is defined as the end effector's error between the actual pose and the nominal one, which is mainly caused by geometric tolerances, such as assembly and has multiple passive joints, all of them can be treated as active ones as in the series robot, if sensors are attached to them during calibration. In this way, the calibration method by the local POE formula used for serial robots can be directly applied to parallel manipulators with passive joints. Specifically, to save the number of sensors used in the revolute joints, the angle of the second revolute joint of each branched chain is computed numerically according to the kinematic constraints of the mechanism. Eventually, a linear-in-parameter error model is synthesized, which greatly improves the efficiency in the coming least squares estimation of parameters. Due to the symmetry of this mechanism, the idle pair of branched chains can be used for cross-validation of the kinematic and error models, thereby ensuring the model accuracy. A simulation study of the calibration algorithm is performed to verify the effectiveness of the algorithm. Last but not least, the pre-processing compensation procedure of the command pose is suggested for this 3T1R parallel manipulator.

Kinematics of a Symmetrical 3T1R Parallel Manipulator
In order to establish a kinematic model for the symmetrical 3T1R parallel manipulator, the basic structure of this parallel manipulator is firstly introduced. Subsequently, its forward kinematics based on the local POE formula is derived. It lays the foundation for the subsequent establishment of the kinematic error model for this parallel manipulator.

Structure of the 3T1R Parallel Manipulator
As shown in Figure 1, it has two identical branched chains, yielding symmetrical structures. Each branched chains has two identical sub-branches, so the whole mechanism can be treated as a "binary-tree" structure. This parallel manipulator is driven by four linear actuators installed on the base, and two linear actuators on one side share one stator. The entire manipulator mechanism looks like the shape of the letter "M", and it can move along the X, Y, and Z axes and rotate about the Z axis.

Collinear linear modules
Moving platform

Dyad kinematics Based on local POE formula
The Dyad branched chain is shown in Figure 2. Based on the traditional local POE formula method of building in Figure 2a, the Dyad kinematics can be obtained as (1): where T j−1,j ∈ SE(3), represents the pose of {O j } relative to {O j−1 } and uses the representation method of the local coordinate system.ŝ j is the corresponding element for T j−1,j in se (3).
Dyad kinematics based on different methods of establishing coordinate systems for the local POE formula. (a) to establish a coordinate system arbitrarily, (b) to establish a coordinate system which leads to reduced number of kinematic parameters.
According to the characteristics of the local POE formula, it is known that the local coordinate system of the POE formula can be assigned arbitrarily, that is the pose of {O j } is arbitrary, so it requires 13 parameters to completely define T j−1,j . Among them, the 13 parameters include: 6 parameters within initial pose T j−1 (0) ∈ SE(3), 6 parameters within joint twistŝ j ∈ se(3), and 1 parameter within joint variable q j ∈ R 1X1 . Compared with the D-H method, which only needs four parameters to define the kinematics, this method of establishing the kinematic model requires too many parameters.
In order to reduce the number of parameters in the kinematic model, a new method of assigning the local coordinate systems is proposed, as shown in Figure 2b. In this method, the origin of the coordinate system {O j−1 } is at the center of the joint, and the direction of the z axis is along the joint line. Furthermore, an additional local coordinate system {O j−1 } is on the same link at the joint j. Now, (1) is written as: where T j−1,j (0) represents the pose transformation of {O j−1 } relative to {O j−1 }. eŝ j q j represents the pose transformation of {O j } relative to {O j−1 }. According to the transformation between Lie groups and Lie algebras, T j−1,j (0) ∈ SE(3), at least onet j ∈ se(3) exists, making T j−1,j (0) ∈ SE(3). Therefore, T j−1,j (0) = eˆt j , and (2) is written as:

Branched Chain Kinematics Based on the Local POE Formula
Based on (2), consider a single branched chain of the parallel manipulator with (n + 1) links, numbered sequentially in the order 0, 1, 2 · · · , n from base coordinate system {0} to tool coordinate system {n}; the forward kinematics is expressed as:

Kinematics of the Parallel Manipulator Branched Chain i Based on the Local POE Formula
The schematic of the 3T1R parallel manipulator is shown in Figure 3a. It has a symmetrical structure, dividing it into two identical parts, I and II. This parallel manipulator contains four linear actuators B 1 , B 2 , B 3 , and B 4 . Each linear actuator to the midpoint of the end effector p is regarded as a branched chain. On side I, the branched chains formed by B 1 and B 2 to the end p point are recorded as i and i, respectively, where i = 1; on side II, the branched chains formed by B 3 and B 4 to the end p point are recorded as i and i, respectively, where i = 2. Any one of the branched chains from B 1 and B 2 and any one from B 3 and B 4 are chosen to build the entire parallel manipulator kinematics, so there are four possible combinations. In this article, B 1 and B 3 are chosen for kinematic calibration and recorded as the i th branched chain, i = 1, 2, which is the structure indicated by the dark solid line in Figure 1. Surely, because this parallel manipulator structure is symmetrical, the remaining two branch chains can be used to carry out the calibration as well for cross-validation. This also helps to enhance the accuracy of the kinematic model.
There are four joint modules within each branched chain i (i = 1, 2), which are one active translation joint and three passive rotation joints. The coordinate systems are established as shown in Figure 1, and the way of naming is referenced in Nomenclature. The joint ij is an active joint for i = 1, 2 and j = 1, and the joint ij is a passive joint for i = 1, 2 and j = 2, 3, 4. The pose of the end effector coordinate system {P} relative to the parallel manipulator base coordinate system {B} is defined as the forward kinematics. The kinematic equation of the mechanism is given as:

Establishing the Kinematic Error Model for this 3T1R Parallel Manipulator
Since this parallel manipulator is a closed-chain mechanism, the kinematic error at the end effector of the parallel manipulator should be the same if they are calculated from individual branched chains.

Establishment of the Kinematic Error Model from a Single Branched Chain
An error model of the 3T1R parallel manipulator, considering the kinematic errors of the single branched chain, is established. Due to the part manufacturing and assembly error, the actual pose of the end effector is different from its nominal value. Since: is the kinematic transformation from coordinate system {i1} to coordinate system {B}. The kinematic Equation (5) of a single branched chain is simplified as: Ift ∈ se(3) for a given T(0) ∈ SE(3), then eˆt = T holds. Therefore, for the initial pose T i(j−1),ij (0), , so (6) is simplified to: From (6), it is known that the modeling method using the local POE formula, whose forward kinematics g i is a function of the initial pose Differentiate the kinematic parameters of (8) to obtain the calibration model, yielding: where δT · T −1 represents the end effector's pose error of this parallel manipulator base coordinate system {B} through the branched chain i, and the motion parameter error is δT(0), δs, δq. The purpose of the kinematic calibration algorithm is to find the appropriate δT(0), δs, δq, achieving the least squares fitting between two sides of (9), which is: From earlier analysis, the kinematic error model consists of 13 parameters from the initial poses, joint twists, and joint variables. In order to simplify the calibration model, we assign the coordinate systems in the way shown in Figure 2b. Since the initial pose coordinate system in the local POE formula can be established at any point on the link, the initial pose T(0) is floating. Therefore, the kinematic error is regarded as only due to the deviation of the initial poses T(0), while the joint twists and joint variables are regarded as being accurate [12]. Therefore, (9) and (10) are simplified to: By using two branched chains for kinematic calibration, the goal of optimization is revised as: From (7), δg i used in (11) and (13) is obtained as: If the parameter error is relative to its local coordinate system, δ(e t ij ) = eˆt ij δt ij , so (14) is further expressed as: δg i = eˆt i1 δt i1 eŝ i1 q i1 · eˆt i2 eŝ i2 q i2 · eˆt i3 eŝ i3 q i3 · eˆt i4 eŝ i4 q i4 · eˆt P + eˆt i1 eŝ i1 q i1 · eˆt i2 δt i2 eŝ i2 q i2 · eˆt i3 eŝ i3 q i3 · eˆt i4 eŝ i4 q i4 · eˆt P +eˆt i1 eŝ i1 q i1 · eˆt i2 eŝ i2 q i2 · eˆt i3 δt i3 eŝ i3 q i3 · eˆt i4 eŝ i4 q i4 · eˆt P + eˆt i1 eŝ i1 q i1 · eˆt i2 eŝ i2 q i2 · eˆt i3 eŝ i3 q i3 · eˆt i4 δt i4 eŝ i4 q i4 · eˆt P +eˆt i1 eŝ i1 q i1 · eˆt i2 eŝ i2 q i2 · eˆt i3 eŝ i3 q i3 · eˆt i4 eŝ i4 q i4 · eˆt P δt P .
By multiplying (15) with g −1 i on both its left and right sides, where g −1 From the adjoint transformation Ad X y = XyX −1 , (16) is simplified to: where: For (17), the left term δg i · g −1 i ∈ se (3) represents the end effector's pose error in the base coordinate system {B}, according to the definition of the logarithm of the matrix on SE (3): where g a i is the pose measured at the end of the branched chain i and g −1 i is the inverse of nominal pose matrix at the end of the branched chain i. From (17) and (18), we get: (19) is written as a linear-in-parameter form: where y i is the pose error of the end effector from the branched chain i, J i is the error Jacobian matrix of the branched chain i, and x i is the kinematic parameter error of the branched chain i, given as:

Establishment of the Overall Kinematic Error Model
The kinematic error of the end effector calculated from the kinematic error models of different branched chains should be identical. This yields: where: In the above error model, it has a total of 60 error parameters to describe the kinematic error. Through the discernibility analysis of the error model deduced from the POE formula, the maximum number of independent motion parameters of a universal non-over-constrained parallel manipulator after eliminating redundant error components is 4r + 2p + 6, where r and p represent equivalent rotation joints and translation joints, respectively [12]. In the 3T1R parallel manipulator, there are six equivalent rotation joints and two translation joints. Therefore, the maximum number of identifiable parameters after eliminating redundant error components is 34. According to (21), the number of errors obtained using the method described in this article is 60, so the error model does have redundancy in terms of the number of parameters.

Method to Reduce the Number of Sensors Used in Passive Joints
In general, since the parallel manipulator contains many passive joints, a sensor will be installed on each passive joint to measure the corresponding joint angle throughout the calibration process. The required measurements for calibration are accomplished by both built-in sensors and external measurement device. The built-in sensors are used for measuring rotation of the passive joints, and a laser tracker is used for measuring the end effector's pose. As the 3T1R parallel manipulator has many passive joints and such a measurement scheme is costly, a recursive method is proposed to estimate one of the joint angles in each branched chain, so that the number of sensors being used is reduced.
This recursive method is based on the coordinates of the end point A i and uses the local POE formula method to establish the branched chain kinematics of the end point A i , such as (22): where p i represents the coordinates of point A i on the branched chain i relative to the base coordinate system {B} and p i represents the coordinates of the point A i on the branched chain i relative to the coordinate system {i3}. By setting: (22) is re-written as: The difference between the nominal and actual distance of two revolute joints at the end effector is given by: where − −− → A 1 A 2 represents the nominal distance of between A 1 and A 2 and d 12 represents the actual distance between A 1 and A 2 . Here, Differentiating (26), we get: With (25), (26), (27), and (28), it yields: where: (29) is written in a recursive form as: where k is the number of iterations. This means, by utilizing the geometric constraints of the mechanism, the proposed recursive method is able to estimate the rotation angle of the second revolute joint of each branched chain. This saves one rotation sensor for each branched chain during calibration.

A Recursive Least Squares Method to Identify the Parameters in the Kinematic Error Model
In order to improve the efficiency, a recursive least squares method is used in the parameter identification process for the kinematic error model. Suppose that when measuring positions of m points, there will be m position error vectors, and the Jacobian matrix will be expanded to m terms, so the kinematic error vector X will remain unchanged. Hence, (31) is obtained as: The solution in the least squares sense for X is given as: (32) can be further simplified by recursion. Once the kinematic error parameter X is determined, the initial pose T i(j−1),ij (0) will be updated by substituting X into the following equation: until the norm X of X approaches zero. Therefore, the pose computed from the calibrated kinematic model approaches the actual pose.

Simulation Results
The kinematic modeling schematic diagram and the establishment of the coordinate system are shown in Figure 3, where l u = 0.7, l d = 0.7, a = 0.3, b = 1.0 (units: m). According to the modeling method using the local POE formula given in Nomenclature, the values of the initial pose and joint twist of adjacent links of two branched chains are given as: In order to simulate the deviation of kinematic parameters from the nominal ones, the errors δt, δs, and δq are applied to the initial poses, joint twists, and joint variables in the theoretical model. Noise is introduced in the measurement of joint variables to test the robustness of the parameter identification method. In this way, an actual model is obtained, where the actual initial pose is: the actual joint twist is: the actual joint variable is: Hence, the kinematic model of the actual branched chain i is: Add deviations δt, δs, and δq to the nominal values to get the actual value of the parallel manipulator simulation, as shown in Table 1. Each joint, regardless of being active or passive, is considered as a joint with one degree of freedom, so the error vector of the spin volume should be orthonormal, that is w i + δw = 1,(w i + δw) T (v i + δv) = 0.  Use the numerical forward kinematics algorithm to obtain the joint displacements and joint angles of 20 different parallel manipulator poses; 2.
Assign errors to kinematic parameters, such as δt, δs, and δq, as shown in Table 1.
The actual joint twist is computed as s a ij = s ij + δs ij ; 5.
The actual joint variable is computed as q a ij = q ij + δq ij ; 6.
The recursive calibration algorithm is used to identify the kinematic errors of the parallel manipulator.
The initial pose after calibration is recorded as T c i(j−1),ij (0), as shown in Table 2. When trace(R) = −1, 1 + 2 cos φ = trace(R), and ω < π, we get log(T) = t. That is, where,ω If φ is small,ω ≈ (R − R T )/2. When trace(R) = −1, log R = (2k + 1)πν. In this case, k is any integer, and ν is a unit eigenvector of log R with an eigenvalue of one.  The kinematic error between the calibrated value T c and the actual value T a is obtained through the above method, as shown in Table 2, including the orientation and position error at the end effector. The orientation error is the angular difference between the actual value of the end effector and the calibrated one. Position error is the norm of the coordinate difference between the actual value of the midpoint of the moving platform and the calibration one. Note that the derivation of kinematic parameters was arbitrarily assumed in this simulation example, and it may not truly reflect the actual situation in the parameter space due to manufacturing and assembly error. However, through the end effector's pose error simulation result as in Figure 4, it is seen that both the orientation and position errors converged, and they were reduced by at least 10 times, thereby verifying the effectiveness of the calibration algorithm. Compared with other methods based on the deviation of the distance between points on the end-effector of parallel manipulators [26,27], our method gave better fitness of the pose error at the central point of the end effector. Since the model was linear-in-parameter, the pose fit with higher matching and better fitting conversion effects. Therefore, in the case of a linear error model, the pose error was used for fitting so that the error reduction only needed one iteration.  To guide the future calibration experiments, the flowchart of the calibration process is given in Figure 5. The meanings of various symbols in the flowchart can be referenced in Nomenclature.

Pre-Processing Compensation
Since some robot controllers cannot directly utilize the calibrated active and passive joint angles to drive the robot, the joint variables were preprocessed to obtain the pose of the end effector that could be directly utilized by the robot.
First, the initial pose T c i(j−1),ij (0) after calibration of each link was obtained through the calibration process described above, as shown in Table 2. Combined with the commanded end pose T co , the calibrated joint angle q c ij was calculated using the inverse kinematics after calibration. Since the manipulator could not directly use q c ij , the following conversion was required. By bringing q c ij into the positive kinematics solution based on the nominal kinematics parameters, the compensated end effector pose T cp could be obtained. At this time, T cp could be directly used for offline programming. The flowchart of this process is shown in Figure 6. The meaning of various symbols is listed in Nomenclature.

Conclusions
In this paper, a kinematic model of a symmetrical 3T1R parallel manipulator was proposed based on the local POE formula. By using local POE formula for kinematic error modeling, the joint twists and joint variables could maintain their nominal values, and the errors were set to be concentrated in the matrices of initial poses. As the end effector's pose errors of the parallel manipulator established through different branched chains were identical, the entire pose error of the parallel manipulator's end effector was formed by equaling the pose error of the end effector calculated in individual branched chains.In addition, the joint angle of the second revolute joint was calculated by a recursive algorithm according to the inherent kinematic constraints. This saved the rotary sensors installed at the second revolute joint of each branched chain. Through simulation, it was shown that the deviation between the actual and calibrated end effector's pose could be reduced by at least 10 times with the proposed kinematic error modeling and parameter identification methods. The procedures to use this error model for profile pre-compensation in the actual testbed were also proposed.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: 3T1R parallel manipulator The parallel manipulator can achieve three degrees of freedom of translation along the X, Y, and Z axes and one degree of freedom of rotation around the Z axis POE formula The product of exponentials formula D-H convention The Denavit-Hartenberg convention

Nomenclature
For the ease of reference, the definition of the symbols used in this article are listed as follows: n The number of passive joints i Represents the parallel manipulator B 1 , B 3 branched chains (i = 1, 2) i Represents the parallel manipulator B 2 , B 4 branched chains (i = 1, 2) j The number of joints on each branched chain (j = 1, 2, 3, 4) g i Forward kinematics of these parallel manipulator branched chains i T i(j−1),ij (0) The initial pose of the coordinate system {ij} relative to the coordinate system {i(j − 1)} T c i(j−1),ij (0) After calibration, the initial pose of the coordinate system {ij} relative to the coordinate system {i(j − 1)} T a i(j−1),ij (0) The actual initial pose of the coordinate system {ij} relative to the coordinate system {i(j − 1)} q ij Standard representation of joint variable; it represent the joint angle or joint displacement q n ij Standard representation of the nominal joint variable q c ij Standard representation of the joint variable after calibration

{B}
Base coordinate system of the parallel manipulator {i0}, {i0} Base coordinate system of these parallel manipulator branched chains i and i {i1}, {i1} The translational motion coordinate system of these parallel manipulator branched chains i and i {i2}, {i2} The rotational motion coordinate system where it is connected to the modules of these parallel manipulator branched chains i and i {i3} The coordinate system of l d link rotates around the D point of these parallel manipulator branched chains i and i

{i4}
The parallel manipulator moving platform rotational motion coordinate system around the A i point {P} Midpoint coordinate system of the end effector g c Forward kinematics after calibration of the parallel manipulator T n The nominal pose of the end effector T a The actual pose of the end effector T co The command pose of the end effector during error compensation T cp The end pose after compensation Ad g i The adjoint transformation of g i , also written as Ad (g i ) , g i ∈ SE(3) The representation of the end error on these parallel manipulator branched chains i in the base coordinate system {B} δt ij The error of t in coordinate system {ij} J i Error Jacobian matrix of these parallel manipulator branched chains i J Error Jacobian matrix of the entire parallel manipulator