INTRODUCTION
Comprehensive studies and works have been made in the introduction of robotic
technology in clinical medicine^{[1]}. Two main areas are addressed
in Paolo et al.^{[1]} where robotic systems have been proposed
to serve as an aid to the doctor in some surgical operations or to serve as
aid systems to people with disabilities. Despite these comprehensive studies,
in fact most of these comprehensive robotic systems are still research devices
and only few of them have been put into practice. Beside this, it is worth mentioning
that most robotic systems in medical applications are of serial manipulators
because most of them are simply derived from industrial robotics. For example,
the 6 degree of freedom robot for proton therapy center^{[2]}, standardknife
and Xrayknife, the assistant robot Tou developed by Univeritat Politecnica
de Catalunya^{[3]}, CT, etc., their servo system are all of CNClike
serial mechanisms.
Recently, some efforts have been directed towards the investigation of alternative manipulator designs based on the concepts of closed kinematic chain due to the following advantages as compared with the traditional open kinematic chain manipulators: more rigidity, higher accuracy due to lack of cantilever like structure and greater dynamic capabilities.
In this study, a novel 3RPS parallel manipulator (Fig. 1), is introduced to an ultrasound therapy machine called ultrasoundknife, which is undergoing at the Huazhong University of Science and Technology right now. The machine serves the applications of HIFU (HighIntensityFocusedUltrasound). The manipulation approach we analyzed is based on an inparallel actuated tripodlike manipulator which has two degrees of orientation freedom and one degree of translational freedom. From the mechanical point of view, such 3 DOF manipulator cannot satisfy the basic clinical requirements. The curing knife must be well immobilized and precisely aligned with the tumor to give full advantage of the dose potential of this facility. To achieve this, the mobile platform is supposed to have two rotational degrees of freedom and three translational degrees of freedom. Since 3RPS manipulator has only 3 DOF, a traditional 2DOF XY worktable may be placed over the mobile platform to fulfill the requirement of DOFs. It can be seen that such hybrid mechanism has 5 DOF and is much more simpler and lighter than the traditional 5 DOF serial manipulators. Thus the cost of such machine could be reduced a lot.
Design requirements and scheme of the machine: The general schematic of the machine is shown in Fig. 2 .The patient is placed on a couch attached to the XY worktableand kept still all the time during the whole curing process. A 3RPS manipulator is installed below the patient.
The machine is of a hybrid 5axis mechanism with three spatial degrees of freedom and two rotational degrees of freedom. It’s worth mentioning that the two translational DOFs in the XY plane, provided by the worktable, are mainly to overcome the shortcoming of the 3RPS manipulator; although it also promises a good working space. For HIFU applications, the required absolute position accuracy must be no more than ±0.5 mm. Large error may be dangerous to the patient.

Fig. 1: 
3RPS parallel manipulator 
Accuracy is important because the treatment plans are usually made beforehand according to the patient coordination.
The design requirements of this machine are similar to those of any traditional radiationcuring machine^{[2]}. The only thing we need to mention is the oblique movement of mobile platform around Zaxis, which carries the sound transducers that varies from 0 to 20°. The maximum working space height is about 200 mm.
It’s obvious that the performance of this machine directly depends on the performance of the parallel mechanism below the worktable. If great execution error exists, calibration or compensation methods must be carried out to make up the difference.
For convenience, we put forward the nominal kinematic parameter values of 3RPS
manipulator under consideration: R = 200 mm, r = 140 mm, minimal working
height of the mobile platform is 771 mm and maximal working height is 971 mm.
As for the radiator focus, it is assumed to be 244 mm right now. The practical
machine may be a little different from the values above, but almost the same.
Literature survey: As a novel parallel manipulator, 3RPS mechanism and its application has already been studied by a few researchers. Lee and Shah^{[4,5]} put forward the inverse, forward kinematic and dynamic analysis of 3RPS manipulators. Physical constraints imposed by the limits of ball joints and the link lengths have also been discussed. At the end, Lee and Shah^{[4]} concluded that ‘typical applications are automated assembly, contour machiningand material handling’.
Their contribution makes 3RPS manipulators available in practical applications^{[68]}. Dujun et al.^{[8]} analyse the end effector’s position error with respect to components of 3RPS mechanism by matrix differentiation method. Qualitative conclusions have been reached, listed as follows: (1) Large position error will be caused with structural parameter deviation exists; (2) Position error depends on the critical position of mobile platform. The closer the mobile platform approaches the edge of the workspace, the larger the position error will be observed; (3) Largest component of position error is along the Z axis; (4) The component of position error along the +X axis is larger than that along theX axis when the mobile platform is moving along X axis; (5) The components of position error along the X axis and Z axis is symmetric with respect to the Y axis, while the position error along the Y axis is linear to the geometric deviations of three spherical joints, if such deviations are all the same. After all, Dujun et al.^{[8]} don’t prompt a calibration method or compensation method so that the actual kinematic parameters could be estimated and then be used to modify the kinematic model residing in the manipulator controller.
There are quite a few references concerning parallel manipulator calibration. While few of them mention hybrid manipulator calibration. Fortunately, for hybrid manipulators, the basis is all the same.

Fig. 2: 
Scheme of the Ultrasoundknife machine 
Wang and Masory^{[9]} and Masory et al.^{[10]} and others discussed the accuracy issue of a Stewart platform. At the end, they state that the platform pose error has been reduced at least one order of magnitude, while the proposed algorithm is computationally expensive and has to be performed on a powerful computer for practical implementation. Vischer and Clavel^{[11]} and Masory et al.^{[10]} take ‘implicit calibration’ method as the standard method to calibrate parallel robots. Implicit calibration means that the residuals for the identification process are defined based on the implicit formulation of the closure equations only.
Alternatively, Zhuang and Liu^{[12]}, show that by installing 2 or more redundant sensors at passive joints on the Stewart platform, the system is able to perform selfcalibration. It is worth mentioning because no extra measuring devices are needed, thus seemed to be more valuable in a practical circumstance. The drawback of this method is that, as stated by the authors, computation is more intensive than some existing methods. Since calibration is usually carried out every year or at certain period of time. Such computation burden is acceptable.
Iurascu and Park^{[13]} study the existing calibration methodsand conclude that the majority of these approaches are formulated in terms of a linearized set of kinematic equations, typically either the inverse kinematics or the loop closure equations. They put forward geometric method as a unified, highlevel framework for closed chain kinematic calibration. Such approach rests on viewing the configuration space of the kinematic chain as an embedded submanifold of an ambient manifoldand formulating error measurements based on the natural metric in this ambient manifold. Unlike many of the previous approaches to calibration that begin with a suitable linearization of the relevant kinematic equations, followed by a recasting into a leastsquares problem, in their approach they directly minimize the geometrically defined objective function.
In this paper, we would like to borrow the idea of selfcalibration and geometric method to carry out our calibration for hybrid 3RPS mechanism.
Error sources: Discussion of error sources is the basis of error modeling for any mechanism. We neglect any unrepeatable kinds of error, such as thermo deflection, backlash in motor box, etc. Such kinds of error depend on specific operating circumstancesand they are so hard to be modeled that it is very difficult or expensive to promote execution accuracy just by modifying the program inside the controller.
For those repeatable errors, we mainly consider the following 3 categories:
• 
Machining errors: Machining error exists in any machined components
that may be installed into the robot because of machining tolerance. Such
error is inevitable but can be limited within a certain range. Once machined,
it could be represented by a constant geometric variable (length, height,
etc.), if we neglect thermo effects. For 3RPS mechanism, the mobile platform
cannot be just the same as the equilateral triangle as we specified. Three
edges may be different. So are the three legs. The leg length offsets may
be different too. 
• 
Assembly errors: Assembly errors include position and angular errors,
which appear during installation. This means the installation of the base
of legs may not be at the right place and at the right angle. 
• 
Joint errors: Usually, high quality of joints is one of the basic requirements
for robot design. As joints are the key components for any robot, intricate
joint errors may be enlarged or/and accumulated, which leads to comparatively
large end effector’s position or orientation error. 
Performance analysis
Parameterization and error modeling: There are 8 components in a
3RPS mechanism: static platform, 6 linkagesand mobile platform.
Mobile platform: As for the mobile platform, the key problem is the position of the center of three spherical joints.
In order to make spherical joints work smoothly, certain kind of backlash should be promised. What’s more, installation introduces larger position errors between each other. In order to combine these errors together and decrease the number of variables, we take the distances between the joints as variables: r_{12}, r_{23}, r_{31},
where
r is nominal radius of the circumcircle of the mobile platform.
On the other hand, r_{12}, r_{23}, r_{31} only defines the shape of the triangleand doesn’t provide information about the position and orientation of the mobile platform. Three vertexes’ coordination should be added:
It’s very obvious that
and r_{12}, r_{23}, r_{31} must
satisfy following constraints, (Fig. 3).

Fig. 3: 
Topreview of 3RPS manipulator 
Static platform: The variables representing static platform are different
from the mobile platform. The positions of three revolute joints are constant
once they have been installed, however, their critical coordinates are not the
nominal ones. In fact, even if the installation of revolute joints is accurate,
assembly error between revolute joints and prismatic joints will make the actual
revolution center deviate from its nominal place. (Fig. 3,
4) nominal revolute center is A, but the actual center is
A’, which depends on both the angular installation errorβ and the
length of the linkage. In order to utilize the nominal value of revolute joints’
positions, the static platform is represented by a triangle whose vertices are:
The height of revolute joints’ center is comparatively easy to control,
thus the error along the Zaxis is very small compared with that along the X
or Y axis, since they are installed on the same plane and the joints are standard
joints. So Δ1_{z}, Δ_{2z}, Δ_{3z} could
be eliminated from the variable set, or to say that they are assumed to be zero.
The variables chosen to represent the base are simplized to be:
Linkage legs: Linkage legs are active actuators in 3RPS mechanism. The end effector’s position and orientation directly depends on the input lengths of links. Linkages are modeled as extensible lines, whose lengths are variable.
Machining errors, assembly errors and joint errors exert their influence on the links simultaneously. Machining errors are very easy to be seen, which results in the offset errors in lengths.
Prismatic joint errors could be ignored, because the prismatic joints are designed to be of ball screws, which are very precise and could be taken as ideal prismatic joints. Only linkage offset errors Δ_{l1}, Δ_{l2}, Δ_{l3} are taken into account. On the other hand, the orientation errors of these links cannot be ignored (Fig. 4). Tiny angular deviation from its nominal place may cause an enlarged position error in the spherical joint center’s position. Such errors mainly originate from the revolute joint’s orientation. In this paper, we would like to model these angular errors by three variables: α_{1}, α_{2}, α_{3} (Fig. 5) .
Link length offsets and its orientation introduces two sets of constraints, which are:
Link length offsets constraints:
Orientation constraints:
Forward Kinematic model with error: Equations above define the practical behavior of 3RPS mechanism. Of course, they are unlinear and very hard to have an analytical result. What is more, the position and orientation of the mobile platform is implied in the 9 parameters of the 3 vertices: x_{1}, y_{1}, z_{1}, x_{2}, y_{2}, z_{2}, x_{3}, y_{3}, z_{3}.
Orientation is described by the equation below:
Furthermore, we care much more about the position of the curing point or the
focus point of the transducers installed on the mobile platform. (Fig.
6), the focus point could be taken as a point attached rigidly to the mobile
platform without knowing its position or its focus length. Certainly, its nominal
position and focus length could be known beforehand. In the mobile platform
frame xoy, the focus point could be assumed to be (x_{f}, y_{f},
z_{f}, )’, where the nominal value is (0, 0, f). In this paper,
f = 244 mm.
where

is the orientation of the mobile platform; 

is the x direction of mobile framework; 

is the y direction of mobile framework; and 

is the focus position. 

Fig. 4: 
Revolute joint model 

Fig. 5: 
Focus above the mobile palatform 

Fig. 6: 
Focus above the mobile plateform 
Performance evaluation without calibration: We assume the size of each
error set and select a random value for evaluation. With enough number of tests
(or simulations), we could know the performance of 3RPS manipulator, thus that
of our hybrid mechanism is clear too.
Limitations on each error set: We have 4 sets of error variables: position
errors of 3 vertices on mobile platform r_{12}, r_{23},
r_{31}; position errors of static platform: Δ_{1x},
Δ_{1y}, Δ_{2x}, Δ_{2y}, Δ_{3x},
Δ_{3y} ; orientation and offset errors of 3 links: α_{1},
α_{2}, α_{3}, Δ_{l1}, Δ_{l2},
Δ_{l3} and position error of focus caused by installation: x_{f},
y_{f}, z_{f}. r_{12}, r_{23},
r_{31} are mainly caused by manufacturing and assembly. With
nominal distance √3r≈ 242.5 mm, assuming manufacturing precision
is of H7, the distance error will be 0.02 mm. Assembly process introduces much
more error in distance between each spherical joint. It is difficult to predict
the assembly range right now, as no prototype machine exists. Here, we assume
r_{12}, r_{23}, r_{31} random in the
range of 1mm.Δ_{1x}, Δ_{1y}, Δ_{2x},
Δ_{2y}, Δ_{3x}, Δ_{3y} are mainly assembly
errors. Nominal distance between each revolute joint is 346.5 mm. Position error
could be of the range 0.02 mm. The practical revolute center may not coincide
with the physical revolute joint center. Such error may be enlarged and assumed
to be as great as 2 mm. α_{1}, α_{2}, α_{3},
is controllable and could be limited inside the range of 0.05°. Δ_{l1},
Δ_{l2}, Δ_{l3} are along the link lengths, which fall
in the range of 1mm. As for the x_{f}, y_{f}, z_{f},
their error limitation is 5 mm.
100 random error sets have been studied. In this way, we could observe the position error bound of the focus point.
The simulation results coincide with Dejun et al.^{[8]}. The maximum focus position error is observed to be almost as large as 12 mm at the height of 971 mm. That is too much from our expectation. It shows that the position accuracy of 3RPS manipulator in fact is practically of the same order of traditional serial manipulators, although theoretically speaking it should be better. Kinematic calibration has to be carried out to find the suitable parameter values to approximate machine’s practical performance. So the first thing we should do is to find the machine’s practical performance. Certain kind of measurement has to be carried out.
Calibration Measurement
Measurement setup: For any measurement process, we have to carry
out measurement setup first. Since the measurement method corresponds to the
calibration method directly, we would like to describe in detail how to measure
and calibrate this hybrid manipulator.
First, select a reference plane, which can be considered as parallel to the base. This plane is called the measuring plane. Corresponding measuring plane height could be tested beforehand according to its nominal position. This height needs not to be very precise. In fact, the nominal height is quite a good choice. We will show that the final calibration result has very little relation with the precision of measuring plane height.
What we care about most is the position of the focus point, but the focus point is a virtual point. Without any help, it’s almost impossible to test the critical position of the focus point. Here, we employ a ballend probe to represent the right place of the focus point. This is a good idea because the focus point itself could actually be described as a box of 2x2x8 mm^{3} in case of touching point area being smaller than 2x2 mm^{3}. During the ultrasoundtransducer installation, we can adjust the focus point just at the ball end (Fig. 7).
Practical measurement procedure: When three leg length inputs are all the same with the end of the probe ball touching the measuring plane, make a mark at the plane. This point could be taken as the origin of the working table.

Fig. 8: 
Simulate measurement process 
Now, with the ball end touching the measuring plane, move the mobile platform. Of course, the touch point would slide on the plane. We can move the working table to ensure that the mark on the measuring plane coincides with the ball end. What we need is the XY offset value of the working table and corresponding three leg lengths. Repeat this operation, we can have a group of data sets, which will be used for parameter identification processing below.
Numerical measurement simulation: To verify this method, we developed a set of program to simulate the measurement process (Fig. 8).
Hybrid inverse/forward kinematics is a little different from that of pure parallel mechanism. The height of the measuring plane should be taken into account beside the mobile platform orientation.
Forward kinematics of this hybrid mechanism has already been put forward. The only difference is the measuring noise that should be taken into account. In this paper, we assume the measuring noise is no more than 0.5 mm.
For inverse kinematics, one of the inputs is the orientation of the mobile
platform. Obviously, the calculated orientation of mobile platform should be
equivalent to that of the input. So we have the following equations to describe
the inverse kinematics of our hybrid mechanism:
There’re only two independent equations in
yet. The last independent equation is of the focus height.
where

is the orientation input, 

is the calculated orientation, 

is the x direction of mobile framework, while 

is the y direction. 
Coordinate origin: Numeric measurement simulation output is the key step in our study of hybrid mechanism calibration. Its validity, reliabilityand stability are the crucial factors. There’re two things worth mentioning: the first one is the origin of mobile framework; the second one is the origin of the working table.
As for the origin of mobile framework, when we analyze 3RPS manipulator’s performance, the origin is the geometric center of the mobile triangle whether it’s an ideal equilateral triangle or not. Because the parameter deviations of mobile platform itself have already been generalized by three edge lengths’ deviationsand have nothing to do with the mobile framework. What’s more, it is still worth being pointed out that the focus position is described in the mobile framework. The focus position deviation is constant but the parameters describing deviation depends directly on the origin of the mobile framework. We choose geometric center as the origin just because it’s very easy to calculate and has no negative effect on the final result.
The origin of the working table is another problem. Working table origin corresponds to the 3RPS manipulator origin in xy direction. We set the three leg lengths to have the same value, thus the mobile platform can only move parallelly up and down. If all the parameters are of their nominal value, the center of the mobile platform will not move in neither x direction nor y directionand this position (in xy plane) is just the origin of the hybrid mechanism. Of course, this point does not need to be reached at the same time as 2D working table reaches its mechanical origin. If the xy position of this point should be decided, we only need to record a certain point on the working table to coincide with this position. And all the transducer outputs at this time could be taken as corresponding to the origin point of the hybrid mechanism. Unfortunately, all the practical parameters are not of their nominal value and there’s no way to test such xy position. Anyway, at measuring stage, we still follow this idea to position the working table origin. In fact, the working table offset has very little effect on the final result in case we follow the same mobile framework. Just like CNC machine’s mechanical origin, even if we don’t know the critical position of machine origin, we still can set up a working coordination for operation. What’s important is the movement of the working table in xy direction, not merely its origin point.
Parameter identification: Now, we have a group of measurement data, which contains the 3 leg length inputs and corresponding working table offsets, which stand for focus offsets in xy direction. The aim of calibration is to find an appropriate set of parameters to satisfy the relationship between the input values and the endeffector outputs, according to a certain set of functions (forward kinematics or inverse kinematics or loop closure equations). As Zhuang and liu^{[12]} mentioned, this is a classical curvefitting or surfacefitting problem. Here, we can assume the 3 leg lengths are the horizontal coordinates and the 2 working table offset values are the vertical coordinates. Our aim is to find the right curve that suits this data set (Fig. 9).
First of all, we should give an objective function, which evaluates the residuals between the idea curve and calculated curv.
Suppose measurement data are written in such format:

Fig. 9: 
Parameter Identification 

Fig. 10: 
Precision analysis at height of measuring plane compared with
practical focus position
(A) Nominal parameters
(B) Calibrated parameters with measuring noise of 0.5 mm
(C) alibrated parameters with measuring noise of 0.2 mm
(D) Calibrated parameters without measuring noise 
Table 1: 
Practical and Identified values for the mechanism parameters 


Fig. 11: 
Precision analysis at the height of highest plane compared
with practical focus position Calibration results acquired at the height
of measuring plane
(A) Nominal parameters
(B) Calibrated parameters with measuring noise of 0.5 mm
(C) Calibrated parameters with measuring noise of 0.2 mm
(D) Calibrated parameters without measuring noise
(The number of checkpoints/measurements is 20.)
(Mobile platform ‘rotation’ around Z axis ranges from 0° to
360° with the interval of 1.1°.)
(The angle between Z axis and ideal focus axial ranges from 0°to 20°
with the interval of 5°.) 

Fig. 12: 
Effect of measurement number and measurement noise on calibration
result 
The objective function is the measurement of the accumulate residual between the ideal curve and the calculated curve. There're three components in the residual including the residual of focus height in addition to those of x and y offsets.
With the help of forward kinematics, we can have the focus position with respect to certain set of parameter values. The x, y focus offset could be compared with the measurement valuesand the z focus offset should always equal measuring plane height. The differences between them are taken as the residual.
For curvefitting problem, the standard solution is to minimize the residuals in a leastsquare sense. It can be solved with the LevenbergMarquardt (LM) algorithm that is implemented in the ‘optimization toolbox’ of Matlab. The LMalgorithm is a mixture between the GaussNewton and the steepest descent algorithm aiming at keeping the quadratic convergence rate of GaussNewton algorithm by avoiding the problem of rank deficiency of the identification Jacobian. As for the initial values for iterative LMalgorithm, the nominal parameters of the hybrid mechanism could be taken and all the noise is set to zero.
Assume the parameter deviation values are the same as we have used in error
analysis, within 27 iterations and about 4 minutes calculation time (PIII 450
MHz), the LMmethod has identified the following parameter set (Table
1) with 20 measurements. Of course, the calculation time is the function
of number of measurements. The larger the number of measurements, the longer
the time needed to complete the identification.
Performance evaluation with calibration: We want to make sure how good our calibration result is.
There’re are three things to be mentioned: the first is that we want to learn how much has been improved in focus point position according to our calibration method; the second is that as measurements are always carried out at the same height, thus we want to know the performance of such result at the other heights; the third is about the appropriate number of measurements.
1. For the first problem, we do the verification at the measuring height, but with much more measurement points and different size of measurement noise. Both the focus position error based on calibrated parameters and the focus position error based on nominal parameters are compared with that based on practical ones. We can see clearly that the measurement noise has an obvious effect on the quality of the calibration result. If no measurement noise exists, the calibration result will be accurate. In a practical circumstance, the measurement could not be larger than 0.2 mm for an experienced operator (Fig. 10).
2. Another question is that: since the calibration is carried out only on the
measuring plane, is it still working for other planes at different heights?
Now, we can look into the most extreme cases, which are the highest and lowest
plane. In our study, the lowest plane is the measuring plane. So we only care
about the highest one whose height is 200 mm greater. Following the same route
and the same calibration results, we have: 50 calibrated values with the measuring
noise of 0.5 mm, 0.3 mm.and 0.2 mm have been tried (Fig. 11).
A. The greatest focus position error takes place at the largest
angle between Z axis and focus with the same mobile platform center height.
In our application, it’s on the 20° rotation along Z axis.
B. The focus position precision has been promoted almost 10 times better
after calibration. The maximum focus position error is observed to be about
0.65mm, no more than 1 mm in 20timesmeasurementcalibration.
C. If no measurement noise exist, the calibrated result will be just the
same as practical ones. And no focus position error could be observed.
3. Another problem is about the measurement number. The calibration result has something to do with the measurement noise. From the theoretical point of view, the true value could be reached if we can do infinite number of measurement. All the measurements noise, could be eliminated by finding the mean value. But considering the practical implementation difficulty and numerical performance, we find that the calibration result has not been improved when measurement exceeds certain number (Fig. 12). Horizontal coordination stands for the number of checkpoints, while the vertical coordination stands for maximum and minimum focus position error observed. We start from 5 points. Measurement noise is of magnitude 0.2 mm.
50 times of measurement seems to be a good choice. The focus point position error is observed to be 0.2748±0.089 mm on the highest planeand 0.1945±0.018 mm on the lowest plane in terms of 0.2 mm measurement error and 20 times of simulation.
The performance evaluation and kinematic calibration has been carried out. Numeric simulation shows that the practical performance of our hybrid manipulator is not as good as theoretical prediction because of machining errors, assembly errors, joint errors and other kind of random errors.
To close the difference between the theoretical prediction and the practical performance, we develop a kind of calibration method to enhance the performance of our hybrid manipulator. The numeric performance is reasonably good and within 50 sets of measurements, we could have a set of parameter values, which approximates the practical mechanism. The focus position error after calibration is mainly decided by the measurement number. The more measurements we have, the better the calibration model we approach. However, the computational expense will increase, as the number of measurement data increases.