# Building a shape model from data

The goal in this tutorial is to learn how to build a Statistical Shape Model from meshes in correspondence. Furthermore, we discuss the importance of rigid alignment while doing so.

The following resources from our online course may provide some helpful context for this tutorial:

• Learning a model from example data (Video)
##### Preparation

As in the previous tutorials, we start by importing some commonly used objects and initializing the system.

import scalismo.geometry._
import scalismo.common._
import scalismo.ui.api._
import scalismo.mesh._
import scalismo.io.{StatisticalModelIO, MeshIO}
import scalismo.statisticalmodel._
import scalismo.registration._
import scalismo.statisticalmodel.dataset._

scalismo.initialize()
implicit val rng = scalismo.utils.Random(42)

val ui = ScalismoUI()


Let us load (and visualize) a set of face meshes based on which we would like to model shape variation:

val dsGroup = ui.createGroup("datasets")

val meshFiles = new java.io.File("datasets/nonAlignedFaces/").listFiles

val (meshes, meshViews) = meshFiles.map(meshFile => {
val meshView = ui.show(dsGroup, mesh, "mesh")
(mesh, meshView) // return a tuple of the mesh and the associated view
}) .unzip // take the tuples apart, to get a sequence of meshes and one of meshViews


You immediately see that the meshes are not aligned. What you cannot see, but is very important for this tutorial, is that the meshes are in correspondence. This means that for every point on one of the face meshes (corner of eye, tip of nose, …), we can identify the corresponding point on other meshes. Corresponding points are identified by the same point id.

Exercise: verify that the meshes are indeed in correspondence by displaying a few corresponding points.

#### Rigidly aligning the data:

In order to study shape variations, we need to eliminate variations due to relative spatial displacement of the shapes (rotation and translation). This can be achieved by selecting one of the meshes as a reference, to which the rest of the datasets are aligned. In this example here, we simply take the first mesh in the list as a reference and align all the others.

val reference = meshes.head
val toAlign : IndexedSeq[TriangleMesh[_3D]] = meshes.tail


Given that our dataset is in correspondence, we can specify a set of point identifiers, to locate corresponding points on the meshes.

val pointIds = IndexedSeq(2214, 6341, 10008, 14129, 8156, 47775)
val refLandmarks = pointIds.map{id => Landmark(s"L_\$id", reference.pointSet.point(PointId(id))) }


After locating the landmark positions on the reference, we iterate on each remaining data item, identify the corresponding landmark points and then rigidly align the mesh to the reference.

val alignedMeshes = toAlign.map { mesh =>
val landmarks = pointIds.map{id => Landmark("L_"+id, mesh.pointSet.point(PointId(id)))}
val rigidTrans = LandmarkRegistration.rigid3DLandmarkRegistration(landmarks, refLandmarks, center = Point(0,0,0))
mesh.transform(rigidTrans)
}


Now, the IndexedSeq of triangle meshes alignedMeshes contains the faces that are aligned to the reference mesh.

Exercise: verify visually that at least the first element of the aligned dataset is indeed aligned to the reference.

### Building a discrete Gaussian process from data

Now that we have a set of meshes, which are in correspondence and aligned to our reference, we can turn the dataset into a set of deformation fields, from which we then build the model:

val defFields = alignedMeshes.map{ m =>
val deformationVectors = reference.pointSet.pointIds.map{ id : PointId =>
m.pointSet.point(id) - reference.pointSet.point(id)
}.toIndexedSeq

DiscreteField[_3D, UnstructuredPointsDomain[_3D], EuclideanVector[_3D]](reference.pointSet, deformationVectors)
}


Learning the shape variations from this deformation fields is done by calling the method createUsingPCA of the DiscreteLowRankGaussianProcess class. Note that the deformation fields need to be interpolated, such that we are sure that they are defined on all the points of the reference mesh.

val interpolator = NearestNeighborInterpolator[_3D, EuclideanVector[_3D]]()
val continuousFields = defFields.map(f => f.interpolate(interpolator) )
val gp = DiscreteLowRankGaussianProcess.createUsingPCA(reference.pointSet, continuousFields)


Exercise: display the mean deformation field of the returned Gaussian Process.

Exercise: sample and display a few deformation fields from this GP.

Exercise: using the GP’s *cov method, evaluate the sample covariance between two close points on the right cheek first, and a point on the nose and one on the cheek second. What does the data tell you?*

By combining this Gaussian process over deformation fields with the reference mesh, we obtain a Statistical Mesh Model:

val model = StatisticalMeshModel(reference, gp.interpolate(interpolator))


Notice that when we visualize this mesh model in Scalismo-ui, it generates a GaussianProcessTransformation and the reference mesh in the Scene Tree of Scalismo-ui. .

val modelGroup = ui.createGroup("model")
val ssmView = ui.show(modelGroup, model, "model")


### An easier way to build a model.

Performing all the operations above every time we want to build a PCA model from a set of files containing meshes in correspondence can be tedious. Therefore, Scalismo provides a more easy to use implementation via the DataCollection data structure.

The DataCollection class in Scalismo allows grouping together a dataset of meshes in correspondence, in order to make collective operations on such sets easier.

We can create a DataCollection by providing a reference mesh, and a sequence of meshes, which are in correspondence with this reference.

val dc = DataCollection.fromMeshSequence(reference, alignedMeshes)._1.get


For each mesh, the data collection automatically computes the transformation, which transforms the reference into the respective mesh:

val item0 :DataItem[_3D] = dc.dataItems(0)
val transform : Transformation[_3D] = item0.transformation


Now that we have our data collection, we can build a shape model as follows:

val modelNonAligned = StatisticalMeshModel.createUsingPCA(dc).get

val modelGroup2 = ui.createGroup("modelGroup2")
ui.show(modelGroup2, modelNonAligned, "nonAligned")


Here again, a PCA is performed based the deformation fields retrieved from the data in correspondence.

Notice that, in this case, we built a model from misaligned meshes in correspondence.

Exercise: sample a few faces from the second model. How does the quality of the obtained shapes compare to the model built from aligned data?

Exercise: using the GUI, change the coefficient of the first principal component of the nonAligned shape model. What is the main shape variation encoded in the model?