\section{Anisotropic total variation}
The method considered will build on the total variation regularization
-method of Section \ref{sec:total_variation}. From anisotropic diffusion
-in Section \ref{sec:anisotropic_diffusion} we borrow the idea of making
+method of Section~\ref{sec:total_variation}. From anisotropic diffusion
+in Section~\ref{sec:anisotropic_diffusion} we borrow the idea of making
the regularization in each point directionally dependent. We introduce
the anisotropic total variation
\begin{equation}
\label{eq:aniso_tv_sqrt}
\end{equation}
for all $u \in C^1(\Omega)$. If $A(x)$ is the identity
-matrix we get the regular total variation found in Definition
-\ref{def:tv}. One problem with the regular total
-variation method is that when reducing the total variation, it will also
-try to reduce the variation over known edges in the image, which can
-lead to contrast loss, especially in thin details. By controlling $A(x)$
-such that $\nabla u(x)$ is weighted down across known edges, we hope to
+matrix we get the regular total variation found in
+Definition~\ref{def:tv}. One problem with the regular total variation
+method is that when reducing the total variation, it will also try to
+reduce the variation over known edges in the image, which can lead to
+contrast loss, especially in thin details. By controlling $A(x)$ such
+that $\nabla u(x)$ is weighted down across known edges, we hope to
retain the regularization properties of the method while reducing some
of the negative effects. If the variation across an edge is ``ignored''
by the energy functional, there is no gain in reducing the height of the
&= \sup_{\eta^T A^{-1} \eta \leq 1} \int_\Omega u \diver \eta \,
dx,
\end{align}
+
where $\xi$ and $\eta = \Ahalf \xi$ are in $C_c^\infty(\Omega,
\mathbb{R}^2)$, the space of $C^1$ vector fields with compact support.
If we let $\norm{\xi}_A = \sqrt{\xi^T A \xi}$ and $\norm{\eta}_A^* =
\subsection{Anisotropy tensor}
\label{sec:anisotropy_tensor}
-\fixme{
+\fixme{%
This tensor can be visualized, and can be used for edge
detection. There is probably a bit to be said about the smoothing we
do. We also discuss whether to use the noisy image (yes, probably),
Edges can be found in many different ways, but as suggested by Weickert
in his book on Anisotropic Diffusion \cite{weickert1998anisotropic}, and
-briefly mentioned in Section \ref{sec:anisotropic_diffusion}, a good
+briefly mentioned in Section~\ref{sec:anisotropic_diffusion}, a good
starting point is the \emph{edge detector} $\nabla \tilde{f}_\sigma$ and
the structure tensor
\begin{equation}
noise and small details.
\begin{figure}
- \centering
+ \centering{}
\includegraphics[width=0.4\textwidth]{fig/edges.png}
- \caption{
+ \caption{%
The largest eigenvalue of the structure tensor is $\abs{\nabla
\tilde{f}_\sigma(x)}^2$ and is as we see a good edge detector.
}
as the edge detector itself. Its eigenvalues will be 0 and $\labs{\nabla
\tilde{f}_\sigma(x)}^2$ with corresponding eigenvectors $v_1$ and $v_2$
perpendicular and parallel to $\nabla \tilde{f}_\sigma(x)$ respectively.
-Figure \ref{fig:edges} shows that the largest eigenvalue of the
+Figure~\ref{fig:edges} shows that the largest eigenvalue of the
structure tensor is a good edge detector.
We are also interested in identifying and being sensitive to features in
\begin{equation}
\begin{aligned}
\sigma_1 &= \left(1 +
- \frac{(\lambda_1 - \lambda_2)^2}{\gamma^2}\right)^{-1}, \\
+ \frac{(\lambda_1 - \lambda_2)^2}{\omega^2}\right)^{-1}, \\
\sigma_2 &= 1.
\end{aligned}
\label{eq:sigma_construction}
\sigma_2 = 1$ and the anisotropic total variation is close to the
regular total variation.
-\fixme{
+\fixme{%
The numerical problems should be discussed somewhere but maybe not
here. The eigenvalues are extracted such that
}
takes values either close to or above $F(u)$. More specifically, for a
sequence $u_k$ converging to $u$, we have $F(u) \leq \liminf_k F(u_k)$.
For a function $f : \mathbb{R} \to \mathbb{R}$ this can be visualized as
-in Figure \ref{fig:lower_semicont}.
+in Figure~\ref{fig:lower_semicont}.
Since our space $L^2(\Omega)$ is of infinite dimensions things become a
little bit problematic here. The problem lies in the fact that a
%\TVA(u) = \sup_{\norm{\xi}_A^* \leq 1} \int_\Omega u \, \diver \xi \,
%dx
\TVA(u) = \sup \left\{\int_\Omega u \, \diver \xi \, dx : \xi \in
- C_c^\infty(\Omega, \mathbb{R}^2), \norm{\xi}_A^* \leq 1 \right\}
+ C_c^\infty(\Omega, \mathbb{R}^2), \norm{\xi}_A^* \leq 1 \right\}
\end{equation}
This is again is a
supremum of weakly continuous functionals. Thus the regularization term
1 & \text{if } t > s + r
\end{cases}
\end{equation}
- visualized in Figure \ref{fig:eta_r} such that its derivative takes the
- form shown in Figure \ref{fig:eta_r_diff}. By composing the function
+ visualized in Figure~\ref{fig:eta_r} such that its derivative takes the
+ form shown in Figure~\ref{fig:eta_r_diff}. By composing the function
$\eta_r$ with our image $u$ and Green's identity we obtain
\begin{equation}
\int_\Omega - \eta_r(u) \diver \xi \, dx
perimeter will be computed using a graph cut machinery.
\fixme{well, if we assumed differentiability, we wouldn't need the sup
-definition of the TV...}
+definition of the TV.}
\section{Cauchy--Crofton formulas}
Intuitively, a long curve will intersect more lines.
We write $\mathcal{L}$ for the set of all straight lines in the plane,
-and parametrize them as shown in Figure \ref{fig:line_param}. Thus a
+and parametrize them as shown in Figure~\ref{fig:line_param}. Thus a
line is parametrized by the angle $\phi \in [0, 2\pi)$ of the normal going to the
origin, and the distance $\rho \in [0, \infty)$ from origin to the line. Sometimes it is
more convenient to consider a unit vector $\nu$ giving the direction of
0 & \sigma_2
\end{pmatrix}.
\end{align}
- As $U$ and $U^T$ correspond to rotations they do not contribute to
- the Jacobian, except for changing the input angle of the operator
- $\Lambda$ such that $J_M(\ell_{\phi, \rho}) = J_{\Sigma^2}(U^T \ell_{\phi,
- \rho})$. Given a line
+ As $U$ and $U^T$ correspond to rotations and our measure
+ $\mathcal{L}$ is invariant under rigid motions, $U$ and $U^T$ do not
+ have direct contributions to the Jacobian. They do however affect the
+ input angle of the operator $\Lambda$ such that $J_M(\ell_{\phi,
+ \rho}) = J_{\Sigma^2}(U^T \ell_{\phi, \rho})$. Given a line
\begin{equation}
\ell_{\phi, \rho} = \begin{pmatrix}
\rho\cdot \cos \phi \\
$\abs{J_{\Sigma^2}(\ell_{\phi, \rho})} = \partial_\phi\theta \cdot
\partial_\rho \eta$. Differentiation yields
\begin{align}
- \partial_\phi \theta &= \frac{ \frac{\sigma_1}{\sigma_2}
+ \partial_\phi \theta &= \frac{\frac{\sigma_1}{\sigma_2}
\sec^2 \phi}{1 + \frac{\sigma_1^2}{\sigma_2^2}
\tan^2 \phi}
= \frac{\sigma_1 \sigma_2}{\sigma_1^2 \sin^2 \phi + \sigma_2^2
= \frac{\sigma_1^2 \sigma_2^2}
{\left(\sigma_1^2 \sin^2 \phi + \sigma_2^2 \cos^2 \phi
\right)^{\sfrac{3}{2}}}
- = \frac{\sigma_1^2 \sigma_2^2
+ = \frac{\sigma_1^2 \sigma_2^2
}{\left(\sigma_1^2 \nu_x^2 + \sigma_2^2
\nu_y^2\right)^{\sfrac{3}{2}}}
- = \frac{\det \Sigma^2}{(\nu^T \cdot \Sigma^2 \cdot
+ = \frac{\det \Sigma^2}{(\nu^T \cdot \Sigma^2 \cdot
\nu)^{\sfrac{3}{2}}}.
\end{equation}
We are interested in the Jacobian of the whole transformation, so
\begin{equation}
\abs{C}_{M_\pi}
= \int_C \abs{\dot{C}(t)}_{M_\pi} \, dt
- = \int_C \sqrt{ \dot{C}(t)^T M_\pi\big( C(t) \big) \dot{C}(t) }
+ = \int_C \sqrt{\dot{C}(t)^T M_\pi\big( C(t) \big) \dot{C}(t)}
\, dt.
\label{eq:piecewise_curve_length}
\end{equation}
$\labs{\dot{C}(t)}_M$. We have assumed bounds on the eigenvalues
of $M(x)$ such that according to the Rayleigh principle
\begin{equation}
- K \geq \lambda_1 = \max_\xi \frac{\xi^T M_\pi(x) \xi}{ \xi^T \xi
+ K \geq \lambda_1 = \max_\xi \frac{\xi^T M_\pi(x) \xi}{\xi^T \xi
}
\implies \xi^T M_\pi(x) \xi \leq K \, \xi^T \xi \quad \forall \xi
\end{equation}
function $w_\pi$ is then bounded such that
\begin{equation}
\sum_{x \in \ell_{\nu, \rho} \cap C} w_\pi(\nu, x)
- \leq
+ \leq
\sum_{x \in \ell_{\nu, \rho} \cap C}
\frac{K^2}{k^{\sfrac{3}{2}}}
=
= g(\ell_{\nu, \rho}).
\end{equation}
This is integrable following from the Euclidean Cauchy--Crofton
- formula in Theorem \ref{thm:euclidean_cauchy_crofton} and the fact
+ formula in Theorem~\ref{thm:euclidean_cauchy_crofton} and the fact
that we assumed $C$ to be of finite length
\begin{equation}
\int_\mathcal{L} g(\ell_{\nu, \rho}) \, d\mathcal{L}(\ell_{\nu,
\end{equation}
\end{proof}
-With the anisotropic coarea formula in Theorem
-\ref{thm:anisotropic_coarea} we have a way to calculate the anisotropic
-total variation by integrating the anisotropic perimeter of each level
-set of the image as introduced in \eqref{eq:functional_per}. In the
-Euclidean case, which here would amount to setting the anisotropy tensor
-$A$ equal to the identity matrix $I$, the perimeter coincides nicely
-with the length of the boundary curve, assuming some regularity for the
-boundary.
+With the anisotropic coarea formula in
+Theorem~\ref{thm:anisotropic_coarea} we have a way to calculate the
+anisotropic total variation by integrating the anisotropic perimeter of
+each level set of the image as introduced in \eqref{eq:functional_per}.
+In the Euclidean case, which here would amount to setting the anisotropy
+tensor $A$ equal to the identity matrix $I$, the perimeter coincides
+nicely with the length of the boundary curve, assuming some regularity
+for the boundary.
In the general case we need
to be more careful. As can be seen in \eqref{eq:perimeter_calc}, the
If $P$ is a 90\textdegree{} rotation matrix we have
\begin{equation}
\begin{aligned}
- \PerA(U; \Omega) &= \int_{\partial U} \sqrt{ \langle
+ \PerA(U; \Omega) &= \int_{\partial U} \sqrt{\langle
\nu_{\partial U}, A(x) \nu_{\partial U} \rangle} \, dt \\
- &= \int_{\partial U} \sqrt{ \langle P \nu_{\partial U}, P A(x)
+ &= \int_{\partial U} \sqrt{\langle P \nu_{\partial U}, P A(x)
P^T \, P \nu_{\partial U} \rangle } \, dt.
%&\stackrel{\mathclap{\gamma \hskip 0.07em = \hskip 0.07em
%\partial U}}{=}
normal $\nu_{\partial U}$ gives us the tangent $\pm \dot{\gamma}$.
Inserting this into \eqref{eq:per_to_length1} we get
\begin{equation}
- \PerA(U; \Omega) = \int_\gamma \sqrt{ \langle \dot{\gamma}, M(x) \,
+ \PerA(U; \Omega) = \int_\gamma \sqrt{\langle \dot{\gamma}, M(x) \,
\dot{\gamma} \rangle }.
\label{eq:per_to_length2}
\end{equation}
of $\Omega$.
Now we make sure that all the assumptions of the Riemannian
-Cauchy--Crofton formula in Theorem \ref{thm:riemannian_cauchy_crofton}
+Cauchy--Crofton formula in Theorem~\ref{thm:riemannian_cauchy_crofton}
are fulfilled so that it can be applied to the curve length integral we
have constructed in \eqref{eq:per_to_length2}.
-The structure tensor is constructed as described in Section
-\ref{sec:anisotropy_tensor}
+The structure tensor is constructed as described in
+Section~\ref{sec:anisotropy_tensor}
\begin{equation}
S_\rho(x) = \left(K_{\rho} * \left( \nabla \tilde{f}_{\sigma}
\otimes \nabla \tilde{f}_{\sigma} \right)\right)(x),
\end{equation}
Because of the convolutions with the Gaussian function, this is a smooth
continuous map from $\bar{\Omega}$ to $\mathbb{R}^{2 \times 2}$. As we
-can see in \eqref{eq:s_eigenvalues}, the eigenvalues
+can see in \eqref{eq:s_eigenvalues}, the eigenvalues
depend continuously on the coefficients of the
elements in
the structure tensor $S_\rho(x)$. The extreme value theorem \fixme{ref}
our anisotropy tensor $A(x)$ is bounded away from zero as
\begin{equation}
\lambda_2 = \left(1
- + \frac{(s_1 - s_2)^2}{\gamma^2}\right)^{-1} \geq \left(1 +
- \frac{s_1^2}{\gamma^2}\right)^{-1} \geq k > 0.
-\end{equation}
-\fixme{bah, now gamma is a curve...}
+ + \frac{(s_1 - s_2)^2}{\omega^2}\right)^{-1} \geq \left(1 +
+ \frac{s_1^2}{\omega^2}\right)^{-1} \geq k > 0.
+\end{equation}
Hence, our metric tensor $M(x) = P A(x) P^T$ and curve length
calculation in
-\eqref{eq:per_to_length2} fulfill all the assumptions of Theorem
-\ref{thm:riemannian_cauchy_crofton}. \fixme{did we argue for
+\eqref{eq:per_to_length2} fulfill all the assumptions of
+Theorem~\ref{thm:riemannian_cauchy_crofton}. \fixme{did we argue for
continuous?} Thus we can apply the Cauchy--Crofton formula of
\eqref{eq:riemannian_cauchy_crofton} to calculate the curve length in
\eqref{eq:per_to_length2} and through the anisotropic coarea formula in
-Theorem \ref{thm:anisotropic_coarea} obtain our final expression for the
+Theorem~\ref{thm:anisotropic_coarea} obtain our final expression for the
continuous anisotropic total variation
\begin{equation}
\TVA(u) = \int_{-\infty}^{\infty} \int_\mathcal{L} \sum_{x \in
- \ell_{\nu, \rho} \cap \gamma_s } \frac{\det M(x)}
+ \ell_{\nu, \rho} \cap C_s } \frac{\det M(x)}
{2 \left( \nu^T \cdot M(x) \cdot \nu \right)^{\sfrac{3}{2}}}
\, d\mathcal{L}(\ell_{\nu, \rho}) \, ds,
\end{equation}
-where $\gamma_s = \partial \{ u > s \}$ is the boundary of the level set
+where $C_s = \partial \{ u > s \}$ is the boundary of the level set
for level $s$. Note that $P$ does not affect the determinant, i.e.\
$\det A = \det PAP^T = \det M$, and from our eigendecomposition in
\fixme{ref} we see that the transformation $PAP^T \to M$ actually
just absorb it into the $\beta$ parameter of
\eqref{eq:first_anisotropic_functional}. The perimeter is then
calculated using a discretized version of the Cauchy--Crofton formula
-introduced in Theorem \ref{thm:riemannian_cauchy_crofton}. Again, we
+introduced in Theorem~\ref{thm:riemannian_cauchy_crofton}. Again, we
stop the sum at $L-2$ since the level set $\{ u > L - 1\}$ is empty and
has zero perimeter.
\begin{figure}
\input{fig/line_disc}
\end{figure}
-By approximating the integral Theorem
-\ref{thm:riemannian_cauchy_crofton} by a discrete sum we obtain
+By approximating the integral
+Theorem~\ref{thm:riemannian_cauchy_crofton} by a discrete sum we obtain
the approximation
\begin{equation}
\begin{aligned}
The set of lines $\mathcal{L}$ has been discretized to the lines
$\mathcal{L}_D$. Note that we are approximating the length of the
\emph{differentiable} curve $C$ in $\Omega$. Being a difference in the
-$\rho$ parameter of our line discretization in Figure
-\ref{fig:line_param}, the difference $\Delta \rho$ represents the
-distance from one line to the next in a line family as shown in Figure
-\ref{fig:line_family}. The difference $\Delta \phi$ is taken to be the
-average of the distance to the two neighboring line families as shown in
-Figure \ref{fig:line_neigh}
+$\rho$ parameter of our line discretization in
+Figure~\ref{fig:line_param}, the difference $\Delta \rho$ represents the
+distance from one line to the next in a line family as shown in
+Figure~\ref{fig:line_family}. The difference $\Delta \phi$ is taken to
+be the average of the distance to the two neighboring line families as
+shown in Figure~\ref{fig:line_neigh}
In the discrete setting our domain $\Omega$ is discretized as a regular
grid $\mathcal{G}$. Our image is then reduced to a function $u :
a line of that family intersecting it.
The set of lines can then be represented by the neighborhood of a pixel as
-shown in Figure \ref{fig:line_neigh}. Extending the edges shown in the
-figure gives all lines going through the point considered. Figure
-\ref{fig:line_family} shows all lines of a given family, i.e.\ lines
-having the same angle parameter $\phi$.
+shown in Figure~\ref{fig:line_neigh}. Extending the edges shown in the
+figure gives all lines going through the point considered.
+Figure~\ref{fig:line_family} shows all lines of a given family, i.e.\
+lines having the same angle parameter $\phi$.
Thus not only have we discretized the set of lines, but each line is
made up of \emph{edges} going from one grid point to the next. We will
``did $e$ cross $C$ or not?'' This amounts to checking whether the
terminals of $e$ lie on different sides of the curve $C$, and the
approximation is exact for zero or one intersection points, but will, as
-we see in Figure \ref{fig:curve_edge}, not be entirely correct when we
+we see in Figure~\ref{fig:curve_edge}, not be entirely correct when we
have more.
\begin{figure}
terminal in $\{ u > \lambda \}$ and the other outside. Thus we rewrite
the sum over $e \cap C$ such that
\begin{equation}
- \PerA(\{u > \lambda\}; \Omega) \approx \sum_{e_{ab}} \abs{
- u^\lambda_a - u^\lambda_b} \frac{\det M(e_{ab}) \norm{e_{ab}}^2
- \, \delta^2 \, \Delta\phi}{2 \left(e_{ab}^T \cdot M(e_{ab}) \cdot
- e_{ab}\right)^{\sfrac{3}{2}}}.
+ \PerA(\{u > \lambda\}; \Omega) \approx \sum_{e_{ab}}
+ \abs{u^\lambda_a - u^\lambda_b} \frac{\det M(e_{ab})
+ \norm{e_{ab}}^2 \, \delta^2 \, \Delta\phi}{2 \left(e_{ab}^T
+ \cdot M(e_{ab}) \cdot e_{ab}\right)^{\sfrac{3}{2}}}.
\label{eq:per_approx2}
\end{equation}
The absolute value $\abs{u^\lambda_a - u^\lambda_b}$ is one if one of
increasingly refined discrete domain $\mathcal{G}$, there exists a choice for
$\mathcal{L}_D$ that leads to a consistent Cauchy--Crofton formula. For
convenience we will use a neighborhood representation of $\mathcal{L}_D$
-similar to the one in Figure \ref{fig:line_neigh}.
+similar to the one in Figure~\ref{fig:line_neigh}.
If we consider the edges $e$ of each family separately, the curve length
approximation in \eqref{eq:cauchy_crofton_approx2} can be written
point.
For each $\phi$ parameter, our discretization in the $\rho$ dimension
-can be regarded as a midpoint rule as shown in Figure
-\ref{fig:line_midpoint}. Thus if $\sup \Delta \rho \to 0$, this part of
-the discretization is fine.
+can be regarded as a midpoint rule as shown in
+Figure~\ref{fig:line_midpoint}. Thus if $\sup \Delta \rho \to 0$, this
+part of the discretization is fine.
\begin{figure}
\input{fig/line_midpoint}
The discretization in the $\phi$ dimension can also be regarded as a
version of the \emph{rectangle method}, although not the midpoint rule.
-As shown in Figure \ref{fig:circ_rule}, the summand is evaluated on the
+As shown in Figure~\ref{fig:circ_rule}, the summand is evaluated on the
endpoint of the partition intervals $[\phi_k, \phi_{k+1}]$ and the
difference is taken to be $\Delta \phi = \phi_{k+1} - \phi_k$. Thus if
$\sup \Delta \phi \to 0$, this discretization is also consistent.
To show that all these properties can be fulfilled, we look at a
particular neighborhood stencil construction.
Consider a square centered around a grid point with side lengths
-$\sqrt{\delta}$ as shown in Figure \ref{fig:square_cons}. As $\delta$
+$\sqrt{\delta}$ as shown in Figure~\ref{fig:square_cons}. As $\delta$
goes to zero, the size of this square will go to zero. Inside this
square we can fit a square of $n^2 = \lfloor 1 / \sqrt{\delta} \rfloor^2$ grid
points. This means that the number of grid points along the outer edge
our neighborhood a grid point having the same angle $\phi$ to the
$x$-axis. This means either including the actual grid point at the
outer edge, or one having the same angle, just closer to the center.
-This construction can be seen in Figure \ref{fig:square_cons} for $n =
+This construction can be seen in Figure~\ref{fig:square_cons} for $n =
5$. The maximal $\Delta \phi$ will then be between the horizontal or
-vertical edge and its neighbors, shown in Figure \ref{fig:square_cons}
+vertical edge and its neighbors, shown in Figure~\ref{fig:square_cons}
as angle $a$. These angles can be calculated to be
\begin{equation}
\sup \Delta \phi = \arctan \frac{1/n}{n/2} = \arctan
\begin{equation}
\norm{e} \leq \sqrt{\delta / 2} \to 0.
\end{equation}
-And finally we know from Lemma \ref{lem:delta_rho} that for each line family
+And finally we know from Lemma~\ref{lem:delta_rho} that for each line family
$\delta^2 = \Delta \rho \norm{e}$ and the fact that $\norm{e} \geq
\delta$. Thus for the inter-line distance $\Delta \rho$ we have
\begin{equation}
summand converges to the integrand, and the differences $\Delta \phi$
and $\Delta \rho$ go to zero. Thus the perimeter approximation in
\eqref{eq:per_approx1} is consistent with the continuous formulation in
-Theorem \ref{thm:riemannian_cauchy_crofton}.
+Theorem~\ref{thm:riemannian_cauchy_crofton}.
Note that as we will work with digital images with fixed resolutions, we
do not really have the chance to refine our discretization. We do
We will construct a graph in such a way that if a variable
$u^\lambda_x$ ends up on the source side of the cut we set $u^\lambda_x
= 0$, and if it ends up on the sink side we set $u^\lambda_x = 1$, as in
-Definition \ref{def:graph_representable}. This is an arbitrary choice,
+Definition~\ref{def:graph_representable}. This is an arbitrary choice,
but still something we have to keep in mind through the rest of the
section.
$u^\lambda_x \in \{0,1\}$. The energy $E_\lambda^x(u^\lambda_x)$
and minimum cut obtaining this configuration is shown. The last
two columns show the capacities of the cut for each of the two
- graph constructions in Figure \ref{fig:norm_subgraph}. We
+ graph constructions in Figure~\ref{fig:norm_subgraph}. We
verify that for each of the two graph constructions, the cut
capacities are equal to the energies, plus a constant.
}
The fidelity term of our energy function in \eqref{eq:total_energy}
simplifies to
\begin{align}
- E_\lambda^x(0) &=
+ E_\lambda^x(0) &=
N_x(\lambda + 1) -
N_x(\lambda)
\label{eq:fidelity_energy0} \\
where $E_\lambda^x(0)$ might be positive or negative depending on
$\lambda$ and the pixel value $v_x$.
-Figure \ref{fig:norm_subgraph} shows how graphs can be constructed to
+Figure~\ref{fig:norm_subgraph} shows how graphs can be constructed to
represent this part of the total energy. The construction differs
-depending on whether $E_\lambda^x(0)$ is positive or negative. Table
-\ref{tab:fid_energy} shows how the cuts correspond to the values of
-$u^\lambda_x$ and we can easily verify that the constructed graph
+depending on whether $E_\lambda^x(0)$ is positive or negative.
+Table~\ref{tab:fid_energy} shows how the cuts correspond to the values
+of $u^\lambda_x$ and we can easily verify that the constructed graph
actually represents the fidelity term in the energy function.
\subsubsection{Regularization term}
For our neighboring relation in \eqref{eq:total_energy} of the form
\begin{equation}
- E^{x,y}(u^\lambda_x, u^\lambda_y) =
+ E^{x,y}(u^\lambda_x, u^\lambda_y) =
w_{xy} \abs{u^\lambda_x - u^\lambda_y}
\end{equation}
we have
\end{aligned}
\label{eq:neigh_energies}
\end{equation}
-and by Theorem \ref{thm:graph_rep_id} our energy function is graph
+and by Theorem~\ref{thm:graph_rep_id} our energy function is graph
representable. In \cite{kolmogorov2002energy}, Kolmogorov and Zabih
presents a way to construct a graph for any graph representable function
-on the form shown in Theorem \ref{thm:graph_rep_id}. Since the
+on the form shown in Theorem~\ref{thm:graph_rep_id}. Since the
energies in \eqref{eq:neigh_energies} are especially simple, the
construction and presentation is simplified.
-Figure \ref{fig:neigh_subgraph} shows two different ways of how a graph
-can be constructed to represent the regularization term. See Table
-\ref{tab:neigh_energy} for an overview of how the two values of
-$u^\lambda_x$ corresponds to cuts in the graph.
+Figure~\ref{fig:neigh_subgraph} shows two different ways of how a graph
+can be constructed to represent the regularization term.
+See Table~\ref{tab:neigh_energy} for an overview of how the two values
+of $u^\lambda_x$ corresponds to cuts in the graph.
\begin{table}[t]
\centering
configuration the corresponding energy and the cut yielding this
configuration is shown. The last two columns show the capacities
of the cut in the two alternative graph constructions shown in
- Figure \ref{fig:neigh_subgraph}. We verify that for each of the
+ Figure~\ref{fig:neigh_subgraph}. We verify that for each of the
two graph constructions, the cut capacities are equal to the
energies, plus a constant.
}
\label{tab:neigh_energy}
\end{table}
-Figure \ref{fig:big_graph} shows a visualization of how the final
+Figure~\ref{fig:big_graph} shows a visualization of how the final
graph might look with all its edges. The source will have a lot of
outgoing edges, one for each pixel, while the sink has one incoming edge
for each pixel. The vertices corresponding to the pixels are only
\label{eq:net_flow_cut}
\end{equation}
Note how this definition differs from the capacity of a cut $c(S,T)$ in
-Definition \ref{def:s_t_cut}. While the capacity of a cut represents how
+Definition~\ref{def:s_t_cut}. While the capacity of a cut represents how
much flow it is maximally possible to send from $S$ to $T$, the net flow
across a cut represents the net amount of flow going across the cut,
counting negatively the flow that goes back from $T$ to $S$.
\end{enumerate}
\end{theorem}
See \cite{cormen2009introduction} for a proof, and remark that because
-of the inequality in \eqref{eq:flow_leq_cut}, the cut in Statement
-\ref{it:max_flow_min_cut} is a minimum cut. The theorem does not tell us
-how to find such a cut, and there are multiple ways. One possibility is
-to take $S$ to be all vertices reachable from the source in the residual
-graph and $T = V - S$.
+of the inequality in \eqref{eq:flow_leq_cut}, the cut in
+Statement~\ref{it:max_flow_min_cut} is a minimum cut. The theorem does
+not tell us how to find such a cut, and there are multiple ways. One
+possibility is to take $S$ to be all vertices reachable from the source
+in the residual graph and $T = V - S$.
-Figure \ref{fig:aug_flow} shows a simple graph which already has five
+Figure~\ref{fig:aug_flow} shows a simple graph which already has five
units flowing from $s$ to $t$. The marked path is a possible augmenting
path, and note that it follows an edge in $E$ in the reverse direction,
made possible by the construction of the residual graph just
\subsection{Ford-Fulkerson}
The Ford-Fulkerson algorithm is the most basic augmenting path
algorithm, which can be extended to more advanced algorithms. It is
-stated in pseudocode in Algorithm \ref{alg:ford_fulkerson}, and the idea
+stated in pseudocode in Algorithm~\ref{alg:ford_fulkerson}, and the idea
is to augment the flow along paths from $s$ to $t$ until it is no longer
possible.
\begin{algorithm}
The push procedure moves excess flow from an active vertex along an edge
$(u, v) \in E_f$ for which $d(u) = d(v) + 1$, i.e.\ to a vertex with a
smaller distance label. We call such edges \emph{admissible}. See
-Algorithm \ref{alg:push} for a pseudocode implementation of the push
+Algorithm~\ref{alg:push} for a pseudocode implementation of the push
operation.
Assuming that $f$ is a valid preflow, it is easy to verify that the
%\State $e(v) \mathrel{+}= f_\text{aug}$
\EndFunction
\end{algorithmic}
-\caption{
+\caption{%
\sf The push procedure of the Push-Relabel algorithm.
}
\label{alg:push}
The relabel procedure is our tool for changing the distance labels of
the vertices. It changes the label of a vertex to the greatest possible
value, which is one more than the lowest label among its neighbors in
-the residual graph. See Algorithm \ref{alg:relabel} for a pseudocode
+the residual graph. See Algorithm~\ref{alg:relabel} for a pseudocode
implementation.
\begin{algorithm}
\begin{algorithmic}
A vertex $u$ can only be successfully relabeled to obtain a new label if
the outgoing edges of $u$ in the residual graph have changed since the
previous relabeling. This is why the push and relabel procedures often
-are combined into a \emph{discharge} procedure as shown in Algorithm
-\ref{alg:discharge}. When it is run on an active vertex $u$, we push as
-much as possible of the excess flow to other vertices before the vertex
-is relabeled.
+are combined into a \emph{discharge} procedure as shown in
+Algorithm~\ref{alg:discharge}. When it is run on an active vertex $u$,
+we push as much as possible of the excess flow to other vertices before
+the vertex is relabeled.
\begin{algorithm}
\begin{algorithmic}
\Function{Discharge}{$u$}
heuristics are the gap and global relabeling heuristics, both aiming to
reduce the total number of relabeling steps.
-The gap relabeling heuristic aims to find a label $k$ as in Theorem
-\ref{thm:cut_identification} such that no vertex has that label.
+The gap relabeling heuristic aims to find a label $k$ as in
+Theorem~\ref{thm:cut_identification} such that no vertex has that label.
From vertices $u$ with $d(u) > k$, there are no unsaturated edges going
to vertices with smaller distance labels, so no more flow can ever find
its way from these vertices to the sink. These vertices are therefore
given the label $N$ and never considered again as they will never become
-active. Algorithm \ref{alg:gap} shows a pseudocode representation of
+active. Algorithm~\ref{alg:gap} shows a pseudocode representation of
what is done once a gap $k$ is found.
\begin{algorithm}
When running the push-relabel algorithm with the gap heuristic, we can
be sure that there will never be a vertex $u$ with label $d(u) = N-1$ at
the end of the algorithm, i.e.\ we know that there will always be a gap
-at label $N-1$. This can be seen using the same reasoning as in Theorem
-\ref{thm:cut_identification}, because if there was a vertex with label
+at label $N-1$. This can be seen using the same reasoning as in
+Theorem~\ref{thm:cut_identification}, because if there was a vertex with label
$N-1$, there would only be $N-3$ vertices possibly having labels in
$\{1, \ldots, N-2\}$, so a gap must exist somewhere in that interval.
When using the gap relabeling heuristic, such a gap can not exist, so we
can conclude that there is no vertex with label $N-1$.
-Using Theorem \ref{thm:cut_identification} we can then conclude that the
+Using Theorem~\ref{thm:cut_identification} we can then conclude that the
sets $S = \{ u \in V : d(u) \geq N\}$ and $T = V - S$ form a minimum cut
of the graph.
\subsection{Parametric push-relabel algorithm}
Now we have an algorithm for finding a minimum $s$-$t$-cut in a graph,
-so let's return to the graph constructed in Section
-\ref{sec:graph_construction}. For every level $\lambda \in \{0,
+so let's return to the graph constructed in
+Section~\ref{sec:graph_construction}. For every level $\lambda \in \{0,
\ldots, L\}$ we want to find a minimum $s$-$t$-cut which gives us the
thresholded image $u^\lambda$. These can then hopefully be stacked
together to form the final image $u$.
\input{fig/norm_evolution}
\end{figure}
-Going back to the graph representations in Figure
-\ref{fig:norm_subgraph} and Figure \ref{fig:neigh_subgraph} we know that
-only the capacity of edges from sub-graphs representing the fidelity
-term depend on our level parameter $\lambda$. From the expression in
-\eqref{eq:fidelity_energy0}, visualized in Figure
-\ref{fig:norm_evolution}, we see that the energy term $E_\lambda^x(0)$
-increases monotonically with increasing $\lambda$ parameter. Let $u, v
-\in V - \{s, t\}$. Since the edges in Figure \ref{fig:norm_subgraph} now
-are the only ones depending on $\lambda$, the following is true for
-\emph{decreasing} values of $\lambda$
+Going back to the graph representations in
+Figure~\ref{fig:norm_subgraph} and Figure~\ref{fig:neigh_subgraph} we
+know that only the capacity of edges from sub-graphs representing the
+fidelity term depend on our level parameter $\lambda$. From the
+expression in \eqref{eq:fidelity_energy0}, visualized in
+Figure~\ref{fig:norm_evolution}, we see that the energy term
+$E_\lambda^x(0)$ increases monotonically with increasing $\lambda$
+parameter. Let $u, v \in V - \{s, t\}$. Since the edges in
+Figure~\ref{fig:norm_subgraph} now are the only ones depending on
+$\lambda$, the following is true for \emph{decreasing} values of
+$\lambda$
\begin{description}
\item[Edges from $s$ to $u$\textmd{:}]
- As seen in Figure \ref{fig:norm_subgraph_neg} the capacity of
+ As seen in Figure~\ref{fig:norm_subgraph_neg} the capacity of
these edges will increase monotonically with decreasing
$\lambda$ parameter.
\item[Edges from $u$ to $v$\textmd{:}]
These edges have no $\lambda$-dependence and will remain
unchanged.
\item[Edges from $v$ to $t$\textmd{:}]
- As seen in Figure \ref{fig:norm_subgraph_pos} the capacity of
+ As seen in Figure~\ref{fig:norm_subgraph_pos} the capacity of
these edges will decrease monotonically with decreasing
$\lambda$ parameter.
\end{description}
\subsubsection{Output image construction}
-We mentioned already in Section \ref{sec:total_energy} that in order to
+We mentioned already in Section~\ref{sec:total_energy} that in order to
be able to construct our output image $u$, the thresholded images
$u^\lambda$ would have to stack one on top of the other as shown in
-Figure \ref{fig:img_decomp}. Because of the reuse of the distance labels
+Figure~\ref{fig:img_decomp}. Because of the reuse of the distance labels
between the iterations of the push-relabel algorithm, we can guarantee
that this is possible.
divide and conquer approach only yields improved performance when using
the $L^2$ norm in the fidelity term. This has to do with the fact that
the fidelity term for the $L^1$ norm only changes once for each pixel,
-as we can see in Figure \ref{fig:norm_evolution}, reducing the amount of
+as we can see in Figure~\ref{fig:norm_evolution}, reducing the amount of
work that has to be done in each iteration of the regular parametric
push-relabel algorithm.
\subsection{Implementation}
-A \cpp\ implementation can be found in Appendix
-\ref{app:c++implementation}. It uses the open computer vision library
-OpenCV \cite{opencv_library} to load and save image files.
+A \cpp\ implementation can be found in
+Appendix~\ref{app:c++implementation}. It uses the open computer vision
+library OpenCV \cite{opencv_library} to load and save image files.
Note that when implementing maximum flow algorithms it is not a good
idea, memory- and performance-wise, to actually construct the residual
projects / master / commitdiff