Joint Center: Functional Method
Up ] Orientation Angles ] Computation of the Orientation Angles ] Joint Angle vs. Orientation Angles ] Angular Velocity vs. Orientation Angles ] Electromagnetic Motion Sensors ] [ Joint Center: Functional Method ] Computation of the Rotation Matrix ] Helical (Screw) Axis ] Motion Sensors: Joint Center ]

HU Gamage & Lasenby (2002)
Nelder-Mead Downhill Simplex Method
Halvorsen, et al. (1999)
Newton's Method
Simple Linear Method
References and Related Literature


The functional method is useful in locating the center of rotation of a ball-and-socket joint or the axis of rotation of a hinge joint. It computes the center/axis of rotation numerically rather than relying on a set of predetermined ratios or prediction equations.

The hip joint is a ball-and-socket joint and has 3 degrees of freedom. In an ideal situation, markers attached to the thigh will form spheres with the hip joint being the center if the pelvis is stationary. Thus, it is possible to numerically compute the center of the spheres based on the marker coordinates.

Four methods (algorithms) are included in this page: 2 closed-form least square algorithms by HU Gamage & Lasenby (2002), and Halvorsen et al. (1999), respectively, the Nelder-Mead downhill simplex method, and the Newton's method.


HU Gamage & Lasenby (2002)

Acknowledgement: The author is grateful to Dr. Nels Madsen of Auburn University for his invaluable contribution to this page. He kindly provided a detailed rundown of the equations.

Fig. 1 shows the geometric relationship between the hip joint and the markers attached to the thigh. Each marker on the thigh will form a sphere with the hip joint being its center. Vector in Fig. 1 is the position of the hip joint center while vector is the position of marker j at frame i.

    [Fig. 1]

The error function can be defined as

,    [1]

where Rj = radius of the sphere formed by marker j. Thus, the cost function for the least square approximation of the hip joint center becomes 

.    [2]

The cost function U in [1] must be minimized to obtain the position of the center of the spheres (hip joint). 

Now, let's define the mean position of marker j shown in Fig. 1 as:

.    [3]

The relative position of the marker to the mean position at frame i can be expressed in terms of the mean position: 



.    [5]

Vector shown in Fig. 1 is the relative position of the joint center to the mean position of marker j:

.    [6]



The cost function shown in [1] can be then rewritten as

.    [8]

Minimization of the cost function U yields



.    [10]

From [9]:

,    [11]

where .

From [10]:

.    [12]

Note that [5] and [11] were used in deriving [12]. From [12]:

.    [13]

Rewrite [13] in the matrix form and divide both sides by n:

,    [14]

where { } = column matrix of the vector. [14] basically has the following form:

,    [15]




.    [17]

Aj = (3 x 3) while Bj = (3 x 1). From [15], [16], and [17]:

.    [18]

[18] is basically in the form of 

,    [19]



Expand [19] and [20] to all markers:

,    [21]

where C = (3 x 3) while D = (3 x 1). [19] is a linear system of equations of the position of the hip joint.

From [21]:

.    [22]

This method was introduced by H. U. Gamage and Lasenby (2002). The advantage of this method is that it can handle multiple markers at once. The basic assumption is that the spheres formed by the markers have a common center (hip joint).


Nelder-Mead Downhill Simplex Method

Another algorithm that can be used in finding the hip joint center is the so-called downhill simplex method following Nelder and Mead. This method performs multidimensional minimization of the given function. In other words, it finds the minimum of the given cost function and the corresponding position of the hip joint center. Since this method does not involve derivatives, the cost function can be written as either




Although [24] looks more complex than [23], it requires less computations since the magnitude of a vector is typically computed from its square. that minimizes U is the position of the hip.

The uniqueness of this method is in the method of systematic perturbation of through iterations in multi dimensions to obtain the minimum of function U. It uses the so-called simplex for this and tetrahedron defined by 4 vertices (points) is the simplex in three- dimension. In other words, four position vectors instead of one will be used in obtaining the position of the hip joint center. If , the initial guess of , is available, the four vertices of the starting simplex can be given as

,    [25]

where = constant, and i, j, & k = unit vectors of the axes. For example, the overall mean position of the markers can be used as the initial guess.

Each iteration consists of the following steps:

  1. Identification of the vertices that provide the largest, second largest, and smallest function value: high, next high, and low point, respectively.
  2. Convergence check. If it is sufficiently converged, read the low point coordinates as the hip joint center and stop the iterations.
  3. Total number of iterations. See if it exceeded the maximum allowed. If yes, stop the iterations.
  4. Determination of the type of simplex manipulation (Fig.2). Fig. 2a shows the simplex at the start of the current iteration. 

        [Fig. 2] (Adapted from Press et al., 1986)

There are 6 different types of simplex manipulations possible:

  1. Reflection: moving the high point through the opposite face to a lower point (Fig. 2b) or

,    [26]

where vr = reflected vertex position, j = the high point, = Kronecker delta, and v = vertex position.

  1. Reflection & expansion (Fig. 2c): extrapolation of the reflected point by a factor of 2, or

,    [27]

where vre = reflected and expanded vertex position.

  1. Contraction (Fig. 2d): one-dimensional contraction of the high point toward the opposite face, or

,    [28]

where vc = contracted vertex position.

  1. Reflection & contraction (Fig. 2e): one-dimensional contraction of the reflected point, or

,    [29]

where vrc = reflected and contracted vertex position.

  1. Multiple contraction (Fig. 2f): multi-dimensional contraction of the vertices toward the low point, or

,    [30]

where vmci = contracted position of vertex i toward the low point (k), and vk = position of the low point.

  1. Reflection and multiple contraction (Fig. 2f): multi-dimensional contractions of the vertices toward the low point after reflection, or

,    [31]

where vrmci = reflected and contracted position of vertex i toward the low point (k).

Through the different simplex manipulations described above, the method approaches to the 'valley floor' and provides minimum of the cost function and location of the hip joint center. See Press et al. (1986) for details of the downhill simplex method and the sample function Amoeba. The iterations can be stopped if the change in the function value is fractionally smaller than a given tolerance:

,    [32]

where j = high point, k = low point, and U = cost function.

It is customary to start the iteration once again w/ the low point obtained from the previous session being one of the starting vertices to make sure that the converged low point is the true solution.


Halvorsen et al. (1999)

This method is based on a simple geometric relationship between the marker positions at two different instants. If the pelvis is stationary, the thigh markers always lie on the spheres they form with the hip joint being the center. In this case the line connecting two different positions of the same marker at two different instants is always perpendicular to the line drawn from the hip joint to their midpoint (Fig.3):

    Fig. 3

.    [33]

The error function is then


Let the cost function U be

,    [35]

which needs to be minimized. Minimization of the cost function U shown above has two meanings: (a) keep the angle between the two vector as close to 90 degrees as possible, and (b) keep the magnitude of the vector drawn from the joint center to the midpoint as small as possible. From [33]:

.    [36]

Covert [36] to the equivalent matrix form:

,    [37]

where {} = column matrix operator. Expand [37] to m markers:

,    [38]


,    [39]


.    [40]

Form [38]:

.    [41]

Theoretically, [41] can be expanded to multiple two-frame combinations to increase the redundancy. One issue in this method is to determine the frames to be used in the computation.


Newton's Method

For a given marker, let the random error function for marker j at frame i be

.    [42]

For n frames and m markers, we will have a nonlinear system of m x n equations. The partial derivatives of function are

,    [43]


.    [44]

Therefore, the Jacobian matrix of the nonlinear system of equations of interest can be written as

.    [45]

The Newton's method involves iterations and each iteration consists of three steps. The first step is to solve the following linear system of equations:

,    [46]

where, k = the current iteration, = the solution from the previous iteration,

,    [47]


.    [48]

In order for [46] to work, the initial guess of rc and Rj's, , must be provided. The mean coordinates of the marker may be used as the initial value.

The second step in each iteration is to update rc and Rj's using the solution obtained in the previous step:

.    [49]

The last step is to check if the solution is sufficiently converged:

,    [50]

where TOL = tolerance. The iteration will stop if the solution is sufficiently converged.


Simple Linear Method

Acknowledgement: This method was suggested by Dr. Yong-San Yoon of Korea Advanced Institute of Science and Technology (KAIST).

From [1]:


[51] can be rewritten as

,    [52]


.    [53]

[52] is basically a linear equation of rc and . Expand [52] to n frames and m markers to obtain

.    [54]

Use the typical least square approach to solve the system of equations. Rj's of the markers can be computed afterward using [53], if necessary.


References and Related Literature

Gander, W., Golub, G. H., & Strebel, R. (1994). Fitting of circles and ellipses: least squares solution. Tech Report, Departement Informatik, ETH Zurich.

Halvorsen, K., Lesser, M., & Lundberg, A. (1999). A new method for estimating the axis of rotation and the center of rotation. J. Biomechanics 32, 1221-1227.

Hiniduma Udugama Gamage, S.S., & Lasenby, J. (2002). New least squares solutions for estimating the average center of rotation and the axis of rotation. J. Biomechanics 35, 87-93.

Press, W.H., Flannery, B.P., Teukolsky, S.A., & Vetterling, W.T. (1986). Numerical recipes: the art of scientific computing. New York, NY: Cambridge University Press.

Thompson, M.S., Dawson, T., Kuiper, J.H., Northmore-Ball, M.D., & Tanner, K.E. (2000). Acetabular morphology and resurfacing design. J. Biomechanics 33, 1645-1653.


Young-Hoo Kwon, 1998-