Next Article in Journal
Gold Nanorod-Assisted Photothermal Therapy and Improvement Strategies
Next Article in Special Issue
How do Paraspinal Muscles Contract during the Schroth Exercise Treatment in Patients with Adolescent Idiopathic Scoliosis (AIS)?
Previous Article in Journal
Review of Biomechanical Studies and Finite Element Modeling of Sternal Closure Using Bio-Active Adhesives
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Turmell-Meter: A Device for Estimating the Subtalar and Talocrural Axes of the Human Ankle Joint by Applying the Product of Exponentials Formula

1
Facultad de Ciencias Básicas e Ingeniería, Universidad de los Llanos, Villavicencio 500002, Colombia
2
Instituto Universitario de Automática e Informática Industrial (Instituto ai2), Universitat Politècnica de València, 46022 Valencia, Spain
*
Author to whom correspondence should be addressed.
These authors contributed equally to this work.
Bioengineering 2022, 9(5), 199; https://doi.org/10.3390/bioengineering9050199
Submission received: 3 February 2022 / Revised: 4 April 2022 / Accepted: 20 April 2022 / Published: 4 May 2022
(This article belongs to the Special Issue Biomechanics-Based Motion Analysis)

Abstract

:
The human ankle is a complex joint, most commonly represented as the talocrural and subtalar axes. It is troublesome to take in vivo measurements of the ankle joint. There are no instruments for patients lying on flat surfaces; employed in outdoor or remote sites. We have developed a “Turmell-meter” to address these issues. It started with the study of ankle anatomy and anthropometry. We also use the product of exponentials’ formula to visualize the movements. We built a prototype using human proportions and statistics. For pose estimation, we used a trilateration method by applying tetrahedral geometry. We computed the axis direction by fitting circles in 3D, plotting the manifold and chart as an ankle joint model. We presented the results of simulations, a prototype comprising 45 parts, specifically designed draw-wire sensors, and electronics. Finally, we tested the device by capturing positions and fitting them into the bi-axial ankle model as a Riemannian manifold. The Turmell-meter is a hardware platform for human ankle joint axes estimation. The measurement accuracy and precision depend on the sensor quality; we address this issue by designing an electronics capture circuit, measuring the real measurement with a Vernier caliper. Then, we adjust the analog voltages and filter the 10-bit digital value. The Technology Readiness Level is 2. The proposed ankle joint model has the properties of a chart in a geometric manifold, and we provided the details.

1. Introduction

Taking in vivo measurements in the human ankle joint is troublesome because the ankle is a complex mechanism [1]. Deviations in the axis increase the pronation or supination moments, causing instability and enhancing injuries risk. In this work, we present a device intended for the study of the human ankle joint (HAJ). Modeling and measuring this lower limb joint is essential in physiology, biomechanics, and rehabilitation (also in humanoid robotic limb development).
Our primary aim is to develop a device for the two axes model estimation of the human ankle joint. Secondary objectives are: it must be non-invasive, compact, energy-efficient, and easy to set up and transport. It should also be compatible with laying positions, such as with the foot in the elevated position. To accomplish the objectives, we followed a plan, first by understanding the ankle movements. Then, we used statistics for dimensional determination. We also use a modern approach, such as the Product of Exponential (POE) formula. We then designed the structure based on embedded non-invasive distance sensors.
Our contribution to the ankle joint axis localization is the holistic development of a specific device. Draw-wire sensors measure distance, are composed of a wire wound around a drum, and are attached to a potentiometer and a spring. They are retractile with constant tension. For bias correction and gain calibration, we designed a capture system. We adjust the voltage to avoid the maximal value of the analog-to-digital conversion. We calibrated each sensor through direct measuring with a Vernier caliper. Then, we measured the voltage and adjusted the offset and gain by a calibration program in Processing (Software). Limitation measurements are by 10-bit analog-to-digital converters and digitally filtered in the acquisition board. Technology Readiness Level (TRL) is 2.
We highlight our approach over traditional methods because we apply the POE formula to the ankle kinematic model. Furthermore, we estimate the ankle axis localization by a geometric approach, solved algebraically. We computed it from the pseudo-inverse application. For the talocrural and subtalar axes estimation, we use circle fitting. As an alternative ankle joint representation, we propose a Riemannian chart. We have limited the scope to the human ankle joint (HAJ) model. There are applications in physical therapy and HAJ mobility diagnosis.
The state of the art in the ankle localization is detailed in [1,2,3,4,5,6,7,8,9,10,11,12,13,14].
There are different HAJ models in the literature; we focus on the two-axes approach. The approach is recommended by the International Society of Biomechanics (ISB) [15], anatomy and biomechanics books [3,16,17,18,19], and simulation software [20]. We found models of the ankle joints in several articles [14,21,22,23,24,25,26]. Contributions to the study of the ankle joint axes are in [2,8,9,27]. The most cited research about the subtalar axis are in [5,7,10,11,12,13]. A literature review of functional representations is in [4].
Draw-wire sensors (DWS) are distance measurement sensors, who use a wire coiled on a drum attached to a potentiometer and a spiral spring that are retractile at constant tension. Similar robotic applications are in [28,29,30], also in linear position tracking [31], and easy robot programming [32]. Inertial measurement units (IMU) were post-processed and complemented with other sensors [33,34,35,36,37]. We shall employ our device for the HAJ bi-axial measurements and for other models as well [38]. BiodexTM and HumacnormTM are manufacturers of general kinetics machines.
We divide the materials and methods section into two subsections: the motion theory and the mechatronics design. In the first section, we study anatomy, statistics, proportions, and anthropometry to understand the functional HAJ movements and standard dimensions. Then we perform the HAJ simulation using the POE formula. Here, we do not include a deep study of infinitesimal kinematics. We intend to design a device for a healthy HAJ with no singularities with a continuous range of movement. We describe the trilateration method to find the platform pose. It is a geometrical method based on tetrahedrons; we avoid numerical solutions that depend on finite derivative terms. The tetrahedron is a well-defined 3D geometrical structure. Solving tetrahedron geometry is the expansion of planar trigonometry. Knowing the sides allows us to find the height of a tetrahedron. We attach the platform to the foot; the sensors are passive elements and do not support or add high tensile forces. We have selected the first seven sensor configurations 3-2-2 (seven sensors) instead of 3-3-3 (nine sensors) or 3-2-1 (five sensors) for hardware limitations, sensor redundancy, and symmetrical design (for both limb use).
The device’s mechatronics design and implementation are in the second subsection. We used Draw-wire sensors to measure the tetrahedron sides. These sensors have a constant tension because they comprise a drum attached to a spiral spring. We limit them to the maximal distance, and the precision depends on the potentiometer and electronics signal conditioning with a high common-mode rejection ratio (CMRR). The calibration process deals with accuracy and precision. First, we made rough adjustments to the acquisition system. Second, the software calibration process makes fine adjustments. Our proposed method avoids numerical errors because it uses geometric formulas. We validate the position through sensor redundancy. We conduct calibration and testing in a healthy patient and represent the HAJ movements as a manifold chart. The complementary source code was uploaded to [39].

2. Materials and Methods

This section is grouped in two main subsections, first the motion theory, and second the mechatronic system. For the first part, we show the simulation using anthropometric values and the POE formula. Using the plots, we estimate the DWS maximal length. Next, we present the device’s geometrical design and the trilateration method. Finally, we compute the axis position by circle fitting and modeling the ankle joint as a Riemannian manifold chart. In the second subsection, we describe the mechanics design and implementation, we used SolidWorks® (2017–2018 Student Edition, Dassault Systèmes, Vélizy-Villacoublay, France), KiCad©(6.0.4, Jean-Pierre Charras and KiCad developers, CERN, Linux Foundation), and FreeCad© (0.19, Jürgen Riegel, Werner Mayer, Yorik van Havre and others) StepUp tools addon.

2.1. Motion Theory

For the simulation with the POE formula, we adapt the data from [40], proportions from [41,42], and statistics from [43].

2.1.1. References Assignation

Figure 1 presents the reference points and the mean distances taken from [40].
A, B, and C are the triangle’s vertices in a platform fixed to the foot, the K, L, and O distances from the most medial and lateral points from the black-filled to the white-filled marker. M 1 and M 2 define the talocrural (TC) axis. We show top-transverse and right-lateral views in Figure 2 with distances Q, W, and w. N 1 and N 2 determine the subtalar (ST) axis.
Table 1 enumerates the mean values of Figure 1 and Figure 2.
In Figure 3, we show the ST and TC axes from several viewpoints. The TC axis refers to the sagittal plane and the ST to the transverse plane.

2.1.2. Anatomical and Geometrical Correspondence

We define the sagittal (lateral) plane as the X-Z plane (perpendicular to the y-axis). The coronal (frontal) plane is the Y-Z plane (x-axis is normal to it); the transverse (axial) plane is the X-Y plane (perpendicular to the z-axis). Figure 4, left, shows this corresponding references.
With this reference frame, we can define the TC axis orientation from a unitary vector in the z-direction. We first rotate it −80° around the x-axis; then we turn it −6° around the z-axis. A unitary vector in the x-axis direction defines the ST axis, rotating 41° about the y-axis, followed by a 23° rotation around the z-axis.
We show the fibula, tibia, talus, calcaneus 3D position, reference points, TC, and ST axes in Figure 4, right.
In this image, A 0 , B 0 , and C 0 are the vertices from the platform fixed to the foot, and P M is the triangle’s center. S 1 , S 2 , and S 3 are fixed to the shank relative to the origin point P 0 . M 1 and M 2 define the TC axis; N 1 and N 2 correspond to the ST axis. We define r 1 and r 2 as the sagittal plane intersection with the TC and ST axes.

2.1.3. Size and Dimensions

This first part help us to determine the HAJ axes direction and orientation for some cases. However, it is difficult to design a device that fits all humans, and we cannot make a device that fits 90 % percentiles; we intend to design a device scalable and adjustable in a defined population group. We also make an effort to design adjustable foot and shank attachments. To do so, we select the device dimensions using the proportions extracted from [41]. The heigh is H, and the proportions we use are: distance from the knee to the foot is 0.285H, the distance from the ankle to the foot is 0.039H and the foot widht is 0.055H and the foot length is 0.152H.
We select the origin of coordinates between the knee and the ankle, d m is the distance from P M to P O . This distance is proportional to the body’s height H . To do so, we define d m as follows:
d m = P 0 P M = 0.285 0.039 2 + 0.039 · H = 0.162 · H .
For the sake of obtaining the prototype dimensions, we use statistics for a specific population. In [43], the mean height H of an adult male is 175 cm; by substituting this value into the equation, the knee-ankle distance is 28.35 cm. The distance d p 12 between points r 1 and r 2 about the TC and ST axes on the sagittal plane is:
d p 12 = r 1 r 2 = Q ,
the projection of the most medial point (MMP) on the sagittal plane is:
P MMP = ( x MMP , 0 , z MMP ) ,
and for the most lateral point is:
P MLP = ( x MLP , 0 , z MLP ) .
The point M 1 p is the projection of M 1 on the sagittal plane; we calculate it from the P and O values.
M 1 p = ( x MMP P , 0 , z MMP O ) ,
also, M 2 p is M 2 ; we estimate the projection from L and K through:
M 2 p = ( x MLP L , 0 , z MLP K ) .
Therefore, the segment M 2 M 1 ¯ has the sagittal projection M 2 p M 1 p ¯ ; it has the same proportional relation R = W / w in respect to M 2 p r 1 ¯ , then:
M 2 M 1 M 2 r 1 = W w = R ,
solving for r 1 gives the following:
r 1 = M 2 M 2 M 1 R .
By knowing the distance Q projected in the sagittal plane and r 1 , the angle 41° we calculate r 2 from:
r 2 = Q cos ( 41 ) , 0 , sin ( 41 ) + r 1 ,
The distance from the origin P O to the plantar surface of the foot is d m , we choose a circumscribed equilateral triangle with vertices A 0 , B 0 , C 0 as the platform base. The coordinates of A 0 are:
A 0 = r p , 0 , d m ,
for B 0 are:
B 0 = r p cos 60 , r p sin 60 , d m ,
and for C 0 :
C 0 = r p cos 60 , r p sin 60 , d m ,
where r p is proportional to H, then:
r p = 2 3 · H .
In summary, we estimate P 0 , r 1 , r 2 ; and the platform’s vertices A 0 , B 0 , and C 0 . They are not arbitrarily selected, on the contrary, we employed anthropometry, statistics, and proportions.

2.1.4. Product of Exponentials Formula

In this section, we employ the PoE formula. We follow the intuitive concept that inter-bone contact surfaces determine HAJ movements. Therefore, we represent these movements as a Special Euclidean group SE(3) in matrix form:
g = R p ^ T 0 1 × 3 1 ,
where R 3 × 3 is the rotation matrix and p ^ T is the translation vector.
For the initial point A 0 :
g A ( 0 ) = I 3 × 3 A ^ 0 0 1 × 3 1 ,
for B 0 :
g B ( 0 ) = I 3 × 3 B ^ 0 0 1 × 3 1 ,
and for C 0
g C ( 0 ) = I 3 × 3 C ^ 0 0 1 × 3 1
We define ω ^ 1 = ( ω x 1 , ω y 1 , ω z 1 ) as a unitary vector for the TC axis direction given by:
ω ^ 1 = M 2 M 1 M 2 M 1 ,
and a directed vector r ^ 1 from P O to r 1 is:
r ^ 1 = r 1 P O ,
then, an orthogonal vector to r ^ 1 and ω ^ 1 is:
ν ^ θ 1 r 2 z = ω ^ 1 × r ^ 1 ,
together, ω ^ 1 and ν ^ θ 1 r 2 z compound the six-dimensional vector ξ ^ 1 :
ξ 1 ^ = v ^ 1 ω ^ 1 .
In the same way, there are correspondent vectors for the TC axis:
ω ^ 2 = N 2 N 1 N 2 N 1 ,
r ^ 2 = r 2 P O ,
ν ^ 2 = ω ^ 2 × r ^ 2 ,
and:
ξ 2 ^ = v ^ 2 ω ^ 2 .
We compute R for each joint i = 1 , 2 from the Rodrigues’ formula:
e ( Ω i θ i ) = I 3 × 3 + Ω sin θ i + Ω 2 1 cos θ i ,
where Ω is the skew symmetric matrix:
Ω = 0 ω z i ω y i ω z i 0 ω x i ω y i ω x i 0 .
The exponential formula is:
e ξ i θ i = e Ω i θ i τ i 0 1 × 3 1 ,
and, τ i is translation vector:
τ i = I 3 × 3 e ω ^ i θ i ω ^ i × ν ^ + ω ^ i ω ^ i T ν ^ i θ i
Points A, B, and C have invariant relative positions, and there are two rotating joints; the PoE formula for A is:
g A = e ξ ^ 1 θ 1 e ξ ^ 2 θ 2 g A ( 0 ) = R p ^ A 0 1 ,
where p ^ A is the instantaneous position vector of A, the PoE for B:
g B = e ξ ^ 1 θ 1 e ξ ^ 2 θ 2 g B ( 0 ) = R p ^ B 0 1 ,
and the PoE for C is:
g C = e ξ ^ 1 θ 1 e ξ ^ 2 θ 2 g C ( 0 ) = R p ^ C 0 1 ,
θ 1 is the TC rotation angle from the zero position, and θ 2 is the ST rotation from the zero position. For the sake of clarity, we show the section of the ankle with the vectors r ^ 1 , ω ^ 1 , ν ^ 1 and r ^ 2 ; also the points A, B, C, and P O in Figure 5.

2.1.5. Forward Kinematics

In this subsection, we show the simulation of the movements of the ankle by using the measurements and the PoE. The code is in SageMath Computer Algebraic System (CAS), which lets us manage symbolic notation, and interactive plotting in a Jupyter notebook. All source was uploaded to Git-Hub [39].
The simulation plot for the platform’s central point is in Figure 6a. We show the points P O , A 0 , B 0 , C 0 , r 1 , r 2 , and the surfaces representing each group of movements. The forward kinematics with θ 1 r a n g e = θ 2 r a n g e = [ 15 , 15 ] and θ 1 = θ 2 = 10 is in Figure 6b. For θ 1 r a n g e = θ 2 r a n g e = [ 10 , 10 ] and θ 1 = θ 2 = 5 is in Figure 6c.
Such a representation lets us compute the ankle joint ROM in all directions. Groups of A, B, C, and PM movements are smooth surfaces or geometric manifolds. They have two DOF, with a limited domain due to the axes ROM.

2.1.6. Geometric Design and Trilateration Method

In the last section, we note that in a healthy ankle, the range of motion from three points in a platform attached to the foot, pertain to a surface without singularities. Moreover, we note that we can trace tetrahedrons from the reference base on the shank to the platform points. Tetrahedrons can be solved by knowing the triangle base, and the sides. We choose the complete tetrahedron with three distance sensors in the point A for symmetry conservation in the case of using the device in the left or right foot. We avoid the use of numerical methods such as the Newton–Raphson (NR), for reducing time of computation. Furthermore, we choose the symmetry and redundancy in the apexes B and C. We realize that by knowing the platform dimensions, two sensors and the apex A coordinates, we can define a plane rotated with respect to the base and solve other tetrahedrons corresponding to the B and C apexes. We also take a holistic approach, we knew that micro-controller systems often have two cores and eight or ten analog to digital converter channels. We used 7 channels, leaving three for temperature, battery level, and voltage input detection.
Finally, based on such considerations, we show a geometric design in the Figure 7a platform center, and in Figure 7b are the vertices.
By considering the distances between the origin and the vertices, we estimate the DWS maximal length in every module.
l m a x = max p A ( θ 1 , θ 2 ) A + r m
Here, l m a x is the maximal possible length from the triangular inequality, p A is the positions group in g A , r m is the module’s radius, and A B is the base point.
The main design requirement is the localization of three points attached to the foot. We estimate the actual position employing a DWS array in a tetrahedral structure to find the apex, which is a platform vertex. In Figure 8 we show the design structure.
P O and P M are the base and platform reference frames. The platform has known dimensions and the number of sensors is seven. First, we compute A p from three distances: l A 1 = A p A 1 , l A 2 = A p A 2 , and l A 3 = A p A 3 . Then, we compute B p and B p apexes after A p employing two DWS. We summarize the method in a flowchart; Figure 9.

2.1.7. Finding the Apex in Tetrahedron A

In this section, we compute the tetrahedron T A with base A = [ A 1 , A 2 , A 3 ] and apex A p . Figure 10 shows the method we use.
In Figure 10 we see that triangles 132 = [ A 1 , A 3 , A p 132 ] and 231 = [ A 2 , A 3 , A p 231 ] are two sides of the tetrahedron T A developed on the base plane.
We compute the A p 132 and A p 231 orthogonal projection on each adjacent side of the module base triangle [ A 1 , A 2 , A 3 ] by tracing a circle centered on A 1 with radius A p A 2 and the circle centered on A 3 with radius A p A 3 ; resulting in A p 132 and A p 131 intersection points. In addition, the circle centered on A 2 with radius A p A 2 intersects the circle centered in A 3 at points A p 231 and A p 232 . The segment from A p 132 to A p 131 intersects the points defined by A p 231 and A p 232 at A pxy . In the case of tetrahedron T A , we determine A pxy = ( Apx , Apy , 0 ) as A p projection on the base plane. It is easy to realize that the height of T A is the absolute value of the A pz coordinate. Then, we can find the distance from A pxy to A 3 as a triangle [ A pxy , A 3 , A p ] side; the other is A pz , and the hypotenuse is the distance l A 3 = A p A 3 , then, A pz is:
A pz = l A 3 2 ( A pxy A 3 ) 2

2.1.8. Tetrahedrons B and C Apexes

In this subsection, we show that, by knowing A p , the point B p needs two sensors to be found. To determine the result of the tetrahedron T B , we consider the base of a triangle B 1 , B 3 , A p in Figure 11a.
We compute the angle α from the XY plane to a normal vector n ^ A p B :
n ^ A p B = ( B 3 A p ) × ( B 3 A p ) ( B 3 A p ) × ( B 3 A p ) ,
and, the angle α is:
α = a c o s ( n ^ A p B · n ^ z ) ,
where n ^ z is the unitary vector normal to the XY plane.
The tetrahedron sides are the lengths l B 1 = B p B 1 , l B 3 = B p B 3 , and d A p B p = B p A p . The rotation axis is in the direction B 1 B 3 . The B p r is B p s rotated α in angle about this axis. In Figure 11b, we show how to find the B p r apex, similarly to that of a tetrahedron T A . Finally, when B p r is found, the contrary rotation about the axis B 1 B 3 gives the B p s .
There are two possible apex values: B p s 1 over, and B p s 2 below of the XY plane. We show the B p r apex below the XY plane in Figure 12.
We use the same method to solve the T C apex. For the correct apex selection, the condition when the side of the platform distance d C p B p is:
d C p B p = B p s C p s .

2.1.9. Procedure for Found Platform Positions

We must fix the shank and the foot to the base and platform. Then we mark the MMP and the MLP. To do so, we design a detachable reference point from module A. Initially, we attach the foot and the shank to the device, and then we mark and record the MMP and MLP; Figure 13 shows the detailed view.
We compute the platform position from the seven sensor lengths. The main steps for capture data series are:
1.
Capture the initial position at horizontal relative position from d r = I M U 2 I M U 1 readings;
2.
Compute jerk j r k = | d r i d r i 1 | ;
3.
Move the foot continuously until jerk crosses zero again.
First, we capture the sensor lengths by activating a button in the computer software. Every time, we compute the absolute difference from IMU2 to the IMU1 readings. If the differences are constants, then there is no platform and base relative movement. We compute the jerk by relative acceleration differentiation. The data capturing process ends when the acceleration change crosses zero. Jerk changes activate the capture of IMU data.
The symbolic equations to find A p , B p , and C p from the captured data, were found by the SageMath CAS. By using the prototype dimensions and the sensor lengths, we compute the platform’s position and orientation. Here, the origin is from the initial DWS lengths l M i 0 , where M is the module A, B, or C; and i is the sensor number i = 1 , 2 , 3 .
After MLP and MMP registering, we attach the apex of module A to the platform, define the sagittal plane perpendicular to the ABC base plane, and intersect point A. By implementing the trilateration method mentioned before, we compute the points A 0 , B 0 , and C 0 .
Figure 14a illustrate the point positions with the device in the initial portable configuration. The apexes’ computation are in Figure 14b–d.

2.1.10. Computing the Axis Position and Direction

From the anthropometric values [40], we put a mean model in the Turmell-Meter. The TC axis will be defined by M 1 0 and M 2 0 . The sagittal plane intersection with the M 1 M 2 ¯ segment is r 1 . For example, the TC axis approximation is computed from most lateral point (MLP), the most medial point (MMP) and L, K, P, and O:
M 1 0 = M L P [ L , 0 , K ] ,
and:
M 2 0 = M M P [ P , 0 , O ] ,
from these values, we solve for r 1 from the plane y = 0 intersection with the line L T C :
L T C = V M 1 = ρ ( M 2 M 1 ) M 2 M 1 ,
where V is a point pertaining to L T C .
The ST axis sagittal intersection r 2 initial point is:
r 2 = r 1 + 2 2 [ Q , 0 , Q ] .
Here, r 1 and r 2 are reference values computed from the previously mentioned anthropometric mean values. Such initial points are for reference, comparison and validation of the trilateration and regression method. The tracked trajectory data set is processed offline. We use the least squares normal vector to the plane, this direction is similar to the circle approximation. From here, we compute the TC axis first, and then the ST axis. To do so, we compute the TC axis position by employing dorsiflexion and plantarflexion, because the TC axis is the most dominant in such movements. The method used is circle fitting in a plane containing the trajectory points. A further model refinement can be made with optimization, and machine learning methods, such as gradient descent and the symbolic product of exponential formula.
First, we found the TC axis orientation ω 1 by registering several trajectories. For each trajectory, we have a list of data points P = [ x , y , z ] , which pertain to a plane:
a x + b y + c z + d = 0 ,
where a, b, and c are the components of a direction vector perpendicular to a plane containing the points. Solving out for z, we have the system:
x 0 y 0 1 x 1 y 1 1 x n 1 y n 1 1 a b d = c z 0 z 1 z n 1
which has the form:
A x = B
there are more equations than unknowns. From linear algebra and least squares we knew that the pseudo inverse is A + = ( A T A ) 1 A T , then a normal vector is:
a b d = ( A T A ) 1 A T B
Now, we compute c by replacing a, b, d in the plane equation, and finally we get n ^ = [ a , b , c ] T . We found the angle between the normal plane and the X-Y plane, after knowing the normal vector by applying the Rodrigues’ formula, v ^ = n ^ × k ^ , with k ^ = [ 0 , 0 , 1 ] T
P r = P cos ( θ ) + ( v ^ × P ) sin ( θ ) + v ^ ( v ^ · P ) ( 1 c o s θ ) .
where θ = arccos n ^ · k ^ n ^ .
After this, we estimate the plane, and rotate all the data points onto the X-Y plane. We search for a circle in the X-Y plane, and rearrange the equation for least squares estimation by using a variable substitution.
( x x c ) 2 + ( y y c ) 2 = r 2 ( 2 x c ) x + ( 2 y c ) y + ( r 2 x c 2 y c 2 ) = x 2 + y 2 c 0 x + c 1 y + c 2 = x 2 + y 2
where c = [ c 0 , c 1 , c 2 ] T with c 0 = 2 x c , c 1 = 2 y c , and c 2 = r 2 x c 2 y c 2 .
By taking the rotated points, P r we have a linear system:
x 0 y 0 1 x 1 y 1 1 x n 1 y n 1 1 c 0 c 1 c 2 = x 0 2 + y 0 2 x 1 2 + y 1 2 x n 1 2 + y n 1 2 .
that has the form:
A c = b
In this system, we have more equations than unknowns, then, we search for the c values that minimize the squared difference b A c 2 .
arg min c R 3 b A c 2 .
We found the center point C p = [ x c , y c ] and radius r by solving:
2 x c = c 0 2 y c = c 1 r 2 x c 2 y c 2 = c 2 .
Finally, we apply a rotation to the center in respect to the original plane. This point pertains to the TC axis. For each trajectory A, B, C, we get three planes, and three centers, the TC line direction is parallel to the planes’ normal vectors. The information is complete by determining the plane orientation.
The ST axis estimation is similar, but employs trajectories from inversion and eversion movements.
This is a basic estimation; by conducting optimization on the product of the exponential formula, we enhance the accuracy of the axis position estimation.

2.1.11. Ankle Joint Movements as a Manifold

In this subsection, we explain how the centers r 1 , r 2 and directions ω 1 , ω 2 define a manifold representing the HAJ movements. The circle center points calculated pertain to the TC and ST axes; they are the initial data to fit the product of the exponential formula. In Figure 15a, we show the complete platform’s center point manifold. It is topologically similar to a torus.
A manifold chart represents the range of motion limits, we show an example of the geodesic as a trajectory on the manifold in Figure 15b; this explains how to map ankle coordinates, and a straight trajectory with initial velocity and no external force action. We have the data necessary for the line intersection with the sagittal plane, the center points, and the direction gives a line:
p ^ l = l ^ 0 + l ^ d ,
where p ^ l is the parametric line, l ^ is a parallel vector to it, l ^ 0 is a known vector in such line, and d R , replacing the parametric equation in the plane equation:
( p ^ l p ^ 0 ) · ( n ^ p ) = 0 ,
where p ^ 0 is a known vector in the plane, and n ^ p is the plane’s normal vector, solving for d, gives:
d = ( p ^ 0 l ^ 0 ) · n ^ p l ^ · n ^ p ,
and replacing in the TC axis line equation:
r 1 = c 1 + ω 1 d ,
where r 1 is the TC axis intersection with the sagittal plane. The point c 1 is the center, and the axis direction ω 1 , both were found by circle fitting. Furthermore, packing in six dimensional Plücker line coordinates, we have:
m ^ 1 = r 1 × ω 1 ,
and the l 1 six dimensional vector is:
l 1 = [ ω 1 x : ω 1 y : ω 1 z : m 1 x : m 1 y : m 1 z ] .
We include those data for the PoE formula simulation and the manifold representation.

2.2. Mechatronic System Design

In this section, we design DWS to measure the lengths of the tetrahedron sides; they are arranged as structural parts. Their maximal length estimation is from the forward kinematics simulation. We design the shank attachment from the dimensions, proportions, and statistical data.

2.2.1. Draw-Wire Sensor

We use flat springs. They are not exposed to a high load against gravity, and are in two or three concurrent groups. In Figure 16, we depict the design, composed of three 3D printed parts, potentiometer, flat spring, bolts, and nuts.
A two-coil winch drives the potentiometer; a flat spring retracts a wire attached to the winch. When we pull the wire, the spring retracts it. The value of each turn is from the nominal value of the potentiometer, R n = 2.2 k Ω , divided into ten turns, that is 220 Ω per turn. The diameter is D = 3.8 cm , the spring could be compressed in four turns. The maximal length is as follows:
l m a x = 4 · D · π
Which is 47.75 cm approximately, this value is greater than l m a x for all groups of movements.

2.2.2. Mechanical Parts

The attachment on the calf has a size according to the simulation. We use the mesh model of a leg to guide the shape of the calf support, as in Figure 17a. We also scale and divide this structure into seven parts for 3D printing. An aluminum tube is the support structure, as in Figure 17b, and a neoprene band attaches the shank to the support with Velcro fabric.
All the DWS modules are in a plate, the A module has three DWS, B and C modules has two DWS, as in Figure 18a. The design of the foot attachment is from standard measurements to adjust the foot’s length and width, as in Figure 18b.

2.2.3. Electronics

Two operational amplifiers in instrumentation configuration are the base block of the acquisition system, as Figure 19 shows. We employ the KiCad software for the circuit design.
The voltage gain in the instrumentation amplifier is:
A v = v o v i = 1 = R 2 R 1 + 2 R 2 R 1 ,
By selecting R 2 = 100 k Ω , R 1 = 1 k Ω , and R G = 5 k Ω , the voltage gain is 141. With 34 mV as voltage input, we get:
v o = A v · v i = 4.794 V
The final acquisition circuit has seven instrumentation amplifiers, with bias and gain trimmers for calibration. We design the printed board circuit as an ™Arduino Mega 2560 Shield, and assemble the components to the board by throw-hole soldering. We feed the circuits with a power system with two 18650 Li-Ion batteries in series, a backup pack, a Battery Management System (BMS); a 5V buck and a 12V boost converters. The Figure 20 shows the schematics. Finally, we add connectors for the MPU, OLED, and Bluetooth modules.

2.2.4. Electronics Casing

We export the KiCad printed circuit design to FreeCAD StepUp to design the case containing all the components, focusing on a compact configuration design. The two main electronic components are the Arduino Mega 2560 and the Orange Pi One single board computer. We place the components, such as the Dual Pole Dual Throw (DPDT) toggle switches, symmetrically on the box sides. Figure 21 shows the main sides and the final assembly of the electronics case.
Each box has attached components to optimize the space. We test every component, and then install the support structure.

2.2.5. Final Mechanical Assembly

The prototype consist of 45 3D printed parts, the union of main components is by an 8 mm steel threaded rod. The sub-assemblies uses M3 bolts and nuts. Figure 22 shows the assembly CAD.

2.2.6. Calibration and Validation Software

Calibration is with the Arduino board connected to the PC, running a calibration program in processing. The basic program reads the IMU measurements and captures readings from the draw-wire sensors through the ADC inputs. The raw data are integer values with signs 2 bytes wide, the two 1-byte registers converted to 2-byte integers. An exponentially weighted moving average (EWMA) algorithm filters the raw signals and sends them to the PC via a serial port. The lengths computed are from the initial values plus the scaled sensor inputs with:
l i M j = d i M j + m i M j s i M j ,
here, l i M j is the length in cm from the i wire to the j module, d i M j is the initial distance, m i M j is the measured digital value, and s i M is the scale factor in digital units per cm.
We present a rendered image with a scaled 175 cm model in Figure 23.

3. Results

We organize this section as follows: first, we show the simulation; second, the final prototype; third, the trilateration and axis orientation; and finally, an ankle manifold representation.

3.1. Simulation Results

In this subsection, we use different values from Table 1 to estimate the work-space and range of motion. First, we show the variation of mean value results, and second the platform position simulation by changing the range of movement and angles.

Changing Statistical Mean Values

Figure 24a shows the complete manifold, taking into account the intervals θ 1 , θ 2 180 , 180 . It also shows the platform’s initial position, the TC axis reference, the initial ST reference, the initial orientation, and a parametric trajectory with equal angle rate variation. In Figure 24b is the attaching point A simulation; Figure 24c,d depicts the simulations of B and C, respectively.
In Figure 25a we show the platform’ central point simulation with variations of 10% below the statistical mean values; Figure 25b shows the simulation changing 10% over the statistical mean values; Figure 26a is the attaching point A simulation adding the 10% mean values; and Figure 26b subtracts 10% of the mean values. Figure 27a,b are the results for the platform attaching point B. We show the results for the attaching point C in Figure 28a,b.
Finally, by changing the range of maximum and minimum angles, an example of the interactive simulation is in Figure 29a,b. We capture the view of the sliders and also show the simulation rendering result.

3.2. Final Prototype

In this section, we describe the results of the TM design, which are the assembled device and calibration. We try several designs and finally the CAD model is in [44]. First, we show images of the connected electronics parts. Second, we assemble the structure and perform calibrations. Third, we probe the device in a healthy patient to validate the prototype adaptability. We print the structural parts using ABS and the draw-wire sensor using PLA; PETG is in the supports and the case.

3.2.1. Printed and Connected Electronics

We place the electronics in each side. In Figure 30, the connections and box sides and charge of the batteries.

3.2.2. Printed and Assembled Structure

We assemble all structural components carefully, putting them together with stainless-steel threaded rods; then we place the draw-wire sensors, the acquisition board, connections, and final structure for calibration. Figure 31 shows the assembly.

3.2.3. Calibration Results

We calibrate the system by using a personal computer. The resulting calibration, and measures of the lengths, are in Figure 32. The lecture is at the initial position, then we compare with the SolidWorks® model measurements and the Vernier caliper real measurements for each DWS. The Table 2 shows the calibration results.
Figure 33a shows the length with a SolidWorks® Measurement tool for module A, sensor 1; the lecture for sensor 2 is in Figure 33b. In Figure 33c, is the sensor 3 length. Table 3 shows the error measured in the real prototype and in SolidWorks®.

3.3. Trilateration Results

In this section, we use the measurements from the sensors to compute trilateration, then we compare them with the simulation results. The foot and shank fit in the adjustable platform and support structure, respectively, as is shown in the initial procedure in Figure 13. By introducing the DWS lengths to the virtual model, we compute the A, B and C coordinates in four consecutive positions. In the Table 4 are the seven sensors lenghts, and in the Table 5, we show the A, B and C coordinates for the four positions. The resulting figures for the first two positions are in Figure 34a–b, and for the latest two positions in Figure 35a,b. We show the base triangles, the points, the sensors, the platform and the circles on the base.

3.4. TC Axis Circle Fitting

The results of circle fitting for trajectories A, B, and C are in Table 6, corresponding to ankle joint plantar/dorsiflexion movements. We show the circle fitting for trajectories A, B, C, and PM in the Figure 36a–d.

3.5. ST Axis Circle Fitting

The results of ST circle fitting for trajectories A, B, C, and PM are in Table 7, corresponding to ankle joint inversion movements. We show the circle fitting for trajectories A, B, C, and PM in the Figure 37a–d.

3.6. Ankle Manifold Representation

In this section, we show the results in the software SageMath Manifolds. We load the model and visualize it as a manifold, we show the axis and the sagittal plane intersection. With the model parameters loaded, r 1 , r 2 , ω 1 , ω 2 , and the origin established in the center of the base modules. We apply the equation:
r ^ 1 = c ¯ 0 + n ^ p · d
where c ¯ 0 is the median center computed from trajectories A, B, and C center fitting, and n ^ p is the median planes’ normal vectors containing the circles. Table 8 shows values for the TC axis in the PM chart. In Table 9, we show the Plucker coordinates for the TC and ST axes.
Finally, Figure 38a shows the ankle manifold, and Figure 38b, the chart representing the range of movement and angle coordinates.

4. Discussion

In this work, we addressed the human ankle joint model from an alternative approach. We used statistical measurements for the development of a new device, specially designed to capture the human ankle joint movements. In animal joints, it is difficult to place encoders and linear sensors to measure the range of movement of complex joints in each internal living tissue reference frame. The product of exponential formulas uses only two frames, and it is useful in this case. Furthermore, in our work, we used a trilateration method for finding the device’s platform position, which is an analytic method. Therefore we avoid numerical approximations that can diverge and reduce rounding errors. We proposed the ankle joint model as a Riemannian manifold. We can define a chart as a subset of such a manifold with angle coordinates for measuring the range of movement. Our presented device is lightweight, non-invasive, and can be used in remote places, on beds, or on the floor. By characterizing the ankle parameters, we can conduct symmetry studies by correlating the left and right ankle joints. We can enhance the device configuration in future versions by replacing the draw-wire sensors used from potentiometers to digital encoders connected by a CAN bus, reducing wiring, space, weight, and energy consumption. We will use the model for the synthesis and reconfiguration of an ankle parallel rehabilitation robot, programmed by symmetrical movements at the opposite ankle. By employing the axis location and the screw theory, forces, and torques, we will study the ankle dynamics by using reciprocal screws to the axis location in a re-configurable platform. The robot will be lightweight because of the use of cable-driven actuators, inspired by antagonistic muscles that work with reciprocal inhibition for energy optimization. The robot will reconfigure the structure, considering the ankle joint as a central mast, and referenced it with MMP and MLP markers.
Figure 39 shows a schematic of the re-configurable approach.
Other applications are, for example, by visualizing the platform trajectories one can explain how the calcaneal Achilles insertion is near to the platform’s A point. The platform’s normal vector changes abruptly near this region, as was depicted in Figure 24b and Figure 26a,b. Furthermore, Riemannian models have different properties. We will explore diagnosis and treatments based on the model and metrics by employing machine learning algorithms. This approach can be applied to other joints in humans and other animals, by designing specialized re-configurable hardware and software. Tracking the parameters in different ages and weight conditions, and comparing the ankle models in healthy and injured people.

5. Conclusions

Computer tomography (CT) and magnetic resonance (MR) images have greater precision and accuracy. Measurements in medical imaging will help us compare the errors (RMS) in the HAJ. In biomechanics, we have not found an ideal model for error comparison. Then, we will compare the error with an accurate measurement. The device has limitations regarding mechanical precision and deformation of its parts. We face up to the error through the electronic design system. The calibration process is imperative for enhancing accuracy.
The calibration process is human-dependent. We read the digital measurement and compare it with caliper measurements directly in the sensor. Then, we register the data in a table to find the equivalence. An electronic board with trimmers avoids saturation, bias, and calibration; a 10 bit ADC and an exponentially weighted moving average (EWMA) filter the noise signals. We have implemented a processing (Software) calibration interface. We avoid adding more specific technical data, such as CMRR, ADC speed, mechanical tolerance, and other issues inherent to the measuring devices.
Digital sensors, communications, and POE function fitting use machine learning techniques.
The ankle is the most commonly injured joint of the lower limb, fundamental to the human body’s balance; it is necessary to measure the range of motion by in vivo methods for patients in lying positions in reduced or remote places. The device’s development considers ankle anatomy and anthropometry. We propose a Riemannian manifold model based on the device’s data readings. Performing simulations enabled us to design the size of the device and the maximal length of the wires. We present a trilateration algorithm, projecting the tetrahedron’s sides on the base plane. The sensors are modular and part of the device’s lightweight and portable structure. The electronic system is modular, replaced by other single-board computers (SBC) and microcontroller unities. We will also use the TM for ankle characterization and diagnosis for rehabilitation robotics, prosthesis, and orthosis design. The prototype is not a finished product (the TRL is 2). The work’s scope is to validate the use of a modern alternative biomechanic representation of the human ankle joint. It is a platform for testing an alternative trilateration method that employs draw-wire sensors (DWS). Such sensors have a constant tension, coiled on a drum attached to a potentiometer, and a flat spiral spring. We also attempted to develop a flexible device design for several foot sizes. We are working on a newer device version with an enhanced attachment system, a more compact design, and digital DWS compatible with a configurable robot. Machine learning and edge computing will assist in disease diagnosis and rehabilitation of patients.

Author Contributions

Conceptualization, J.V.-R., Á.V. and Ó.A.-V.; methodology, Á.V. and Ó.A.-V.; software, J.V.-R.; validation, Á.V. and Ó.A.-V.; formal analysis, J.V.-R.; investigation, J.V.-R.; resources, Á.V. and Ó.A.-V.; data curation, Á.V. and Ó.A.-V.; writing—original draft preparation, J.V.-R.; writing—review and editing, Á.V. and Ó.A.-V.; visualization, J.V.-R.; supervision, Á.V.; project administration, Á.V.; funding acquisition, J.V.-R., Á.V. and Ó.A.-V. All authors have read and agreed to the published version of the manuscript.

Funding

This research was partially funded by Colciencias-Colfuturo PhD Scholarships Program Educational Credit Forgivable grant number 568, and by Vicerrectorado de Investigación de la Universitat Politècnica de València (PAID-11-21).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Informed consent was obtained from all subjects involved in the study.

Data Availability Statement

The code and CAD electronics and mechanical designs are available.

Acknowledgments

The authors thank the Colfuturo Colciencias Collaboration for supporting this work, as well as the Universitat Politècnica de València and the Universidad de los Llanos.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

The following abbreviations are used in this manuscript:
HAJHuman Ankle Joint
ISBInternational Society of Biomechanics
DWSDraw-Wire Sensors
IMUInertial Measurement Units
PoEProduct of Exponentials
DoFDegrees of Freedom
RoMRange of Motion
BMSBattery Management System

References

  1. Krähenbühl, N.; Horn-Lang, T.; Hintermann, B.; Knupp, M. The subtalar joint. EFORT Open Rev. 2017, 2, 309–316. [Google Scholar] [CrossRef] [PubMed]
  2. Nichols, J.A.; Roach, K.E.; Fiorentino, N.M.; Anderson, A.E. Predicting tibiotalar and subtalar joint angles from skin-marker data with dual-fluoroscopy as a reference standard. Gait Posture 2016, 49, 136–143. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  3. Xie, S.S. Kinematic and Computational Model of Human Ankle. In Advanced Robotics for Medical Rehabilitation: Current State of the Art and Recent Advances; Xie, S.S., Ed.; Springer Tracts in Advanced Robotics; Springer International Publishing: Cham, Switzerland, 2016; pp. 185–221. [Google Scholar] [CrossRef]
  4. Jastifer, J.R.; Gustafson, P.A. The subtalar joint: Biomechanics and functional representations in the literature. Foot 2014, 24, 203–209. [Google Scholar] [CrossRef]
  5. Van Alsenoy, K.; De Schepper, J.; Santos, D.; Vereecke, E.; D’Août, K. The Subtalar Joint Axis Palpation Technique: Part 1—Validating a Clinical Mechanical Model. J. Am. Podiatr. Med. Assoc. 2014, 104, 238–246. [Google Scholar] [CrossRef] [PubMed]
  6. Van Alsenoy, K.K.; D’Août, K.; Vereecke, E.E.; De Schepper, J.; Santos, D. The Subtalar Joint Axis Palpation Technique: Part 2: Reliability and Validity Results Using Cadaver Feet. J. Am. Podiatr. Med. Assoc. 2014, 104, 365–374. [Google Scholar] [CrossRef]
  7. De Schepper, J.; Van Alsenoy, K.; Rijckaert, J.; De Mits, S.; Lootens, T.; Roosen, P. Intratest reliability in determining the subtalar joint axis using the palpation technique described by K. Kirby. J. Am. Podiatr. Med. Assoc. 2012, 102, 122–129. [Google Scholar]
  8. Parr, W.C.H.; Chatterjee, H.J.; Soligo, C. Calculating the axes of rotation for the subtalar and talocrural joints using 3D bone reconstructions. J. Biomech. 2012, 45, 1103–1107. [Google Scholar] [CrossRef] [Green Version]
  9. Leitch, J.; Stebbins, J.; Zavatsky, A.B. Subject-specific axes of the ankle joint complex. J. Biomech. 2010, 43, 2923–2928. [Google Scholar] [CrossRef]
  10. Lewis, G.S.; Cohen, T.L.; Seisler, A.R.; Kirby, K.A.; Sheehan, F.T.; Piazza, S.J. In vivo tests of an improved method for functional location of the subtalar joint axis. J. Biomech. 2009, 42, 146–151. [Google Scholar] [CrossRef]
  11. Lewis, G.S.; Kirby, K.A.; Piazza, S.J. Determination of subtalar joint axis location by restriction of talocrural joint motion. Gait Posture 2007, 25, 63–69. [Google Scholar] [CrossRef]
  12. Spooner, S.; Kirby, K. The subtalar joint axis locator: A preliminary report. J. Am. Podiatr. Med. Assoc. 2006, 96, 212–219. [Google Scholar] [CrossRef] [PubMed]
  13. Kirby, K.A. Subtalar Joint Axis Location and Rotational Equilibrium Theory of Foot Function. J. Am. Podiatr. Med. Assoc. 2001, 91, 465–487. [Google Scholar] [CrossRef] [PubMed]
  14. Dul, J.; Johnson, G.E. A kinematic model of the human ankle. J. Biomed. Eng. 1985, 7, 137–143. [Google Scholar] [CrossRef]
  15. Wu, G.; Siegler, S.; Allard, P.; Kirtley, C.; Leardini, A.; Rosenbaum, D.; Whittle, M.; D’Lima, D.D.; Cristofolini, L.; Witte, H.; et al. ISB recommendation on definitions of joint coordinate system of various joints for the reporting of human joint motion—Part I: Ankle, hip, and spine. J. Biomech. 2002, 35, 543–548. [Google Scholar] [CrossRef]
  16. Mann, R.A. Biomechanics of the Ankle. In Joint Surgery Up to Date; Hirohata, K., Kurosaka, M., Cooke, T.D.V., Eds.; Springer: Tokyo, Japan, 1989; pp. 73–81. [Google Scholar] [CrossRef]
  17. Winter, D.A. Biomechanics and Motor Control of Human Movement; John Wiley & Sons: Hoboken, NJ, USA, 2009. [Google Scholar]
  18. Dawe, E.J.C.; Davis, J. (vi) Anatomy and biomechanics of the foot and ankle. Orthop. Trauma 2011, 25, 279–286. [Google Scholar] [CrossRef]
  19. Coughlin, M.J.; Saltzman, C.L.; Mann, R.A. Mann’s Surgery of the Foot and Ankle E-Book: Expert Consult-Online; Elsevier Health Sciences: Amsterdam, The Netherlands, 2013. [Google Scholar]
  20. Delp, S.; Loan, J.; Hoy, M.; Zajac, F.; Topp, E.; Rosen, J. An interactive graphics-based model of the lower extremity to study orthopaedic surgical procedures. IEEE Trans. Biomed. Eng. 1990, 37, 757–767. [Google Scholar] [CrossRef]
  21. Donatelli, R. Normal Biomechanics of the Foot and Ankle. J. Orthop. Sport. Phys. Ther. 1985, 7, 91–95. [Google Scholar] [CrossRef]
  22. Bähler, A. The biomechanics of the foot. Clin. Prosthetics Orthot. 1986, 10, 8–14. [Google Scholar]
  23. Lundberg, A.; Svensson, O.K.; Németh, G.; Selvik, G. The axis of rotation of the ankle joint. J. Bone Jt. Surg. Br. Vol. 1989, 71, 94–99. [Google Scholar] [CrossRef]
  24. Singh, A.K.; Starkweather, K.D.; Hollister, A.M.; Jatana, S.; Lupichuk, A.G. Kinematics of the Ankle: A Hinge Axis Model. Foot Ankle Int. 1992, 13, 439–446. [Google Scholar] [CrossRef]
  25. Leardini, A.; O’Connor, J.; Catani, F.; Giannini, S. A geometric model of the human ankle joint. J. Biomech. 1999, 32, 585–591. [Google Scholar] [CrossRef]
  26. Brockett, C.L.; Chapman, G.J. Biomechanics of the ankle. Orthop. Trauma 2016, 30, 232–238. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  27. Bruening, D.; Richards, J. Skiing-Skating: Optimal ankle axis position for articulated boots. Sport. Biomech. 2005, 4, 215–225. [Google Scholar] [CrossRef] [PubMed]
  28. Andrade-Cetto, J.; Thomas, F. A Wire-Based Active Tracker. IEEE Trans. Robot. 2008, 24, 642–651. [Google Scholar] [CrossRef] [Green Version]
  29. Thomas, F.; Ottaviano, E.; Ros, L.; Ceccarelli, M. Coordinate-free formulation of a 3-2-1 wire-based tracking device using Cayley-Menger determinants. In Proceedings of the 2003 IEEE International Conference on Robotics and Automation (Cat. No.03CH37422), Taipei, Taiwan, 14–19 September 2003; Volume 1, p. 355. [Google Scholar] [CrossRef] [Green Version]
  30. Thomas, F.; Ros, L. Revisiting trilateration for robot localization. IEEE Trans. Robot. 2005, 21, 93–101. [Google Scholar] [CrossRef] [Green Version]
  31. Salleh, S.; Rahmat, M.F.; Othman, S.M.; Abidin, H.Z. Application of draw wire sensor in the tracking control of an electro hydraulic actuator system. J. Teknol. 2015, 73, 51–57. [Google Scholar] [CrossRef] [Green Version]
  32. Jiafan, Z.; Jinsong, L.; Liwei, Q.; Dandan, Z. Kinematic analysis of a 6-DOF wire-based tracking device and control strategy for its application in robot easy programming. In Proceedings of the 2009 IEEE International Conference on Robotics and Biomimetics (ROBIO), Guilin, China, 9–23 December 2009; pp. 1591–1596. [Google Scholar] [CrossRef]
  33. Bulling, A.; Blanke, U.; Schiele, B. A Tutorial on Human Activity Recognition Using Body-worn Inertial Sensors. ACM Comput. Surv. 2014, 46, 1–33. [Google Scholar] [CrossRef]
  34. Chermak, L.; Aouf, N.; Richardson, M.A.; Visentin, G. Real-Time Smart and Standalone Vision/IMU Navigation Sensor; Springer Science and Business Media LLC.: Berlin/Heidelberg, Germany, 2016. [Google Scholar] [CrossRef] [Green Version]
  35. Ong, Z.C.; Noroozi, S. Development of an Economic Wireless Human Motion Analysis Device for Quantitative Assessment of Human Body Joint; Elsevier BV: Amsterdam, The Netherlands, 2018. [Google Scholar] [CrossRef] [Green Version]
  36. Porciuncula, F.; Roto, A.V.; Kumar, D.; Davis, I.; Roy, S.; Walsh, C.J.; Awad, L.N. Wearable Movement Sensors for Rehabilitation: A Focused Review of Technological and Clinical Advances. PM&R 2018, 10, S220–S232. [Google Scholar] [CrossRef] [Green Version]
  37. Wahyudi, W.; Listiyana, M.S.; Sudjadi, S.; Ngatelan, N. Tracking Object based on GPS and IMU Sensor. In Proceedings of the 5th International Conference on Information Technology, Computer, and Electrical Engineering (ICITACEE), Semarang, Indonesia, 26–28 September 2018; pp. 214–218. [Google Scholar] [CrossRef]
  38. Gregorio, R.D.; Parenti-Castelli, V.; O Connor, J.J.; Leardini, A. Mathematical models of passive motion at the human ankle joint by equivalent spatial parallel mechanisms. Med. Biol. Eng. Comput. 2007, 45, 305–313. [Google Scholar] [CrossRef]
  39. Vargas-Riaño, J.H. Turmell-Meter. 2022. Available online: https://github.com/juliohvr/Turmell-Meter (accessed on 2 February 2022).
  40. Isman, R.E.; Inman, V.T.; Poor, P.M. Anthropometric studies of the human foot and ankle. Bull. Prosthet. Res. 1969, 11, 97–129. [Google Scholar]
  41. Drillis, R.; Contini, R.; Bluestein, M. Body Segment Parameters; New York University, School of Engineering and Science: New York, NY, USA, 1966. [Google Scholar]
  42. Hebbelinck, M.; Ross, W.D. Kinanthropometry and biomechanics. In Biomechanics IV; Springer: Berlin/Heidelberg, Germany, 1974; pp. 535–552. [Google Scholar]
  43. Fryar, C.D.; Carroll, M.D.; Gu, Q.; Afful, J.; Ogden, C.L. Anthropometric Reference Data for Children and Adults: United States, 2015–2018; Vital & Health Statistics. Series 3 Analytical and Epidemiological Studies; CDC: Atlanta, GA, USA, 2021; pp. 1–44. [Google Scholar]
  44. Vargas Riaño, J.; Valera, Á.; Agudelo Varela, O. Turmell-Metre|3D CAD Model Library|GrabCAD. Available online: https://grabcad.com/library/turmell-metre-1 (accessed on 2 February 2022).
Figure 1. Reference points from anthropometric values K, L, O, and P.
Figure 1. Reference points from anthropometric values K, L, O, and P.
Bioengineering 09 00199 g001
Figure 2. Q, W, and w distances from lateral and transverse views.
Figure 2. Q, W, and w distances from lateral and transverse views.
Bioengineering 09 00199 g002
Figure 3. Mean relative position of the ST and TC axis.
Figure 3. Mean relative position of the ST and TC axis.
Bioengineering 09 00199 g003
Figure 4. Planes, axes, and points of corresponding references.
Figure 4. Planes, axes, and points of corresponding references.
Bioengineering 09 00199 g004
Figure 5. Vectors and points on the sagittal plane.
Figure 5. Vectors and points on the sagittal plane.
Bioengineering 09 00199 g005
Figure 6. Forward kinematics for (a) initial position and (b) θ 1 r a n g e = θ 2 r a n g e = [ 15 , 15 ] , θ 1 = θ 2 = 10 and (c) θ 1 r a n g e = θ 2 r a n g e = [ 10 , 10 ] , θ 1 = θ 2 = 5 .
Figure 6. Forward kinematics for (a) initial position and (b) θ 1 r a n g e = θ 2 r a n g e = [ 15 , 15 ] , θ 1 = θ 2 = 10 and (c) θ 1 r a n g e = θ 2 r a n g e = [ 10 , 10 ] , θ 1 = θ 2 = 5 .
Bioengineering 09 00199 g006
Figure 7. Geometric design: (a) is the platform center, base, and r 1 , r 2 ; and (b) platform vertices with talocrural and subtalar axis.
Figure 7. Geometric design: (a) is the platform center, base, and r 1 , r 2 ; and (b) platform vertices with talocrural and subtalar axis.
Bioengineering 09 00199 g007
Figure 8. Geometric design of the DWS arrays.
Figure 8. Geometric design of the DWS arrays.
Bioengineering 09 00199 g008
Figure 9. Tetrahedron trilateration flowchart.
Figure 9. Tetrahedron trilateration flowchart.
Bioengineering 09 00199 g009
Figure 10. Finding the apex A p .
Figure 10. Finding the apex A p .
Bioengineering 09 00199 g010
Figure 11. Rotation of α angle about the axis B 1 B 3 : (a) original tetrahedron, T B (b) rotated tetrahedron.
Figure 11. Rotation of α angle about the axis B 1 B 3 : (a) original tetrahedron, T B (b) rotated tetrahedron.
Bioengineering 09 00199 g011
Figure 12. Finding the apex B pr .
Figure 12. Finding the apex B pr .
Bioengineering 09 00199 g012
Figure 13. Adjusting the foot, the shank and the Most Medial Point reference.
Figure 13. Adjusting the foot, the shank and the Most Medial Point reference.
Bioengineering 09 00199 g013
Figure 14. Computed positions from sensor lengths at portable configuration: (a) the rest position, (b) apex A, (c) apex B, and (d) apex C.
Figure 14. Computed positions from sensor lengths at portable configuration: (a) the rest position, (b) apex A, (c) apex B, and (d) apex C.
Bioengineering 09 00199 g014
Figure 15. Simulation of the platform central point with variations in the mean statistical values: (a) platform’s center point manifold, (b) manifold chart and a geodesic.
Figure 15. Simulation of the platform central point with variations in the mean statistical values: (a) platform’s center point manifold, (b) manifold chart and a geodesic.
Bioengineering 09 00199 g015
Figure 16. Draw-wire sensor design.
Figure 16. Draw-wire sensor design.
Bioengineering 09 00199 g016
Figure 17. Mechanical attachment: (a) calf support and (b) aluminum tube structure.
Figure 17. Mechanical attachment: (a) calf support and (b) aluminum tube structure.
Bioengineering 09 00199 g017
Figure 18. Base and platform: (a) DWS modules support (b) platform with foot’s size adjustment.
Figure 18. Base and platform: (a) DWS modules support (b) platform with foot’s size adjustment.
Bioengineering 09 00199 g018
Figure 19. Two Op. Amp. instrumentation amplifier.
Figure 19. Two Op. Amp. instrumentation amplifier.
Bioengineering 09 00199 g019
Figure 20. Power system with backup, BMS, boost, and buck converters.
Figure 20. Power system with backup, BMS, boost, and buck converters.
Bioengineering 09 00199 g020
Figure 21. Modular electronics casing.
Figure 21. Modular electronics casing.
Bioengineering 09 00199 g021
Figure 22. Complete prototype.
Figure 22. Complete prototype.
Bioengineering 09 00199 g022
Figure 23. Rendered image with a 175 cm height patient.
Figure 23. Rendered image with a 175 cm height patient.
Bioengineering 09 00199 g023
Figure 24. Simulation of all points: (a) platform’s central point, (b) attachment a, (c) attachment b, and (d) attachment c.
Figure 24. Simulation of all points: (a) platform’s central point, (b) attachment a, (c) attachment b, and (d) attachment c.
Bioengineering 09 00199 g024
Figure 25. Simulation of the platform central point with variations in the mean statistical values: (a) 10% below, and (b) 10% over.
Figure 25. Simulation of the platform central point with variations in the mean statistical values: (a) 10% below, and (b) 10% over.
Bioengineering 09 00199 g025
Figure 26. Simulation of the platform’s attaching point A: (a) mean values plus 10%, (b) mean values minus 10%.
Figure 26. Simulation of the platform’s attaching point A: (a) mean values plus 10%, (b) mean values minus 10%.
Bioengineering 09 00199 g026
Figure 27. Attaching point B simulation: (a) adding 10% to the statistic mean values, (b) subtracting 10%.
Figure 27. Attaching point B simulation: (a) adding 10% to the statistic mean values, (b) subtracting 10%.
Bioengineering 09 00199 g027
Figure 28. Simulation results for C: (a) mean values plus 10%, (b) mean values minus 10%.
Figure 28. Simulation results for C: (a) mean values plus 10%, (b) mean values minus 10%.
Bioengineering 09 00199 g028
Figure 29. Interactive simulation example: (a) sliders, (b) rendering.
Figure 29. Interactive simulation example: (a) sliders, (b) rendering.
Bioengineering 09 00199 g029
Figure 30. Connections and electronics.
Figure 30. Connections and electronics.
Bioengineering 09 00199 g030
Figure 31. Assembled structure.
Figure 31. Assembled structure.
Bioengineering 09 00199 g031
Figure 32. Processing calibration interface.
Figure 32. Processing calibration interface.
Bioengineering 09 00199 g032
Figure 33. Measuring in SolidWorks (2017–2018 Student Edition, Dassault Systèmes, Vélizy-Villacoublay, France)®: (a) sensor 1, (b) sensor 2, (c) sensor 3.
Figure 33. Measuring in SolidWorks (2017–2018 Student Edition, Dassault Systèmes, Vélizy-Villacoublay, France)®: (a) sensor 1, (b) sensor 2, (c) sensor 3.
Bioengineering 09 00199 g033
Figure 34. First two trilateration results: (a) position 1, (b) position 2.
Figure 34. First two trilateration results: (a) position 1, (b) position 2.
Bioengineering 09 00199 g034
Figure 35. Latest two trilateration results: (a) position 3, (b) position 4.
Figure 35. Latest two trilateration results: (a) position 3, (b) position 4.
Bioengineering 09 00199 g035
Figure 36. TC axis circle fitting: (a) trajectory A, (b) trajectory B, (c) trajectory C, (d) trajectory PM.
Figure 36. TC axis circle fitting: (a) trajectory A, (b) trajectory B, (c) trajectory C, (d) trajectory PM.
Bioengineering 09 00199 g036
Figure 37. ST axis circle fitting: (a) trajectory A, (b) trajectory B, (c) trajectory C, (d) trajectory PM.
Figure 37. ST axis circle fitting: (a) trajectory A, (b) trajectory B, (c) trajectory C, (d) trajectory PM.
Bioengineering 09 00199 g037
Figure 38. Ankle joint manifold. (a) Manifold for PM, (b) chart with ankle axis coordinates.
Figure 38. Ankle joint manifold. (a) Manifold for PM, (b) chart with ankle axis coordinates.
Bioengineering 09 00199 g038
Figure 39. Re-configurable cable-driven robot concept.
Figure 39. Re-configurable cable-driven robot concept.
Bioengineering 09 00199 g039
Table 1. Mean values of anthropometric measurements.
Table 1. Mean values of anthropometric measurements.
VariableK (cm)L (cm)O (cm)P (cm)Q (cm)R = W/w
Mean1.2 cm1.1 cm1.6 cm1.0 cm0.5 cm0.54 cm
Table 2. Calibration results with digital measurements and real measurements.
Table 2. Calibration results with digital measurements and real measurements.
Measurementsl1M1l2M1l3M1l1M2l2M2l1M3l2M3
BCD value239330246265177252242
Vernier Caliper, cm8.0 cm5.3 cm6.9 cm13.0 cm8.4 cm7.8 cm11.5 cm
Table 3. Error compared with SolidWorks® measurements.
Table 3. Error compared with SolidWorks® measurements.
Measurementsl1M1l2M1l3M1
Measured distance7.622 cm5.33 cm6.384 cm
Error in cm0.38 cm−0.030 cm0.52 cm
Table 4. Sensor measurements in four different positions.
Table 4. Sensor measurements in four different positions.
Positionsl1M1l2M1l3M1l1M2l2M2l1M3l2M3
Pos1., cm11.0 cm12.6 cm12.5 cm14.8 cm10.8 cm15.2 cm11.9 cm
Pos2., cm10.2 cm11.7 cm11.6 cm15.2 cm11.3 cm15.5 cm12.2 cm
Pos3., cm9.40 cm10.8 cm10.8 cm15.6 cm11.7 cm15.8 cm12.5 cm
Pos4., cm8.56 cm9.89 cm9.95 cm16.0 cm12.2 cm16.0 cm12.7 cm
Table 5. A, B and C coordinates computed from the four positions.
Table 5. A, B and C coordinates computed from the four positions.
PositionsABC
Pos1., cm(−11.7, −1.06, −11.0) cm(6.11, −9.77, −8.76) cm(5.54, 8.81, −9.35) cm
Pos2., cm(−12.1, −0.93, −10.2) cm(5.62, −9.92, −9.37) cm(4.83, 9.46, −10.1) cm
Pos3., cm(−12.4, −0.65, −9.39) cm(5.27, −9.79, −9.68) cm(4.94, 9.03, −10.2) cm
Pos4., cm(−12.7, −0.48, −8.53) cm(4.68, −10.0, −10.3) cm(3.54, 10.7, −11.1) cm
Table 6. TC axis circle fitting.
Table 6. TC axis circle fitting.
TrajectoryCenterDirectionRadius
A(0.08649, 2.138, −6.712) cm(−0.089, −0.95, 0.31)7.666
B(0.5713, 5.531, −7.824) cm(−0.089, −0.95, 0.31)5.246 cm
C(−0.2442, −2.669, −5.315) cm(−0.089, −0.95, 0.31)7.206 cm
PM(0.1552, 1.642, −6.683) cm(−0.089, −0.95, 0.31)5.375 cm
Table 7. ST axis circle fitting.
Table 7. ST axis circle fitting.
TrajectoryCenterDirectionRadius
A(4.444, 1.825, −9.008) cm(−0.75, −0.28, 0.60)2.428 cm
B(1.757, 0.6768, −6.925) cm(−0.75, −0.28, 0.60)6.567 cm
C(0.1578, 0.1819, −5.807) cm((−0.75, −0.28, 0.60)6.935 cm
PM(2.087, 0.8882, −7.281) cm(−0.75, −0.28, 0.60)3.875
Table 8. Axis estimation data.
Table 8. Axis estimation data.
AxisMedian CenterMedian Normalr ω
TC(1.92, 0.783, −7.10) cm(−0.750, −0.280, 0.600)(−0.174, 0.000, −5.43) cm(−0.750, −0.280, 0.600)
ST(0.121, 1.89, −6.70) cm(−0.0890, −0.950, 0.310)(−0.0562, 0.000, −6.08) cm(−0.0890, −0.950, 0.310)
Table 9. Plucker line coordinates.
Table 9. Plucker line coordinates.
AxisPlucker Line Coordinates
TC[−0.750: −0.280: 0.600: −1.52: 4.17: 0.0487]
ST[−0.0890: −0.950: 0.310: −5.78: 0.559: 0.0534]
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Agudelo-Varela, Ó.; Vargas-Riaño, J.; Valera, Á. Turmell-Meter: A Device for Estimating the Subtalar and Talocrural Axes of the Human Ankle Joint by Applying the Product of Exponentials Formula. Bioengineering 2022, 9, 199. https://doi.org/10.3390/bioengineering9050199

AMA Style

Agudelo-Varela Ó, Vargas-Riaño J, Valera Á. Turmell-Meter: A Device for Estimating the Subtalar and Talocrural Axes of the Human Ankle Joint by Applying the Product of Exponentials Formula. Bioengineering. 2022; 9(5):199. https://doi.org/10.3390/bioengineering9050199

Chicago/Turabian Style

Agudelo-Varela, Óscar, Julio Vargas-Riaño, and Ángel Valera. 2022. "Turmell-Meter: A Device for Estimating the Subtalar and Talocrural Axes of the Human Ankle Joint by Applying the Product of Exponentials Formula" Bioengineering 9, no. 5: 199. https://doi.org/10.3390/bioengineering9050199

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop