The Turmell-Meter: A Device for Ankle Joint Axis Estimation Appliyng Product of Exponentials

The human ankle is a complex joint, most commonly represented as talocrural and 1 subtalar axes. It is difficult to locate and take in vivo measurements of the ankle joint. There are 2 no instruments for patients lying on a bed or the floor; that can be used in outdoor or remote 3 sites. We have developed a "Turmell-meter" to address these issues. We started with the study of 4 ankle anatomy and anthropometry, then we used the product of exponentials’ formula to visualize 5 the movements. Furthermore, we built a prototype using human proportions and statistics. For 6 pose estimation, we used a trilateration method by applying tetrahedral geometry. Additionally, 7 we computed the axis direction by fitting 3D circles, plotting the manifold and chart as an ankle 8 joint model. We presented the results of simulations, a prototype comprising 45 parts, specifically 9 designed draw-wire sensors, and electronics. Finally, we tested the device by capturing positions 10 and fitting them into the bi-axial ankle model as a Riemannian manifold. The Turmell-meter is 11 intended to be a hardware platform for human ankle joint axis estimation, it is adjustable and 12 has an easy setup. The proposed model has the properties of a chart in a geometric manifold, we 13 provided the details 14


Introduction
In this work, we present a device intended for the study of the human ankle joint 20 (HAJ). Modeling and measuring this lower limb joint is essential in physiology, biome- 21 chanics, and rehabilitation (also in humanoid robotic limb development). 22 23 The HAJ is fundamental for human locomotion. And the ankle joint sprain is 24 the most common lower limb injury in sports [1], football [2], basketball [3], high school 25 sports [4], military academy [5], service [6] and physical training [7]. It causes chronic 26 ankle instability [8], and is costly for the healthcare system [9]. Treatment and healing 27 require measuring the range of motion. 28 29 Similar to other characteristics, the HAJ model is unique for each person. Individual 30 variations and anthropometric measurements depend on gender, age, and phenotype. 31 There are few types of equipment for in vivo ankle joint measurements on reduced spaces 32 such as beds or for patients laying on the floor in remote places. 33 34 85 For the simulation with the PoE formula, we adapt the data from [82], proportions 86 from [83,84], and statistics from [85]. 87 2.1.1. References Assignation 88 Figure 1 presents the reference points and the mean distances taken from [82].    Table 1. Mean values of anthropometric measurements .

97
We define the sagittal (lateral) plane as the X-Z plane (perpendicular to the y-98 axis). The coronal (frontal) plane is the Y-Z plane (x-axis is normal to it); the transverse 99 (axial) plane is the X-Y plane (perpendicular to the z-axis). Figure 4, left, shows this 100 corresponding references.  With this reference frame, we can define the TC axis orientation from a unitary 102 vector in the z-direction. We first rotate it -80°around the x-axis; then we turn it -6°103 around the z-axis. A unitary vector in the x-axis direction defines the ST axis, rotating 104 41°about the y-axis, followed by a 23°rotation around the z-axis. 105 We show the fibula, tibia, talus, calcaneus 3D position, reference points, TC, and ST 106 axes in Figure 4, right. 107 In this image, A 0 , B 0 , and C 0 are the vertices from the platform fixed to the foot, 108 and P M is the triangle's center. S 1 , S 2 , and S 3 are fixed to the shank relative to the origin 109 point P 0 . M 1 and M 2 define the TC axis; N 1 and N 2 correspond to the ST axis. We define 110 r 1 and r 2 as the sagittal plane intersection with the TC and ST axes. 111 112 We estimate the device dimensions from anthropometric proportions in [83] and 113 use the segments proportions shown in Fig. 5  We select the origin of coordinates between the knee and the ankle, d m is the 115 distance from P M to P O . This distance is proportional to the body's height H. To do so, 116 we define d m as follows:

Size and Dimensions
and, according to [85], the mean height H of an adult male is 175 cm; by substituting 118 this value into the equation, the knee-ankle distance is 28.35 cm. The distance d p12 119 between points r 1 and r 2 about the TC and ST axes on the sagittal plane is: the projection of the most medial point (MMP) on the sagittal plane is: and for the most lateral point is: also, M 2p is M 2 ; we estimate the projection from L and K through: Therefore, the segment M 2 M 1 has the sagittal projection M 2p M 1p ; it has the same 125 proportional relation R = W/w in respect to M 2p r 1 , then: solving for r 1 gives the following: By knowing the distance Q projected in the sagittal plane and r 1 , the angle 41°we The distance from the origin P O to the plantar surface of the foot is d m , we choose 130 a circumscribed equilateral triangle with vertices A 0 , B 0 , C 0 as the platform base. The 131 coordinates of A0 are for B 0 are: and for C 0 : Where r p is proportional to H, then: In summary, we estimate P 0 , r 1 , r 2 ; and the platform's vertices A 0 , B 0 , and C 0 . They 133 aren't arbitrary selected, on the contrary, we employed anthropometry, statistics, and 134 proportions. In this section, we employ the PoE formula. We follow the intuitive concept that 137 inter-bone contact surfaces determine HAJ movements. Therefore, we represent these 138 movements as a Special Euclidean group SE(3) in matrix form: where R 3×3 is the rotation matrix andp T is the translation vector.

140
For the initial point A 0 : for B 0 : and for C 0 142 Preprints (www.preprints.org) | NOT PEER-REVIEWED | Posted: 7 February 2022 doi:10.20944/preprints202202.0080.v1 We defineω 1 = (ω x1 , ω y1 , ω z1 ) as a unitary vector for the TC axis direction given 143 by: and a directed vectorr 1 from P O to r 1 is: then, an orthogonal vector tor 1 andω 1 is: together,ω 1 andν θ 1 r 2z compound the six-dimensional vectorξ 1 : In the same way, there are correspondent vectors for the TC axis: and: We compute R for each joint i = 1, 2 from the Rodrigues' formula: where Ω is the skew symmetric matrix: The exponential formula is: and, τ i is translation vector: Points A, B, and C have invariant relative positions, and there are two rotating 154 joints; the PoE formula for A is: and the PoE for C is: θ 1 is the TC rotation angle from the zero position, and θ 2 is the ST rotation from the 158 zero position. For the sake of clarity, we show the section of the ankle with the vectors 159r 1 ,ω 1 ,ν 1 andr 2 ; also the points A, B, C and P O in Figure 6.  The simulation plot for the platform's central point is in Figure 7(a). We show the 167 points P O , A 0 , B 0 , C 0 , r 1 , r 2 , and the surfaces representing each group of movements.

168
The forward kinematics with θ 1range = θ 2range = [− 15,15] and θ 1 = θ 2 = 10 • is in Figure   169 7(b). For θ 1range = θ 2range = [−10, 10] and θ 1 Figure 7 Here, l max is the maximal possible length from the triangular inequality, p A is the 179 positions group in g A , r m is the module's radius, and A B is the base point. 180 The main design requirement is the localization of three points attached to the foot. 181 We estimate the actual position employing an array of DWS in a tetrahedral structure to 182 find the apex, which is a platform vertex. In Figure 9 we show the design structure.   Figure 10 shows the method we use. In Figure 10  tetrahedron T A , we determine A pxy = (Apx, Apy, 0) as A p projection on the base plane.

199
It is easy to realize that the height of T A is the absolute value of the A pz coordinate. Then, 200 we can find the distance from A pxy to A 3 as a triangle △[Apxy, A 3 , A p ] side ; the other is 201 A pz , and the hypotenuse is the distance l A3 = ∥A p − A 3 ∥, then, A pz is:

203
In this subsection, we show that, by knowing A p , the point B p needs two sensors to 204 be found. To determine the result of the tetrahedron T B , we consider the Figure 11(a) the 205 base of a triangle △ B 1 , B 3 , A p . 206 We compute the angle α from the XY plane to a normal vectorn ApB : and, the angle α is: wheren z is the unitary vector normal to the XY plane. that of a tetrahedron T A . Finally, when Bpr is found, the contrary rotation about the axis 213 B1 − B3 gives the Bps.

214
There are two possible apex values: Bps1 over, and Bps2 below of the XY plane. We

215
show the Bpr 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 217 condition when the side of the platform distance d CpBp is:

219
In this section, we design DWS to measure the lengths of the tetrahedrons sides; 220 they are arranged as structural parts. Their maximal length estimation is from the 221 forward kinematics simulation. We design the shank attachment from the dimensions, 222 proportions and statistical data. We use flat springs. They aren't exposed to a high load against gravity, and are in 225 two or three concurrent groups. In Figure 13 we depict the design, composed of three 226 3D printed parts, potentiometer, flat spring, bolts, and nuts. All the DWS modules are in a plate, the A module has three DWS, B and C modules 241 has two DWS as in Figure 15(a). The design of the foot attachment is from standard 242 measurements to adjust the foot's length and width, as in Figure 15 The voltage gain in the instrumentation amplifier is: By selecting R 2 = 100kΩ, R 1 = 1kΩ, and R G = 5kΩ, the voltage gain is 141. With 249 34mV as voltage input, we get: The final acquisition circuit has seven instrumentation amplifiers, with bias and 251 gain trimmers for calibration. We design the printed board circuit as an ™Arduino Mega  Each box has attached components to optimize the space. We test every component, 265 and then install the support structure. The prototype consist of 45 3D printed parts, the union of main components is by an 268 8 mm steel threaded rod. The sub-assemblies uses M3 bolts and nuts. Figure 19 shows 269 the assembly CAD.
here, l iMj is the length in cm from the i wire to the j module, d iMj is the initial 279 distance, m iMj is the measured digital value, and s iM is the scale factor in digital units 280 per cm. 281 We present a rendered image with a scaled 175 cm model in the Figure 20.

299
The symbolic equations to find A p , B p and C p from the captured data, were found 300 by the SageMath CAS. By using the prototype dimensions, and the sensor lengths, we  After MLP and MMP registering, we attach the apex of the module A to the platform.

304
And define the sagittal plane perpendicular to the ABC base plane and intersecting the 305 A point. By implementing the trilateration method before mentioned, we compute the 306 points A 0 , B 0 , and C 0 . is r 1 . We estimate the initial TC axis approximation with: and from these values, we solve for r 1 from the intersection of the line: and the plane y = 0.

316
The ST axis sagittal intersection r 2 initial point is: This is a statistical initial estimation, we use for comparison with the orientation 318 captured by doing circle fitting to A, B and C trajectories. First, we found the TC axis or in all n data is in a n × 3 matrix, solving for z we have which has the form: there are more equations, the pseudo inverse is A + = (A T A) −1 A T , and the normal 324 vector is: By using the Rodrigues' formula we found the rotated points to the XY plane by where θ = arccos n·k ∥n∥ .

328
The circle in XY plane can be rearranged by where c = [c 0 , c 1 , c 2 ] T with c 0 = 2x c , c 1 = 2y c , and c 2 = r 2 − x 2 c − y 2 c .

330
By taking the rotated points, P r we have a linear system . .
has the form Again, we have more equations than unknowns, then, we search for the c values 332 that minimize the squared difference ∥b − Ac∥ 2 .
We found the center point C p = [x c , y c ] and radius r by solving: Finally, we apply a rotation to the center in respect to the original plane. This point 335 pertain to the TC axis, for each trajectory A, B, C we get three planes, and three centers, wherep is the parametric line point,l 0 is a known point in the line, and d ∈ R, 353 replacing the parametric equation in the plane equation: wherep 0 is a known point in the plane, andn p is the plane's normal vector. 355 solving for d, gives: and replacing in the TC axis line equation: where r 1 is the TC axis intersection with the sagittal plane. The point c 1 is the 358 center, and the axis direction ω 1 , both were found by circle fitting. Also, packing in six 359 dimensional Plücker line coordinates, we have: and the l 1 six dimensional vector is: We include those data for the PoE formula simulation and the manifold representa-362 tion.

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

386
In this section, we describe the results of the TM design, which are the assembled 387 device and calibration. We try several designs and finally the CAD model is in [86]. First, 388 we show images of the connected electronics parts. Second, we assemble the structure 389 and perform calibration. Third, the device calibration results and probe the device in a 390 healthy patient to validate the prototype adaptability. We print the structural parts using 391 ABS, the draw-wire sensor using PLA; PETG is in the supports and the case.  393 We place the electronics in each side. In Figure 30, then connect the box sides; and 394 charge the batteries.

396
We assemble all structural components carefully, and putting together with stainless-397 steel threaded rods; then we put the draw-wire sensors, the acquisition board, connec-398 tions, and final structure for calibration. The Figure 31 shows the assembly.  The Table 2 shows the calibration results. The Figure 33(a) shows the length with SolidWorks Measurement tool for module 405 A, sensor 1; the lecture for sensor 2 is in Figure 33(b). In Figure 33(c) is the sensor 3 406 length.

407
The Table 3 shows the error measured in real prototype and in SolidWorks.

409
In this section we use the measurements from the sensors to compute trilateration,

417
In this section, we show the results of model fitting to solve the axis orientation.

418
First, we show the TC axis, and then the ST axis. The results of circle fitting for trajectories A, B, and C are in the The The results of ST circle fitting for trajectories A, B, C, and PM are in the Table 6, 426 corresponding to ankle joint inversion movements. With the model parameters loaded, r 1 , r 2 , ω 1 , ω 2 , and the origin established in the center 433 of the base modules. We apply the equation: where c 0 is the median center computed from trajectories A, B and C center fitting, 435 andn p is the median planes normal vectors containing the circles. The Table 7 shows 436 values for the TC axis in PM chart. In the Table 8, we show the Plucker coordinates for the TC and ST axis.  Finally, the Figure 37(a) shows the ankle manifold, and the Figure 37        The following abbreviations are used in this manuscript: In this section, we start with the ankle description, which presents a complex 513 movement. First, we study the shank, ankle, and foot bones. Then, we analyzed the 514 ankle movements based on the anatomic spatial and functional representation.

516
We start with an understanding of inter-bone contact surfaces when studying ankle 517 movements. In Fig. A1, we identify the names of the bones of the left and right feet. In Fig. A2, we use the right-hand rotation convention and present the movements 519 systematically. Also, we organize those movements into two rows, corresponding to 520 pronation and supination. In addition, we show the hindfoot and midfoot are the most 521 involved segments in ankle movements.  The most accepted approach for modeling the ankle joint is the biaxial movement.

524
It results from the interaction of several bones, such as the fibula, tibia, talus, calcaneus, 525 navicular, cuboid, and three cuneiform bones. As shown in Fig. A3, the first axis 526 corresponds to the rotation from the talus regarding the tibia-fibula fixed joint.

Listing 1: Sagemath Forward Kinematics
For visualization and interactive view in ®Acrobat Reader, we exported the symbolic code to Asymptote, the listing B2 can be executed in http://asymptote.ualberta.ca/