Calibration procedure and biomechanical validation of an universal six degree-of-freedom robotic system for hip joint testing

Abstract

Purpose

Among various test methods for different human joints, the use of robot systems has attracted major interest and inherits the potential to become a gold standard in biomechanical testing in the future. A key issue associated with those robot-based platforms is the accurate definition of parameters, e.g., tool center point (TCP), length of tool or anatomical trajectories of movements. These must be precisely correlated to the physiological parameters of the examined joint and its corresponding bones. Exemplified for the human hip joint, we are creating an accurate calibration procedure for a universal testing platform by using a six degree-of-freedom (6 DOF) robot and optical tracking system for recognition of anatomical movements of the bone samples.

Methods

A six degree-of-freedom robot (TX 200, Stäubli) has been installed and configured. The physiological range of motion of the hip joint composed of a femur and a hemipelvis was recorded with an optical 3D movement and deformation analysis system (ARAMIS, GOM GmbH). The recorded measurements were processed by automatic transformation procedure (created in Delphi software) and evaluated in 3D CAD system.

Results

The physiological ranges of motion were reproduced for all degrees of freedom with the six degree-of-freedom robot in adequate accuracy. With the establishment of a special calibration procedure by using a combination of different coordinate systems, we were able to achieve a standard deviation of the TCP depending of the axis between 0.3 and 0.9 mm and for the length of tool between + 0.67 and − 0.40 mm (3D CAD processing) resp. + 0.72 mm to − 0.13 mm (Delphi transformation). The accuracy between the manual and robotic movement of the hip shows an average deviation between − 0.36 and + 3.44 mm for the points on the movement trajectories.

Conclusion

A six degree-of-freedom robot is appropriate to reproduce the physiological range of motion of the hip joint. The described calibration procedure is universal and can be used for hip joint biomechanical tests allowing to apply clinically relevant forces and investigate testing stability of reconstructive osteosynthesis implant/endoprosthetic fixations, regardless of the length of the femur, size of the femoral head and acetabulum or whether the entire pelvis or only the hemipelvis will be used.

Introduction

Due to the upright posture and bipedal gait that have resulted during the human evolution process, the development of osteoarthritis in the hip joint can be considered as a causative consequence and is nowadays generally considered to be a one of the most frequent civilization diseases [1, 2]. In Germany, nearly 18% of adults aged 18 years and older report to suffer from osteoarthritis in the last 12 months, with the prevalence being higher in women (21.8%) than in men (13.9%). The percentage of osteoarthritis patients increases significantly with age, among those aged 65 and older, almost half of the female population (48.1%) and a third of men (31.2%) are affected [3]. In order to develop new treatment strategies, tailored to the individual underlying pathology, a profound understanding of the biomechanics is essential.

In the literature, there are various approaches studying biomechanics and kinematics of the hip. Bergmann et al. as well as other authors have impressively shown that the dynamic loads carried by the human hip joint exceed the weight of the body almost by four times during walking and may even reach eightfold increased values during jumping or stumbling [4, 5]. Some authors reproduce physiological movements and force loads in the joint [6,7,8], while others focus on researching different implant materials, for example, to determine the corrosion and emergence of metal particles released by friction of implant parts [9, 10]. However, there is still no consistent methodology in the literature for biomechanical testing of joints, especially those associated to the hip joint.

The first biomechanical simulators were developed in the 50’s of the last century [11]. Until their appearance, various types of test machines (usually standard strength testing machines) were utilized using standard test methods [12,13,14,15,16]. For the improvement of the reproducibility and documentation of the precision, some authors reported the use of different 3D measuring devices allowing the reconstruction of simple movements in the hip joint [17,18,19,20]. Others used more complex measurement systems equipped with additional elements to emulate or measure muscle response [21, 22]. The common challenge of these solutions is the limited ability to reproduce complex joint movements, the restricted possibility of active control and the difficulty of reconfiguration of the device for other types of joints, movements or loading forces [17, 20, 23,24,25].

As van Arkel et al. described in their study [16], many research laboratories have developed new methods to study hip joint biomechanics based on various techniques, such as digital image correlation [26, 27], real-time contact-pressure measurement [28, 29], optical tracking motion analysis [30, 31], 3D digital reconstructions combining CT scans and motion tracking [32, 33], custom built rigs in servo-hydraulic actuators/materials testing devices [16, 26, 34, 35] or 6 DOF robots with actuators [36, 37]. Interesting advanced systems were also developed by Colbrunn et al. [36] and Al-Haifi et al. [38]. They constructed a three-dimensional testing system based on a six degrees of freedom parallel robot platform (Steward platform). Colbrunn et al. added a rotary stage to their system as a seventh degree of freedom [36]. The use of six degrees of freedom robots in biomechanical studies has a relatively short history [39, 40]. Due to their kinematic abilities, computer control and high positioning repeatability, they are very useful for precisely reproducing anatomical movements in joints [39, 41,42,43] and therefore are becoming the new gold standard [39, 44]. There are already studies that have been answered clinical relevant questions by using a 6 DOF robots in the knee joint [45, 46].

While most of these biomechanical testing devices focused on limited range of motion under constant cyclic loads [26, 34, 47], only few of them allowed variable loads while simultaneously reproducing more complex physiological movements [16, 18]. This establishes a strong demand for a uniform and standardized test setup (universal joint kinemator), in order to study complex biomechanical testing approaches for different motion and kinematics of the hip joint. For precise control of the robot platforms, a clear definition of a coordinate system in relation to the robot and the joint is mandatory. This is a precondition for exact 3D analysis and measurement as well as transmission of variable forces [42].

The purpose of this study was therefore to establish a standardized calibration procedure for a hip joint testing platform, which is also transferable to other joints (s.c. ‘universal joint kinemator’). For this purpose, we have used quantitatively assessed parameters, i.e., range of motion (trajectory of the motion), joint rotation center (defined for robots as tool center point) and length of bone specimen. Precise determination of these parameters is an essential precondition for development and setup of a 6 DOF robots, potentially allowing to perform subsequent biomechanical experiments under reproducible, dynamic real-time conditions with physiologically measured loadings.

Methods

Six degree-of-freedom robotic system for testing hip joint

For the establishment of an ‘universal joint kinemator,’ a six degree-of-freedom robotic system (TX 200, Stäubli, Switzerland) was installed (Fig. 1), as it has been used for investigations of the knee and hip [41, 42]. With this robot, a maximum dynamic load of 1500 N and maximum static load of 9000 N can be applied, with a maximum range of 219.4 cm and a position repeatability (ISO 9283) of ± 0.06 mm [48]. In front of the robot, bidirectional (2 DOF—movements possible in X and Y direction) rail fastenings and a steel mounting plate with grid holes were installed to allow maximum variability of the object attachment. A force torque (FT) sensor (FTN-Omega160-IP-R15-NETB SI-1000–120, Schunk GmbH & Co. KG) with ± 1000 N maximum force measurement range in X and Y direction, ± 2500 N in Z direction and torque ± 120 Nm in Rx, Ry and Rz direction (measurement resolution of forces: 0.25 N, and torques: 0.025 Nmm for Rx, Ry and 0.015 Nmm for Rz direction) was additionally mounted on the robot flange to measure loads (Fig. 1).

Fig. 1

General view of the testing stand and the main elements with coordinate systems (left), presentation of the ‘length of tool’ (LOT) in relation to the ‘center of robot flange’ and the ‘tool center point’ (TCP) (right)

Samples

Artificial hemipelvis-femur-bone samples with normal anatomical proportions (right hip, foam cortical shell, SKU: 1172, medium size, Sawbones, Malmö, Sweden) were used for the testing of the robot movement. The distal femurs were removed 5 cm from the joint line with an oscillating saw to ensure a secure fixation in the holding adapter and were placed into bone cement (Dental Plaster Typ 4, Excalibur, water/plaster ratio 22:100, Siladent, Dr. Böhme & Schöps GmbH) in an exact position according to the mechanical axis of the femur in relation to the center of rotation of the hip joint. In terms of alignment, we have been oriented to the International Society of Biomechanics (ISB) recommendations for definition of joint coordinate systems [49]. The hemipelvis-femur was then installed on the robot testing setup. The distal femur (cemented in holder) was attached to the force torque sensor and the hemipelvis to the rail fastening system (Fig. 1). Subsequently, the fixation of the joint was removed, so that the joint could be moved manually for calibration of future movements by the robot.

A hemipelvis with native hip joint and complete femur from a male (81 years) specimen was used for the final investigation of the accuracy between manual and robotic movements (positive evaluation for the use of human specimens for biomechanical testing: EK 51,012,019).

Photo-optical and kinematic calibration

Photo-optical measurements were carried out by a high-resolution 3D-stereo photo-optical system (ARAMIS ASX, GOM GmbH, Braunschweig, Germany) with a frequency up to 10 Hz, a measurement volume of 1100 × 800 × 500 mm and an accuracy of 0.001 mm to evaluate the axes and finally the precision of the movements. The system was aligned and positioned at 90° to the sagittal plane (y-axis) of the hemipelvis-femur according to the recommended orientation conditions [49].

Establishment of the coordinate systems

The spatial understanding of these coordinate systems is essential, since a conversion of the coordinates between each other is necessary for the transformation of the motion between the investigated specimen and the robot. The structure of the ‘universal joint kinemator’ has four main coordinate systems defined in the Cartesian space related to the basic elements (Fig. 1):

• Robot world coordinate system—associated with the robot base,

• Robot flange coordinate system—associated with the robot flange (sixth axis),

• ARAMIS coordinate system—associated with the optical ARAMIS measuring system,

• Bone sample coordinate system—associated with the testing bone sample.

The Robot World coordinate system of the TX 200 Stäubli is located on the first moveable element (rotation between the robot base and the robot shoulder) of the robot (Fig. 1) [48]. Additionally connected with the Robot World coordinate system is the so-called tool center point (TCP) specifying the tool tip position (point around which the robot arm moves) of the robot tool and corresponds to the joint rotation center (JRC). ‘Tool Length’ is related to the distance between the TCP point (JRC) and center of robot flange (Fig. 1), which corresponds to the most distal movement point of the robot. Practice shows that most often the TCP point does not coincide with the Z-axis of robot flange but is spatially shifted in relation to it. The ARAMIS coordinate system is associated with a set of three groups of optical markers (two on the left and one on the right ARAMIS origin element) that form the reference system (Fig. 1). The bone sample coordinate system is related to the bone sample mounted to the bidirectional rail system (Fig. 1).

Determination of motion parameters of manual movements (stage I)

The aim of stage I is to determinate motion parameters of the specimen plane (bone sample coordinate system) by using standardized hemipelvis-femur-bone samples (Sawbones, Malmö, Sweden). For this purpose, the trajectories of physiological anatomical movements (range of motion—ROM) were evaluated and used to define the TCP. The femur was moved manually without use of external forces (preloads) and the range of motion was recorded with the ARAMIS system (Fig. 2). During the test, four types of physiological movements were measured: flexion–extension, abduction–adduction, external and circular movement of femur. The collected data was obtained from the measurement system in data sets of X-, Y-, Z coordinates and angles of rotations α, β and γ. To obtain correct data, a specially configured ‘calibration adapter’ (aluminum element) with 3D/6 DOF definition and precisely marked points (seven markers of the ARAMIS system) was manufactured (Fig. 4, bottom) and mounted on the top of the femur holder (Fig. 2).

Fig. 2

Calibration Stage I—Evaluation of the manual movements of hemipelvis femur bone sample (Sawbone) with mounted calibration adapter: flexion–extension (1), (1ꞌ), abduction–adduction (2), (2ꞌ), visualization of femur trajectory for all movements including circular movements (3), (3ꞌ) images recorded by ARAMIS system

To evaluate the accuracy of the setup for the bone sample coordinate system, the tests were repeated in six different positions with an offset of 50 mm from the initial position, while maintaining the full range of motion of the robot (Fig. 3). This results in a change in the position and tilt angle of the Aramis system, while maintaining the correct measurement range of the system in relation to the robot and the examined bone sample.

Fig. 3

The varying settings of the robot test platform during the experiment with six different mounting positions of the bone sample (base plate with holes raster 50 × 50 mm)

Determination of motion parameters of the robot (stage II)

For second stage calibration process, the ‘Calibration adapter’ (see stage I) was mounted to the robot flange (Fig. 4, (2)). As a primary step, the robot's YZ plane was defined and measured, by moving the robot to five different localizations in space (points: P1–P5) documented optically by ARAMIS system (Fig. 4, (1)). As a variable different rotations of the robot flange around, these points were used, while all points were localized on constant X value (rotation performed sequentially at individual points P1-P11 in a sequence of rotations of the robot flange around the X-, Y- and Z-axis).

Fig. 4

Stage II calibration of the robot: calibration points P1–P5 localized on plane YZ (1), calibration points P6–P11 localized on plane XZ (1ꞌ), examples of robot flange manipulations—two different positions (2) and (2ꞌ), with mounted calibration adapter (images recorded by ARAMIS system)

In a secondary step, the XZ plane was defined as a vertical axis localized at the intersection of the YZ axes and through the point P1 (from step 1). The procedure was repeated analogously to the previous one generating the points P6–P11 (Fig. 4, (1’)).

Evaluation of tool center point and length of tool (stage III)

The main objective of stage III was to obtain fundamental parameters for the definition of the position of the tool center point (TCP) and the length of tool (LOT) in the robot world coordinate system and therefore consequently in the coordinate system of the robot control panel. For evaluation of the calibration process, the following parameters were determined:

• Distance between calibration points (P1–P5 in YZ and P6–P11 in XZ robot plane) in comparison with the theoretical set of coordinates in the robot control unit, defining the error of measuring and reconstructing the position of the robot flange during calibration procedure

• Reconstructed Robot World origin (point xR, yR, zR = 0, 0, 0),

• Maximum deviation of the center of the Robot World origin,

• Main localization of the Robot World origin (defined as an average value of all calibrations points) related to the ARAMIS origin.

Establishment of a numerical transformation for robot controlling (Delphi transformation)

The data obtained from the measurements using the ARAMIS optical system were automatically processed (by geometric transformations) in an automated transformation procedure developed for this purpose (in the Delphi computer programming language). The automatic transformation procedure allows the transformation of the input data from the calibration to obtain the following data: conversion of ARAMIS coordinates into Robot world coordinates, length of tool (length of bone sample, LOT), movements trajectory of the bone sample around the joint, TCP defined in Robot world coordinates and the error of the obtained results (Figs. 8 and 10).

For the analysis of the calibration measurement, digital measurement data are given for the different robot positions i = 1,…,n, in robot world coordinates $x_{R,i} ;y_{R,i} ;z_{R,i}$ as well as in ARAMIS coordinates $x_{A,i} ;y_{A,i} ;z_{A,i}$ (Fig. 5).

Fig. 5

Bases of ARAMIS—and robot world coordinate system

We aimed to determine a transformation matrix to convert the coordinates in any robot position from one coordinate system into the other one. This transformation was carried out in two steps, separated into translation shift and rotation shift of both coordinate systems. The translation shift of both coordinate systems is performed quite simply by calculation of the coordinate differences of both measured coordinate sets with the first data set i = 1 of the calibration measurement, where all coordinates are shifted to a Zero-point. With the robot control values of the first date set of the calibration measurement $x_{R,1} ;y_{R,1} ;z_{R,1}$, it follows the zeroed robot coordinates and the zeroed ARAMIS coordinates.

$$\begin{array}{*{20}c} {\tilde{x}_{R,i} = x_{R,i} - x_{R,1} } & {\tilde{x}_{A,i} = x_{A,i} - x_{A,1} } \\ \end{array}$$

$$\begin{array}{*{20}c} {\tilde{y}_{R,i} = y_{R,i} - y_{R,1} } & {\tilde{y}_{A,i} = y_{A,i} - y_{A,1} } \\ \end{array}$$

$$\begin{array}{*{20}c} {\tilde{z}_{R,i} = z_{R,i} - z_{R,1} } & {\tilde{z}_{A,i} = z_{A,i} - z_{A,1} } \\ \end{array}$$

The rotation matrix is calculated via matrix operations. For that purpose, the available measurement data are collected in two 3 × n-matrices, where n is the total number of the calibration measurement data sets.

Both zeroed measurement data matrices are

$$\begin{array}{*{20}c} {\underline {R} = \left[ {\begin{array}{*{20}c} {\tilde{x}_{R,1} } & {\tilde{x}_{R,2} } & \cdots & {\tilde{x}_{R,n} } \\ {\tilde{y}_{R,1} } & {\tilde{y}_{R,2} } & \cdots & {\tilde{y}_{R,n} } \\ {\tilde{z}_{R,1} } & {\tilde{z}_{R,2} } & \cdots & {\tilde{z}_{R,n} } \\ \end{array} } \right]} & {\underline {A} = \left[ {\begin{array}{*{20}c} {\tilde{x}_{A,1} } & {\tilde{x}_{A,2} } & \cdots & {\tilde{x}_{A,n} } \\ {\tilde{y}_{A,1} } & {\tilde{y}_{A,2} } & \cdots & {\tilde{y}_{A,n} } \\ {\tilde{z}_{A,1} } & {\tilde{z}_{A,2} } & \cdots & {\tilde{z}_{A,n} } \\ \end{array} } \right]} \\ \end{array}$$

The relation between both measurement data sets is defined by the—initially unknown—transformation matrix $\underline {T}_{A \to R}$, where index A → R indicates the transformation direction from ARAMIS coordinates into robot coordinates:

$$\underline {R} = { }\underline {T}_{A \to R} \cdot \underline {A}$$

For solution of this equation to the transformation matrix, at first the transposed matrix $\underline {A}^{T}$ needs to be multiplied from the right side:

$$\underline {R} \cdot \underline {A}^{T} = { }\underline {T}_{A \to R} \cdot \underline {A} \cdot \underline {A}^{T}$$

Because the resulting 3 × 3-matrix $\underline {A} \cdot \underline {A}^{T}$ is general invertible, it can be continued with

$$\underline {R} \cdot \underline {A}^{T} \cdot \left( {\underline {A} \cdot \underline {A}^{T} } \right)^{ - 1} = { }\underline {T}_{A \to R} \cdot \underline {A} \cdot \underline {A}^{T} \cdot \left( {\underline {A} \cdot \underline {A}^{T} } \right)^{ - 1}$$

Afterward because of $\underline {A} \cdot \underline {A}^{T} \cdot \left( {\underline {A} \cdot \underline {A}^{T} } \right)^{ - 1} = \underline {E}$, the equation for the calculation of the transformation matrix from ARAMIS coordinates into robot coordinates follows as

$$\underline {T}_{A \to R} = \underline {R} \cdot \underline {A}^{T} \cdot \left( {\underline {A} \cdot \underline {A}^{T} } \right)^{ - 1}$$

The transformation matrix in the opposite direction from robot coordinates into ARAMIS coordinates is the inverse matrix

$$\underline {T}_{R \to A} = \left( {\underline {T}_{A \to R} } \right)^{ - 1}$$

Determination of sphere equation of the tested hip joint

The hip joint measurements delivered measurement data sets $x_{i} ;y_{i} ;z_{i}$ of possible positions of the center point of the robot adapter plate. If the trajectory of the measured data is similar to a spheroidal shape, only the parameters of the related sphere—middle point $x_{m} ;y_{m} ;z_{m}$ and radius $R$ of the sphere (Fig. 6) need to be given to the control system of the robot.

Fig. 6

Sphere parameters x_m, y_m, z_m, R and one hip joint measurement data set ${x}_{i},{y}_{i},{z}_{i}$

If the trajectory of the measured data is not similar to a spheroidal shape, only the trajectory of an ‘allowed coordinate trace’ can be given to the control system of the robot. The proof of the likelihood of the hip joint measurement data to a spheroidal shape is carried out by the following procedure, which can be performed using ARAMIS coordinates as well as in robot coordinates since the conversion between both coordinate system is possible with the procedure as described above. In general, a sphere can be described by equation

$$\left( {x - x_{m} } \right)^{2} + \left( {y - y_{m} } \right)^{2} + \left( {z - z_{m} } \right)^{2} = R^{2}$$

which can be expanded to

$x^{2} - 2 \cdot x \cdot x_{m} + x_{m}^{2} + y^{2} - 2 \cdot y \cdot y_{m} + y_{m}^{2} + z^{2} - 2 \cdot z \cdot z_{m} + z_{m}^{2} = R^{2}$.

The sorting of all terms with the parameters of the sphere $x_{m}$ $y_{m}$ $z_{m}$ and R to the left side of the equation and all other terms to the right side of the equation leads to

$x_{m}^{2} + y_{m}^{2} + z_{m}^{2} - R^{2} - 2 \cdot x_{m} \cdot x - 2 \cdot y_{m} \cdot y - 2 \cdot z_{m} \cdot z = - \left( {x^{2} + y^{2} + z^{2} } \right)$.

The factors related to 1, x, y, and z on the left side of the equation can be substituted to

$$p_{1} = x_{m}^{2} + y_{m}^{2} + z_{m}^{2} - R^{2}$$

$$p_{2} = - 2 \cdot x_{m}$$

$$p_{3} = - 2 \cdot y_{m}$$

$p_{4} = - 2 \cdot z_{m}$.

With these substitutions, it follows

$$p_{1} + p_{2} \cdot x + p_{3} \cdot y + p_{4} \cdot z = - \left( {x^{2} + y^{2} + z^{2} } \right)$$

The four unknown factors $p_{1}$ $p_{2}$ $p_{3}$ and $p_{4}$ require at least four equations. These four equations can be obtained by a sufficient repetition number of the hip joint measurement. If measurement data sets of more than four hip joint measurements are available, these measurements can be written as matrices

$$\left[ {\begin{array}{*{20}c} 1 & {x_{1} } & {y_{1} } & {z_{1} } \\ 1 & {x_{2} } & {y_{2} } & {z_{2} } \\ \vdots & \vdots & \vdots & \vdots \\ 1 & {x_{i} } & {y_{i} } & {z_{i} } \\ \vdots & \vdots & \vdots & \vdots \\ 1 & {x_{n} } & {y_{n} } & {z_{n} } \\ \end{array} } \right] \cdot \left[ {\begin{array}{*{20}c} {p_{1} } \\ {p_{2} } \\ {p_{3} } \\ {p_{4} } \\ \end{array} } \right] = \left[ {\begin{array}{*{20}c} { - \left( {x_{1}^{2} + y_{1}^{2} + z_{1}^{2} } \right)} \\ { - \left( {x_{2}^{2} + y_{2}^{2} + z_{2}^{2} } \right)} \\ \vdots \\ { - \left( {x_{i}^{2} + y_{i}^{2} + z_{i}^{2} } \right)} \\ \vdots \\ { - \left( {x_{n}^{2} + y_{n}^{2} + z_{n}^{2} } \right)} \\ \end{array} } \right]$$

, respectively, $\underline {K} \cdot \underline {p} = \underline {R}$.

The searched parameter vector $\underline {p}^{T} = \left[ {p_{1} ,p_{2} ,p_{3} ,p_{4} } \right]$ can be calculated by use of the pseudo-inverse matrix via.

$$\underline {p} = \left( {\underline {K}^{T} \underline{ \cdot K} } \right)^{ - 1} \cdot \underline {K}^{T} \cdot \underline {R} .$$

The parameters of the sphere equation $x_{m}$, $y_{m} ,$ $z_{m}$ and R result from the backward substitutions

$$x_{m} = - p_{2} /2$$

$$y_{m} = - p_{3} /2$$

$$z_{m} = - p_{4} /2$$

$$R = \sqrt {x_{m}^{2} + y_{m}^{2} + z_{m}^{2} - p_{1} } .$$

While this matrix operation is very fast and effective, the analysis of the measurement data sets from the hip joint measurements did not lead to useful results. A possible reason might be a partially very unfavorable condition number of the linear equation system as well as the disturbance by some inaccurate data sets of the hip joint measurement.

To overcome this problem, the calculation of the sphere equation was separated into all possible four-point combinations of the entire hip joint measurement. That means, if the entire hip joint measurement consists of $n$ single measurements, all these single measurements will be combined into all their possible four-point combinations. For example, if only five single measurements exist, the possible four-point combination is:

• 1. Measurement 1, measurement 2, measurement 3, measurement 4

• 2. Measurement 1, measurement 2, measurement 3, measurement 5

• 3. Measurement 1, measurement 2, measurement 4, measurement 5

• 4. Measurement 1, measurement 3, measurement 4, measurement 5

In general, if $n$ single measurements exist, the total number of possible four-point-combinations is factorial of $n$ divided by factorial of four:

$$n!/4! = n!/24 .$$

For each combination, the four equations for the calculation of the sphere equation were solved

$$p_{1} + p_{2} \cdot x_{1} + p_{3} \cdot y_{1} + p_{4} \cdot z_{1} = - \left( {x_{1}^{2} + y_{1}^{2} + z_{1}^{2} } \right)$$

$$p_{1} + p_{2} \cdot x_{2} + p_{3} \cdot y_{2} + p_{4} \cdot z_{2} = - \left( {x_{2}^{2} + y_{2}^{2} + z_{2}^{2} } \right)$$

$$p_{1} + p_{2} \cdot x_{3} + p_{3} \cdot y_{3} + p_{4} \cdot z_{3} = - \left( {x_{3}^{2} + y_{3}^{2} + z_{3}^{2} } \right)$$

$p_{1} + p_{2} \cdot x_{4} + p_{3} \cdot y_{4} + p_{4} \cdot z_{4} = - \left( {x_{4}^{2} + y_{4}^{2} + z_{4}^{2} } \right)$.

Substitution and backward substitution of $p_{1}$, $p_{2}$, $p_{3}$ and $p_{4}$ was carried out as described before. In result, the parameters of the sphere equation have been calculated for each four-point combination of the hip joint measurement (Fig. 7, upper left). As can be seen, some single four-point combinations led to huge errors in the calculation of the sphere equation parameters. To exclude them from the further analysis, the mean value and the standard deviation of all parameters were calculated. All results were removed, which exceeded a defined threshold of plus/minus standard deviation from the mean value. The result of the first application of this threshold rule is depicted in Fig. 7 (upper right). While the result of the first application of the threshold rule was still not sufficient, the same rule was applied subsequently for a third (Fig. 7, bottom left) and fourth time, see Fig. 7 (bottom right).

Fig. 7

Parameters of the sphere equation from four-point combination

After this repeated application of the threshold rule, analysis results close to the median of all remaining results. These results were used to calculate the arithmetic average of equation $x_{m}$, $y_{m}$, $z_{m}$ and $R$ which are considered as the resulting parameters of the sphere equation.

Furthermore, the standard deviation of all remaining single values of $x_{m}$, $y_{m}$, $z_{m}$ and $R$ was calculated. If these standard deviations are small, the trajectory of the hip joint measurement is similar to a sphere shape. If these standard deviations are large, the trajectory of the hip joint measurement is not similar to a sphere shape. The decision, which standard deviation can be evaluated as ‘small’ or ‘large,’ depends on the specified tolerances for the following hip joint test.

In addition, 3D CAD illustration (Rhinoceros 3D, Robert McNeel & Associates) was used to visualize the input data and, respectively, the evaluated results. Subsequently, the evaluated data of the bone sample calibration procedures were processed, to identify the TCP of the bone sample (artificial hemipelvis-femur).

Evaluation of accuracy between manual (calibration) and robotic (automatic reproduction) movement:

Based on the proposed methodology, a motion reproduction for human hip joint specimens was performed. Therefore, a human hemipelvis cadaver was mounted on the test platform by using special holders. Primarily, the trajectories of the physiological movements of the hip were evaluated (stage I). The femur was moved manually (without external force) and the physiological range of motion was recorded with the ARAMIS system (Fig. 8). The motion was recorded with a frequency of 10 Hz, which allows for a large number of measurement points with good reproduction of the pathway.

Fig. 8

Documentation of manual flexion–extension of the human hip joint: images (1), (2), (3) present three selected positions from the movement

Secondary, the robot's calibration movements were performed according to the aforementioned procedure. The robot's movements are always the same and are determined in the robot's control program according to calibration points P1–P11 (stage II).

In the next step, the developed numerical transformation (Delphi Transformation) was used to convert the coordinates of the manual movement of the human hip specimen (recorded in optical ARAMIS coordinates system), into the Robot world coordinate system, which is necessary to control the robot's movement adequately (stage III).

Three data sets were used as input for the Delphi transformation (Fig. 9):

• (1) Coordinates from the robot control program including positions of the robot for the calibration procedure—points P1–P11 (data in robot world coordination system)

• (2) Coordinates of the calibration plate attached to the robot flange measured during robot calibration movements done by the ARAMIS optical measurement system (points recorded in the ARAMIS coordinate system)

• (3) Coordinates of the markers of the calibration plate attached to the femur bone sample holder obtained during manual movement (points recorded in the ARAMIS coordinate system).

Fig. 9

Block diagram of the Delphi Transformation: (1) input data of robot calibration positions P1-P11 from robot control program, (2) input data of robot positions during calibration movement from ARAMIS system, (3) input data of hip specimen manual movements from ARAMIS system, (a) output data of TCP, (b) output data of length of tool (LOT), (c) output data of coordinates points (trajectories) transformed into robot world coordinate system

As the results of Delphi Transformation, the following output is obtained (Fig. 9, Additional file 1):

• (a) TCP which corresponds to the JRC,

• (b) LOT which corresponds to the length of the femur bone with the holder,

• (c) Trajectories (set of coordinates) of the femur cadaver movement transformed into robot world coordinate system.

The reconstruction of the anatomical movements can be achieved using data from the Delphi Transformation in two ways:

Method 1 uses a designated TCP and LOT. This method is reserved solemnly for joints that can be defined as ball joints (i.e., the hip joint). In the presented investigation, the coordinates of the TCP were determined as xR = 1009.45, yR = − 124.54, zR = 59.39 (coordinates written in the Robot world coordinate system), while the LOT was evaluated to be 383.11 mm.

Method 2 uses data in the form of a trajectory of motion of the end of the femur bone (farthest tracker from the JRC), without using the JRC and without using definition of tool length (TCP set as value 0 mm + length of force torque sensor). This method is more universal and can be applied to various types of joints, e.g., the elbow joint. In the described investigation, a series of coordinates (xR, yR, zR) were obtained for more than 230 points showing the movement of flexion–extension, abduction–adduction and circular movement of the femur.

The output data obtained from the Delphi Transformation (either according to method 1 or 2) is entered into the robot control system. The robot is then positioned by the operator at a point coinciding with a single point in the motion trajectory. Next, the attachment of the femur holder to the robot flange is performed. Now the physiological movements in automatic mode are feasible (Fig. 10).

Fig. 10

Physiological movements of the human hip specimen reproduced by the robot for flexion–extension: images (1), (2), (3) present three selected positions from the movement

To evaluate the quality of the robot's movement reproduction, measurements of the robot's automatic movements were made and analyzed using the same optical ARAMIS system. All measurement procedures with the ARAMIS system were performed on the same day of the investigation, with equal conditions and in the same position according to the robot. For better understanding and visualization of the obtained results, the trajectories of movements during manual calibration and robotic motion were compiled and analyzed in a three-dimensional CAD system.

Results

Determination of the accuracy (repeatability) of designation of the reference point of the world robot coordinate system

In the present study, reconstruction of Robot World origin was created individually for each of the six measurement positions (Fig. 3). For this purpose, the xA, yA, zA reference point coordinates measured by the ARAMIS optical system for all six bone sample mounting positions were compared with each other. The average value of deviation around the theoretical position of Robot world origin xR, yR, zR = 0, 0, 0 was  ± 0.3 mm. The maximum values did not exceed ± 0.7 mm (Table 1).

Table 1 Average values of deviation for reconstruction of the Robot World origin—values based on all six measurement positions

Accuracy of calibration procedure of robot parameters TCP and LOT

Tool center point—TCP (hip rotation center)

Since the fixed holder of the bone samples was shifted by six various positions on the mounting plate (Fig. 3), the values obtained for the TCP also reflect these different spatial positions. For easier interpretation of the results, the data were normalized according to the shift values resulting from the 50 × 50 mm hole matrix and the fixing positions of the bone sample holder (Table 2).

Table 2 TCP position normalized to Position 0,0,0 in World Robot coordination system for six positions of bone samples (all values in mm)

In addition, the data with the TCP coordinates were compiled and transformed with the values of the main setting in the robot control program. The deviations were adjusted to a mean value of deviation at zero. The obtained results represent the deviation values of TCP position as defined in the robot control program (Table 3).

Table 3 Deviation of the TCP position referred to the settings in robot control program (all values in mm)

Length of tool (LOT)

The last parameter analyzed was the determination of the tool length (Table 4) which defines the distance between the TCP and the center of robot flange after attachment of the femur holder to the robot (Fig. 1).

Table 4 Results of the tool length (LOT) determination for bone sample (all values in mm)

The average value for the LOT measured with 3D CAD was approximately 366.47 mm (including the length of the femur bone sample together with the femur holder) with median the range of deviation not exceeding 0.04 mm. The extreme values of deviations are about maximum/minimum deviation of + 0.67 mm and − 0.4 mm, respectively. Similarly, values could also be determined for the Delphi transformation, where the LOT was approx. 367.24 mm with mean deviation of 0.48 and the maximum/minimum deviation of 0.72 mm and 0.13 mm.

Comparative investigation of the accuracy between manual and robotic movement with hip specimens

In the final stage of the investigation, the quality of the motion trajectory was compared between manual and robotic movement (Figs. 8, 10 and 11, Table 5). For each movement, the deviations of the manual trajectory curve points xA, yA, zA (recorded by Aramis during calibration) vs. the curve points reproduced by the robot in the xR, yR, zR axis directions were evaluated. In each case, the average, maximum and minimum deviation of the points of the compared curves were determined.

Fig. 11

Comparison of motion trajectories recorded for flexion–extension movement: manual motion of the bone sample recorded during calibration (red curve), reconstruction of the motion of the bone sample made by the robot (blue curve)

Table 5 Results of the trajectory curve coordinates deviation recorded for flexion–extension motion

The results obtained lead to the conclusion that the highest deviation of the movement trajectory can be evaluated for the x-axis. The average deviation of the points on the trajectory curves varies between 3.44 and 0.17 mm, with the lowest deviation in the z-axis and the highest in x-axis direction. The maximum deviation with 3.66 mm was also in the x-axis, while the lowest (0.24 mm) was evaluated in the y-axis. The average value for the Euclidean distance in 3D space was obtained at 3.59 mm (value determined for all X, Y, Z coordinates).

Discussion

Biomechanical and kinematic understanding of complex joints (i.e., the hip joint) is essential for further comprehension of pathologies, injuries and surgical procedures. The main limitation of common biomechanical investigations is the mostly uniplanar force or tensile application [50,51,52,53]. Therefore, 3D motion and force input is necessary for measuring complex movement patterns and kinematics [4, 5, 54]. The use of a 6 DOF robots guarantees the possibility of free motion in 3D-space while maintaining high positioning repeatability [41, 44, 55, 56]. Additionally, the opportunity of developing appropriate control algorithms allows reproducibility of the complex 3D movements of the test sample in accordance with the anatomical structure of the examined joint and bone. This will lead to a better understanding of the biomechanics and kinematics of the hip joint and, in the long term, may help to investigate endoprostheses biomechanically and make them safer for patients.

The main focus of the current study was to develop an accurate and fast calibration method for this testing approach by using a robotic platform without the need for detailed measurements of the geometrical joint parameters (femoral head, acetabulum) and characteristic points of specimens (to determine the axes defined by the ISB). Therefore, an algorithm was developed to establish a universal applicable biomechanical platform using a numerical transformation. Within the presented calibration, measurements are obtained using an optical measurement system and data processing takes about 10 min to generate the coordinates (trajectory) for robot movements.

We were able to show that by combining different coordinate systems, a standard deviation of the TCP between 0.3 and 0.9 mm depending on the axis could be achieved. The maximum deviation of the LOT was evaluated between + 0.67 and − 0.4 mm (3D CAD processing) and 0.72–0.13 mm (Delphi transformation). The accuracy between the manual and robotic movement of the hip shows an average deviation between − 0.36 and + 3.44 mm for the points on the movement trajectories. The presented error represents the values as a total error for the entire ‘universal kinemator’ system. The total error results from adding up error of control and positioning of the robot, measurement error of the optical measurement (ARAMIS) system, the error in defining individual planes and axes of the robot and bone samples, the error of joint imperfections and error in defining the calibration element and its mounting.

Goldsmith et al. examined the reproducibility by using the hip joint model and was able to demonstrate that a 6 DOF robotic system allows to identify hip passive paths (without external loads) in a highly repeatable manner (median RMS error of < 0.1 mm and < 0.4 degrees) [41, 42]. Moreover, the robotically simulated clinical exams were found to be consistent and repeatable (rotational RMS error ≤ 0.8 degrees) in determining hip ranges of motion. After converting the degree values to our specified LOT, the deviation seems to be similar in the underlying study (1 degree of error results in a position error of about 3.5 mm at the end of the femur). However, the indication of the maximum deviation appears to be superior, since maximum divergent robot trajectory can cause possible peak forces with potential damage to the specimen under investigation. Compared to Goldsmith et al. in the aforementioned method, there is no need to measure the femoral head or hip acetabulum for defining the position of the JRC relative to the bone landmarks [41]. In addition, there is no need for intraarticular preparation and no risk of damaging the hip joint and subsequently altering its biomechanical parameters, which may have negative effects on TCP when a best-fit sphere is used to describe the locations recorded on the surface of the femoral head with the intact capsule [41].

In order to ensure correct fixation and load on the examined hip joint, it is necessary to properly define the individual anatomical and mechanical axes of the tested sample, as well as define the hip joint center [49]. For determination of HJC, special procedures and protocols must be used [57]. The calibration procedure presented in the current study allows to define basic parameters such as HJC (robot definition: TCP), geometric dimensions of the bone sample mounted in the handles (robot definition: LOT) without risk of accidental damage to the joint/bone sample or other elements of the test set up. The advantage of our proposed method in comparison with others is that the movement can be mimicked by two methods, using the TCP, the length of tool (based on an approximation of the shape of the joint to a sphere) or using the recorded physiological movement trajectory of the femur [41, 55]. Moreover, in most setups a compressive force (axially and medially 10 N) is required to ensure contact between the femoral head and the articular surface of the acetabulum [42]. In our methodology, this is not necessary (measurements are taken without external forces) since the entire joint is intact and provides full anatomical mobility.

While there are no consistent calibration procedures in terms of the establishment of biomechanical platforms using a 6 DOF, the current study demonstrates a transparent and universal applicable calibration pathway for further research approaches. The proposed innovative trajectory-based approach is a universal solution independent of the shape of the joint. As a result, the method can be used to accurately represent joint function, regardless of the fact that the femoral head is naturally conchoidal and slightly aspherical in shape [58]. This can be particularly important when performing tests on pathological joints (i.e., after trauma). This approach also has the advantage of making the system robust to inaccuracies in the mounting of the specimens in the grips relative to the planes and the axes defined using the ISB coordinate system [49].

A limitation of our study is the small number of bone samples examined and thus associated failure to account for the variables of age and gender. The recorded inaccuracies in the realization of the movements may be due to some systematic error in the whole system, including errors in movement registration, robot calibration, sample mounting, data processing and generation of the control program for the robot. In addition, the effect of changing the length of the test sample (LOT) to the error of trajectory reconstruction was not defined. In the underlying study, the length of the femur including the handle was close to 367 mm, which may also affect the deviation values. These values can be different in comparison with tests presented by other researches where the length of tool is much shorter [55]. However, with the current method, it is possible to perform complex movements in different axes simultaneously—for example, rotational movement around the JRC.

Conclusion

A six degree-of-freedom robot is appropriate to reproduce the physiological range of motion in hip joint. The described procedure of calibration is universal and can be used for hip joint biomechanical tests, regardless of the length of the femur, size of the femoral head and acetabulum or whether the entire pelvis or hemipelvis is used.

Clinical implications of this universal calibration procedure may include further robot studies that are aimed to apply physiological loading forces [4, 54]. The use and translation of these in vivo measured data, indicating force transmission and loads in the joints during activities of daily living, allows much more realistic biomechanical stability testing of endoprosthetic implants and osteosynthesis [54]. Pursuing this approach may add much more value to the currently used and commonly available test procedures. These and other potential applications of the procedure for other kinematic and biomechanical joint examinations will be the subject of further investigations.