Embedding shapes with Green’s functions for global shape matching

ARTICLE IN PRESS

JID: CAG

[m5G;August 10, 2017;5:38]

Computers & Graphics xxx (2017) xxx–xxx

Contents lists available at ScienceDirect

Computers & Graphics journal homepage: www.elsevier.com/locate/cag

Technical Section

Embedding shapes with Green’s functions for global shape matching Oliver Burghard∗, Alexander Dieckmann, Reinhard Klein

Q1

Universität Bonn, Friedrich-Ebert-Allee 144, Bonn 53113, Germany

a r t i c l e

i n f o

Article history: Received 13 December 2016 Revised 14 June 2017 Accepted 21 June 2017 Available online xxx Keywords: Shape understanding Non-isometric shape matching Functional maps

a b s t r a c t We present a novel approach for the calculation of dense correspondences between non-isometric shapes. Our work builds on the well known functional map framework and investigates a novel embedding for the alignment of shapes. We therefore identify points with their Green’s functions of the Laplace–Beltrami operator, and hence, embed shapes into their own function space. In our embedding the L2 distances are known as the biharmonic distances, so that our embedding preserves the intrinsic distances on the shape. In the novel embedding each point-to-point map between two shapes becomes and can be represented as an affine map. Functional constraints and novel conformal constraints can be used to guide the matching process.

1

1. Introduction

2

Finding correspondences between two or more shapes is an important sub-task for a variety of applications, in which information has to be transferred or correlated between shapes. For example, local deformations can be transferred for shape editing [1–3], and correlations between corresponding regions can be exploited to compress dynamic meshes [4,5] and to create generative shape models [6–8]. Finding correspondences is especially interesting between nonisometric shapes, see Fig. 1 for some examples. Many previous approaches tailored to register isometric shapes fail in this case. Extrinsic non-rigid ICP [7,9] and variants [1,10–12] suffer from unreliable correspondences on extrinsic distances and from difficulties in solving the non-linear deformation models. The Blended Intrinsic Maps method [13,14] replaces the extrinsic metric by an intrinsic one and delivers good registration results by assuming the resulting maps to be locally conformal. A problem of BIM is that it is not clear how to incorporate a priori constraints which might be necessary to guide the method to the correct map out of the multiple reasonable ones (Fig. 2). Furthermore, the stitching of local maps leads to inconsistencies at their boundaries. Another group of approaches embeds shapes into a high dimensional space, where L2 distances approximate intrinsic distances. Although most of them allow to incorporate additional constraints, many share the major drawback that their embedding requires a non-linear alignment.

3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25

∗

Corresponding author. E-mail address: [email protected] (O. Burghard).

© 2017 Elsevier Ltd. All rights reserved.

Functional maps [15] overcome this problem by constructing an embedding, in which shapes can be aligned with a linear deformation. Unfortunately L2 distances of delta-distributions, that are typically used to embed points, only approximate intrinsic distances between intrinsically close points, see Fig. 4. This is especially important when only a few implicit constraints are available, such as when matching non-isometric shapes. In contrast to functional maps [15] we identify points with their Green’s functions of the Laplace–Beltrami operator. In this embedding the L2 distances are the well known biharmonic distances [16], which are an intrinsic distance metric on the shape. They are invariant to isometric shape deformations so that pose deformations have little influence on the matching process. We calculate correspondences by aligning these embeddings with an affine deformation, which can be computed reliably and efficiently. There is a linear relation between the Green’s alignment and the (pullback) functional map [15], so that we can incorporate functional constraints and operator commutativity into our setting. Last but not least, we can include additional constraints on the alignment, which require the resulting map to be close to conformal. The main contributions of our paper are (a) a novel embedding of shapes in the functional map framework by identifying points with their Green’s functions, (b) combining our novel embedding with functional constraints and (c) including conformality constraints into functional shape matching. The paper is organized as follows: Section 2 discussed the related work. Section 3 introduces the alignment of shapes with Green’s functions and relates it to the functional maps framework. We further motivate the novel embedding by a comparing Green’s functions and delta-distributions in Section 4. Then Section 5 shows how to utilize functional constraints and

http://dx.doi.org/10.1016/j.cag.2017.06.004 0097-8493/© 2017 Elsevier Ltd. All rights reserved.

Please cite this article as: O. Burghard et al., Embedding shapes with Green’s functions for global shape matching, Computers & Graphics (2017), http://dx.doi.org/10.1016/j.cag.2017.06.004

26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56

JID: CAG 2

ARTICLE IN PRESS

[m5G;August 10, 2017;5:38]

O. Burghard et al. / Computers & Graphics xxx (2017) xxx–xxx

Fig. 1. A variety of results gathered with our method. For each class (fourlegged, teddy, humans, birds) there is a single source (black contour) and a variety of target shapes. Using the sparse correspondences depicted by small spheres we calculate a dense map from the source to each target shape, which we then use to transfer a color field. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

2.1. ICP

Fig. 2. Failure case of Blended Intrinsic Maps, that are difficult to resolve without predefined constraints.

Fig. 3. (Left) Biharmonic distances db (p, x) from a fixed point p and (right) its Green’s function gp and three other Green’s functions.

Fig. 4. Biharmonic distances and L2 distances on delta-distributions for different numbers of eigenvectors.

57 58 59

conformality in the matching process. After discussing the discretization in Section 6 we describe a shape matching algorithm in Section 7, which is evaluated in Section 8.

60

2. Related work

61

Estimating correspondences between different shapes is a challenging task that has been addressed intensively in literature. In this section, we only provide a brief overview on directly related works and kindly refer the interested reader to the recent surveys [17,18].

62 63 64 65

Initially the problem of shape matching appeared in the context of registering sequential point cloud scans of a static scene. This led to the development of rigid ICP algorithms [19], which alternate between detecting corresponding points and rigidly aligning shapes. Due to the local optimization, these techniques depend strongly on the initial correspondences and on heuristics to prune novel correspondences. A variety of methods extend the original ICP metaphor to match deformed shapes by allowing non-rigid deformations in R3 for the alignment [1,7,9–12], A common shortcoming of these methods is the detection of corresponding points based on extrinsic instead of intrinsic distances. For deformable shapes extrinsic distances can be small even for intrinsically distant points. As a consequence these methods typically require a large number of point-to-point constraints to begin with and utilize sophisticated heuristics to prune novel correspondences. The Blended Intrinsic Maps method [13] obtains good results by concatenating and blending multiple conformal maps into a single global map. However, it cannot incorporate user constraints, which are sometimes necessary to solve ambiguities (e.g. Fig. 2). Furthermore, at the boundaries of the local maps the results often exhibit discontinuities. Additionally the method is difficult to generalize to point-clouds or shapes of genus other than zero. Other methods map shapes by parameterizing them on a common domain and then aligning their parameterizations so that either an intrinsic measure of stretch from source to target becomes minimal [21–25] or so that the integrated stretch along a sequence of deformations from the source onto the target [26–28] becomes minimal. These methods deliver continuous maps of high quality, but are often computational demanding and their application on non-simple topologies is non-trivial. Yet another class of methods uses an ICP-like alignment after embedding shapes into a high dimensional space where L2 distances approximate intrinsic ones. Shapes have been embedded with the eigenvectors of an affinity matrix [29,30], with an embedding approximating geodesic distances [31], with an embedding based on electrostatic repulsion [32] and with delta-functions [33]. All of these methods use non-linear maps to align the embeddings. Slightly different are the methods [34,35], where the alignment of shapes is avoided by directly embedding one shape into the other by minimizing a non-linear functional.

67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106

2.2. Functional maps

107

The remarkably successful functional maps framework was introduced in [15]. To the best of our knowledge this paper was the first to fully exploit the fact that a linear alignment of a

108

Please cite this article as: O. Burghard et al., Embedding shapes with Green’s functions for global shape matching, Computers & Graphics (2017), http://dx.doi.org/10.1016/j.cag.2017.06.004

66

109 110

Q2

ARTICLE IN PRESS

JID: CAG

[m5G;August 10, 2017;5:38]

O. Burghard et al. / Computers & Graphics xxx (2017) xxx–xxx

111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155

functional embedding is sufficient to represent arbitrary non-linear alignments in R3 . This significantly simplifies the matching process. Additionally, the authors demonstrate the usefulness of functional constraints, such as matching labeled regions. Particularly for isometric shapes, where many functional constraints are available [36,37] and the alignment is a rotation commuting with the Laplacian, superior results have been achieved. Furthermore the resulting maps can be optimized with an ICP-like alignment algorithm after embedding the shapes with delta-distributions. Typically when matching non-isometric shapes there are few a priori constraints available and the ICP-like alignment becomes especially important. In this case a drawback of their embedding emerges, namely the L2 distances on delta-distributions, which are a critical ingredient for a ICP-like method, are not well-defined. As we describe in Section 4, the distances are only defined after projecting the delta-distributions onto the first eigenvectors of the Laplace operator and the distance metric depends strongly on the number of eigenvectors that were used. The more eigenvectors are used the more localized these distances become. In our work we therefore provide a novel embedding for shape matching that respects intrinsic distances. Functional maps initiated a series of publications, such as improving the extraction of correspondences [33,38,39], utilizing lowrank assumptions on the functional map [40,41] (which still requires multiple initial constraints typically not available for nonisometric shapes) or investigating the matching of shape collections [38,42]. Previous work represented points as Green’s functions [16] and as the related Global Point Signatures [43]. The invariance of the Dirichlet energy to conformal deformations was prominently explored in [44]. Yet to the best of our knowledge, we are the first to compute correspondences by aligning shapes represented with Green’s functions and the first to use functional correspondences to approximate the functional map representation of a conformal map. 3. Embedding with Green’s functions We start our exposition by describing the embedding of a shape M using Green’s functions. Let L2 (M ) = { f : M → R |  2 M f (x ) dx < ∞} be the set of square integrable, real-valued functions on M and let δ p be the delta-distribution at a point p, i.e. δ p , f  = f ( p) ∀ f ∈ L2 (M ). Let further  : L2 (M ) → L2 (M ) be the Laplace–Beltrami operator and its spectral decomposition have the eigenvectors φi ∈ L2 (M ) and the eigenvalues λi ∈ R (λ1 ≤ λ2 ≤ . . . ). Then the Green’s function gp of the Laplace operator  at point p is the solution of the equation:

g p = δ p 156 157 158 159 160 161 162 163 164

(1)

This equation has a solution if and only if δ p is orthogonal to the null-space of . For a simple exposition, we assume that 0 = λ1 < λ2 = 0, which is the case for a compact, simply-connected shape. Hence the null-space of  is the subspace of constant functions, which is spanned by the first eigenvector φ 1 . We further write  for the orthogonal projection on the complement of this nullspace, i.e. ( f ) := f − φ1 φ1 , f . Therefore in our context we define the Green’s function of the Laplace operator of point p ∈ M as the solution to the slightly different equation

g p = (δ p ) = δ p − φ1 φ1 , δ p  165 166

(2)

which now always has a solution that can be written with the pseudo-inverse + as

g p = + δ p =

∞  i=2

167

 g p , φ1  = 0 ,

φi

φi ( p ) . λi

We state a few properties, for the upcoming discussion:

(3)

3

Proposition 1.

168

(a) The mapping from a surface M onto its Green’s functions gM p is injective, thus an embedding. 2 (b) The functions {gM p + φ1 | p ∈ M} form a basis of L (M ). (c) The functions {gM | p ∈ M } are linearly dependent and span a p subspace of L2 (M ) of co-dimension 1(see the Appendix). Several Green’s functions can be seen in Fig. 3. The further points are located on the surface, the more their Green’s functions differ. L2 distances on the Green’s functions are called biharmonic distances [16] and the left of Fig. 3 shows the distance fields of several points. 3.1. Matching shapes





∀p ∈ M

G



 φ1M = 0

G˜

gM p

+φ

M 1



=

gN T ( p)

+φ

N 1

∀ p ∈ M.

one can infer T. Yet not every map G is induced by some point-wise map T. If G is induced by some point-wise map T, then it aligns the Green’s functions:

  ∀ p ∈ M ∃q p ∈ N : G gM = gN p qp

174 175 176 177 178

180 181 182 183 184 185 186 187 188 189 190

191 192 193 194 195 196 197 198 199

200 201 202 203 204 205 206

(6)

On the other hand, if G fulfills the last equation, we can reconstruct T via T : M → N , p → q p . The above results are summarized in the following proposition: Proposition 2. For each point-wise map T, Eq. (4) defines a unique affine map G, that maps the Green’s functions of M onto the Green’s functions on N . Such an alignment of Green’s functions supports Eq. (6), which can be used to restore T from G. Therefore there is an oneto-one relation of point-wise maps T and the Green’s alignment G supporting Eq. (6). 3.2. Pullback functional maps

207 208 209 210 211 212 213 214 215 216

The affine map G is related to the functional maps introduced in [15], which we shortly introduce. Each point-wise map T : M → N induces a pullback functional map F, that pulls function values from N onto M:

f → f ◦ T

(7)

Please cite this article as: O. Burghard et al., Embedding shapes with Green’s functions for global shape matching, Computers & Graphics (2017), http://dx.doi.org/10.1016/j.cag.2017.06.004

173

(5)

  Using G˜ one can write G( f ) = B( f ) + t as t = G˜ φ1M − φ1N and B( f ) = G˜ ( f ) −  f, φ1M (φ1 + 2t ). G is unique, because it can be used to write G˜ using G˜ ( f ) = B( f ) +  f, φ1M (φ1 + 2t ) and from G˜

F : L2 ( N ) → L2 ( M )

171 172

(4)

G is well-defined as due to Proposition 1a different points have different Green’s functions. Actually there is one unique affine map G which fulfills Eq. (4). It can be written as G( f ) = B( f ) + t where B : L2 (M ) → L2 (N ) is a linear map (i.e. B( f + λg) = B( f ) + λB(g) ∀λ ∈ R∀ f, g ∈ L2 (M )) and t ∈ L2 (N ). Choosing G as an affine map simplifies solving for G and allows the inclusion of least squares constraints. To see why choosing G as an affine map is possible, note that due to Proposition 1b there is a unique linear map G˜ with



170

179

In the following we utilize the properties of Green’s functions for the construction of a mapping T : M → N between the shapes M and N . As Green’s functions define a distance field on the surface, we can recover a point from its Green’s function simply by determining the surface point whose Green’s function is most similar. Therefore, Green’s functions represent the intrinsic location of a point on the shape. We can therefore solve for a map T by solving for a map G : L2 (M ) → L2 (N ) which aligns the Green’s functions of M onto the Green’s functions of N and only later calculate T from G. For clarification we add superscripts to the involved quantities on M and N . Then such an alignment G has to fulfill:

G gM = gN p T ( p)

169

217 218 219 220

ARTICLE IN PRESS

JID: CAG 4

221 222 223 224 225 226

[m5G;August 10, 2017;5:38]

O. Burghard et al. / Computers & Graphics xxx (2017) xxx–xxx

Both F as well as G map functions between M and N and are thus “functional maps”. For the sake of a clear notation we refer to T as a point-wise map, to F as the corresponding pullback functional map and to G as the corresponding alignment of Green’s functions. Because F preserves function values at corresponding points, it aligns (dual) delta-distributions ( · T denotes transposed): N δ M p , F ( f )M = F ( f )( p) = f (T ( p)) = δT ( p) , f N

⇒ 227 228 229

N F T δM p = δT ( p)

∀p ∈ M

(8)

Note: The original publication on functional maps [15, Section 6.1] assumes that functional maps defined by Eq. (7) align delta-distributions: N δM p = F δT ( p) ∀ p ∈ M

230 231 232 233 234 235 236 237 238

This is true only for area-preserving map, which fulfills = FT . Interestingly the difference of Eqs. (8) and (9) seems to have had limited effect in previous works. One reason might be that for any functional map F defined by Eq. (7) there is another map F  = F −T that aligns delta-distributions. The difference of F and F depends on the change in the area-form and is often small. Proposition 3. The alignment of Green’s functions G can be written in terms of the pullback functional map F (AM is the area of M, see the Appendix) as T G ( f ) = + N F M ( f ) +

239

241 242 243 244 245 246 247 248 249

1



+N F T φ1M

AM

+ T MB

N +





AM φ

M T 1 t

N +

AM M N ,T φ φ . AN 1 1

251 252 253 254 255 256

(11)

Each constraint on F has an equivalent constraint on G and vice versa. The final correspondences do not depend on whether one optimizes for F or G. They depend on the involved constraints to solve for F and G and on the method to extract correspondences from F and G. The point-wise maps T induce only a fraction of all possible functional maps. Not only does any induced functional map align dual delta-distributions (Eq. (8)), but also any functional map F which does so according to Eq. (12) is induced by a point-wise map T. Thus, the alignment of dual delta-functions by F and the alignment of Green’s functions by G are equivalent conditions. N ∀ p ∈ M ∃q p ∈ N : F T δ M p = δq p

250

(10)

and F can be written in terms of G( f ) = B( f ) + t as

F = 240

(9) F −1

(12)

Proposition 4. Let F be a pullback functional map and G an alignment of Green’s functions which are related by Eq. (10) or (11), then the following statements are equivalent: (i) G aligns Green’s functions by Eq. (6). (ii) F aligns dual delta-distributions by Eq. (12). (iii) There is a point-wise map T inducing F by Eq. (7) or G by Eq. (4).

257

4. Green’s functions vs. delta-distributions

258

Next we discuss the differences between our novel embedding and the original embedding with delta-distributions and why it matters for non-isometric shape matching. Aligning shapes with an ICP-like algorithm is a minimization of the non-linear ICP energy functional:

259 260 261 262



E[G; T ] = 263 264 265

M

 M  N 2 G g p − g T ( p ) d p 2

Fig. 5. Matching two shapes from the “fourlegged” dataset from point-to-point constraints using either delta-distributions or Green’s functions.

2 delta-distributions are not square integrable (δ M p ∈ L (M )) and their L2 distances are therefore not well-defined. To replace the Green’s functions in Eq. (13) with delta-distributions they first have to be projected onto a finite subspace. For example, the authors of [15, Section 6.1] embed the shapes with delta-distributions projected onto the first k eigenvectors of the Laplace operator:

p ∈ M → δ p(k ) =

k 

φi φi ( x ) ∈ L 2 ( M )

267 268 269 270 271

(14)

i=1

Unfortunately, the emerging distances approximate intrinsic distances only for a small number of eigenvectors (e.g. k < 10), while for large k the L2 distances do no more approximate intrinsic distances, but only discriminate points:

k (k )  k→∞ 0 δ p − δq(k) 2 = (φi (x ) − φi (y ) )2 −−−→ i=1

x=y ∞ else

272 273 274 275

(15)

Eventually also the Green’s functions have to be approximated by a finite basis. In contrast to L2 distances on delta-distributions the biharmonic distances are well approximated with the first few eigenvectors of the Laplace operator. Approximations of both distance fields with different numbers of eigenvectors are shown in Fig. 4. For 10 0 0 eigenvectors delta-distributions only discriminate the query point in agreement with Eq. (15). The effect of the different embeddings on the alignment process is shown in Fig. 5, where two shapes are matched with an ICP-like method that will be described in Section 7 (no conformality, α = 0). Using a small number of eigenvectors (k = 5) the results using delta-distributions and Green’s functions are similar, while for 20 eigenvectors the alignment with delta-distributions becomes discontinuous, even after iterating. This is especially the case, when only few a priori correspondences are known, which is typical for non-isometric matching. 5. Functional constraints and conformality

276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291

292

(13)

This functional depends on the L2 distances in the embedding. For Green’s functions these distances are known as biharmonic distances—a well-defined intrinsic distance metric [16]. In contrast,

Next we utilize the linear relation between F and G to transfer functional constraints and operator commutativity into our setting. Afterwards we propose novel functional constraints for conformal maps.

Please cite this article as: O. Burghard et al., Embedding shapes with Green’s functions for global shape matching, Computers & Graphics (2017), http://dx.doi.org/10.1016/j.cag.2017.06.004

266

293 294 295 296

ARTICLE IN PRESS

JID: CAG

[m5G;August 10, 2017;5:38]

O. Burghard et al. / Computers & Graphics xxx (2017) xxx–xxx

297

5.1. Functional constraints

298

Instead of localizing correspondences at single points, it is often more appropriate to determine corresponding regions, or equivalently to require indicator functions to match. This is an instance of a functional constraint, where the map T is known to pull a function fN ∈ L2 (N ) back onto another function hM ∈ L2 (M ):

299 300 301 302

F f N = hM

5

(16)

306

Functional constraints emerge in other applications as well. For example isometric shapes have the same heat- and wave-kernelsignatures [36,37]. As functional constraints are linear constraints in F, they can be written as a linear constraints in G using Eq. (11).

Fig. 6. Smoothed versions of the Stanford bunny colored by the Green’s function of a vertex at the nose. To illustrate the preservation of Green’s functions, the Green’s functions from the deformed meshes were mapped onto the original mesh in the lower row.

307

5.2. Operator commutativity

−T  ( f ) + √1 + F −T φ N is the inverted where G+ ( f ) = + N MF M 1

342

Green’s alignment.

343

303 304 305

308 309 310 311 312 313 314 315

If we have a functional operator on each shape and these operators have the same effect on equivalent functions, then these operators commute with the functional map F. In this case applying the first operator followed by a projection onto the other shape has the same effect as first projecting onto the other shape followed by applying the second operator there. A typical example for operator commutativity are intrinsic symmetries or the Laplace–Beltrami operator for isometric shapes where we have:

M F = F  N 317 318

5.3. Conformal maps

319

Additionally we introduce novel constraints on F and G by assuming that T approximates a conformal map. A conformal map is a locally angle preserving, continuous map. The class of conformal maps includes the isometric maps and is general enough to match typical shapes of equivalent topology. At the same time a conformal map is already determined by a few known correspondences [14]. Furthermore, conformal maps were successfully used in previous work to match near-isometric shapes [13,14]. Conformality does not only restrict the map T, but also the functional map F and the Green’s alignment G. For two functions f, h ∈ L2 (N ) the Dirichlet energy is defined by

320 321 322 323 324 325 326 327 328 329

ED [ f, h] = ∇N f, ∇N hN = N f, hN 330 331 332

and measures how much their gradients agree. A point-wise map is conformal if and only if its functional representation preserves the Dirichlet energy [44]:

M F ( f ), F (h )M = N f, hN 333

or equivalently

F T M = N F −1 . 334 335 336 337

Proposition 5. Let T : M → N be a conformal, point-wise map with the corresponding functional map F and the corresponding Green’s alignment G, then F adheres to Eq. (18). If F fulfills a functional constraint F fN = hM , then it also fulfills Eq. (19) and if G fulfills a pointto-point correspondence G gM = gN p q , then it also fulfills Eq. (20). Area preserving maps (F T = F −1 ) have similar constraints, whose investigation we leave for future work.

344 345 346 347 348 349 350

(17)

Using Eq. (11) this linear constraint on F becomes a linear constraint on G.

316

AN

(18)

Eq. (18) is an instance of a non-linear Procrustes problem, /2 namely that 1M/2 F −1 is orthogonal [45,46]. We propose to use N the current functional and point-to-point constraints to transfer Eq. (18) into linear constraints instead. Each functional constraint F fN = hM is equivalent to fN = F −1 hM , which is combined with Eq. (18) into

5.4. Conformal maps and Green’s functions

351

It is interesting how Green’s functions change under a conformal map. Let T : M → N be a conformal map with the corresponding functional map F and the corresponding Green’s alignment G, then the Green’s function of a point p ∈ M and the Green’s function of the mapped point image T ( p) ∈ N differ by (see the Appendix):

352

1 N N N −1 M T M F −1 gM g p  + + p − g T ( p ) = φ1  φ 1 , F N F φ1  AM = φ1N

∞  i=2



F

−1



λM i

1i

+

1

∞ 

AM

i=2

 



φiM

(F )1i λNi

353 354 355 356 357

(21) 358

(22)

The values F1i and F −1 represent the mean of the mapped 1i eigenvectors on N , which typically are small if the area form changed little, in which case we expect little difference between the original and the mapped Green’s function. Fig. 6 shows the Stanford bunny and multiple increasingly smoothed versions of it. During the smoothing a conformal map between the meshes was obtained [47]. Projecting the Green’s functions from the smoothed meshes onto the original mesh shows that Green’s functions of a single point indeed changed little. Apart from the observation above our work does not exploit the connection of conformal maps and Green’s functions further, but leaves this interesting issue open for future work.

359 360 361 362 363 364 365 366 367 368 369 370

6. Discretization

371 372

340

and each point-to-point constraint G gM = gN p q is equivalent to

To represent functions and with them the Green’s alignment G a basis of the function spaces L2 (M ) and L2 (N ) is required. We use the first k eigenvectors of the Laplace operator [φ1 , φ2 , . . . , φk ], which is a common choice in literature [15]. Practical evaluation shows that biharmonic distances can be well approximated with as little as five eigenvectors (Fig. 4).

341

+ gN gM , which we combined with Eq. (18) into p =G T ( p)

6.1. Spectral decomposition

378

Next we discuss the calculation of the eigenvalues λi and eigenvectors φ i of a shape’s Laplace–Beltrami operator . Let functions

380

338 339

(+N F T M )( fM ) = N ( fN ), 





(19)







1

(+M BT N ) gNT ( p) + 

AN

+M F −T φ1N = gM p

(20)

Please cite this article as: O. Burghard et al., Embedding shapes with Green’s functions for global shape matching, Computers & Graphics (2017), http://dx.doi.org/10.1016/j.cag.2017.06.004

373 374 375 376 377

379

ARTICLE IN PRESS

JID: CAG 6

381 382 383 384 385 386 387

O. Burghard et al. / Computers & Graphics xxx (2017) xxx–xxx

on the surface be represented by vectors, whose coefficients are the function values at the mesh vertices or point-cloud points. Let W be the matrix representation of the scalar product on functions and L be the matrix representation of the Laplace–Beltrami operator, which is typically written as L = W −1C, where C is the Dirichlet energy on functions. The eigenvectors and eigenvalues of L are defined by the following generalized eigenvalue problem:

C φi = λiW φi so that 388 389 390 391

φiT W φ j = δi j and λ1 ≤ λ2 ≤ . . .

For triangle meshes we utilize the well-known Cotan Laplacian [48,49]. We estimate the vertex areas as a third of the sum of the adjacent triangle areas and set W to be a diagonal matrix of the estimated vertex areas. Furthermore

 1 (Cu )i = (cotαi j + cotβi j )(ui − u j ) 2

(23)

j∈N (i )

392 393 394 395 396 397 398

defines the Dirichlet energy, where N (i ) is the 1-neighborhood of the vertex i and α ij , β ij are the two angles opposing the edge from vertex i to vertex j. For point-clouds (e.g. Fig. 11a) we use a variant of the Laplace operator from Belkin and Liu [50,51]. They define L so that for a small time t heat-diffusion matches the euclidean one. With the point positions pi ∈ R3 the Dirichlet energy is:

 pi − p j 2 1  (Cu )i = exp − 2 4t 4π t



( ui − u j )

411

W is again the diagonal matrix of the estimated point areas. We estimate the area of a point as π r2 /3, where r is the average euclidean distance to its six nearest points. This simple heuristic worked well in our experiments and allows us to avoid the more complicated area estimations of Belkin and Liu [50,51]. We choose t so that exp(−d¯2 /4/t ) = 1/10, where d¯ is the average euclidean distance of all the points to their 10 closest neighbors. Next we sparsify C by removing small elements as follows. First we mark all coefficients of C larger than 1/(10 · 4π t2 ). Then in each row we mark the 10 largest coefficients (excluding the diagonal). For symmetry we mark Cij if Cji is marked. Lastly we set all unmarked coefficients to 0 and update the diagonal entries of C so that C (1, . . . , 1 )T = (0, . . . , 0 )T .

412

7. Alignment algorithm

413 415

Next we propose a concrete method to calculate a Green’s alignment of two shapes and therefore a point-wise map T from a few known point-to-point or functional constraints.

416

7.1. Alignment energy

399 401 402 403 404 405 406 407 408 409 410

Fig. 7. Correspondences before and after iterating (20 eigenvectors, α = 0).

7.2. Initial solving

427

Typically there are so few a priori point-to-point and functional constraints, that in the first iteration the energy is underconstrained. We therefore add further regularization constraints. Translation depends on the area scale and is typically small, so that we assume t = t˜ = 0. Furthermore the distortion of G should be as little as possible, so that we add the regularizer B2F ( = 10−6 ). In our experiments this simple choice lead to consistently good results and clearly outperformed other possible terms, such as F 22 or M F − F N 22 . In conclusion the initial alignment G for the a priori point-to-point C0 and functional constraints D0 is the unique minimum of:

428

ECp0 [B, 0, 0] + EDf 0 [B, 0, 0] + B2F ;

417 418 419 420 421 422 423

We describe the affine alignment with the variables B ∈ R(k−1 )×(k−1 ) ; t, t˜ ∈ Rk−1 and we define a quadratic energy, where we incorporate point-to-point and functional constraints in a least squares sense. The parameter α ∈ [0, 1] encodes whether we assume conformality. For the point-to-point correspondences C = (( p1 , q1 ), . . . , ( pn , qn )) ∈ (M × N )n the matching energy ECp [B; t; t˜] is

1 n 424 425

n  i=1

 N 2 N M 2 BgM ˜ pi + t − gqi 2 + α M F gqi + t − g pi 2

i=1

426

(25)

constraints D = (( f1 , h1 ), . . . , ( fm , hm )) ∈ and2 for the2 functional m f L (M ) × L (N ) the matching energy ED [B; t] is m  1  F hNi − fiM 2 + α +N F T M fiM − N hNi 2 2 2 m

where F, B and t are related by Eqs. (10) and (11).

431 432 433 434 435 436 437 438

(27) 439

Once we have an alignment G we map the Green’s functions of M onto the Green’s functions of N and determine novel point-topoint correspondences with a nearest neighbor search:

440

ρM :M → N , ρN :N → M,

(26)

N q → arg miny∈M B gM y + t − gq 2

It is sufficient to calculate correspondences from and onto a subset of both shapes SM and SN , which were initially calculated using farthest point sampling in the Green’s embedding, i.e. using biharmonic distances. For efficient nearest neighbor queries we use k-d trees with the “sliding mid-point rule” [52], which adapts to the low intrinsic dimensionality of the data. In summary, the current point-to-point correspondences C1 are inferred from G by:



p, ρM ( p)

441 442

N p → arg minx∈N B gM p + t − gx 2



 

| p ∈ SM ∪



ρN ( q ) , q | q ∈ S N



7.4. Iterative alignment

445 446 447 448 449

450

We refine the alignment by alternating between solving for the alignment G using the current point-to-point constraints and calculating novel point-to-point constraints C1 from the current alignment. From the second iteration on solving for G is over-constraint and no further regularization is required. G is then defined as the minimum of:

ECp0 [B, t, t˜] + EDf 0 [B, t] + ECp1 [B, t, t˜]

443 444

(28)

451 452 453 454 455 456

(29)

In our experiments we use 20 iterations, which was enough to converge. The effect of iterating is shown in Fig. 7 and the entire alignment procedure is shown in Algorithm 1.

458

8. Evaluation

460

An important application of our method is matching nonisometric shapes, which we evaluate on the Shrec dataset [53]. The dataset contains several shape classes, in which we create maps and evaluate their quality using the dense inter-class correspondences from [11] as ground truth.

Please cite this article as: O. Burghard et al., Embedding shapes with Green’s functions for global shape matching, Computers & Graphics (2017), http://dx.doi.org/10.1016/j.cag.2017.06.004

430

7.3. Updating correspondences

C1 = 414

t = t˜ = 0

429

(24)

j

400

[m5G;August 10, 2017;5:38]

457 459

461 462 463 464 465

ARTICLE IN PRESS

JID: CAG

[m5G;August 10, 2017;5:38]

O. Burghard et al. / Computers & Graphics xxx (2017) xxx–xxx

Algorithm 1 Matching two shapes. 1: 2: 3: 4: 5: 6:

7:

Input: shapes M, N ; initial point-to-point C0 and functional D0 constraints; α ∈ [0, 1] SM /SN ← Farthest_Point_Sampling(M/N , 200 ) B, t ← solve Eq. (27) for an initial alignment using C0 , D0 , α , = 10−6 for i = 1…20 do C1 ← Matches(B, t, SM , SN ) B, t ← solve Eq. (29) for a novel alignment using C0 , D0 , C1 , α return Matches(B, t, M, N )

 Sub routines: function Matches ( B, t, M, N):  Eq. (28)   N return ( p, q ) | p ∈ M, q = arg miny∈N B gM 9: p + t − gy 8:

2  N ∪ ( p, q ) | q ∈ N, p = arg minx∈M B gM x + t − gq 2 

function Farthest_Point_Sampling(M, n) S ← { p} with a random point p on M for i = 2 . . .n do  db biharmonic distances 12: S ← S ∪ arg maxx∈M miny∈S db (x, y ) 13: 10:

7

Table 1 Overview of constraints, embedding w.r.t. F, FT , G and GT . Constraints on F and GT

Functional constraints; embedding with primal delta-distributions

Constraints on FT and G

Conformal functional constraints; embedding with dual delta-distributions; embedding with Green’s functions

of-freedom of the optimization, as without these requirements any (continuous) map is a valid solution. The degrees-of-freedom are reduced by either decreasing the number of eigenvectors k or by increasing the weight α of the conformality constraints. The influence of the degrees-of-freedom on the optimization is shown in Fig. 8d. For three initial point-to-point constraints results improve when degrees-of-freedom are reduced, while for six initial pointto-point constraints the results improve when degrees-of-freedom are increased.

505 506 507 508 509 510 511 512 513

11:

14:

466

return S

8.1. Initial results

474

The first results in Figs. 1, 5 and 7 show that good results can be obtained from few sparse point-to-point correspondences. As discussed in Section 4 and shown in Fig. 5, the well-defined distances on the Green’s functions result in a better alignment than distances on delta-distributions. The novel embedding gives smoother maps already in the first iteration, which are further smoothed by the iterative ICP-like alignment. Figs. 5 and 7 show the effect of iterative alignment depicting maps without iterating and after 20 iterations.

475

8.2. Solving with point-to-point constraints

476

We compare our method to functional maps (FM), which aligns delta-distributions, and Blended Intrinsic Maps (BIM), which is a state-of-the-art method for automatic shape matching and does not use predefined constraints. More precisely, when matching with functional maps, we use Eq. (9) to iteratively align the deltadistributions and extract correspondences via nearest neighbor search. We omit the additional rigid ICP alignment described in [15] as it requires isometric shapes and in principle can be applied to our method as well. In each class of the Shrec dataset we build between 30 random maps. Then we build the point-to-point constraints by geodesic farthest point sampling on the source shape and mapping these points onto the target with the ground truth map. To evaluate the maps we equally distributed 200 points on the source, mapped them onto the target and measured the deviation from the ground-truth. Fig. 8a shows the results of the three methods with either three or six point-to-point correspondences and Fig. 8b gives details for three classes. Independent of the number of constraints our method outperforms the functional maps method. The number of eigenvectors used for the calculations influences the results. Fig. 8c shows a steady improvement as the number of eigenvectors is increased (six point-to-point constraints, α = 0). In principle maps are not restricted when represented with Green’s functions. For every map T there is a Green’s alignment G and for a finite set of constraints there is an infinite number of possible maps. Yet in our case there are two additional requirements for the solution, namely the alignment must be represented with k eigenvectors and for α > 0 the solution must fulfill the conformality conditions. These two requirements can be seen as the degrees-

467 468 469 470 471 472 473

477 478 479 480 481 482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 504

8.3. Solving with functional constraints

514

We labeled three regions on a representative shape of the classes ‘birds’ and ‘fourlegged’ and transfer these labels onto all shapes of the same class using the ground truth correspondences [11]. Each label consists of four colors, so that we can build four functional constraints from the indicator functions. Fig. 9 shows the labels on the bird class as well as the matching results. Here we intentionally have chosen an example where BIM fails to demonstrate the usefulness of predefined constraints. Fig. 10 depicts a quantitative evaluation of the matching process and further results are depicted in the additional material. Our method differs from the FM method in two aspects, namely in the novel embedding and in the addition of conformality constraints. We also considered a third method, that uses dual deltadistributions (to constrain and extract point-to-point correspondences) and conformal functional constraints. This third method differs from each of the other two formulations in only one aspect. The results show that both the novel embedding and the conformality constraints improve the results. While a variety of different combinations of embeddings and constraints are possible only few achieved good results in our experiments. For example, the embedding with dual-delta distributions (F T δ p = δq ) does not work well with functional constraints without conformality (F f = g). The reasons might be that the first is a constraint on F, while the second is a constraint on FT and mixing constraints on F and FT does not work well (see Table 1).

515

8.4. Point clouds and topological changes Due to its underlying simplicity, our method is rather general. It works with point clouds (Fig. 11a) and complicated topologies (Fig. 11b). Fig. 11c shows an example of matching in the presence of severe topological differences, where meshes were sewed together at self-intersections. Apparently our method can be seen as a global intrinsic alignment. Where the matching succeeded, most errors were localized around areas where the topology has changed.

517 518 519 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548

8.5. Isometric matching

549

The embedding with Green’s functions is most useful when there are few known constraints. If the functional constraints already determine the functional map, there is little difference between embedding with Green’s functions and delta-distributions. For example, when matching isometric shapes the heat kernel signatures and wave kernel signatures [36,37] provide a magnitude

550

Please cite this article as: O. Burghard et al., Embedding shapes with Green’s functions for global shape matching, Computers & Graphics (2017), http://dx.doi.org/10.1016/j.cag.2017.06.004

516

551 552 553 554 555

JID: CAG 8

ARTICLE IN PRESS

[m5G;August 10, 2017;5:38]

O. Burghard et al. / Computers & Graphics xxx (2017) xxx–xxx

Fig. 8. Quantitative comparison of approximating point-to-point correspondences on the Shrec dataset [53] using ground truth correspondences of [11] (20 iterations, no functional constraints, no conformality α = 0, sub-sampling 20 0 0 points).

Fig. 9. Labels on the source (b) and target (c) were used to map colors from the source (a) with the various methods (d, e). We use the four functional constraints depicted here and in more detail in the additional material. No point-to-point constraints were used (20 eigenvectors, with conformality α = 1, sub-sampling 20 0 0 points).

556 557 558 559 560 561

of functional constraints. Fig. 12 shows a quantitative evaluation of matching several classes from the TOSCA [54] dataset using heat kernel signatures and a single point-to-point constraint due to the intrinsic symmetry (see also additional material). Here the results of our novel embedding are indeed similar to the original functional maps. Note the results of both algorithms can be further

Fig. 10. Quantitative evaluation of matching with functional correspondences (parameters as in Fig. 9).

improved by applying rigid ICP on the delta-distributions as described in [15]. This was not done here as our focus is on nearisometries and as it improves both methods in the same way.

564

8.6. Blending shapes and limitations

565

A good demonstration of the quality of the correspondences is their utilization to linearly blend the source triangulation onto the target shape, which results a novel triangulation of the target

566

Please cite this article as: O. Burghard et al., Embedding shapes with Green’s functions for global shape matching, Computers & Graphics (2017), http://dx.doi.org/10.1016/j.cag.2017.06.004

562 563

567 568

JID: CAG

ARTICLE IN PRESS

[m5G;August 10, 2017;5:38]

O. Burghard et al. / Computers & Graphics xxx (2017) xxx–xxx

9

Fig. 11. Demonstrating the generality of the approach: matching point-clouds (a), shapes of higher genus (b) and shapes with severe differences in their topology and inner metric (c).

Fig. 13. Visualization of the correspondences by a linear interpolation of the source triangulation (λ = 0) onto the target shape (λ = 1). Please see the additional material for further results (20 iterations, 6 point-to-point constraints, no functional constraints, no conformality α = 0, sub-sampling 20 0 0 points).

Fig. 12. Quantitative results of matching isometric shapes with heat-kernelsignatures. As described in the text, the alignment of Green’s functions and deltadistributions yield similar results if there are sufficient constraints, as it is the case for isometric shapes when using heat-kernel-signatures (20 iterations, 20 eigenvectors, 1 point-to-point correspondence, 40 HKS functional constraints, conformality α = 1, sub-sampling 20 0 0 points).

585

shape. This new triangulation depends on the correspondences and shows their quality. Typically the following errors occur: (1) If the correspondence map is not surjective, i.e. not all target vertices have a corresponding source vertex, then these vertices are not contained in the novel triangulation at all (e.g. see missing vertices on the teddy ear). (2) If neighboring source vertices are mapped differently, then this leads to intersecting edges, outstanding long edges, and edges that are not located on the target shape. Fig. 13 shows such a blending using our results. While in principle our results are of good quality, there are at least two reasons for the occurring misalignments. One is that Green’s functions are very similar in proximity of thin extrusions (e.g. horse legs) and another is that only a local but not a global optimization was used (e.g. teddy arm). Our results can be further improved using a method for affine point cloud alignment such as Coherent Point Drift [55], which is also shown in Fig. 13 and in much more detail in the additional material.

586

9. Concluding remarks

569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584

587 588 589 590 591 592 593 594

We expanded on the understanding of the intrinsic alignment of shapes by reformulating the alignment in an euclidean space with a well-defined metric. Our novel embedding preserves the important advantages of the functional framework, namely that it can be aligned with an linear map, that it is invariant to shape deformations preserving the intrinsic metric and that it can incorporate functional constraints. Using the Green’s embedding has proven to be especially useful for matching non-isometric shapes,

where typically only a few correspondences are known. Additionally to the best of our knowledge, we are the first to include conformality as functional constraints into the matching process. Our evaluation shows that to match non-isometric shapes the novel embedding is superior to the previously utilized deltadistributions. Due to its simplicity and generality our method works on shapes of higher genus, on point clouds and to some degree even in the presence of severe changes in the topology and the inner metric. It therefore might as well serve as a basis for further development of techniques for shape matching. Uncited Reference

596 597 598 599 600 601 602 603 604 605

[20]

606

Appendix

607

2 M Proposition 1a: Injectivity of gM · : M → L (M ), p → g p : Let M M 2 p, q ∈ M with g p = gq and f ∈ L (M ), then 0 =  f, (gM p − M M gM q ) =  f, δ p − δq  = f ( p) − f (q ) ⇒ p = q. Proposition 1b: Let β ∈ L2 (M ) be the coefficients of a linear  combination h(y ) = q β (q )(φ1 + gq (y )) dq ∈ L2 (M ), then

φi , h = φi , β φ1 , h =



 q

√

A 1/λi

i=1 : i≥1

β (q )φ1 , φ1 + gq  dq =

 q

√

β (q ) dq = Aφ1 , β

608 609 610 611 612

613

“i ≥ 2 :

λi φi , h = φi , h = φi , h =     φi , β (q )(φ1 + gq ) dq = φi , β (q )δq dq = q q 

  φi ,  β (q )δq dq = φi , (β ) = φi , β 

q

“Independent”: Let h be a linear combination as defined above with h(y ) = 0, then also β = 0.

Please cite this article as: O. Burghard et al., Embedding shapes with Green’s functions for global shape matching, Computers & Graphics (2017), http://dx.doi.org/10.1016/j.cag.2017.06.004

595

614 615

ARTICLE IN PRESS

JID: CAG 10

616 617 618 619 620 621 622 623 624 625 626 627

O. Burghard et al. / Computers & Graphics xxx (2017) xxx–xxx

“Spans L2 (M )”: Because the φ i form a basis of L2 (M ) it is sufficient to represent the φ i with linear combinations of the gM p . Choosing β = φi results in h = μφi for some μ = 0.  Proposition 1c: Let h(y ) = q β (q )gq (y ) dq ∈ L2 (M ) then analog to above φ1 , h = 0 and φi , h = 1/λi φi , β for i ≥ 2. Therefore  the set is not independent ( q gq (y ) dq = 0) and has co-dimension 1 because φ 1 is not contained in the span. Proposition 3, G(+ Mδp ) =

Eq.

(10):

T + N F δ p

From



AM φ1M ( p)

2

=1

we

get

+ T M + N T M M + + N F φ1 φ1 ( p) = N F δ p = N δT ( p) . T M N M Proposition 3, Eq. (11): F δ p = φ1 (T ( p)) + N B+ Mδp + M ) = φ N (T ( p)) +  δ N N N t = φ1N (T ( p)) + N G(+ δ = δ N T ( p) M p 1 T ( p) N N −1 gM  +  F −1 gM and Eq. (22): F −1 gM N p = φ1 φ1 , F p p + −1 + M N F −1 gM M δ p = +N F T M +M δ M p = N N F p =

1

+ N + T M +N F T [Id + ( − Id )]δ M p = N δT ( p) + N F φ1 

AM

628

Supplementary material

629 630

Supplementary material associated with this article can be found, in the online version, at 10.1016/j.cag.2017.06.004.

631

References

632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 Q3656 657 658 659 660 661 662 663 664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684

[m5G;August 10, 2017;5:38]

[1] Sumner RW, Popovic´ J. Deformation transfer for triangle meshes. ACM Trans Graph 2004;23:399–405. ACM. [2] Sumner RW, Zwicker M, Gotsman C, Popovic´ J. Mesh-based inverse kinematics. ACM Trans Graph 2005;24:488–95. ACM. [3] Kilian M, Mitra NJ, Pottmann H. Geometric modeling in shape space. ACM Trans Graph 2007;26(3):64. [4] Sattler M, Sarlette R, Klein R. Simple and efficient compression of animation sequences. In: Proceedings of the 2005 ACM SIGGRAPH/Eurographics symposium on computer animation. ACM; 2005. p. 209–17. [5] Váša L, Marras S, Hormann K, Brunnett G. Compressing dynamic meshes with geometric Laplacians. Comput Graph Forum 2014;33:145–54. Wiley Online Library. [6] Blanz V, Vetter T. A morphable model for the synthesis of 3D faces. In: Proceedings of the 26th annual conference on computer graphics and interactive techniques. ACM Press/Addison-Wesley Publishing Co.; 1999. p. 187–94. [7] Allen B, Curless B, Popovic´ Z. The space of human body shapes: reconstruction and parameterization from range scans. ACM Trans Graph 2003;22:587– 94. ACM. [8] Hasler N, Stoll C, Sunkel M, Rosenhahn B, Seidel H-P. A statistical model of human pose and body shape. Comput Graph Forum 2009;28:337–46. Wiley Online Library. [9] Amberg B, Romdhani S, Vetter T. Optimal step nonrigid ICP algorithms for surface registration. In: IEEE conference on computer vision and pattern recognition, 2007 (CVPR’07). IEEE; 2007. p. 1–8. [10] Yeh I, Lin C-H, Sorkine O, Lee T-Y, et al. Template-based 3D modelfitting using dual-domain relaxation. IEEE Trans Vis Comput Graph 2011;17(8):1178–90. [11] Burghard O., Berner A., Wand M., Mitra N., Seidel H.-P., Klein R.. Compact part-based shape spaces for dense correspondences. arXiv preprint arxiv: 131175352013. [12] Yoshiyasu Y, Ma W-C, Yoshida E, Kanehiro F. As-conformal-as-possible surface registration. Comput Graph Forum 2014;33:257–67. Wiley Online Library. [13] Kim VG, Lipman Y, Funkhouser T. Blended intrinsic maps. ACM Trans Graph 2011;30:79. ACM. [14] Lipman Y, Funkhouser T. Möbius voting for surface correspondence. ACM Trans Graph 2009;28:72. ACM. [15] Ovsjanikov M, Ben-Chen M, Solomon J, Butscher A, Guibas L. Functional maps: a flexible representation of maps between shapes. ACM Trans Graph 2012;31(4):30. [16] Lipman Y, Rustamov RM, Funkhouser TA. Biharmonic distance. ACM Trans Graph 2010;29(3):27. [17] Van Kaick O, Zhang H, Hamarneh G, Cohen-Or D. A survey on shape correspondence. Comput Graph Forum 2011;30:1681–707. Wiley Online Library. [18] Tam GK, Cheng Z-Q, Lai Y-K, Langbein FC, Liu Y, Marshall D, et al. Registration of 3D point clouds and meshes: a survey from rigid to nonrigid. IEEE Trans Vis Comput Graph 2013;19(7):1199–217. [19] Rusinkiewicz S, Levoy M. Efficient variants of the ICP algorithm. In: Third international conference on 3-D digital imaging and modeling, 2001. Proceedings. IEEE; 2001. p. 145–52. [20] Huang Q-X, Adams B, Wicke M, Guibas LJ. Non-rigid registration under isometric deformations. Comput Graph Forum 2008;27(5):1449–57. doi:10.1111/j. 1467-8659.2008.01285.x. [21] Kraevoy V, Sheffer A. Cross-parameterization and compatible remeshing of 3D models. ACM Trans Graph 2004;23:861–9. ACM.

[22] Aigerman N, Poranne R, Lipman Y. Lifted bijections for low distortion surface mappings. ACM Trans Graph 2014;33(4):69. [23] Aigerman N, Poranne R, Lipman Y. Seamless surface mappings. ACM Trans Graph 2015;34(4):72. [24] Aigerman N, Lipman Y. Orbifold Tutte embeddings.. ACM Trans Graph 2015;34(6):190. [25] Aigerman N, Lipman Y. Hyperbolic orbifold Tutte embeddings. ACM Trans Graph 2016;35(6):217. [26] Kurtek S, Srivastava A, Klassen E, Laga H. Landmark-guided elastic shape analysis of spherically-parameterized surfaces. Comput Graph Forum (Eurograph) 2013;32:429–38. Blackwell Publishing Ltd. [27] Kurtek S, Klassen E, Gore JC, Ding Z, Srivastava A. Elastic geodesic paths in shape space of parameterized surfaces. IEEE Trans Pattern Anal Mach Intell 2012;34(9):1717–30. [28] Laga H., Xie Q., Jermyn I.H., Srivastava A.. Numerical inversion of SRNF maps for elastic shape analysis of genus-zero surfaces. arXiv preprint arxiv: 1610045312016. [29] Jain V, Zhang H. Robust 3D shape correspondence in the spectral domain. In: IEEE international conference on shape modeling and applications, 2006 (SMI 2006). IEEE; 2006. 19–19. [30] Sahillioglu Y, Yemez Y. Minimum-distortion isometric shape correspondence using em algorithm. IEEE Trans Pattern Anal Mach Intell 2012;34(11):2203–15. [31] Aflalo Y, Kimmel R. Spectral multidimensional scaling. Proc Natl Acad Sci 2013;110(45):18052–7. [32] Boscaini D, Girdziušas R, Bronstein MM. Coulomb shapes: using electrostatic forces for deformation-invariant shape representation. In: Proceedings of the 7th eurographics workshop on 3D object retrieval. Eurographics Association; 2014. p. 9–15. [33] Rodolà E., Moeller M., Cremers D.. Point-wise map recovery and refinement from functional correspondence. arXiv preprint arxiv: 1506056032015. [34] Bronstein AM, Bronstein MM, Kimmel R. Generalized multidimensional scaling: a framework for isometry-invariant partial surface matching. Proc Natl Acad Sci 2006;103(5):1168–72. [35] Aflalo Y, Dubrovina A, Kimmel R. Spectral generalized multi-dimensional scaling. CoRR 2013. arXiv:1311.2187. [36] Sun J, Ovsjanikov M, Guibas L. A concise and provably informative multi-scale signature based on heat diffusion. Comput Graph Forum 2009;28:1383–92. Wiley Online Library. [37] Aubry M, Schlickewei U, Cremers D. The wave kernel signature: a quantum mechanical approach to shape analysis. In: 2011 IEEE international conference on computer vision workshops (ICCV workshops). IEEE; 2011. p. 1626–33. [38] Nguyen A, Ben-Chen M, Welnicka K, Ye Y, Guibas L. An optimization approach to improving collections of shape maps. Comput Graph Forum 2011;30:1481– 91. Wiley Online Library. [39] Corman E, Ovsjanikov M, Chambolle A. Continuous matching via vector field flow. Comput Graph Forum 2015;34:129–39. Wiley Online Library. [40] Kovnatsky A, Bronstein MM, Bronstein AM, Glashoff K, Kimmel R. Coupled quasi-harmonic bases. Comput Graph Forum 2013;32:439–48. Wiley Online Library. [41] Kovnatsky A, Bronstein MM, Bresson X, Vandergheynst P. Functional correspondence by matrix completion. In: Proceedings of the IEEE conference on computer vision and pattern recognition; 2015. p. 905–14. [42] Huang Q, Wang F, Guibas L. Functional map networks for analyzing and exploring large shape collections. ACM Trans Graph 2014;33(4):36. [43] Rustamov RM. Laplace–Beltrami eigenfunctions for deformation invariant shape representation. In: Proceedings of the fifth eurographics symposium on geometry processing. Eurographics Association; 2007. p. 225–33. [44] Rustamov RM, Ovsjanikov M, Azencot O, Ben-Chen M, Chazal F, Guibas L. Mapbased exploration of intrinsic shape differences and variability. ACM Trans Graph 2013;32(4):72. [45] Gower JC, Dijksterhuis GB. Procrustes problems, 3. Oxford: Oxford University Press; 2004. [46] Viklands T. Algorithms for the weighted orthogonal procrustes problem and other least squares problems, Sweden: Umea University; 2006. Ph.D. thesis. [47] Crane K, Pinkall U, Schröder P. Spin transformations of discrete surfaces. ACM Trans Graph 2011;30:104. ACM. [48] Pinkall U, Polthier K. Computing discrete minimal surfaces and their conjugates. Exp Math 1993;2(1):15–36. [49] Meyer M, Desbrun M, Schröder P, Barr AH. Discrete differential-geometry operators for triangulated 2-manifolds. Visualization and mathematics III. Springer; 2003. p. 35–57. [50] Belkin M, Sun J, Wang Y. Constructing Laplace operator from point clouds in r d. In: Proceedings of the twentieth annual ACM-SIAM symposium on discrete algorithms. Society for Industrial and Applied Mathematics; 2009. p. 1031–40. [51] Liu Y, Prabhakaran B, Guo X. Point-based manifold harmonics. IEEE Trans Vis Comput Graph 2012;18(10):1693–703. [52] Maneewongvatana S, Mount DM. It’s okay to be skinny, if your friends are fat. In: Center for geometric computing 4th annual workshop on computational geometry, 2; 1999. p. 1–8. [53] Giorgi D., Biasotti S., Paraboschi L.. Shrec: shape retrieval contest: watertight models track. http://watertight.ge.imati.cnr.it, 2007. [Online]. [54] Bronstein AM, Bronstein MM, Kimmel R. Numerical geometry of non-rigid shapes. Springer Science & Business Media; 2008. [55] Myronenko A, Song X, Carreira-Perpinán MA, et al. Non-rigid point set registration: coherent point drift. Adv Neural Inf Process Syst 20 07;19:10 09.

Please cite this article as: O. Burghard et al., Embedding shapes with Green’s functions for global shape matching, Computers & Graphics (2017), http://dx.doi.org/10.1016/j.cag.2017.06.004

685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745 746 747 748 749 750 751 752 753 754 755 756 757 758 759 760 761 762 763 764 765 766 767 768 769