FLASH: Fast landmark aligned spherical harmonic map

This code computes landmark-aligned optimized harmonic maps of genus-0 closed triangle meshes using the FLASH method in [1].
Any comments and suggestions are welcome. 

If you use this code in your own work, please cite the following papers:
[1] P. T. Choi, K. C. Lam, and L. M. Lui, 
    "FLASH: Fast Landmark Aligned Spherical Harmonic Parameterization for Genus-0 Closed Brain Surfaces."
    SIAM Journal on Imaging Sciences, vol. 8, no. 1, pp. 67-94, 2015.

Copyright (c) 2013-2022, Gary Pui-Tung Choi
https://math.mit.edu/~ptchoi/

===============================================================
Usage (see demo.m):

[map,s1,s2] = flash(v1,f1,lm1,v2,f2,lm2,lambda,bigtri1,bigtri2)

Input:
v1: mx3 matrix for vertices of the source mesh
f1: px3 matrix for faces of the source mesh
lm1: kx1 cell array with each entry containing a row vector of indices of the landmark curves on the source mesh
v2: nx3 matrix for vertices of the target mesh
f2: qx3 matrix for vertices of the target mesh
lm2: kx1 cell array with each entry containing a row vector of indices of the landmark curves on the target mesh
lambda: a non-negative landmark matching factor (a larger value gives a smaller landmark error)
(Optional) bigtri1: ID of the face to be "punctured" for the spherical map of the source surface
(Optional) bigtri2: ID of the face to be "punctured" for the spherical map of the target surface
If bigtri1 and bigtri 2 are not specified, the algorithm will look for the most regular triangle on each mesh. 
Choosing a triangle too close to the landmark curves may affect the computation. 

Output:
map: mx3 matrix for vertices of the landmark aligned spherical harmonic map
s1: mx3 matrix for vertices of the spherical conformal parameterization of the source mesh
s2: nx3 matrix for vertices of the spherical conformal parameterization of the target mesh