+
+void solveWeights(BlendshapeModel *b, std::map<VertIdx, glm::vec3> manips) {
+
+ int numBlends = b->blendshapes.size();
+ int numManips = manips.size();
+
+ Eigen::VectorXf manipOffset(3 * numManips); // 3 * numManips x 1
+ {
+ int i = 0;
+ for (auto &m : manips) {
+ auto origPos = b->neutral[m.first.second];
+ auto manipOffsetV = m.second - origPos;
+ for (int j = 0; j < 3; j++)
+ manipOffset(i * 3 + j) = manipOffsetV[j];
+ i++;
+ }
+ }
+
+ Eigen::MatrixXf Bbar(3 * numManips, numBlends); // 3 * numManips x numBlends
+ for (int i = 0; i < numBlends; i++) {
+ int j = 0;
+ for (auto &m : manips) {
+ for (int k = 0; k < 3; k++)
+ Bbar(j * 3 + k, i) =
+ b->blendshapes[i].deltas[m.first.second][k];
+ j++;
+ }
+ }
+
+ Eigen::VectorXf prevWeights(numBlends); // numBlends x 1
+ for (int i = 0; i < numBlends; i++)
+ prevWeights(i) = b->blendshapes[i].weight;
+
+ float alpha = .1, mu = 1;
+
+ Eigen::MatrixXf LHS(3 * numManips + numBlends * 2, numBlends);
+ LHS << Bbar, alpha * Eigen::MatrixXf::Identity(numBlends, numBlends),
+ mu * Eigen::MatrixXf::Identity(numBlends, numBlends);
+
+ Eigen::VectorXf RHS(3 * numManips + numBlends * 2);
+ RHS << manipOffset, alpha * prevWeights, Eigen::VectorXf::Zero(numBlends);
+
+ Eigen::VectorXf w =
+ LHS.jacobiSvd(Eigen::ComputeThinU | Eigen::ComputeFullV).solve(RHS);
+
+ for (int i = 0; i < numBlends; i++)
+ b->blendshapes[i].weight = w[i];
+}
+
+void stepBlendshapeAnim(BlendshapeModel *b) {
+ for (int i = 0; i < b->blendshapes.size(); i++) {
+ b->blendshapes[i].weight = b->animation[b->curFrame][i];
+ }
+ b->curFrame++;
+ b->curFrame %= b->animation.size();
+}