Models for Shape

Snakes/Active Contour Models

  • Starts with a circular configuration and adapts to the shape.

  • Useful when there is little prior information about the object.

Shape Templates

  • Starts with the correct shape of the object.

  • Want to estimate the parameters of the transformation that maps the shape onto the current image.

  • Cannot be fit in closed form because we do not know which edge points in the image correspond to each landmark point in the model.

    • Iterative closest point algorithm alternatively matches points in the image to the landmark points and computes the best closed form transformation.

Statistical Shape Models/Active Shape Models/Point Distribution Models

  • Describe the variation within a class of objects.

  • Capable of adapting to novel individual shapes from that class.

  • Generalized procrustes analysis alternatively optimizes the mean shape and the parameters of the transformations that best maps the observed data points to this mean.

Non-Gaussian Statistical Shape Models

  • Gaussian process latent variable model (GPLVM) extends the PPCA model so that the hidden variables are transformed through a fixed nonlinearity before being weighted by the basis functions.

Exercise 17.1

\(\mathbf{C}_2\) is a singular matrix.

[1]:
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt

def conicMatrix(v, xlim=[-2, 2], ylim=[-2, 2], xlab='x', ylab='y', title='', col='k'):
    x = np.linspace(*xlim, 200)
    y = np.linspace(*ylim, 200)
    x, y = np.meshgrid(x, y)

    _ = (v[0] * x**2 + v[1] * x * y + v[2] * y**2 + v[3] * x + v[4] * y + v[5])
    plt.contour(x, y, _, [0], colors=col)
    plt.xlim(*xlim)
    plt.ylim(*ylim)
    plt.xlabel(xlab)
    plt.ylabel(ylab)
    plt.title(title)
    plt.axis('equal')
    plt.show()

conicMatrix([3, 0, 2, 0, 0, -1], xlim=[-1, 1], ylim=[-1, 1], title='Ellipse', col='red')
conicMatrix([0, 0, 0, 1 * 2, 0, -2], title='Vertical Line', col='green')
conicMatrix([-1, 0, 0, 0, 1 * 2, 0], xlim=[-4, 4], ylim=[-1, 9], title='Parabola', col='blue')
<Figure size 640x480 with 1 Axes>
<Figure size 640x480 with 1 Axes>
<Figure size 640x480 with 1 Axes>

Exercise 17.2

The distance transform computes for each pixel the distance to the nearest background pixel whose value is zero. To compute the distance to the nearest non-zero element of the original binary image, apply the same algorithm to the negated binary image.

[2]:
import numpy

def city_block_distance(grid):
    dt = numpy.zeros_like(grid)
    max_distance = grid.shape[0] * grid.shape[1]

    #top-left corner of matrix in numpy
    for i in range(grid.shape[0]):
        for j in range(grid.shape[1]):
            #already at boundary
            if 0 == grid[i,j]:
                continue

            #not at boundary, so initially assume maximum distance to boundary
            dt[i,j] = max_distance
            #examine west and north neighbors
            for k, l in [(i - 1, j), (i, j - 1)]:
                if 0 <= k < grid.shape[0] and 0 <= l < grid.shape[1]:
                    dt[i,j] = min(1 + dt[k,l], dt[i,j])
                else:
                    #values outside grid are zero
                    dt[i,j] = min(numpy.inf, dt[i,j])
    print('Result from forward pass:\n{}'.format(dt))

    #bottom-right corner of matrix in numpy
    for i in range(grid.shape[0] - 1, -1, -1):
        for j in range(grid.shape[1] - 1, -1, -1):
            #already at boundary
            if 0 == dt[i,j]:
                continue

            #examine east and south neighbors
            for k, l in [(i + 1, j), (i, j + 1)]:
                if 0 <= k < grid.shape[0] and 0 <= l < grid.shape[1]:
                    dt[i,j] = min(1 + dt[k,l], dt[i,j])
                else:
                    #values outside grid are zero
                    dt[i,j] = min(numpy.inf, dt[i,j])
    print('Result from backward pass:\n{}'.format(dt))

grid = numpy.asarray([[0, 0, 0, 0, 1, 0, 0],
                      [0, 0, 1, 1, 1, 0, 0],
                      [0, 1, 1, 1, 1, 1, 0],
                      [0, 1, 1, 1, 1, 1, 0],
                      [0, 1, 1, 1, 0, 0, 0],
                      [0, 0, 1, 0, 0, 0, 0],
                      [0, 0, 0, 0, 0, 0, 0]])
city_block_distance(1 - grid)
Result from forward pass:
[[49 49 49 49  0  1  2]
 [49 49  0  0  0  1  2]
 [49  0  0  0  0  0  1]
 [49  0  0  0  0  0  1]
 [49  0  0  0  1  1  2]
 [49  1  0  1  2  2  3]
 [49  2  1  2  3  3  4]]
Result from backward pass:
[[3 2 1 1 0 1 2]
 [2 1 0 0 0 1 2]
 [1 0 0 0 0 0 1]
 [1 0 0 0 0 0 1]
 [1 0 0 0 1 1 2]
 [2 1 0 1 2 2 3]
 [3 2 1 2 3 3 4]]

Exercise 17.3

\[\begin{split}\text{curve}[\mathbf{w}, n] &= -\left( \mathbf{w}_{n - 1} - 2 \mathbf{w}_n + \mathbf{w}_{n + 1} \right)^\top \left( \mathbf{w}_{n - 1} - 2 \mathbf{w}_n + \mathbf{w}_{n + 1} \right)\\ \text{curve}[\mathbf{w}, 2] &= -\left( \mathbf{w}_1 - 2 \mathbf{w}_2 + \mathbf{w}_3 \right)^\top \left( \mathbf{w}_1 - 2 \mathbf{w}_2 + \mathbf{w}_3 \right)\\ &= \begin{bmatrix} 100 - 2x + 200 & 100 - 2y + 300 \end{bmatrix} \begin{bmatrix} 100 - 2x + 200\\ 100 - 2y + 300 \end{bmatrix}\\ &= (300 - 2x)^2 + (400 - 2y)^2\end{split}\]

By inspection, the point at the minimum is \(\mathbf{w}_2 = \begin{bmatrix} 150 & 200 \end{bmatrix}^\top\).

Exercise 17.4

Since the image is a single color, there are no edge pixels so the distance transform will return \(\infty\) for each pixel. This reduces the likelihood in (17.7) to a constant, so only the prior information will affect the snakes model. The prior coefficients can be set to make the points spread out (or collapse) smoothly (or irregularly).

Exercise 17.5

Suppose the center of the shape could be computed from the edge information. (17.5) could use the averaged distance to the center from each landmark, which is essentially automating (17.8).

Exercise 17.6

Recall that \(\mathbf{w} \in \mathbb{R}^{2N}\) is a compound vector containing all the of the \(x\)- and \(y\)-positions of the landmark points in an image.

We could take a maximum likelihood approach and optimize (17.23) with respect to \(\mathbf{h}\). This is equivalent to a least square fitting as shown in (4.14).

Exercise 17.7

Given \(\mathbf{W} = \mathbf{U} \mathbf{L} \mathbf{V}^\top\) where \(\mathbf{U}, \mathbf{V}\) are orthogonal matrices and \(\mathbf{L}\) is a diagonal matrix,

\[\begin{split}\mathbf{W} \mathbf{W}^\top = \mathbf{U} \mathbf{L} \mathbf{V}^\top \left(\mathbf{U} \mathbf{L} \mathbf{V}^\top\right)^\top = \mathbf{U} \mathbf{L} \mathbf{V}^\top \mathbf{V} \mathbf{L} \mathbf{U}^\top = \mathbf{U} \mathbf{L}^2 \mathbf{U}^\top\\\\ \mathbf{W}^\top \mathbf{W} = \left(\mathbf{U} \mathbf{L} \mathbf{V}^\top\right)^\top \mathbf{U} \mathbf{L} \mathbf{V}^\top = \mathbf{V} \mathbf{L} \mathbf{U}^\top \mathbf{U} \mathbf{L} \mathbf{V}^\top = \mathbf{V} \mathbf{L}^2 \mathbf{V}^\top.\end{split}\]

Exercise 17.8

Let \(\boldsymbol{\theta} = \left\{ \boldsymbol{\mu}, \boldsymbol{\Phi}, \sigma^2 \right\}\).

(a)

\[\begin{split}\DeclareMathOperator{\NormDist}{Norm} q_i(\mathbf{h}_i) &= Pr(\mathbf{h}_i \mid \mathbf{w}_i, \boldsymbol{\theta})\\ &= \frac{ Pr(\mathbf{w}_i \mid \mathbf{h}_i, \boldsymbol{\theta}) Pr(\mathbf{h}_i) }{ Pr(\mathbf{w}_i \mid \boldsymbol{\theta}) } & \quad & \text{(7.50)}\\ &= \frac{ \NormDist_{\mathbf{w}_i}\left[ \boldsymbol{\mu} + \boldsymbol{\Phi} \mathbf{h}_i, \sigma^2 \mathbf{I}_{2N} \right] \NormDist_{\mathbf{h}_i}\left[ \boldsymbol{0}, \mathbf{I}_k \right] }{ \NormDist_{\mathbf{w}_i}\left[ \boldsymbol{\mu}, \boldsymbol{\Phi} \boldsymbol{\Phi}^\top + \sigma^2 \mathbf{I}_{2N} \right] } & \quad & \text{(17.23), (17.24), and (17.25)}\\ &= \NormDist_{\mathbf{h}_i} \left[ \hat{\boldsymbol{\mu}}, \hat{\boldsymbol{\Sigma}} \right] & \quad & \text{(a.1)}\end{split}\]

(b)

\[\begin{split}\DeclareMathOperator*{\argmax}{argmax} \DeclareMathOperator{\E}{\mathrm{E}} \DeclareMathOperator{\tr}{\mathrm{tr}} \boldsymbol{\theta}^{[t]} &= \argmax_\boldsymbol{\theta} \left[ \sum_{i = 1}^I \int q_i^{[t]}(\mathbf{h}_i) \log Pr(\mathbf{w}_i, \mathbf{h}_i \mid \boldsymbol{\theta}) d\mathbf{h}_i \right] & \quad & \text{(7.51)}\\ &= \argmax_\boldsymbol{\theta} \sum_{i = 1}^I \E\left[ \log Pr(\mathbf{w}_i \mid \mathbf{h}_i, \boldsymbol{\theta}) \right] & \quad & \text{(7.36)}\\ &= \argmax_\boldsymbol{\theta} \sum_{i = 1}^I -\frac{1}{2} \E\left[ 2N \log(2 \pi) + \log \left\vert \sigma^2 \mathbf{I}_{2N} \right\vert + \left( \mathbf{w}_i - \boldsymbol{\mu} - \boldsymbol{\Phi} \mathbf{h}_i \right)^\top \left( \sigma^2 \mathbf{I}_{2N} \right)^{-1} \left( \mathbf{w}_i - \boldsymbol{\mu} - \boldsymbol{\Phi} \mathbf{h}_i \right) \right] & \quad & \text{(7.37)}\\ &= \argmax_\boldsymbol{\theta} \sum_{i = 1}^I \E\left[ -N \log 2 \pi - 2N \log \sigma - \frac{1}{2} \sigma^{-2} \left( \mathbf{w}_i - \boldsymbol{\mu} - \boldsymbol{\Phi} \mathbf{h}_i \right)^\top \left( \mathbf{w}_i - \boldsymbol{\mu} - \boldsymbol{\Phi} \mathbf{h}_i \right) \right]\\ &= \argmax_\boldsymbol{\theta} \sum_{i = 1}^I \E\left[ -N \log 2 \pi - 2N \log \sigma - \frac{1}{2} \sigma^{-2} \left( \mathbf{w}_i^\top \mathbf{w}_i + \boldsymbol{\mu}^\top \boldsymbol{\mu} + \mathbf{h}_i^\top \boldsymbol{\Phi}^\top \boldsymbol{\Phi} \mathbf{h}_i + 2 \boldsymbol{\mu}^\top \boldsymbol{\Phi} \mathbf{h}_i - 2 \mathbf{w}_i^\top \boldsymbol{\mu} - 2 \mathbf{w}_i^\top \boldsymbol{\Phi} \mathbf{h}_i \right) \right]\\ &= \argmax_\boldsymbol{\theta} L(\boldsymbol{\theta})\end{split}\]

The following relations are useful for (b.1), (b.2), (b.3):

\[\begin{split}\E\left[ \mathbf{h}_i \right] &= \hat{\boldsymbol{\mu}}\\\\ \E\left[ \mathbf{h}_i \mathbf{h}_i^\top \right] &= \hat{\boldsymbol{\Sigma}} + \E\left[ \mathbf{h}_i \right] \E\left[ \mathbf{h}_i \right]^\top = \hat{\boldsymbol{\Sigma}} + \hat{\boldsymbol{\mu}} \hat{\boldsymbol{\mu}}^\top.\end{split}\]

The SVD solution yields the same mean but different principal components.

(a.1)

Substitute in the corresponding values into Exercise 7.9 to give

\[\begin{split}\hat{\boldsymbol{\Sigma}} &= \left( \boldsymbol{\Phi}^\top \boldsymbol{\Sigma}^{-1} \boldsymbol{\Phi} + \mathbf{I}_k \right)^{-1}\\ &= \left( \sigma^{-2} \boldsymbol{\Phi}^\top \boldsymbol{\Phi} + \mathbf{I}_k \right)^{-1}\end{split}\]

and

\[\begin{split}\hat{\boldsymbol{\mu}} &= \hat{\boldsymbol{\Sigma}} \boldsymbol{\Phi}^\top \boldsymbol{\Sigma}^{-1} (\mathbf{w}_i - \boldsymbol{\mu})\\ &= \sigma^{-2} \left( \sigma^{-2} \boldsymbol{\Phi}^\top \boldsymbol{\Phi} + \mathbf{I}_k \right)^{-1} \boldsymbol{\Phi}^\top (\mathbf{w}_i - \boldsymbol{\mu})\\ &= \left( \boldsymbol{\Phi}^\top \boldsymbol{\Phi} + \sigma^2 \mathbf{I}_k \right)^{-1} \boldsymbol{\Phi}^\top (\mathbf{w}_i - \boldsymbol{\mu}).\end{split}\]

(b.1)

\[\begin{split}\frac{\partial L}{\partial \boldsymbol{\mu}} &= \sum_{i = 1}^I \E\left[ -\frac{1}{2} \sigma^{-2} \left( 2 \boldsymbol{\mu} + 2 \boldsymbol{\Phi} \mathbf{h}_i - 2 \mathbf{w}_i \right) \right] & \quad & \text{(C.27), (C.28), (C.33)}\\ 0 &= \sum_{i = 1}^I -\boldsymbol{\mu} - \boldsymbol{\Phi} \E[\mathbf{h}_i] + \mathbf{w}_i & \quad & \E[] \text{ is a linear operator}\\ &= \sum_{i = 1}^I -\boldsymbol{\mu} - \boldsymbol{\Phi} \left[ \left( \boldsymbol{\Phi}^\top \boldsymbol{\Phi} + \sigma^2 \mathbf{I}_k \right)^{-1} \boldsymbol{\Phi}^\top \left( \mathbf{w}_i - \boldsymbol{\mu} \right) \right] + \mathbf{w}_i & \quad & \E[\mathbf{h}_i] = \hat{\boldsymbol{\mu}}\\ &= \sum_{i = 1}^I \left( \mathbf{I}_{2N} - \boldsymbol{\Phi} \left( \boldsymbol{\Phi}^\top \boldsymbol{\Phi} + \sigma^2 \mathbf{I}_k \right)^{-1} \boldsymbol{\Phi}^\top \right) \mathbf{w}_i - \left( \mathbf{I}_{2N} - \boldsymbol{\Phi} \left( \boldsymbol{\Phi}^\top \boldsymbol{\Phi} + \sigma^2 \mathbf{I}_k \right)^{-1} \boldsymbol{\Phi}^\top \right) \boldsymbol{\mu}\\ &= \sum_{i = 1}^I \left( \mathbf{I}_{2N} + \sigma^{-2} \boldsymbol{\Phi} \boldsymbol{\Phi}^\top \right)^{-1} \mathbf{w}_i - \left( \mathbf{I}_{2N} + \sigma^{-2} \boldsymbol{\Phi} \boldsymbol{\Phi}^\top \right)^{-1} \boldsymbol{\mu} & \quad & \text{Sherman-Morrison-Woodbury formula}\\ \boldsymbol{\mu} &= \frac{1}{I} \sum_{i = 1}^I \mathbf{w}_i\end{split}\]

(b.2)

\[\begin{split}\frac{\partial L}{\partial \boldsymbol{\Phi}} &= \sum_{i = 1}^I \E\left[ -\frac{1}{2} \sigma^{-2} \left( 2 \boldsymbol{\Phi} \mathbf{h}_i \mathbf{h}_i^\top + 2 \boldsymbol{\mu} \mathbf{h}_i^\top - 2 \mathbf{w}_i \mathbf{h}_i^\top \right) \right] & \quad & \text{(C.31), (C.29)}\\ 0 &= \sum_{i = 1}^I -\boldsymbol{\Phi} \E\left[ \mathbf{h}_i \mathbf{h}_i^\top \right] - \boldsymbol{\mu} \E\left[ \mathbf{h}_i^\top \right] + \mathbf{w}_i \E\left[ \mathbf{h}_i^\top \right] & \quad & \E[] \text{ is a linear operator}\\ \boldsymbol{\Phi} &= \left( \sum_{i = 1}^I (\mathbf{w}_i - \boldsymbol{\mu}) \E\left[ \mathbf{h}_i^\top \right] \right) \left( \sum_{i = 1}^I \E\left[ \mathbf{h}_i \mathbf{h}_i^\top \right] \right)^{-1}\end{split}\]

(b.3)

\[\begin{split}\frac{\partial L}{\partial \sigma} &= \sum_{i = 1}^I \E\left[ -2N\sigma^{-1} + \sigma^{-3} \left( \mathbf{w}_i^\top \mathbf{w}_i + \boldsymbol{\mu}^\top \boldsymbol{\mu} + \mathbf{h}_i^\top \boldsymbol{\Phi}^\top \boldsymbol{\Phi} \mathbf{h}_i + 2 \boldsymbol{\mu}^\top \boldsymbol{\Phi} \mathbf{h}_i - 2 \mathbf{w}_i^\top \boldsymbol{\mu} - 2 \mathbf{w}_i^\top \boldsymbol{\Phi} \mathbf{h}_i \right) \right]\\ \sigma^2 &= \frac{1}{2NI} \sum_{i = 1}^I (\mathbf{w}_i - \boldsymbol{\mu})^\top (\mathbf{w}_i - \boldsymbol{\mu}) + \E\left[ \mathbf{h}_i^\top \boldsymbol{\Phi}^\top \boldsymbol{\Phi} \mathbf{h}_i \right] - \E\left[ \left( \mathbf{w}_i - \boldsymbol{\mu} \right)^\top \boldsymbol{\Phi} \mathbf{h}_i \right] & \quad & \E[] \text{ is a linear operator}\\ &= \frac{1}{2NI} \sum_{i = 1}^I (\mathbf{w}_i - \boldsymbol{\mu})^\top (\mathbf{w}_i - \boldsymbol{\mu}) + \tr\left( \boldsymbol{\Phi}^\top \boldsymbol{\Phi} \hat{\boldsymbol{\Sigma}} \right) + \hat{\boldsymbol{\mu}}^\top \boldsymbol{\Phi}^\top \boldsymbol{\Phi} \hat{\boldsymbol{\mu}} - \left( \mathbf{w}_i - \boldsymbol{\mu} \right)^\top \boldsymbol{\Phi} \hat{\boldsymbol{\mu}} & \quad & \text{Matrix Cookbook (328), (330)}\end{split}\]

Exercise 17.9

\[\begin{split}L &= \sum_{n = 1}^N -\sigma^{-2} \left( \mathbf{y}_n - \mathbf{A} \mathbf{w}_n - \mathbf{b} \right)^\top \left( \mathbf{y}_n - \mathbf{A} \mathbf{w}_n - \mathbf{b} \right) - \log \mathbf{h}^\top \mathbf{h}\\ &= \sum_{n = 1}^N -\sigma^{-2} \left( \mathbf{y}_n^\top \mathbf{y}_n - \mathbf{y}_n^\top \mathbf{A} \mathbf{w}_n - \mathbf{y}_n^\top \mathbf{b} - \mathbf{w}_n^\top \mathbf{A}^\top \mathbf{y}_n + \mathbf{w}_n^\top \mathbf{A}^\top \mathbf{A} \mathbf{w}_n + \mathbf{w}_n^\top \mathbf{A}^\top \mathbf{b} - \mathbf{b}^\top \mathbf{y}_n + \mathbf{b}^\top \mathbf{A} \mathbf{w}_n + \mathbf{b}^\top \mathbf{b} \right) - \log \mathbf{h}^\top \mathbf{h}\\ &= \sum_{n = 1}^N -\sigma^{-2} \left( \mathbf{y}_n^\top \mathbf{y}_n - 2 \mathbf{y}_n^\top \mathbf{A} \mathbf{w}_n + \mathbf{w}_n^\top \mathbf{A}^\top \mathbf{A} \mathbf{w}_n - 2 \mathbf{b}^\top \mathbf{y}_n + 2 \mathbf{b}^\top \mathbf{A} \mathbf{w}_n + \mathbf{b}^\top \mathbf{b} \right) - \log \mathbf{h}^\top \mathbf{h}\\ &= \sum_{n = 1}^N -\sigma^{-2} \left[ \left( \mathbf{y}_n - \mathbf{b} \right)^\top \left( \mathbf{y}_n - \mathbf{b} \right) + \mathbf{w}_n^\top \mathbf{A}^\top \mathbf{A} \mathbf{w}_n - 2 \left( \mathbf{y}_n - \mathbf{b} \right)^\top \mathbf{A} \mathbf{w}_n \right] - \log \mathbf{h}^\top \mathbf{h}\end{split}\]

where \(\mathbf{w}_n = \boldsymbol{\mu}_n + \boldsymbol{\Phi}_n \mathbf{h}\) and

\[\begin{split}\frac{\partial L}{\partial \mathbf{h}} &= \sum_{n = 1}^N -2 \sigma^{-2} \left( \boldsymbol{\mu}_n^\top \mathbf{A}^\top \mathbf{A} \boldsymbol{\Phi}_n + \mathbf{h}^\top \boldsymbol{\Phi}_n^\top \mathbf{A}^\top \mathbf{A} \boldsymbol{\Phi}_n \right) + 2 \sigma^{-2} \left( \mathbf{y}_n - \mathbf{b} \right)^\top \mathbf{A} \boldsymbol{\Phi}_n - \frac{2}{\mathbf{h}^\top \mathbf{h}} \mathbf{h}^\top & \quad & \text{(C.28), (C.33)}\\ \sigma^2 \mathbf{h}^\top + \sum_{n = 1}^N \mathbf{h}^\top \boldsymbol{\Phi}_n^\top \mathbf{A}^\top \mathbf{A} \boldsymbol{\Phi}_n &= \sum_{n = 1}^N \left( \mathbf{y}_n - \mathbf{A} \boldsymbol{\mu}_n - \mathbf{b} \right)^\top \mathbf{A} \boldsymbol{\Phi}_n & \quad & \left\Vert \mathbf{h} \right\Vert_2 = 1\\ \mathbf{h} &= \left( \sigma^2 \mathbf{I} + \sum_{n = 1}^N \boldsymbol{\Phi}_n^\top \mathbf{A}^\top \mathbf{A} \boldsymbol{\Phi}_n \right)^{-1} \sum_{n = 1}^N \boldsymbol{\Phi}_n^\top \mathbf{A}^\top \left( \mathbf{y}_n - \mathbf{A} \boldsymbol{\mu}_n - \mathbf{b} \right).\end{split}\]

Exercise 17.10

The activation function (9.66) can use \(\mathbf{h}\) to learn the coefficients of the gender classification model (9.65).

The ICP approach to (17.31) already handles the lack of landmark points.

Exercise 17.11

Each point in the point distribution could have mass-spring constraints, which can be used to iteratively redistribute the points towards the hidden portion of the shape.

Exercise 17.12

[GH+96] provides more details on Modeling Complex Data Densities.

[TB99] elaborates on the book’s coverage of probabilistic principal component analysis by showing how to train a mixture of them.

[PM00] extends the book’s analysis of t-distributions to mixture models.

[dRF03][KD04] are alternative nonlinear shape models that should be read after the foregoing materials.

Exercise 17.13

Rewrite \(\mathbf{x}_\cdot \in \mathbb{R}^2\) as

\[\begin{split}\mathbf{x}_\cdot &= \alpha \mathbf{a}_\cdot + \beta \mathbf{b}_\cdot + \gamma \mathbf{c}_\cdot\\ &= \alpha \mathbf{a}_\cdot + \beta \mathbf{b}_\cdot + (1 - \alpha - \beta) \mathbf{c}_\cdot & \quad & \alpha + \beta + \gamma = 1\\ &= \alpha (\mathbf{a}_\cdot - \mathbf{c}_\cdot) + \beta (\mathbf{b}_\cdot - \mathbf{c}_\cdot) + \mathbf{c}_\cdot\\ \mathbf{x}_\cdot - \mathbf{c}_\cdot &= \begin{bmatrix} \mathbf{a}_\cdot - \mathbf{c}_\cdot & \mathbf{b}_\cdot - \mathbf{c}_\cdot \end{bmatrix} \begin{bmatrix} \alpha\\ \beta \end{bmatrix}\end{split}\]

When \(\mathbf{x}_\cdot \in \mathbb{R}^3\) then

\[\begin{split}\begin{bmatrix} \mathbf{x}_\cdot\\ 1 \end{bmatrix} &= \begin{bmatrix} \mathbf{a}_\cdot & \mathbf{b}_\cdot & \mathbf{c}_\cdot\\ 1 & 1 & 1 \end{bmatrix} \begin{bmatrix} \alpha\\ \beta\\ \gamma \end{bmatrix}.\end{split}\]

Assuming that the original and transformed positions have the same barycentric coordinates, warping the whole image is equivalent to

  1. computing each position’s barycentric coordinates,

  2. warping each triangle,

  3. solving for the position relative to the warped triangle.

Exercise 17.14

(i)

\[\begin{split}0 = \tilde{\mathbf{w}}^\top \begin{bmatrix} \mathbf{A} & \mathbf{b}\\ \mathbf{b}^\top & c \end{bmatrix} \tilde{\mathbf{w}} &= \begin{bmatrix} \tilde{\mathbf{x}}^\top & s \end{bmatrix} \begin{bmatrix} \mathbf{A} & \mathbf{b}\\ \mathbf{b}^\top & c \end{bmatrix} \begin{bmatrix} \tilde{\mathbf{x}}\\ s \end{bmatrix}\\ &= \begin{bmatrix} \tilde{\mathbf{x}}^\top & s \end{bmatrix} \begin{bmatrix} \mathbf{A} \tilde{\mathbf{x}} + \mathbf{b} s\\ \mathbf{b}^\top \tilde{\mathbf{x}} + cs \end{bmatrix}\\ &= \tilde{\mathbf{x}}^\top \mathbf{A} \tilde{\mathbf{x}} + 2 \tilde{\mathbf{x}}^\top \mathbf{b} s + cs^2\\ &= \gamma + \beta s + \alpha s^2\end{split}\]

(ii)

There exists a real unique solution for \(s\) when the discriminant is zero:

\[\begin{split}0 = \beta^2 - 4 \alpha \gamma &= \left( 2 \tilde{\mathbf{x}}^\top \mathbf{b} \right)^2 - 4c \left( \tilde{\mathbf{x}}^\top \mathbf{A} \tilde{\mathbf{x}} \right)\\ &= \tilde{\mathbf{x}}^\top \mathbf{b} \mathbf{b}^\top \tilde{\mathbf{x}} - c \tilde{\mathbf{x}}^\top \mathbf{A} \tilde{\mathbf{x}}\\ &= \tilde{\mathbf{x}}^\top \left( \mathbf{b} \mathbf{b}^\top - c \mathbf{A} \right) \tilde{\mathbf{x}}\\ &= \tilde{\mathbf{x}}^\top \mathbf{C} \tilde{\mathbf{x}}\end{split}\]

If the camera has intrinsic matrix \(\boldsymbol{\Lambda}\), then \(\mathbf{w} = \boldsymbol{\Lambda}^{-1} \mathbf{x} = \tilde{\mathbf{x}}\) and \(\mathbf{C} = \boldsymbol{\Lambda}^{-\top} \left( \mathbf{b} \mathbf{b}^\top - c \mathbf{A} \right) \boldsymbol{\Lambda}^{-1}\).