Realtime inverse kinematics for the upper limb: a modelbased algorithm using segment orientations
 Bence J. Borbély^{1}Email author and
 Péter Szolgay^{1}
DOI: 10.1186/s129380160291x
© The Author(s) 2017
Received: 16 August 2016
Accepted: 30 November 2016
Published: 17 January 2017
Abstract
Background
Model based analysis of human upper limb movements has key importance in understanding the motor control processes of our nervous system. Various simulation software packages have been developed over the years to perform model based analysis. These packages provide computationally intensive—and therefore offline—solutions to calculate the anatomical joint angles from motion captured raw measurement data (also referred as inverse kinematics). In addition, recent developments in inertial motion sensing technology show that it may replace large, immobile and expensive optical systems with small, mobile and cheaper solutions in cases when a laboratoryfree measurement setup is needed. The objective of the presented work is to extend the workflow of measurement and analysis of human arm movements with an algorithm that allows accurate and realtime estimation of anatomical joint angles for a widely used OpenSim upper limb kinematic model when inertial sensors are used for movement recording.
Methods
The internal structure of the selected upper limb model is analyzed and used as the underlying platform for the development of the proposed algorithm. Based on this structure, a prototype marker set is constructed that facilitates the reconstruction of modelbased joint angles using orientation data directly available from inertial measurement systems. The mathematical formulation of the reconstruction algorithm is presented along with the validation of the algorithm on various platforms, including embedded environments.
Results
Execution performance tables of the proposed algorithm show significant improvement on all tested platforms. Compared to OpenSim’s Inverse Kinematics tool 50–15,000x speedup is achieved while maintaining numerical accuracy.
Conclusions
The proposed algorithm is capable of realtime reconstruction of standardized anatomical joint angles even in embedded environments, establishing a new way for complex applications to take advantage of accurate and fast modelbased inverse kinematics calculations.
Keywords
Inverse kinematics (IK) Inertial measurement unit (IMU) Realtime Upper limb OpenSim Embedded systems WearableBackground
Quantitative movement analysis is a key concept in understanding processes of the human movement system. Evolved, high precision measurement devices have advanced research activity in movement rehabilitation [1–4], performance analysis of athletes [5] and general understanding of the motor system [6–9] during the last decades by making movement pattern comparison possible. This advancement was further accelerated by modelbased analysis approaches that enabled explicit characterization of the studied movement patterns [10–15].
The most widely used measurement systems apply lineofsight (LoS) methods (optical or ultrasoundbased) that require a fixed markersensor structure. In these cases passive (optical) or active (ultrasound) markers are placed on anatomically relevant locations of the studied subject. Having the markers in place, measurements have to be performed in a specialized laboratory environment where locations and orientations of the sensor elements (i.e. cameras or microphones) are known and invariant (at least across trials). This means that even the spatial positions of the markers can be determined with good accuracy—especially with optical systems—the possible range of motion will always be constrained by the actual measurement volume covered by the sensors of the system. Although this property is not an issue for many movement analysis scenarios, there are cases when a measurement method allowing unconstrained free space movement would be more beneficial (e.g. various outdoor activities or ergonomic assessment of work environments, to name a few).
Advancements in the field of inertial sensor technology have given rise to new development directions in laboratoryfree movement analysis methods. The main difference between LoS and inertial systems is the recorded modality: while LoS methods determine the spatial locations of markers based on planar position (optical) or timing (ultrasound) information, inertial sensors give their orientation in space by measuring physical quantities acting on them directly. These quantities are linear acceleration and angular velocity in most cases while they are supplemented with magnetic field measurements in more complete setups. To obtain orientation from raw inertial measurements, various sensor fusion algorithms have been developed utilizing Kalmanfilters [16, 17], gradient descent methods [18], complementary filters [19, 20] and other techniques [21], most of them being capable for realtime operation in embedded systems. In addition, recent evolution of chipscale inertial sensors based on MEMS technology further widened the possibilities of wearable measurement device development by making the core sensing elements available for better integration. While there is no gold standard among fusion algorithms and sensor chips as compromises have to be made in aspects of accuracy, system complexity and computational demand of the fusion algorithm, it can be stated that inertial sensor technology is taking a more and more growing part in human movement measurements (a good example for this progress is Xsens’ product portfolio).
In addition to accurate measurement, proper evaluation of the recorded motion is an other key building block of human movement analysis. Considering the kinematic assessment of upper limb movements by using inertial sensors, reliable methods have been proposed to estimate joint kinematics in custom local coordinate systems, also dealing with the problem of functional calibration and measurement error compensation [22–25]. While these methods may work well in movement analysis scenarios focusing only on the spacial kinematics of the upper limb, the goal of the present study is to extend the applicability of a more complex freely available upper limb model [26] with inertial measurement and realtime joint angle reconstruction capabilities. The chosen model includes 15 degrees of freedom and 50 muscle compartments and in addition to movement kinematics it enables the evaluation of muscletendon lengths, moment arms, muscle forces and joint moments in an anatomically reasonable setup (also conforming to the ISB recommendation [27] as those presented in [23] and [24]). Our motivation of choosing this approach instead of using direct joint kinematics estimates from inertial sensor data lies in the strong belief that the utilization of this upper limb model leads to better understanding of the additional inherent properties of arm movements (e.g. muscle activation patterns) compared to a pure kinematic investigation.
By using the model from [26], two main differences from the direct approach should be considered: (1) because the model was developed mainly for optical motion capture technology, joint angle reconstruction is less coupled with raw measurement data in the sense that it is performed as an inverse task based on locations of arbitrarily placed virtual markers (instead of relying on sensor orientations directly) and (2) as a consequence of the inverse kinematics approach, the model raises more computational demand for joint angle reconstruction. While the first point contributes to a possibly wider adoption of the model (i.e. measurement data of any reasonable source—being it optical, ultrasound based, inertial, etc.—only needs to be expressed in terms of virtual marker locations), the second can be considered as a drawback in certain situations when realtime joint angle reconstruction would be beneficial.
To analyze modeldefined parameters of recorded movements, various software packages have been developed for biomechanical analysis over the past decades. Some of them are tightly integrated into their corresponding measurement system’s ecosystem with only kinematic reconstruction capabilities [28–30] while others are independent of a particular hardware setup (SIMM [31] and OpenSim [32]), even including analysis options for muscle activations and movement dynamics if the corresponding models are available. While these tools may address the joint angle reconstruction problem involving inverse kinematics in different ways resulting in varying performance, OpenSim [32] was chosen as the reference modelbased analysis software for the purposes of the current study, because (as to the best of the authors’ knowledge) it is the only free and open source (thus widely accessible) option in the biomechanical modeling field.
Because the upper limb model used in this study was originally developed for the SIMM software, an other reason for using OpenSim was it’s ability to import and handle SIMM models natively, so combined with the model it provides a complete package for open source arm movement analysis. For kinematic evaluation, OpenSim uses the “standard” offline measurementscalinginverse kinematics pipeline where the actual biomechanical model (single limb to full body) is fitted to measurement data. During this process, positions of virtual markers placed on specific model segments are fitted to experimentally recorded marker positions of the subject with the same arrangement. Scaling is important to generate subjectspecific model instances while inverse kinematics (IK) is performed to extract modeldefined anatomical joint angles that produced the movement. OpenSim uses a text based structured XML model format that contains all information needed for the biomechanical description of the human body [bodies, kinematic constraints and forces (i.e. muscles)] that are accessible through API calls, too.
Complex measurement and analysis of upper limb movements including kinematics and muscle activities is an exciting and growing subfield of human movement analysis [33–37] that promises better understanding of control patterns during specific movements, and as an example benefit may—on the longer term—advance control techniques currently applied to arm and hand prostheses. This process however needs tighter integration of kinematic measurement and reconstruction (from raw data to anatomical joint angles) as the time and computational overhead of the offline measurementscalinginverse kinematics scheme gives a bottleneck in applications where realtime analysis of the control patterns with respect to the actual kinematics would be beneficial.
As a specific example, let us consider a supervised learning based muscle activity classification framework like the one proposed in [38] that is built around the increasingly popular concept of deep neural networks [39]. The goal of this setup is to assess the classification capability of various network architectures in estimating upper limb movement kinematics based only on muscle activity recordings. The large amount of labeled data that is needed to train these networks (even iteratively after initial deployment) should be presented by a method that performs muscle activity data labeling based on the actual (anatomically relevant) kinematic state of the limb automatically as the measurements occur. By using the approach proposed in this paper, labeled measurement data could be collected in realtime and simultaneously from multiple recording devices and a high level of automation might be reached in the process of measurement, network training and deployment.
Purpose of the study
The main goal of this study is to extend the measurement and analysis workflow of human arm movements with a method that allows accurate and realtime calculation of anatomical joint angles for a widely used upper limb model in cases when inertial sensors are used for movement recording. For this purpose a custom kinematic algorithm is introduced that utilizes orientation information of arm segments to perform joint angle reconstruction. Accuracy and execution times of the proposed algorithm are validated against the most widely available biomechanic simulation software’s inverse kinematics algorithm on various platforms.
Methods
Upper limb model
To analyze arm kinematics with OpenSim the most complete model available was chosen known as the Stanford VA Upper Limb Model [26]. It is freely available as part of the Simtk project [40] in SIMM model format [41] that can be imported directly into OpenSim. The model is based on experimental data, includes 15 degrees of freedom and 50 muscle compartments and enables the evaluation of kinematics, muscletendon lengths, moment arms, muscle forces and joint moments in an anatomically reasonable setup. After importing, the structure of the model follows OpenSim’s convention including bodies connected with joints, rotational and translational kinematic constraints and forces defining muscle paths and attributes (for details, see Additional file 1). The 15 degrees of freedom define the kinematics of the shoulder (3), elbow (2), wrist (2), index finger (4) and thumb (4). As the current work focuses on kinematics of the shoulder, elbow and wrist joints only, any muscles and kinematics of the index finger and the thumb will not be taken into account further in this study. The seven degrees of freedom that define the kinematic state of the whole arm excluding the fingers are elevation plane, thoracohumeral (elevation) angle and axial rotation for the shoulder, elbow flexion and forearm rotation for the elbow and deviation and flexion for the wrist.
The model represents the upper limb as a linked kinematic chain of bodies, each having a parent body, a location in the parent’s frame and a joint describing the possible relative motion of the child with respect to the parent. The threedimensional posture of the arm is generated by consecutive rotations of bodies determined by the actual angle values (joint coordinates) in proximal to distal order. As the movement of the shoulder girdle (clavicle, scapula and humerus) is complex and cannot be measured directly in most cases, the model implements regression equations that vary only with the thoracohumeral angle to determine the position of the shoulder joint with respect to the thorax.
Orientation from joint coordinates
 1.
The shaft of the humerus is parallel to the vertical axis of the thorax.
 2.
In case of shoulder elevation, the humerus moves in the plane of shoulder abduction.
 3.
In case of elbow flexion, the forearm moves in the sagittal plane.
 4.
The hand is in the sagittal plane.
 5.
The third metacarpal bone in aligned with the long axis of the forearm.
Markers, scaling and inverse kinematics
Prototype markers
To enable the utilization of the upper limb model with inertial measurements, a prototype marker set was defined (see Additional file 2). For this purpose, orthonormal bases were formed for each anatomical joint (shoulder, elbow and wrist) and markers were placed at specific locations in these bases to reflect the actual compound rotations among the respective degrees of freedom (for the corresponding mathematical definitions, see Appendix 1.
Orthonormal bases
Shoulder The three independent model axes for the shoulder (defined in model bodies humphant, humphant1 and humerus that have the same position), collectively denoted as \(\mathbf {B}_{sh\_orig}\), were good candidates to form a basis because they are unit length vectors (like all axes in the model) and almost orthogonal to each other (pairwise deviations from right angle are 0.00064°, 0.0029° and 0.0002°). For proper operation of the proposed algorithm however, these axes were orthogonalized using QR decomposition (see Appendix 2) to prevent error accumulation during the calculations. This resulted in the orthonormal basis \(\mathbf {B}_{sh\_orth}\).
As a result, rotations in the shoulder can be expressed as elemental rotations of \(\mathbf {B}_{sh\_orth}\) with acceptable angle errors due to the pairwise deviations between the original and new basis vectors after orthogonalization (0.000019°, 0.000655° and 0.002925°, respectively).
Elbow As relative orientation of the two rotation axes in the elbow is not close enough to orthogonal and the axes are defined in different parent bodies (\(\mathbf {r}_\mathtt{{el\_flex}}\) \(\rightarrow \) ulna and \(\mathbf {r}_\mathtt{{pro\_sup}}\) \(\rightarrow \) radius), the orthonormal basis \(\mathbf {B}_\mathtt{{pro\_sup}}\) and the rotation matrix \(\mathbf {R}_{{\mathbf {r}_\mathtt{{el\_flex}}}}^{{\mathbf {B}_\mathtt{{pro\_sup}}}}({\theta _\mathtt{{el\_flex}}})\) were formed to properly express the compound rotation as the product of an axisangle and an elementary rotation about the main axis in \(\mathbf {B}_\mathtt{{pro\_sup}}\). Again, some angle errors are expected while calculating the elbow flexion angle in this basis because \(\mathbf {r}_\mathtt{{el\_flex}}\) is threated as it would belong to the radius body of the model.
Wrist Rotations in the wrist are the most complex among the three anatomical joints. Effects of the two active joint coordinates (deviation and flexion) are distributed among two model bodies (lunate and capitate), each having two nonorthogonal rotations (\(\mathbf {r}_\mathtt{{dev}}\), \(\mathbf {r}_\mathtt{{flex}}\) \(\rightarrow \) lunate and \(\mathbf {r}_\mathtt{{pdr1}}\), \(\mathbf {r}_\mathtt{{pdr3}}\) \(\rightarrow \) capitate) depending on both joint coordinates. To deal with the complexity of this structure, the orthonormal basis \(\mathbf {B}_\mathtt{{pdr3}}\) and rotation axes \(\mathbf {r}_\mathtt{{dev}}^{{\mathbf {B}_\mathtt{{pdr3}}}}\), \(\mathbf {r}_\mathtt{{flex}}^{{\mathbf {B}_\mathtt{{pdr3}}}}\) and \(\mathbf {r}_\mathtt{{pdr1}}^{{\mathbf {B}_\mathtt{{pdr3}}}}\) were constructed to prepare the calculation of \(\theta _\mathtt{{dev}}\) and \(\theta _\mathtt{{flex}}\). Using this approach, \(\mathbf {r}_\mathtt{{dev}}\) and \(\mathbf {r}_\mathtt{{flex}}\) are threated as if they would belong to the capitate body of the model.
Marker placement

Four markers were placed into each orthonormal basis having one at the origin of the actual basis ([0 0 0] in its parent body) and one in each axis of the basis.

The markers were named using the convention PMx_[SHELWR]_[OXYZ] where PM refers to prototype marker, x is the serial number of the basis in which the marker is located (1–3), [SHELWR] refers to the anatomical joint in which the marker is located and [OXYZ] refers to the marker’s location within the actual basis (origin or any of the axes). For example the name of the wrist’s origin marker is PM3_WR_O.
Algorithm description
Shoulder
Although the formulations in (2a) and (2c) could be susceptible to modulo \(\pi \) and sign errors in general, the allowed angle ranges defined in the model (\(\theta _\mathtt{{sh\_elv}}\): \(0^\circ \rightarrow 180^\circ \), \(\theta _\mathtt{{sh\_rot}}\): \(90^\circ \rightarrow 20^\circ \)) keep these equations safe to use until the experimental data does not force the model outside of these ranges.
Elbow
As in the case of the shoulder, (4a) and (4b) should be used with care because of modulo \(\pi \) and sign errors, but again having sufficient joint angle limits in the model (\(\theta _\mathtt{{el\_flex}}\): \(0^\circ \rightarrow 130^\circ \), \(\theta _\mathtt{{pro\_sup}}\): \(90^\circ \rightarrow 90^\circ \)) application of these formulas is safe until experimental data does not force the model outside of these ranges.
Wrist
 1.
\(\left(  {\theta _\mathtt{{flex}}} + \eta \right) \) defines a baseline with constant negative slope for the two possible solutions \(F \left( {\theta _\mathtt{{flex}}}, 1 \right) \) and \(F \left( {\theta _\mathtt{{flex}}}, 1 \right) \).
 2.
Because of the definition of the atan2(y,x) function, the value of \(\text {atan2} \left( \sqrt{\xi  c^2},~c \right) \) will always be positive if \(\sqrt{\xi  c^2}\) is real (i.e. \(c^2 \le \xi \)). This implies that the two solutions to \(F \left( {\theta _\mathtt{{flex}}}, \sigma \right) \) do not cross the baseline but remain “below” (\(\sigma = 1\)) and “above” (\(\sigma = 1\)) of it for all values of \(\theta _\mathtt{{flex}}\).
 1.
The condition in (23) is always met.
 2.
\(\vartheta ^{(k)}_{1,2} ~ ((k = 1) \vee (k = 2))\) are separate real numbers if there is a singularity region in the actual wrist configuration for \(c_k\). In this case \({\theta _\mathtt{{flex}}} = \vartheta ^{(k)}_{1}\) and \({\theta _\mathtt{{flex}}} = \vartheta ^{(k)}_{2}\) indicate singularity border locations directly.
 3.
\(\vartheta ^{(k)}_{1,2} ~ ((k = 1) \vee (k = 2))\) are complex conjugates if there is no singularity region in the actual wrist configuration for \(c_k\). In this case \({\theta _\mathtt{{flex}}} = \text {Re} \left( \vartheta ^{(k)}_{1} \right) = \text {Re} \left( \vartheta ^{(k)}_{2} \right) \) indicates the location where the values of \(F \left( {\theta _\mathtt{{flex}}}, 1 \right) \) and \(F \left( {\theta _\mathtt{{flex}}}, 1 \right) \) are closest to (\(k = 1\)) and furthest from (\(k = 2\)) each other.
 4.
\(\text {Re} \left( \vartheta ^{(2)}_{1,2} \right) \) always remain outside the model defined range of \(\theta _\mathtt{{flex}}\).
 5.
\({\theta _\mathtt{{flex}}} = \eta \) is the “gluing point” of \(F \left( {\theta _\mathtt{{flex}}}, 1 \right) \) and \(F \left( {\theta _\mathtt{{flex}}}, 1 \right) \), meaning that the singularity region for \(c_1\) starts to develop from this location, driving \(F \left( {\theta _\mathtt{{flex}}}, \sigma \right) \) to “stick” to the baseline.
 6.
If there is a singularity region for \(c_1\), \(\text {Re} \left( \vartheta ^{(1)}_{1} \right) \) remains always smaller than \(\mu \) where \(F \left( \mu , \sigma \right) = 0\).
 7.
In cases when the singularity region starts to develop (i.e. \(\left \text {Im} \left( \vartheta ^{(k)}_{1} \right) \right \) is sufficiently small but not zero), two separate roots may appear, but only one being valid.
 8.
\(F \left( {\theta _\mathtt{{flex}}}, \sigma \right) \) will have a valid root at \({\theta _\mathtt{{flex}}} = \mu \) if and only if \(\sigma = \text {sign} \left( \mu  \eta \right) \).
Although the computational demand of wrist angle calculations is higher than of the shoulder and the elbow, the algorithm has still higher overall time efficiency than the optimization approach used by OpenSim’s inverse kinematics tool, as it is shown in the "Results" section.
Algorithm validation

VM1_TH : Thorax marker at the upper end of the sternum.

VM2_AC : Acromioclavicular joint of the shoulder girdle.

VM3_EL_IN : Medial epicondyle of the humerus.

VM4_EL_OUT : Lateral epicondyle of the humerus.

VM5_WR_IN : Distal head of the radius.

VM6_WR_OUT : Distal head of the ulna.

VM7_MC2 : Distal head of the second metacarpal bone.

VM8_MC5 : Distal head of the fifth metacarpal bone.
Simulated movement patterns
Inverse kinematics with OpenSim
To speed up the validation process, OpenSim (v3.3) was compiled from source on a Supermicro server having two Intel^{®} Xeon^{®} E52695 v3 CPUs (with a total of 56 execution threads) and 64 GB RAM, running Ubuntu Server 14.04.2 LTS operating system. Although the inverse kinematics (IK) algorithm in the used OpenSim version do not utilize multicore architectures natively, each IK task can be divided into separate subtasks that can run in parallel thanks to the applied optimization method (there is no data dependency between time frames). To utilize this property, a pipeline was developed using MATLAB and Bash to prepare VMx marker data and the required files for OpenSim and manage file transfers, multithreaded IK execution, results collection and evaluation. One important step before performing the IK calculation in OpenSim is subjectspecific scaling of the used model and relative weighting of the markers. As only simulated data were used in the current study on an unmodified upper limb model, the scaled model file was identical to the original file during IK execution, while all marker weights were equal.
Algorithm implementation
The prototype of the proposed algorithm was implemented in MATLAB and tested with the simulated PMx marker trajectories. Calculation of (6) was performed using MATLAB’s builtin fzero() function. Based on the MATLAB version, the algorithm was implemented in ANSI C to target practical applications. In this case (6) was solved with Brent’s root finding algorithm from [44]. Furthermore, compilation options were included to assess the effects of different data precisions (float or double) on the accuracy and execution time of the algorithm. This was not an option with MATLAB because fzero() cannot be used with float input.
To address possible accuracy problems arising from the lower precision of float data, an additional test case with a simple output continuity check for wrist angles was included, namely when the absolute difference between two successive \(\theta _\mathtt{{flex\_c}}\) values is larger than 5°, the actual \(\theta _\mathtt{{flex\_c}}\) will be the previous \(\theta _\mathtt{{flex\_c}}\) + 0.5°. This modified version of the algorithm is denoted with mod. suffix among the results.
Evaluation platforms
 1.
STM32F407VG  ARM CortexM4 core with single precision floating point unit (FPU), up to 168 MHz core clock, 1 MB Flash memory, 192 KB SRAM.
 2.
STM32F746NG  ARM CortexM7 core with single precision FPU and L1cache, up to 216 MHz core clock, 1 MB Flash memory, 320 KB SRAM.
Performance metrics
To evaluate the overall performance of the algorithm compared to OpenSim’s IK method, accuracy and execution times were analyzed in all cases (OpenSim, MATLAB and C implementations). To assess accuracy, RMS values were computed for the differences between the calculated and simulated joint coordinate trajectories for each trial. Means and standard deviations of these RMS values were then calculated across trials for each platform and precision (where this was applicable).
Running times of OpenSim’s IK evaluation were calculated as a sum of subtask execution times from the IK log output directly. Algorithm execution times were measured by the tic and toc methods in MATLAB, the clock() function from the \(\texttt {<time.h>}\) library for the C implementation on PC and onchip hardware timers clocked at 1 MHz for both MCUs.
Data exclusion from OpenSim trials
Although inverse kinematics in OpenSim was calculated using an unmodified and unscaled model in each trial, there were cases when large step errors occurred at seemingly random locations in the IK output (independently of subtask borders mentioned in "Inverse kinematics with OpenSim" section). This phenomenon may be handled by marker placement adjustment or error checking in measurement data in general. As IK input was strictly controlled by using simulated trajectories and the markers remained intact in the model between trials, further troubleshooting would have been needed to find a solution to this issue. Because the main emphasis of the study is the proposed algorithm and not OpenSim’s internal workings and IK troubleshooting, all OpenSim trials were excluded from final accuracy assessment where any of the resulted joint coordinate RMS errors exceeded 5° to not bias the results with incorrect data. As a result, only 59 trials out of 100 were used to calculate the accuracy of OpenSim’s IK algorithm. This however did not have any effect on the other measurements, so MATLAB and all C results were calculated across 100 trials.
Results
Accuracy
RMS errors
Test environment  \(\theta _\mathtt{{elv}}\) (°)  \(\theta _\mathtt{{sh\_elv}}\) (°)  \(\theta _\mathtt{{sh\_rot}}\) (°)  \(\theta _\mathtt{{el\_flex}}\) (°)  \(\theta _\mathtt{{pro\_sup}}\) (°)  \(\theta _\mathtt{{dev\_c}}\) (°)  \(\theta _\mathtt{{flex\_c}}\) (°) 

OpenSim  0.0429 ± 0.0339  0.0192 ± 0.0053  0.1472 ± 0.0760  0.0764 ± 0.0288  0.6365 ± 0.1701  0.9198 ± 0.2477  2.2916 ± 1.1142 
MATLAB  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0014 ± 0.0007  0.0005 ± 0.0002  0.0008 ± 0.0006  0.0023 ± 0.0041  0.0049 ± 0.0087 
PC double  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0014 ± 0.0007  0.0005 ± 0.0002  0.0008 ± 0.0006  0.0025 ± 0.0051  0.0053 ± 0.0107 
PC float  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0016 ± 0.0013  0.0005 ± 0.0002  0.0008 ± 0.0006  0.4193 ± 0.8995  1.1148 ± 2.4730 
PC float mod.  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0016 ± 0.0013  0.0005 ± 0.0002  0.0008 ± 0.0006  0.0045 ± 0.0092  0.0097 ± 0.0195 
ARM M4 double  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0014 ± 0.0007  0.0005 ± 0.0002  0.0008 ± 0.0006  0.0025 ± 0.0051  0.0053 ± 0.0107 
ARM M4 soft float  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0016 ± 0.0013  0.0005 ± 0.0002  0.0008 ± 0.0006  0.4193 ± 0.8995  1.1147 ± 2.4730 
ARM M4 hard float  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0016 ± 0.0013  0.0005 ± 0.0002  0.0008 ± 0.0006  0.4095 ± 0.9051  1.0944 ± 2.4840 
ARM M4 soft float mod.  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0016 ± 0.0013  0.0005 ± 0.0002  0.0008 ± 0.0006  0.0078 ± 0.0128  0.0170 ± 0.0276 
ARM M4 hard float mod.  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0016 ± 0.0013  0.0005 ± 0.0002  0.0008 ± 0.0006  0.0077 ± 0.0128  0.0167 ± 0.0275 
ARM M7 double  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0014 ± 0.0007  0.0005 ± 0.0002  0.0008 ± 0.0006  0.0025 ± 0.0051  0.0053 ± 0.0107 
ARM M7 soft float  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0016 ± 0.0013  0.0005 ± 0.0002  0.0008 ± 0.0006  0.4193 ± 0.8995  1.1147 ± 2.4730 
ARM M7 hard float  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0016 ± 0.0013  0.0005 ± 0.0002  0.0008 ± 0.0006  0.4095 ± 0.9051  1.0944 ± 2.4840 
ARM M7 soft float mod.  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0016 ± 0.0013  0.0005 ± 0.0002  0.0008 ± 0.0006  0.0078 ± 0.0128  0.0170 ± 0.0276 
ARM M7 hard float mod.  0.0028 ± 0.0003  0.0006 ± 0.0002  0.0016 ± 0.0013  0.0005 ± 0.0002  0.0008 ± 0.0006  0.0077 ± 0.0128  0.0167 ± 0.0275 
Regarding OpenSim it can be seen that errors for each joint coordinate are larger than those provided by our algorithm. The reason for this can lie in the optimization approach of OpenSim that in fact contains hardcoded convergence (0.0001) and iteration (1000) limits. However these limits prevent OpenSim’s IK algorithm to match the simulated movement patterns perfectly, they provide a practical solution to the \(accuracy \leftrightarrow running~time\) tradeoff for the software’s general usage.
MATLAB and C implementations of the proposed algorithm performed equally well for shoulder and elbow angles independent of the used data precision (double / float). This could occur because of the relatively low number of operations needed by these joint coordinates shown in equation groups () and () that prevented considerable error accumulation due to the lower precision of float. Regarding wrist angles a clear distinction can be made between double and float (MATLAB uses double as default). The two main reasons for this phenomenon are (1) the significantly larger computational demand of \(\theta _\mathtt{{dev\_c}}\) and \(\theta _\mathtt{{flex\_c}}\) involving iterative processes that can lead to precision error accumulations and (2) rounding error based mismatch in the root finding process involved in the calculation of \(\theta _\mathtt{{flex}}\) in rare cases when two roots are present in (6). A deeper analysis among the trialwise results revealed that the second reason was more significant as roughly 70% of the trials ended up in no more than 0.1° maximum error with float precision. The rest of the trials contained 1–5 “wrong” samples showing 15°–20° impulselike errors while the remaining samples within the trial did not have this problem. Investigation of the erroneous samples revealed that indeed a wrong root for (6) was found in these cases. To deal with this issue, an output continuity checking step was implemented for float precision in cases denoted with the mod. suffix. This turned out to be a simple yet effective solution to the problem as the corresponding results show the disappearance of the impulselike errors.
Execution time
Execution times
Test environment  Execution time per iteration (ms)  Speedup wrt. OpenSim 

OpenSim  145.0532 ± 10.0669  1x 
MATLAB  2.3656 ± 0.6689  61x 
PC double  0.0111 ± 0.0013  13011x 
PC float  0.0088 ± 0.0008  16416x 
PC float mod.  0.0097 ± 0.0013  14982x 
ARM M4 double  4.8777 ± 0.3554  30x 
ARM M4 soft float  2.7327 ± 0.0928  53x 
ARM M4 hard float  0.9713 ± 0.0214  149x 
ARM M4 soft float mod.  2.7394 ± 0.0930  53x 
ARM M4 hard float mod.  0.9740 ± 0.0216  149x 
ARM M7 double  2.3124 ± 0.1704  63x 
ARM M7 soft float  1.4293 ± 0.0504  101x 
ARM M7 hard float  0.4462 ± 0.0117  325x 
ARM M7 soft float mod.  1.4296 ± 0.0505  101x 
ARM M7 hard float mod.  0.4478 ± 0.0115  324x 
Measurement results show that the optimization approach of OpenSim performed the calculation of a single iteration in 145 ms on average. Because of the application specific nature of the proposed algorithm, its running times considering different implementations (MATLAB/C), data precisions (double/float) and platforms (PC/ARM CortexM) all showed a significant increase in execution performance compared to OpenSim, the worst result being about 5 ms on average for a single iteration.
As expected, the C implementation is more than two orders of magnitude faster than the MATLAB version on the PC, yielding execution times per iteration about 10 μs with all precision variants (double, float and float mod.). Opposed to this, running times on embedded platforms showed more scattered results. The difference between double and float is more expressed in these cases while application of the FPU accelerates float computations even further (hard float entries in Table 2). Regarding the modified algorithm variant it can be seen that even the extra continuity check adds some amount to the execution time per iteration, the possibility to use float precision brings more speed advantage, especially with the FPU enabled. These findings are true for both tested MCUs with the observation that ARM’s M7 architecture is about twice as fast as M4 when running the presented algorithm with the same core clock.
Discussion
Evaluation results of the tested algorithms show that each approach provides proper accuracy for most common arm movement analysis scenarios. One important aspect however is that while OpenSim provides a useful general tool for biomechanical analysis including fields beyond inverse kinematics (e.g. inverse and forward dynamics), the calculation of joint angles from the actual experimental data is rather demanding computationally. As the output of this step gives the basis for all other analysis methods in the software, the amount of time needed for the overall processing pipeline highly depends on the efficiency of this algorithm. As Table 2 shows, the average amount of time needed for OpenSim’s IK algorithm to perform a single iteration would allow about 7 Hz operation that falls behind the generally accepted practice in human movement recording of at least 50 Hz. This property excludes OpenSim from tight integration with systems requiring realtime movement kinematics, however that is not the software’s original target application anyway (up to version 3.3 at least).
Considering the algorithm proposed in the study Tables 1 and 2 show a significant improvement in performance in both accuracy and execution time when compared to OpenSim’s IK method. The main reason for this difference is the algorithm’s application specific nature with the utilization of both the internal structure of the used upper limb model and inertial sensing of movement to determine limb segment orientations directly. As the MATLAB version showed proper accuracy and sufficiently short execution time on the PC, implementation of the algorithm in ANSI C was reasonable to assess its “real” performance without the overhead of a general prototyping tool that MATLAB essentially has. Because accuracy results are the same or very similar across specific variants of the C implementation (i.e. using double/float precision), only execution time differences are discussed later in the text.
Running times of the algorithm’s C implementation showed more than four orders of magnitude speedup on the tested Intel^{®} Core^{®} i5540M processor compared to OpenSim’s IK algorithm on a more recent and higher performance server CPU with Xeon^{®} architecture, yielding about 10 μs execution time per iteration for all variants. However this is an impressive improvement, running the algorithm on PC would still pose problems from practical aspects of possible applications (e.g. total size and mobility of the measurement system or communication overhead between the measuring and processing device), so the real benefit of this speed increase lies in the “spare” performance that opened the way to testing the algorithm in resource constrained embedded environments. Evaluation of the proposed method on high performance MCUs showed that all implementation variants that provided good accuracy (double and [soft/hard] float mod.) had acceptable execution times on both architectures (M4 and M7) for realtime operation, considering 100 Hz as sufficient sampling frequency for human movement analysis. Based on these results, the specific implementation variant should be chosen taking the overall design requirements of the actual practical application into account (i.e. wearable measurement devices like the one presented in [45]) as in most cases the algorithm should fit into a system containing other computationally demanding processes (like sensor fusion algorithms) with power consumption being a critical part of the design for example.
An other practical advantage of the described algorithm is that it enables subjectindependent joint angle reconstruction during the measurements. This means that by taking advantage of the offsetindependent nature of orientation sensing, no scaling is required for the proper calculation of inverse kinematics (opposed to OpenSim) as long as the IMUs are capable to produce good approximations of limb segment orientations.
It needs to be emphasized however that the application specific nature of the algorithm and its dependency on the used upper limb model induce some practical considerations, because having a method that works in a strictly controlled simulation environment does not guarantee its applicability in a real situation. A fundamental thing to consider is whether data provided by real sensors reflect arm segment orientations needed by the algorithm accurately. As this was a key requirement from the beginning of algorithm design, the prototype markers were defined in local bases of the joints that can be directly expressed in terms of sensor orientations (for details, see Additional files 1, 2). As an additional benefit, the definition of an anatomical calibration procedure—often needed when inertial sensors are used for human movement recording—can be avoided as the proposed algorithm does not use segment length information for joint angle reconstruction. What cannot be avoided however is the accurate estimate of sensor orientations, as the whole process depends highly on the precision of this step. Although there is no ultimate solution to the problem of inertial sensor fusion yet, there are continuous algorithmic efforts to reach higher accuracy and reliability (for different highlighted approaches, see "Background" section). But even in cases when the sensors provide accurate orientation information of the measured limb, care must be taken when determining the limb’s reference orientation based on the measurements. The reason for this is mainly intersubject variability in the sense that even the model defines the reference posture clearly, it cannot be assumed that any actual subject will reproduce the same posture very accurately that can lead to offset errors during the measurement. Furthermore, the assumption was made during algorithm development that the measured movement always remains within the valid joint angle ranges defined in the model. As long as this assumption holds (as in the case of simulated movement patterns presented in this study), the algorithm should not have problems with proper joint angle reconstruction. However, if outliers are present in the experimental data (e.g. reference posture errors, inaccuracies in the measurement or the sensor fusion algorithm or extreme anatomical ranges of a subject) undefined output states can occur. This may be handled with a simple saturation technique on the algorithm level but rather should be prevented by applying proper experimental design and calibration methods. In a practical setup this involves proper sensor placement and various steps before the measurements including zero motion offset compensation, hard and soft iron error compensation in the magnetometer and determining relative sensor orientations with respect to the measured segments [46, 47] for example.
Conclusions
With keeping the upper mentioned considerations in mind the proposed algorithm is capable for realtime reconstruction of standardized anatomical joint angles even in embedded environments, opening the way to complex applications requiring accurate and fast calculation of modelbased movement kinematics. Although the presented algorithm is special to the selected upper limb model, the introduced approach by strategically placing the prototype markers can further be applied to other biomechanical models in the future. As a result, the proposed method brings the possibility to widen the application areas of OpenSim with complex models and making its overall analysis pipeline more efficient by accelerating the calculation of inverse kinematics and providing the possibility to perform this step even on the measurement device in cases when accurate inertial movement sensing is applicable.
Abbreviations
 LoS:

line of sight
 MEMS:

microelectromechanical system
 IMU:

inertial measurement unit
 API:

application programming interface
 IK:

inverse kinematics
 XML:

extensible markup language
 ANSI:

American National Standards Institute
 MCU:

microcontroller unit
 FPU:

floating point unit
 PC:

personal computer
 RMS:

root mean square
Declarations
Authors' contributions
Both authors prepared the study and interpreted the results. In addition, BJB developed the mathematical background, designed and implemented the software components required for the validation of the proposed algorithm and drafted the manuscript. Both authors revised the manuscript and gave their final approval for publication. Both authors read and approved the final manuscript.
Acknowledgements
The authors would like to sincerely thank Endre László for his technical assistance and for proof reading the manuscript.
Competing interests
The authors declare that they have no competing interests.
Availability of data and materials
All data generated or analyzed during this study are included in this published article and its supplementary information files.
Funding
This work was financially supported by the Central Funding Program of Pázmány Péter Catholic University under project number KAP150551.1ITK.
Open AccessThis article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
Authors’ Affiliations
References
 Zheng H, Black ND, Harris ND. Positionsensing technologies for movement analysis in stroke rehabilitation. Med Biol Eng Comput. 2005;43(4):413–20. doi:10.1007/BF02344720.View ArticleGoogle Scholar
 Wu CY, Lin KC, Chen HC, Chen IH, Hong WH. Effects of modified constraintinduced movement therapy on movement kinematics and daily function in patients with stroke: a kinematic study of motor control mechanisms. Neurorehabil Neural Repair. 2007;21(5):460–6. doi:10.1177/1545968307303411.View ArticleGoogle Scholar
 Zhou H, Hu H. Human motion tracking for rehabilitation—a survey. Biomed Signal Process Control. 2008;3(1):1–18. doi:10.1016/j.bspc.2007.09.001.View ArticleGoogle Scholar
 Stephenson JL, Lamontagne A, De Serres SJ. The coordination of upper and lower limb movements during gait in healthy and stroke individuals. Gait Posture. 2009;29(1):11–6. doi:10.1016/j.gaitpost.2008.05.013.View ArticleGoogle Scholar
 O’Donoghue P. Research methods for sports performance analysis. Abingdon: Routledge; 2010.Google Scholar
 Yang N, Zhang M, Huang C, Jin D. Motion quality evaluation of upper limb targetreaching movements. Med Eng Phys. 2002;24(2):115–20. doi:10.1016/S13504533(01)001217.View ArticleGoogle Scholar
 Vandenberghe A, Levin O, De Schutter J, Swinnen S, Jonkers I. Threedimensional reaching tasks: effect of reaching height and width on upper limb kinematics and muscle activity. Gait Posture. 2010;32(4):500–7. doi:10.1016/j.gaitpost.2010.07.009.View ArticleGoogle Scholar
 Borbély BJ, Straube A, Eggert T. Motor synergies during manual tracking differ between familiar and unfamiliar trajectories. Exp Brain Res. 2013;232(3):1–13. doi:10.1007/s0022101338010.Google Scholar
 Vignais N, Marin F. Analysis of the musculoskeletal system of the hand and forearm during a cylinder grasping task. Int J Ind Ergon. 2014;44(4):535–43. doi:10.1016/j.ergon.2014.03.006.View ArticleGoogle Scholar
 Song D, Lan N, Loeb GE, Gordon J. Modelbased sensorimotor integration for multijoint control: development of a virtual arm model. Ann Biomed Eng. 2008;36(6):1033–48. doi:10.1007/s1043900894618.View ArticleGoogle Scholar
 Park MS, Chung CY, Lee SH, Choi IH, Cho TJ, Yoo WJ, Park BSMY, Lee KM. Effects of distal hamstring lengthening on sagittal motion in patients with diplegia. Hamstring length and its clinical use. Gait Posture. 2009;30(4):487–91. doi:10.1016/j.gaitpost.2009.07.115.View ArticleGoogle Scholar
 Arnold EM, Ward SR, Lieber RL, Delp SL. A model of the lower limb for analysis of human movement. Ann Biomed Eng. 2010;38(2):269–79. doi:10.1007/s1043900998525.View ArticleGoogle Scholar
 Veeger DHEJ. “ What if”: the use of biomechanical models for understanding and treating upper extremity musculoskeletal disorders. Man Ther. 2011;16(1):48–50. doi:10.1016/j.math.2010.09.004.View ArticleGoogle Scholar
 Gustus A, Stillfried G, Visser J. Human hand modelling: kinematics, dynamics, applications. Biol Cybern. 2012;106:741–55. doi:10.1007/s0042201205324.MathSciNetView ArticleMATHGoogle Scholar
 Bolsterlee B, Veeger DHEJ, Chadwick EK. Clinical applications of musculoskeletal modelling for the shoulder and upper limb. Med Biol Eng Comput. 2013;51(9):953–63. doi:10.1007/s1151701310995.View ArticleGoogle Scholar
 Luinge HJ, Veltink PH. Measuring orientation of human body segments using miniature gyroscopes and accelerometers. Med Biol Eng Comput. 2005;43(2):273–82. doi:10.1007/BF02345966.View ArticleGoogle Scholar
 Sabatelli S, Galgani M, Fanucci L, Rocchi A. A double stage Kalman filter for sensor fusion and orientation tracking in 9D IMU. In: Sensors applications symposium (SAS). New York: IEEE; 2012. p. 1–5.
 Madgwick SOH, Harrison AJL, Vaidyanathan A. Estimation of IMU and MARG orientation using a gradient descent algorithm. IEEE Int Conf Rehabil Robot. 2011;2011:5975346. doi:10.1109/ICORR.2011.5975346.Google Scholar
 Mahony R, Hamel T, Pflimlin JM. Nonlinear complementary filters on the special orthogonal group. IEEE Trans Autom Control. 2008;53(5):1203–18. doi:10.1109/TAC.2008.923738.MathSciNetView ArticleGoogle Scholar
 Tian Y, Wei H, Tan J. An adaptivegain complementary filter for realtime human motion tracking with MARG sensors in freeliving environments. IEEE Trans Neural Syst Rehabil Eng. 2013;21(2):254–64. doi:10.1109/TNSRE.2012.2205706.View ArticleGoogle Scholar
 Olivares A, Górriz JM, Ramírez J, Olivares G. Accurate human limb angle measurement: sensor fusion through Kalman, least mean squares and recursive leastsquares adaptive filtering. Meas Sci Technol. 2011;22(2):25801.View ArticleGoogle Scholar
 Cutti AG, Giovanardi A, Rocchi L, Davalli A, Sacchetti R. Ambulatory measurement of shoulder and elbow kinematics through inertial and magnetic sensors. Med Biol Eng Comput. 2008;46(2):169–78. doi:10.1007/s1151700702965.View ArticleGoogle Scholar
 Kontaxis A, Cutti AG, Johnson GR, Veeger HEJ. A framework for the definition of standardized protocols for measuring upperextremity kinematics. Clin Biomech. 2009;24(3):246–53. doi:10.1016/j.clinbiomech.2008.12.009.View ArticleGoogle Scholar
 de Vries WHK, Veeger HEJ, Cutti AG, Baten C, van der Helm FCT. Functionally interpretable local coordinate systems for the upper extremity using inertial & magnetic measurement systems. J Biomech. 2010;43(10):1983–8. doi:10.1016/j.jbiomech.2010.03.007.View ArticleGoogle Scholar
 Parel I, Cutti AG, Fiumana G, Porcellini G, Verni G, Accardo AP. Ambulatory measurement of the scapulohumeral rhythm: intra and interrater reliability of a protocol based on inertial and magnetic sensors. Gait Posture. 2012;35:636–40. doi:10.3233/9781607500803164.View ArticleGoogle Scholar
 Holzbaur KRS, Murray WM, Delp SL. A model of the upper extremity for simulating musculoskeletal surgery and analyzing neuromuscular control. Ann Biomed Eng. 2005;33(6):829–40. doi:10.1007/s1043900533207.View ArticleGoogle Scholar
 Wu G, Van Der Helm FCT, Veeger HEJ, Makhsous M, Van Roy P, Anglin C, Nagels J, Karduna AR, McQuade K, Wang X, Werner FW, Buchholz B. ISB recommendation on definitions of joint coordinate systems of various joints for the reporting of human joint motion—Part II: Shoulder, elbow, wrist and hand. J Biomech. 2005;38(5):981–92. doi:10.1016/j.jbiomech.2004.05.042. [arXiv:111].View ArticleGoogle Scholar
 Vicon Oxford Foot Model. https://www.vicon.com/products/software/oxfordfootmodel. Accessed 7 Nov 2016.
 Motive: body software for OptiTrack. http://www.optitrack.com/products/motive/body/. Accessed 7 Nov 2016.
 van den Bogert AJ, Geijtenbeek T, EvenZohar O, Steenbrink F, Hardin EC. A realtime system for biomechanical analysis of human movement and muscle function. Med Biol Eng Comput. 2013;51(10):1069–77. doi:10.1007/s115170131076z.View ArticleGoogle Scholar
 Delp SL, Loan JP, Hoy MG, Zajac FE, Topp EL, Rosen JM. An interactive graphicsbased model of the lower extremity to study orthopaedic surgical procedures. IEEE Trans Biomed Eng. 1990;37(8):757–67. doi:10.1109/10.102791.View ArticleGoogle Scholar
 Delp SL, Anderson FC, Arnold AS, Loan P, Habib A, John CT, Guendelman E, Thelen DG. OpenSim: opensource software to create and analyze dynamic simulations of movement. IEEE Trans Biomed Eng. 2007;54(11):1940–50. doi:10.1109/TBME.2007.901024.View ArticleGoogle Scholar
 Muceli S, Farina D. Simultaneous and proportional estimation of hand kinematics from EMG during mirrored movements at multiple degreesoffreedom. IEEE Trans Neural Syst Rehabil Eng. 2012;20(3):371–8. doi:10.1109/TNSRE.2011.2178039.View ArticleGoogle Scholar
 Jiang N, VestNielsen JL, Muceli S, Farina D. EMGbased simultaneous and proportional estimation of wrist/hand dynamics in unilateral transradial amputees. J Neuroeng Rehabil. 2012;9(1):42. doi:10.1186/17430003942.View ArticleGoogle Scholar
 Jiang N, Muceli S, Graimann B, Farina D. Effect of arm position on the prediction of kinematics from EMG in amputees. Med Biol Eng Comput. 2013;51(1–2):143–51. doi:10.1007/s1151701209794.View ArticleGoogle Scholar
 Borbély BJ, Szolgay P. Estimating the instantaneous wrist flexion angle from multichannel surface EMG of forearm muscles. In: 2013 IEEE biomedical circuits and systems conference, BioCAS. New York: IEEE; 2013. p. 77–80.
 Blana D, Kyriacou T, Lambrecht JM, Chadwick EK. Feasibility of using combined EMG and kinematic signals for prosthesis control: a simulation study using a virtual reality environment. J Electromyogr Kinesiol. 2015;29:21–7. doi:10.1016/j.jelekin.2015.06.010.View ArticleGoogle Scholar
 Borbély BJ, Szolgay P. A system concept for emg classification from measurement to deployment. In: 2016 15th international workshop on cellular nanoscale networks and their applications (CNNA). 2016. p. 121–122.
 Schmidhuber J. Deep learning in neural networks: an overview. Neural Netw. 2015;61:85–117. doi:10.1016/j.neunet.2014.09.003. [arXiv:1404.7828].View ArticleGoogle Scholar
 Holzbaur KRS, Murray WM, Delp SL. Upper extremity kinematic model, Simtk resource. https://simtk.org/home/upextmodel. Accessed 6 Jul 2016.
 Delp SL, Loan P, Krystyne B. SIMM 7.0 for windows user’s manual. 2013. http://www.musculographics.com/download/SIMM7.0UserGuide.pdf. Accessed 6 Jul 2016.
 Hicks J. OpenSim documentations: how scaling works. http://simtkconfluence.stanford.edu:8080/display/OpenSim/How+Scaling+Works. Accessed 6 Jul 2016.
 Hicks J. OpenSim documentations: how inverse kinematics works. http://simtkconfluence.stanford.edu:8080/display/OpenSim/How+Inverse+Kinematics+Works. Accessed 6 Jul 2016.
 Press WH, Teukolsky SA, Vetterling WT, Flannery BP. Numerical recipes in the art of scientific computing. 2nd ed. New York: Cambridge University Press; 1992. p. 359–62.MATHGoogle Scholar
 Borbély BJ, Tihanyi A, Szolgay P. A measurement system for wrist movements in biomedical applications. In: 2015 European conference on circuit theory and design (ECCTD). New York: IEEE. p. 1–4.
 Bonnet S, Bassompierre C, Godin C, Lesecq S, Barraud A. Calibration methods for inertial and magnetic sensors. Sens Actuators A Phys. 2009;156(2):302–11. doi:10.1016/j.sna.2009.10.008.View ArticleGoogle Scholar
 Vanegas M, Stirling L. Characterization of inertial measurement unit placement on the human body upon repeated donnings. In: 2015 IEEE 12th international conference on wearable and implantable body sensor networks (BSN). New York: IEEE; 2015. p. 1–6.
 Piovan G, Bullo F. On coordinatefree rotation decomposition: Euler angles about arbitrary axes. IEEE Trans Robot. 2012;28(3):728–33. doi:10.1109/TRO.2012.2184951.View ArticleGoogle Scholar