An Elastic Energy Minimization Framework for Mean Surface Calculation

TL;DR

This paper introduces a novel elastic energy minimization framework for calculating mean surfaces in 3D image analysis, emphasizing reparameterization equations and centroid calculation within a specific mathematical space.

Contribution

It presents new equations for surface reparameterization and centroid calculation using RPSN, enhancing surface averaging methods for 3D image delineations.

Findings

RPSN space is invariant under reparameterization and has a meaningful L2 metric.

Provides reparameterization equations for both RPSN and SRNF representations.

Addresses theoretical aspects of displacement removal in RPSN for surface analysis.

Abstract

As the continuation of the contour mean calculation - designed for averaging the manual delineations of 3D layer stack images - in this paper, the most important equations: a) the reparameterization equations to determine the minimizing diffeomorphism and b) the proper centroid calculation for the surface mean calculation are presented. The chosen representation space: escaled Position by Square root Normal (RPSN) is a real valued vector space, invariant under the action of the reparameterization group and the imposed L2 metric (used to define the distance function) has well defined meaning: the sum of the central second moments of the coordinate functions. For comparision purpose, the reparameterization equations for elastic surface matching, using the Square Root Normal Function (SRNF) are also provided. The reparameterization equations for these cases have formal similarity, albeit the targeted applications differ: SRNF representation suitable for shape analysis purpose whereas RPSN is more fit for the cases where all contextual information - including the relative translation between the constituent surfaces - are to be retained (but the sake of theoretical completeness, the possibility of the consistent relative displacement removal in the RPSN case is also addressed).

Equations

$$d_{\mathbf{P}}^2=\iint {\left({\mathbf{P}}_1-{\mathbf{P}}_2\right)}^{2}dudv$$

$$d_{\mathbf{Q}}^2=\iint {\left({\mathbf{Q}}_1-{\mathbf{Q}}_2\right)}^{2}dudv$$

$$E\left(\mu ,\nu \right)$$

$$\left[J_{\left(\mu ,\nu \right)}\right]$$

$$\left|J_{\left(\mu ,\nu \right)}\right|$$

$$\left[\begin{array}[]{c}\nabla_{\mu}E\\ \nabla_{\nu}E\end{array}\right]=\sqrt{\left|\mathbf{M}\right|\left|\mathbf{N}\right|}\left[J_{\left(\mu,\nu\right)}\right]^{-1}\left[\begin{array}[]{c}\mathbf{m}_{u}\cdot\mathbf{n}-\mathbf{m}\cdot\mathbf{n}_{u}+\frac{1}{2}\mathbf{m}\cdot\mathbf{n}\left(\Gamma_{u}^{\mathbf{R}}-\Gamma_{u}^{\mathbf{S}}\right)\\ \mathbf{m}_{v}\cdot\mathbf{n}-\mathbf{m}\cdot\mathbf{n}_{v}+\frac{1}{2}\mathbf{m}\cdot\mathbf{n}\left(\Gamma_{v}^{\mathbf{R}}-\Gamma_{v}^{\mathbf{S}}\right)\end{array}\right],$$

$${\mathbf{m}}_u\cdot \mathbf{n}-\mathbf{m}\cdot {\mathbf{n}}_u+\frac{1}{2}\mathbf{m}\cdot \mathbf{n}\left({\Gamma }_u^{\mathbf{R}}-{\Gamma }_u^{\mathbf{S}}\right)$$

$${\mathbf{m}}_v\cdot \mathbf{n}-\mathbf{m}\cdot {\mathbf{n}}_v+\frac{1}{2}\mathbf{m}\cdot \mathbf{n}\left({\Gamma }_v^{\mathbf{R}}-{\Gamma }_v^{\mathbf{S}}\right)$$

$$F\left(\mu ,\nu \right)$$

$$\left[\begin{array}[]{c}\nabla_{\mu}F\\ \nabla_{\nu}F\end{array}\right]=\sqrt{\left|\mathbf{M}\right|\left|\mathbf{N}\right|}\left[J_{\left(\mu,\nu\right)}\right]^{-1}\left[\begin{array}[]{c}\mathbf{R}_{u}\cdot\mathbf{S}-\mathbf{R}\cdot\mathbf{S}_{u}+\frac{1}{2}\mathbf{R}\cdot\mathbf{S}\left(\Gamma_{u}^{\mathbf{R}}-\Gamma_{u}^{\mathbf{S}}\right)\\ \mathbf{R}_{v}\cdot\mathbf{S}-\mathbf{R}\cdot\mathbf{S}_{v}+\frac{1}{2}\mathbf{R}\cdot\mathbf{S}\left(\Gamma_{v}^{\mathbf{R}}-\Gamma_{v}^{\mathbf{S}}\right)\end{array}\right]\,.$$

$${\mathbf{R}}_u\cdot \mathbf{S}-\mathbf{R}\cdot {\mathbf{S}}_u+\frac{1}{2}\mathbf{R}\cdot \mathbf{S}\left({\Gamma }_u^{\mathbf{R}}-{\Gamma }_u^{\mathbf{S}}\right)$$

$${\mathbf{R}}_v\cdot \mathbf{S}-\mathbf{R}\cdot {\mathbf{S}}_v+\frac{1}{2}\mathbf{R}\cdot \mathbf{S}\left({\Gamma }_v^{\mathbf{R}}-{\Gamma }_v^{\mathbf{S}}\right)$$

$$\underset{{\mu }_i,{\nu }_i}{\min }\stackrel{[}{i}=0]N-1\sum \iint {\left(\mathbf{Q}-{\mathbf{Q}}_{\mathbf{i}}\right)}^{2}dudv,$$

$${\mathbf{Q}}_{\mathbf{i}}={\mathbf{S}}_{\mathbf{i}}\left({\mu }_i,{\nu }_i\right)\sqrt{\left|J_{\left({\mu }_i,{\nu }_i\right)}\right|\left|{{\mathbf{S}}_{\mathbf{i}}}_{{\mu }_i}\left({\mu }_i,{\nu }_i\right)\times {{\mathbf{S}}_{\mathbf{i}}}_{v_i}\left({\mu }_i,{\nu }_i\right)\right|}$$

$$={\mathbf{S}}_{\mathbf{i}}\left(u,v\right)\sqrt{\left|{{\mathbf{S}}_{\mathbf{i}}}_u\left(u,v\right)\times {{\mathbf{S}}_{\mathbf{i}}}_v\left(u,v\right)\right|}$$

$$\underset{{\mu }_i,{\nu }_i;\delta \mathbf{D}}{\min }\stackrel{[}{i}=0]N-1\sum \iint {\left[\left(\mathbf{S}-\delta \mathbf{D}\right)\sqrt{\left|\mathbf{N}\right|}-\left({\mathbf{S}}_{\mathbf{i}}-\delta \mathbf{D}\right)\sqrt{\left|{\mathbf{N}}_{\mathbf{i}}\right|}\right]}^{2}dudv,$$

$$\frac{\partial }{\partial \delta \mathbf{D}}\stackrel{[}{i}=0]N-1\sum \iint {\left[\left(\mathbf{S}-\delta \mathbf{D}\right)\sqrt{\left|\mathbf{N}\right|}-\left({\mathbf{S}}_{\mathbf{i}}-\delta \mathbf{D}\right)\sqrt{\left|{\mathbf{N}}_{\mathbf{i}}\right|}\right]}^{2}dudv\doteq 0,$$

$$\delta \mathbf{D}=\frac{\stackrel{[}{i}=0]N-1\sum \iint \left({\mathbf{S}}_{\mathbf{i}}\sqrt{\left|{\mathbf{N}}_{\mathbf{i}}\right|}-\mathbf{S}\sqrt{\left|\mathbf{N}\right|}\right)\sqrt{\left|{\mathbf{N}}_{\mathbf{i}}\right|}dudv}{\stackrel{[}{i}=0]N-1\sum \iint {\left(\sqrt{\left|{\mathbf{N}}_{\mathbf{i}}\right|}-\sqrt{\left|\mathbf{N}\right|}\right)}^{2}dudv}=0.$$

$$\underset{\delta {\mathbf{D}}_{\mathbf{i}}}{\min }\stackrel{[}{i}=0]N-1\sum \iint {\left[\mathbf{S}\sqrt{\left|\mathbf{N}\right|}-\left({\mathbf{S}}_{\mathbf{i}}-\delta {\mathbf{D}}_{\mathbf{i}}\right)\sqrt{\left|{\mathbf{N}}_{\mathbf{i}}\right|}\right]}^{2}dudv,$$

$$\frac{\partial }{\partial \delta \mathbf{D}}\stackrel{[}{i}=0]N-1\sum \iint {\left[\mathbf{S}\sqrt{\left|\mathbf{N}\right|}-\left({\mathbf{S}}_{\mathbf{i}}-\delta {\mathbf{D}}_{\mathbf{i}}\right)\sqrt{\left|{\mathbf{N}}_{\mathbf{i}}\right|}\right]}^{2}dudv\doteq 0,\mathbf{i}=\mathrm{0...}N-1$$

$$\delta {\mathbf{D}}_{\mathbf{i}}=\frac{\iint \left({\mathbf{S}}_{\mathbf{i}}\sqrt{\left|{\mathbf{N}}_{\mathbf{i}}\right|}-\mathbf{S}\sqrt{\left|\mathbf{N}\right|}\right)\sqrt{\left|{\mathbf{N}}_{\mathbf{i}}\right|}dudv}{\iint \left|{\mathbf{N}}_{\mathbf{i}}\right|dudv}=0.$$

$${\Gamma }_{ik}^l={\mathbf{S}}^l\cdot {\mathbf{S}}_{ik}$$

$$i,k,l\in \left\{u,v\right\},$$

$${\mathbf{S}}^u$$

$${\mathbf{S}}^v$$

$$\frac{\partial \ln \left|\mathbf{N}\right|}{\partial u}=\frac{1}{\left|\mathbf{N}\right|}\mathbf{n}\cdot \left({\mathbf{S}}_{uu}\times {\mathbf{S}}_v+{\mathbf{S}}_u\times {\mathbf{S}}_{uv}\right)$$

$$={\mathbf{S}}_{uu}\cdot {\mathbf{S}}^u+{\mathbf{S}}_{uv}\cdot {\mathbf{S}}^v$$

$$={\Gamma }_{uu}^u+{\Gamma }_{vu}^v$$

$$\frac{\partial \ln \left|\mathbf{N}\right|}{\partial v}=\frac{1}{\left|\mathbf{N}\right|}\mathbf{n}\cdot \left({\mathbf{S}}_{uv}\times {\mathbf{S}}_v+{\mathbf{S}}_u\times {\mathbf{S}}_{vv}\right)$$

$$={\mathbf{S}}_{uv}\cdot {\mathbf{S}}^u+{\mathbf{S}}_{vv}\cdot {\mathbf{S}}^v$$

$$={\Gamma }_{uv}^u+{\Gamma }_{vv}^v$$

$${\mathbf{n}}_i\cdot {\mathbf{S}}_k=-\mathbf{n}\cdot {\mathbf{S}}_{ik}$$

$$i,k\in \left\{u,v\right\}.$$

$$\mathbf{I}$$

$${\mathbf{P}}_u$$

$${\mathbf{P}}_v$$

Full text

1 Inroduction

Object delineation is an important annotation step to create training data set for the supervised machine learning methods designed for object segmentation. In cases wherever the object boundaries are not definite (blurred, ambiguous), delineations performed by experts often do not agree. A plausible approach to create meaningful annotation samples is to accept the mean of many recommendations excluding some outliers. This approach requires well defined, meaningful metrics on the space of contours (2D) or surfaces (3D). A specifically designed contour representation: the Rescaled Position by Square Velocity (RPSV): $\mathbf{q}\left(t\right)=\mathbf{r}\left(t\right)\sqrt{\dot{\mathbf{r}}\left(t\right)}$ and a complete mathematical framework for the contour averaging (and interpolation) problem were proposed and examined in the paper [5]. RPSV can be considered as the mixing of two: the Kendall’s ’landmark points’ [4] and the Square Root Velocity Function (SRVF) [2][3][6] representations. The description of the contours by preselected position vector set (the landmark points) is the approach of the early shape analysis techniques. The drawback of this representation is that the results (e.g. mean object) is dependent on the predetermined sampling strategy of the landmark points (corresponding to fixed parameterization). The continuous SRVF representation provides the necessary freedom in the form of the optimal reparameterization of the contours. The results are much much more intuitive. On the other hand, since the description of the contours is velocity-vector based, the relative displacement between the constituent contours cannot be retrieved (insomuch as lost by derivation) - an important information for delineation statistics. RPSV contour representation was designed to keep all contextual information (including the relative displacenents of the constituent contours) and covariant description.

Recently, a new representation: the Square Root Normal Function (SRNF), as the 3D generalization of the SRVF was introduced in [1] for elastic surface analysis. In this paper, the 3D generalization

• the ’Rescaled Position by Square root Normal’ (RPSN) - of the RPSV is introduced along with the reparameterization and the proper centroid equations. Only these equations are presented, because the theoretical results as the possibility of the pairwise determination of the optimal parameterization system $\gamma_{i}\left(t\right)$, $\mathbf{q}_{i}\left(t\right)\rightarrow\mathbf{q}_{i}\left(\gamma_{i}\left(t\right)\right)$, $i=0\ldots n-1$ wrt an arbitrarily chosen reference contour (e.g. $\gamma_{0}\left(t\right)\equiv t$) or the consistency of the mean formula $\mathbf{q}\left(t\right)=\frac{1}{n}\underset{i}{\sum}\mathbf{q}_{i}\left(\gamma_{i}\left(t\right)\right)$ with the reparameterization equations, can be repeated and formally transferred to 3D. The reparameterization equations are compared with the corresponding equations of the SRNF. For the sake of theoretical completeness, the possibility of the consistent relative displacement removal (hence the generalization of the landmark point based approach to the continuous case) is also shown.

The structure of the paper is the following. Section 2 summarizes the SRNF and RPSN frameworks. Section 3 is introduces to the reparameterization equations as variational problems, section 4 is dedicated to the proper centroid calculation (for RPSN) with reference to the possibility of relative displacement removal between the constituent surfaces (4.1). Section 5 concludes the paper with discussion. Appendices contain the important derivations.

2 Elastic surface analysis frameworks

We consider the primary space of the smooth surfaces embedded in $\mathbb{R}^{3}$ as the space of two-parameter coordinate function-triplets $X\left(u,v\right),\,Y\left(u,v\right),\,Z\left(u,v\right)$, or - equivalently - two-parameter position vectors $\mathbf{S}\left(u,v\right)=X\mathbf{i}+Y\mathbf{j}+Z\mathbf{k}$ wrt some basis $\mathbf{i},\,\mathbf{j},\,\mathbf{k}$ (known as Gaussian descriprion). The important derived quantities associated with this description are: a) the local (covariant basis): $\mathbf{S}_{u}=\frac{\partial\mathbf{S}}{\partial u}$, $\mathbf{S}_{v}=\frac{\partial\mathbf{S}}{\partial v}$, b) the normal vector of the surface $\mathbf{N}=\mathbf{S}_{u}\times\mathbf{S}_{v}$, c) its length $\left|\mathbf{N}\right|$ - which is also the square root of the determinant of the metric tensor with components $g_{ik}=\mathbf{S}_{i}\times\mathbf{S}_{k},\,i,k\in\left[u,v\right]$, $g=\det\left[g_{ik}\right]$,* i.e.* $\sqrt{g}=\left|\mathbf{N}\right|$ and d) the unit (paremeterization independent) normal vector $\mathbf{n}=\frac{\mathbf{N}}{\left|\mathbf{N}\right|}$.

Both representations SRNF: $\mathbf{P}\doteq\mathbf{n}\sqrt{\left|\mathbf{S}_{u}\times\mathbf{S}_{v}\right|}$ and RPSN: $\mathbf{Q}\doteq\mathbf{S}\sqrt{\left|\mathbf{S}_{u}\times\mathbf{S}_{v}\right|}$ can be considered as the “change of coordinates” in the space of surfaces, albeit these mappings do not necessarily lead bijections between the original and transformed coordinates. In the SRNF case, any constant translation $\mathbf{d}$ between two surfaces $\mathbf{S}$ and $\tilde{\mathbf{S}}=\mathbf{S}+\mathbf{d}$ is removed by derivation, so that in SRNF ’coordinates’ the surfaces are determined only up to an arbitrary translation, naturally representing the quotient space of the surfaces wrt translation - a useful property for the shape analysis.

These representetions have the distinguishing property that is the distance functions defined between two points in the space of contours, using the respective $\mathbb{L}^{2}$ norms:

[TABLE]

and

[TABLE]

are invariant under the action of the reparameterization group $\varGamma$ (in fact invariant quantities under the product group of the rotation and reparameterization as well), because both represent parameterization independent quantities: the surface area (SRNF) and the second central moment (RPSN) of the surfaces. Among the consequences, a unique distance over the equivalence class of the reparameterization group $\varGamma$ can be defined. Also, for multiple surface problems (e.g. averaging) any of the constituent surfaces can be chosen as reference surface (wrt which the best reparameterizations of the other constituents are to be determined). The transformation of the representations under the effect of an arbitrary element $u\rightarrow\mu\left(u,v\right)$, $\nu=\nu\left(u,v\right)$, $\left(\mu,\nu\right)\in\varGamma$ of the reparameterization group are $\mathbf{P}\rightarrow\mathbf{P}\sqrt{\left|J_{\left(\mu,\nu\right)}\right|}$ and $\mathbf{Q}\rightarrow\mathbf{Q}\sqrt{\left|J_{\left(\mu,\nu\right)}\right|}$, where the $\left|J_{\left(\mu,\nu\right)}\right|$ is the determinant of the Jacobian \left[J_{\left(\mu,\nu\right)}\right]=\left[\begin{array}[]{cc}\frac{\partial\mu}{\partial u}&\frac{\partial\mu}{\partial v}\\ \frac{\partial\nu}{\partial u}&\frac{\partial\nu}{\partial v}\end{array}\right] of the reparameterization. Reparameterization acting in the representation space, require the updation of both the lengths and the directions of the points of $\mathbf{P}$ or $\mathbf{Q}$. Alternatively however, the distance minimization can be formulated directly in the surface space (where only the direction of the points need to be updated), the approach pursued in this paper.

3 Optimal reparameterization as variational problem

Let $\mathbf{R}\left(u,v\right)$ with normal vector $\mathbf{M}=\mathbf{R}_{u}\times\mathbf{R}_{v}$ and unit normal $\mathbf{m}=\frac{\mathbf{M}}{\left|\mathbf{M}\right|}$ be the ’reference’ surface and $\mathbf{S}\left(\mu,\nu\right)$ with normal vectors $\mathbf{N}_{\left(\mu,\nu\right)}=\mathbf{S}_{\mu}\times\mathbf{S}_{\nu}$, $\mathbf{n}=\frac{\mathbf{N}}{\left|\mathbf{N}\right|}$ another surface. Let the $\mu\left(u,v\right)$, $\nu\left(u,v\right)$ function duplet is the element of the reparameterization group $\varGamma$. We need the surface $\mathbf{S}\left(u,v\right)=\mathbf{S}\left(\mu\left(u,v\right),\nu\left(u,v\right)\right)$ (with normal $\mathbf{N}=\mathbf{S}_{u}\times\mathbf{S}_{v}=\left|J_{\left(\mu,\nu\right)}\right|\mathbf{N}_{\left(\mu,\nu\right)}$) to be optimally parameterized to the reference surface wrt the distance function (1). Then the minimum distances between surfaces can be formulated as variational problems using the ’direct surface coordinates’ in the space of surfaces instead of their representations. In the case (1) the functional to be minimized is

[TABLE]

where $\left|J_{\left(\mu,\nu\right)}\right|$ stands for the determinant of the Jacobian of the reparameterization:

[TABLE]

i.e. the transformation between the normal vectors $\mathbf{N}=\mathbf{S}_{u}\times\mathbf{S}_{v}=\left|J_{\left(\mu,\nu\right)}\right|\mathbf{N}_{\left(\mu,\nu\right)}$ is used. The Euler-Lagrange equations associated with the minimization problem (3) (see Appendix A) are

[TABLE]

where $\Gamma_{w}^{\mathbf{R}}=\frac{\partial\ln\left|\mathbf{M}\right|}{\partial w}$, $\Gamma_{w}^{\mathbf{S}}=\frac{\partial\ln\left|\mathbf{N}\right|}{\partial w}$, $w\in\left[u,v\right]$ are the Christoffel divergences of the surfaces $\mathbf{R}$ and $\mathbf{S}$ respectively. At around the identity diffeomorphism these equations are simplified to the

[TABLE]

ones that need to be solved for the new point positions on the fix-shaped surface $\mathbf{S}$.

The associated functional to the problem (2) is

[TABLE]

and the associated Euler-Lagrange equations (see Appendix B) are:

[TABLE]

At around the identity diffeomorphism these equations are simplified to

[TABLE]

As in the two dimensional (contour) case, it can be shown that for the mean surface problem:

[TABLE]

the solution is $\mathbf{Q}\left(u,v\right)=\frac{1}{N}\stackrel{{\scriptstyle[}}{{i}}=0]{N-1}{\sum}\mathbf{Q_{i}}\left(u,v\right)$, where the constituents $\left(\mathbf{Q_{i}}\right)$ are all optimally parameterized wrt a freely chosen reference surface (say $\mathbf{Q_{0}}$), i.e. $\mathbf{Q_{i}}\left(u,v\right)=\mathbf{Q_{i}}\left(\mu_{i}\left(u,v\right),\nu_{i}\left(u,v\right)\right)$, $i=1...N-1$, where $\mu_{i}\left(u,v\right),\>\nu_{i}\left(u,v\right)$ are the optimal diffeomorphisms - the solution of the Euler-Lagrange equation systems (9) - determined pairwise between $\mathbf{Q_{0}}=\mathbf{R}\left(u,v\right)\sqrt{\left|\mathbf{R}_{u}\times\mathbf{R}_{v}\right|}$ and

[TABLE]

since $\mathbf{S_{i}}\left(u,v\right)=\mathbf{S_{i}}\left(\mu_{i}\left(u,v\right),\nu_{i}\left(u,v\right)\right)$.

4 Proper centroid

In the RPSN case, the Euler-Lagrange equations for the optimal reparameterization system can be solved for surfaces (that are given as position vectors $\mathbf{S}=X\mathbf{i}+Y\mathbf{j}+Z\mathbf{k}$) wrt any point in the space designated as origin, but the result is dependent on the choice of the origin. However there is an optimal choice of the origin compatible with the minimization problem (10). The new energy comprising the position of this optimal cetroid - wrt the current origin - can be formulated as double-minimization problem:

[TABLE]

where the new representations (wrt the displaced origin) are $\mathbf{\mathbf{Q}^{\delta D}=\left(\mathbf{S}-\delta\mathbf{D}\right)}\sqrt{\left|\mathbf{N}\right|}=\frac{1}{N}\stackrel{{\scriptstyle[}}{{i}}=0]{N-1}{\sum}\mathbf{Q_{i}^{\delta D}}$, $\mathbf{\frac{1}{N}\stackrel{{\scriptstyle[}}{{i}}=0]{N-1}{\sum}\mathbf{Q_{i}^{\delta D}}=\left(\mathbf{S_{i}}-\delta\mathbf{D}\right)}\frac{1}{N}\stackrel{{\scriptstyle[}}{{i}}=0]{N-1}{\sum}\mathbf{S_{i}}\sqrt{\left|\mathbf{N_{i}}\right|}$ (the reference contour $\mathbf{\mathbf{R}=S_{0}}$, $\mathbf{N}=\mathbf{S}_{u}\times\mathbf{S}_{v}$, $\mathbf{N_{i}}=\mathbf{S_{i}}_{u}\times\mathbf{S_{i}}_{v}$). Notation $\delta\mathbf{D}$ indicates that the new origin is to be calculated wrt fixed set of optimally (precedingly) determined reparameterization. The minimization wrt the position of the new origin leads to simple extreme value problem:

[TABLE]

the solution is

[TABLE]

The optimal set of diffeomorphism and the position of the proper centroid are interrelated: the new centroid position involves a new set of optimal reparameterization functions $\mu_{i}\left(u,v\right),\nu_{i}\left(u,v\right)$, $i=1...N-1$, which in turn moves the optimal centroid further away from its ad hoc initial position. Therefore the double optimization problem (11) can be solved iteratively.

4.1 Removing translations

The RPSN representation is intentionally designed to retain the relative displacement information between the constituent surfaces, nevertheless it could be used for shape analysis purpose by appropriately removing the relative displacements between the constituents. Assume the double optimization problem (11) is solved . Then the minimization problem:

[TABLE]

can be used to remove the displacement updating the original positions of the constituents $\mathbf{S_{i}}$ to$\mathbf{S_{i}}-\delta\mathbf{D_{i}}$. From the extreme value problems:

[TABLE]

and the solutions are

[TABLE]

(note the denominator is the surface area of the $i$-th constituent surface.

As the optimal set of diffeomorphism and the position of the proper centroid are interrelated, so the displacements calculated by (16). Iterative solution is possible. Since the quantities that need to be calculated for the proper centroid (13) and the displacement removal equations (16) are largely overlapped, $\delta\mathbf{D}$ and all $\delta\mathbf{D_{i}}$, $\mathbf{i}=0...N-1$ should expediently be calculated in the same step. In this case, however the problem of the convergence remains an open question.

5 Discussion

In this paper an elastic surface mean determination method was presented. The mean surface is calculated from a set of surfaces in a way that all visible information (relative placement, rotation, scale) are retained. At the same time - borrowing the idea from the state of the art shape analysis methods - the parameterization of the surfaces is relaxed. The chosen contour representation (RPSN) and the imposed $\mathbb{L}^{2}$ metric forms a Hilbert space of the contour representations. The metric is chosen to be invariant wrt the reparameterization, the distance function based on it has well defined meaning, the (sum of) the second moment of the surfaces. The mean surface calculation is performed in the quotient space space of surfaces modulo reparameterization group and could be formulated as a double optimization problem: a variational for the system of the optimal parameterization and an extreme value problem for the proper centroid identification.

The work is the direct extension of the 2D contour mean calculation using the representation RPSV (Rescaled Position by Square Velocity) inroduced in the paper [5].

Appendices

In the appendices, the optimal reparameterization for the SRNF (Appendix A) and RPSN (Appendix B) are derived. The common notations and basic formulae used in these appendices are introduced below.

The notations wrt the original parameterization $\left(u,v\right)$ are as follows. The covariant basis vectors of the surface $\mathbf{S}=\mathbf{S}\left(u,v\right)$ are denoted by $\mathbf{S}_{u}=\frac{\partial\mathbf{S}}{\partial u}$, $\mathbf{S}_{v}=\frac{\partial\mathbf{S}}{\partial v}$. The contravariant basis vectors $\mathbf{S}^{u}$, $\mathbf{S}^{v}$ are defined such that their dot (scalar) products with the coordinate basis vectors are $\mathbf{S}^{i}\cdot\mathbf{S}_{k}=\delta_{k}^{i}$, $i,k\in\left\{u,v\right\}$. The normal vector of the surface is denoted by $\mathbf{N}=\mathbf{S}_{u}\times\mathbf{S}_{v}$. $\left|\mathbf{N}\right|$ stands for its length, the unit (hence paremeterization-independent) normal vector is denoted by $\mathbf{n}$ $\left(\mathbf{n}=\frac{\mathbf{N}}{\left|\mathbf{N}\right|}\right)$. The direct (dyadic) product of two vectors $\mathbf{u}$, $\mathbf{v}$ is defined such that its scalar products (contractions) with the vector $\mathbf{w}$ become $\left(\mathbf{uv}\right)\cdot\mathbf{w}=\left(\mathbf{v}\cdot\mathbf{w}\right)\mathbf{u}$ and $\mathbf{w}\cdot\left(\mathbf{uv}\right)=\left(\mathbf{u}\cdot\mathbf{w}\right)\mathbf{v}$. Metric and inverse metric tensor components are defined as $g_{ik}=\mathbf{S}_{i}\cdot\mathbf{S}_{k}$ and $g^{ik}=\mathbf{S}^{i}\cdot\mathbf{S}^{k}$, $i,k\in\left\{u,v\right\}$ respectively. The determinant of the metric tensor is denoted by $g$ $\left(g=\det\left[g_{ik}\right]\right)$, its square root (the ’metric’) $\sqrt{g}=\left|\mathbf{N}\right|$ is used to define surface area element as $dS=\sqrt{g}dudv$. Christoffel symbols (connection components) for embedded surfaces can be defined as

[TABLE]

where vectors $\mathbf{S}_{ik}$ are the second partial derivative of the position vector $\mathbf{S}$ . It can be seen by simple substitution that

[TABLE]

Using (17), (18), the partial derivatives of the logarithm of the metric $\left|\mathbf{N}\right|$ become the Christoffel ’divergences’:

[TABLE]

From identity $\mathbf{n}\cdot\mathbf{S}_{k}\equiv 0$

[TABLE]

Any vector $\mathbf{w}$ can be decomposed in the local basis $\mathbf{S}_{u},\,\mathbf{S}_{v},\,\mathbf{n}$ as $\mathbf{w}=\left(\mathbf{w}\cdot\mathbf{S}^{u}\right)\mathbf{S}_{u}+\left(\mathbf{w}\cdot\mathbf{S}^{v}\right)\mathbf{S}_{v}+\left(\mathbf{w}\cdot\mathbf{n}\right)\mathbf{n}$ or alternativeli in the contravariant basis as $\mathbf{w}=\left(\mathbf{w}\cdot\mathbf{S}_{u}\right)\mathbf{S}^{u}+\left(\mathbf{w}\cdot\mathbf{S}_{v}\right)\mathbf{S}^{v}+\left(\mathbf{w}\cdot\mathbf{n}\right)\mathbf{n}$, it follows that the decomposition of the identity tensor ($\mathbf{I}:\,\mathbf{I}\cdot\mathbf{w}=\mathbf{w}\cdot\mathbf{I}\equiv\mathbf{w}$) can be given in two ways:

[TABLE]

Note that the partial derivatives of the unit normal vector $\mathbf{n}_{u}$ and $\mathbf{n}_{v}$ are the elements of the tangent space hence can be decomposed such that $\mathbf{n}_{i}=\left(\mathbf{n}_{i}\cdot\mathbf{S}^{u}\right)\mathbf{S}_{u}+\left(\mathbf{n}_{i}\cdot\mathbf{S}^{v}\right)\mathbf{S}_{v}$.

The partial derivatives of the SRNF representation $\mathbf{P}=\mathbf{n}\sqrt{\left|\mathbf{S}_{u}\times\mathbf{S}_{v}\right|}$:

[TABLE]

similarly

[TABLE]

The partial derivatives of the RPSN representation $\mathbf{Q}=\mathbf{S}\sqrt{\left|\mathbf{S}_{u}\times\mathbf{S}_{v}\right|}$:

[TABLE]

similarly

[TABLE]

The transformation of the normal $\mathbf{N}=\mathbf{S}_{u}\times\mathbf{S}_{v}\rightarrow\mathbf{N}_{\left(\mu,\nu\right)}=\mathbf{S}_{\mu}\times\mathbf{S}_{\nu}$ can be expressed as

[TABLE]

where

[TABLE]

is the determinant of the Jacobian

[TABLE]

of the reparameterization. The reparameterization $\left(\mu\left(u,v\right),\nu\left(u,v\right)\right)$ considered feasible iff its determinant (27) is not negative for any values of $\left(u,v\right)$ and zero only in isolated points. The inverse of the Jacobian is:

[TABLE]

Using (Appendices) and (Appendices), the transformation of the Christoffel divergences become:

[TABLE]

[TABLE]

Note that one can apply the transformation rule of the Christoffel symbols directly to deduce results (Appendices), (Appendices).

Appendix A

The optimal parameterization for the SRNF problem between surfaces $\mathbf{R}\left(u,v\right)$ (reference) and $\mathbf{S}\left(u,v\right)=\mathbf{S}\left(\mu\left(u,v\right),\nu\left(u,v\right)\right)$ can be formulated as

[TABLE]

$\left(\left|\mathbf{S}_{u}\times\mathbf{S}_{v}\right|=\left|J_{\left(\mu,\nu\right)}\right|\left|\mathbf{S}_{\mu}\times\mathbf{S}_{\nu}\right|,\>\left|J_{\left(\mu,\nu\right)}\right|=\frac{\partial\mu}{\partial u}\frac{\partial\nu}{\partial v}-\frac{\partial\mu}{\partial v}\frac{\partial\nu}{\partial u}\right)$. Since the square of the terms in the parentheses $\iint\left|\mathbf{R}_{u}\times\mathbf{R}_{v}\right|dudv$, $\iint\left|\mathbf{S}_{u}\times\mathbf{S}_{v}\right|dudv$ express the surface area of the surfaces $\mathbf{R}$ and $\mathbf{S}$, problem (38) is equivalent with

[TABLE]

Its Lagrangian is

[TABLE]

$\left(\mathbf{M}=\mathbf{R}_{u}\times\mathbf{R}_{v},\>\mathbf{m}=\frac{\mathbf{M}}{\left|\mathbf{M}\right|}\right)$. From the latter

[TABLE]

where notation $\left|\mathbf{N}\right|$ is exclusively used to denote the normal vector of the surface wrt parameterization $\left(u,v\right)$, i.e. $\mathbf{N}=\mathbf{S}_{u}\times\mathbf{S}_{v}=\left|J_{\left(\mu,\nu\right)}\right|\mathbf{S}_{\mu}\times\mathbf{S}_{\nu}$ (see Eq. (Appendices)). Using the identity (18) and the definitions of the Christoffel symbols (17), (41) can be rearranged as:

[TABLE]

and similarly

[TABLE]

where the

[TABLE]

are the Christoffel divergences , i.e. the relative change of the metric $\left|\mathbf{N}\right|$ wrt the variables $\mu,\,\nu$. Using (35) (36), equations (42) expressed with the original coordinates $\left(u,v\right)$ can be rearranged as

[TABLE]

and simiarly:

[TABLE]

Calculation of the further terms:

[TABLE]

The explicit dependencies from $\left(\mu,\nu\right)$ can be removed using (26):

[TABLE]

Getting rid of the explicit independency from $\left(\mu,\nu\right)$ in (49) greatly simplify the following calculations (using (Appendices) and (23)):

[TABLE]

hence the ’$\mu$ component’ of the Euler-Lagrange equation becomes

[TABLE]

With similar calculation, the ’$\nu$ component’ of the Euler-Lagrange equation is

[TABLE]

Summarizing these components, the Euler-Lagrange equation associated with the reparameterization problem (38) using the inverse of the Jacobian of the reparameterization (35) is

[TABLE]

Appendix B

The optimal parameterization for the RPSN problem between surfaces $\mathbf{R}\left(u,v\right)$ (reference) and $\mathbf{S}\left(u,v\right)=\mathbf{S}\left(\mu\left(u,v\right),\nu\left(u,v\right)\right)$ can be formulated as

[TABLE]

$\left(\left|\mathbf{S}_{u}\times\mathbf{S}_{v}\right|=\left|J_{\left(\mu,\nu\right)}\right|\left|\mathbf{S}_{\mu}\times\mathbf{S}_{\nu}\right|,\>\left|J_{\left(\mu,\nu\right)}\right|=\frac{\partial\mu}{\partial u}\frac{\partial\nu}{\partial v}-\frac{\partial\mu}{\partial v}\frac{\partial\nu}{\partial u}\right)$. Similarly to the SRNF case, the minimization of (56) is equivalent with

[TABLE]

where the Lagrangian is

[TABLE]

$\left(\mathbf{M}=\mathbf{R}_{u}\times\mathbf{R}_{v},\>\mathbf{m}=\frac{\mathbf{M}}{\left|\mathbf{M}\right|}\right)$. From the latter

[TABLE]

$\left(\mathbf{N}=\mathbf{S}_{u}\times\mathbf{S}_{v},\>\mathbf{n}=\frac{\mathbf{N}}{\left|\mathbf{N}\right|}\right)$ and similarly

[TABLE]

where Christoffel divergences $\Gamma_{\mu}^{\mathbf{S}}$, $\Gamma_{\nu}^{\mathbf{S}}$ are defined in (44). Using (35) (36), equations (59) expressed with the original parameters $\left(u,v\right)$ can be rearranged as:

[TABLE]

Calculation of the further terms with the removal of the explicit $\left(\mu,\nu\right)$ dependencies:

[TABLE]

Then using (Appendices) and (25)):

[TABLE]

hence the ’$\mu$ component’ of the Euler-Lagrange equation becomes

[TABLE]

With similar calculation, the ’$\nu$ component’ of the Euler-Lagrange equation is

[TABLE]

Summarizing these components, the Euler-Lagrange equation associated with the reparameterization problem (38) using the inverse of the Jacobian of the reparameterization (35) is

[TABLE]