2.1.1. References Assignation

Figure 1 presents the reference points and the mean distances taken from [40].

**Figure 1.** Reference points from anthropometric values K, L, O, and P.

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<sup>1</sup> and M<sup>2</sup> define the talocrural (TC) axis. We show top-transverse and rightlateral views in Figure 2 with distances Q, W, and w. N<sup>1</sup> and N<sup>2</sup> determine the subtalar (ST) axis.

Table 1 enumerates the mean values of Figures 1 and 2.

**Table 1.** Mean values of anthropometric measurements.


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.

**Figure 3.** Mean relative position of the ST and TC axis.

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.

**Figure 4.** Planes, axes, and points of 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, A0, B0, and C<sup>0</sup> are the vertices from the platform fixed to the foot, and P<sup>M</sup> is the triangle's center. S1, S2, and S<sup>3</sup> are fixed to the shank relative to the origin point P0. M<sup>1</sup> and M<sup>2</sup> define the TC axis; N<sup>1</sup> and N<sup>2</sup> correspond to the ST axis. We define r<sup>1</sup> and r<sup>2</sup> 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<sup>m</sup> is the distance from P<sup>M</sup> to PO. This distance is proportional to the body's height H. To do so, we define d<sup>m</sup> as follows:

$$\mathbf{d}\_{\mathbf{m}} = \left\| \mathbf{P}\_0 - \mathbf{P}\_\mathbf{M} \right\| = \left[ \frac{0.285 - 0.039}{2} + 0.039 \right] \cdot \mathbf{H} = 0.162 \cdot \mathbf{H}.\tag{1}$$

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 dp12 between points r<sup>1</sup> and r<sup>2</sup> about the TC and ST axes on the sagittal plane is:

$$\mathbf{d}\_{\mathbf{p}12} = \|\mathbf{r}\_1 - \mathbf{r}\_2\| = \mathbf{Q} \tag{2}$$

the projection of the most medial point (MMP) on the sagittal plane is:

$$\mathbf{P\_{MMP}} = (\mathbf{x\_{MMP}}, \mathbf{0}, \mathbf{z\_{MMP}})\_{\prime} \tag{3}$$

and for the most lateral point is:

$$\mathbf{P\_{MLP}} = (\mathbf{x\_{MLP}}, \mathbf{0}, \mathbf{z\_{MLP}}).\tag{4}$$

The point M1p is the projection of M<sup>1</sup> on the sagittal plane; we calculate it from the P and O values.

$$\mathbf{M}\_{\rm lp} = (\mathbf{x}\_{\rm MMP} - \mathbf{P}\_{\prime}\mathbf{0}\_{\prime}\mathbf{z}\_{\rm MMP} - \mathbf{O})\_{\prime} \tag{5}$$

also, M2p is M2; we estimate the projection from L and K through:

$$\mathbf{M}\_{\mathbf{2p}} = (\mathbf{x}\_{\text{MLP}} - \mathbf{L}\_{\prime}\mathbf{0}, \mathbf{z}\_{\text{MLP}} - \mathbf{K}).\tag{6}$$

Therefore, the segment M2M<sup>1</sup> has the sagittal projection M2pM1p; it has the same proportional relation R = W/w in respect to M2pr1, then:

$$\frac{\mathbf{M}\_2 - \mathbf{M}\_1}{\mathbf{M}\_2 - \mathbf{r}\_1} = \frac{\mathbf{W}}{\mathbf{w}} = \mathbf{R} \tag{7}$$

solving for r<sup>1</sup> gives the following:

$$\mathbf{r}\_1 = \mathbf{M}\_2 - \frac{\mathbf{M}\_2 - \mathbf{M}\_1}{\mathbf{R}}.\tag{8}$$

By knowing the distance Q projected in the sagittal plane and r1, the angle 41° we calculate r<sup>2</sup> from:

$$\mathbf{r\_2} = \mathbf{Q}[\cos(41^\circ), \mathbf{0}, \sin(41^\circ)] + \mathbf{r\_1} \tag{9}$$

The distance from the origin P<sup>O</sup> to the plantar surface of the foot is dm, we choose a circumscribed equilateral triangle with vertices A0, B0, C<sup>0</sup> as the platform base. The coordinates of A<sup>0</sup> are:

$$\mathbf{A}\_0 = (\mathbf{r}\_\mathbf{p}, \mathbf{0}, -\mathbf{d}\_\mathbf{m}),\tag{10}$$

for B<sup>0</sup> are:

$$\mathbf{B}\_0 = \begin{pmatrix} \mathbf{r}\_\mathrm{P} \cos \mathbf{60}^\circ \text{ } \mathbf{r}\_\mathrm{P} \sin \mathbf{60}^\circ \text{ } -\mathbf{d}\_\mathrm{m} \end{pmatrix} , \tag{11}$$

and for C0:

$$\mathbf{C}\_{0} = \left(\mathbf{r}\_{\rm P}\cos - 60^{\circ}, \mathbf{r}\_{\rm P}\sin - 60^{\circ}, -\mathbf{d}\_{\rm m}\right) \tag{12}$$

where r<sup>p</sup> is proportional to H, then:

$$\mathbf{r\_p} = \frac{2}{5} \cdot \mathbf{H}.\tag{13}$$

In summary, we estimate P0, r1, r2; and the platform's vertices A0, B0, and C0. 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:

$$\mathfrak{g} = \left[ \begin{array}{cc} \mathbf{R} & \hat{\mathfrak{p}}\_T \\ \mathbf{0}\_{1 \times 3} & 1 \end{array} \right] \text{.} \tag{14}$$

where R3×<sup>3</sup> is the rotation matrix and *p*ˆ*<sup>T</sup>* is the translation vector.

For the initial point A0:

$$\mathfrak{g}\_A(\mathbf{0}) = \begin{bmatrix} \mathbf{I}\_{3 \times 3} & \mathbf{\hat{A}}\_0 \\ \mathbf{0}\_{1 \times 3} & \mathbf{1} \end{bmatrix} \tag{15}$$

for B0:

$$\mathfrak{g}\_{\mathcal{B}}(0) = \begin{bmatrix} \mathbf{I}\_{3 \times 3} & \mathbf{\hat{B}}\_{0} \\ \mathbf{0}\_{1 \times 3} & \mathbf{1} \end{bmatrix}' \tag{16}$$

and for C<sup>0</sup>

$$\mathfrak{g}\_{\mathbb{C}}(\mathbf{0}) = \begin{bmatrix} \mathbf{I}\_{3 \times 3} & \mathbf{\hat{C}}\_{0} \\ \mathbf{0}\_{1 \times 3} & \mathbf{1} \end{bmatrix} \tag{17}$$

We define *ω*ˆ <sup>1</sup> = (*ωx*1, *ωy*1, *ωz*1) as a unitary vector for the TC axis direction given by:

$$
\hat{\omega}\_1 = \frac{\mathbf{M}\_2 - \mathbf{M}\_1}{\|\mathbf{M}\_2 - \mathbf{M}\_1\|} \,\tag{18}
$$

and a directed vector rˆ<sup>1</sup> from P<sup>O</sup> to r<sup>1</sup> is:

$$
\mathbf{\dot{r}}\_1 = \mathbf{r}\_1 - \mathbf{P}\_{\mathbf{O}\prime} \tag{19}
$$

then, an orthogonal vector to rˆ<sup>1</sup> and *ω*ˆ <sup>1</sup> is:

$$
\psi\_{\theta\_1 r\_{2\_z}} = -\hat{\omega}\_1 \times \mathbf{f}\_{1\prime} \tag{20}
$$

together, *ω*ˆ <sup>1</sup> and *ν*ˆ*θ*<sup>1</sup> *r*2*z* compound the six-dimensional vector ˆ*ξ*1:

$$
\hat{\xi}\_1 = \begin{pmatrix} \ \hat{\sigma}\_1 \\ \ \hat{\omega}\_1 \end{pmatrix}. \tag{21}
$$

In the same way, there are correspondent vectors for the TC axis:

$$
\hat{\omega}\_2 = \frac{\mathbf{N}\_2 - \mathbf{N}\_1}{\|\mathbf{N}\_2 - \mathbf{N}\_1\|} \,\tag{22}
$$

$$
\mathbf{\dot{r}}\_2 = \mathbf{r}\_2 - \mathbf{P}\_{\bullet \prime} \tag{23}
$$

$$
\mathfrak{d}\_2 = -\mathfrak{d}\_2 \times \mathfrak{k}\_2.\tag{24}
$$

and:

$$
\xi\_2 = \left(\begin{array}{c}\mathfrak{e}\_2\\\mathfrak{a}\_2\end{array}\right).\tag{25}
$$

We compute R for each joint *i* = 1, 2 from the Rodrigues' formula:

$$e^{(\Omega\_{\bar{l}}\theta\_{\bar{l}})} = \mathbf{I}\_{3 \times 3} + \Omega \sin \theta\_{\bar{l}} + \Omega^2 (1 - \cos \theta\_{\bar{l}}) \, \tag{26}$$

where Ω is the skew symmetric matrix:

$$
\Omega = \begin{bmatrix}
0 & -\omega\_{zi} & \omega\_{yi} \\
\omega\_{zi} & 0 & -\omega\_{xi} \\
\end{bmatrix}.
\tag{27}
$$

The exponential formula is:

$$e^{\mathfrak{z}\_i \theta\_i} = \begin{bmatrix} e^{\Omega\_i \theta\_i} & \mathfrak{r}\_i \\ \mathfrak{d}\_{1 \times 3} & 1 \end{bmatrix} \prime \tag{28}$$

and, *τ<sup>i</sup>* is translation vector:

$$\mathbf{r}\_{i} = \left(\mathbf{I}\_{3 \times 3} - e^{\hat{\omega}\_{i}\theta\_{i}}\right)\hat{\omega}\_{i} \times \boldsymbol{\vartheta} + \hat{\omega}\_{i}\hat{\omega}\_{i}^{T}\boldsymbol{\vartheta}\_{i}\boldsymbol{\theta}\_{i} \tag{29}$$

Points A, B, and C have invariant relative positions, and there are two rotating joints; the PoE formula for A is:

$$\mathfrak{g}\_A = e^{\oint\_{\mathbb{S}} \theta\_1} e^{\oint\_{\mathbb{S}} \theta\_2} \mathfrak{g}\_A(0) = \begin{bmatrix} \mathbb{R} & \mathfrak{f}\_A \\ 0 & 1 \end{bmatrix} \, \tag{30}$$

where *p* ˆ *<sup>A</sup>* is the instantaneous position vector of A, the PoE for B:

$$\mathfrak{g}\_{B} = e^{\oint\_{1} \theta\_{1}} e^{\oint\_{2} \theta\_{2}} \mathfrak{g}\_{B}(0) = \begin{bmatrix} \mathbb{R} & \hat{\boldsymbol{\mu}}\_{B} \\ 0 & 1 \end{bmatrix} \tag{31}$$

and the PoE for C is:

$$\mathfrak{g}\_{\mathbb{C}} = e^{\oint\_{\mathbb{S}} \theta\_1} e^{\oint\_{\mathbb{S}} \theta\_2} \mathfrak{g}\_{\mathbb{C}}(0) = \begin{bmatrix} \mathbb{R} & \hat{p}\_{\mathbb{C}} \\ 0 & 1 \end{bmatrix} \tag{32}$$

*θ*<sup>1</sup> is the TC rotation angle from the zero position, and *θ*<sup>2</sup> 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, *ω* ˆ <sup>1</sup>, *ν* ˆ<sup>1</sup> and r ˆ<sup>2</sup> ; also the points A, B, C, and P<sup>O</sup> in Figure 5.

**Figure 5.** Vectors and points on the sagittal plane.

#### 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 PO, A0, B0, C0, r1, r2, and the surfaces representing each group of movements. The forward kinematics with *θ*1*range* = *θ*2*range* = [−15 ◦ , 15 ◦ ] and *θ*<sup>1</sup> = *θ*<sup>2</sup> = 10 ◦ is in Figure 6b. For *θ*1*range* = *θ*2*range* = [−10 ◦ , 10 ◦ ] and *θ*<sup>1</sup> = *θ*<sup>2</sup> = 5 ◦ is in Figure 6c.

**Figure 6.** Forward kinematics for (**a**) initial position and (**b**) *θ*1*range* = *θ*2*range* = [−15◦ , 15◦ ], *θ*<sup>1</sup> = *θ*<sup>2</sup> = 10◦ and (**c**) *θ*1*range* = *θ*2*range* = [−10◦ , 10◦ ], *θ*<sup>1</sup> = *θ*<sup>2</sup> = 5 ◦ .

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.

**Figure 7.** Geometric design: (**a**) is the platform center, base, and *r*<sup>1</sup> ,*r*2; and (**b**) platform vertices with talocrural and subtalar axis.

By considering the distances between the origin and the vertices, we estimate the DWS maximal length in every module.

$$l\_{\max} = \max[\|p\_A(\theta\_1, \theta\_2) - A\| + \mathbf{r}\_{\mathbf{m}}] \tag{33}$$

Here, *lmax* is the maximal possible length from the triangular inequality, *p<sup>A</sup>* is the positions group in g*A*, *r<sup>m</sup>* is the module's radius, and A<sup>B</sup> 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.

**Figure 8.** Geometric design of the DWS arrays.

P<sup>O</sup> and P<sup>M</sup> are the base and platform reference frames. The platform has known dimensions and the number of sensors is seven. First, we compute A<sup>p</sup> from three distances: *lA*<sup>1</sup> = kA<sup>p</sup> − A1k, *lA*<sup>2</sup> = kA<sup>p</sup> − A2k, and *lA*<sup>3</sup> = kA<sup>p</sup> − A3k. Then, we compute B<sup>p</sup> and B<sup>p</sup> apexes after A<sup>p</sup> employing two DWS. We summarize the method in a flowchart; Figure 9.

**Figure 9.** Tetrahedron trilateration flowchart.

2.1.7. Finding the Apex in Tetrahedron A

In this section, we compute the tetrahedron T<sup>A</sup> with base △*<sup>A</sup>* = [A1, A2, A3] and apex Ap. Figure 10 shows the method we use.

**Figure 10.** Finding the apex Ap.

In Figure 10 we see that triangles △<sup>132</sup> = [A1, A3, Ap132] and △<sup>231</sup> = [A2, A3, Ap231] are two sides of the tetrahedron T<sup>A</sup> developed on the base plane.

We compute the Ap132 and Ap231 orthogonal projection on each adjacent side of the module base triangle △[A1, A2, A3] by tracing a circle centered on A<sup>1</sup> with radius kA<sup>p</sup> − A2k and the circle centered on A<sup>3</sup> with radius kA<sup>p</sup> − A3k; resulting in Ap132 and Ap131 intersection points. In addition, the circle centered on A<sup>2</sup> with radius kA<sup>p</sup> − A2k intersects the circle centered in A<sup>3</sup> at points Ap231 and Ap232. The segment from Ap132 to Ap131 intersects the points defined by Ap231 and Ap232 at Apxy. In the case of tetrahedron TA, we determine Apxy = (Apx, Apy, 0) as A<sup>p</sup> projection on the base plane. It is easy to realize that the height of T<sup>A</sup> is the absolute value of the Apz coordinate. Then, we can find the distance from Apxy to A<sup>3</sup> as a triangle △[Apxy, A3, Ap] side; the other is Apz, and the hypotenuse is the distance *lA*<sup>3</sup> = kA<sup>p</sup> − A3k, then, Apz is:

$$\mathbf{A\_{pz}} = \sqrt{l\_{A3}^2 - (\mathbf{A\_{pxy}} - \mathbf{A\_3})^2} \tag{34}$$
