\newcommand{\lnorm}[1]{\lVert #1 \rVert}
\newcommand{\diver}{\mathop{\mathrm{div}}\nolimits}
\newcommand{\TV}{\mathrm{TV}}
-\newcommand{\TVA}{\mathrm{TV}_A}
+\newcommand{\TVA}{\mathrm{TV}_{\!A}}
\newcommand{\BV}{\mathrm{BV}}
\newcommand{\Per}{\mathrm{Per}}
\newcommand{\PerA}{\mathrm{Per}_A}
\chapter{Continuous formulation}
-From the introduction we see that there are many different approaches to
+In the previous chapter we saw that there are many different approaches to
the image restoration problem, all with their own strengths and
weaknesses. The method considered in this thesis is an anisotropic total
-variation formulation, and the hope is that we keep the strengths of the
+variation formulation, and the aim is to keep the strengths of the
anisotropic diffusion and total variation methods, while eliminating
-some of their respective weaknesses.
+some of their weaknesses.
This chapter will be devoted to the continuous formulation of the
method. We will look at the functional we want to minimize and its
anisotropic coarea formula, the anisotropic total variation is rewritten
as an integral of the perimeter of all the level sets of the image.
-Later, the anisotropic Cauchy--Crofton formula is introduced to make it
+Following that, the anisotropic Cauchy--Crofton formula is introduced to make it
feasible to calculate the perimeter of these level sets. All of this
leads up to the discretization of our functional in the next chapter.
before.
Note that $u(x)$ and $A(x)$ are always dependent on the position in the
-image $x$, but we will sometimes drop the $x$, when no misunderstandings
-are possible.
+image $x$, but we will sometimes drop the $x$, when the intended meaning
+is clear.
As we will not always be working with differentiable images, we extend
the definition of the total variation functional. Being symmetric
\end{aligned}
\end{equation}
where $\xi$ and $\eta = \Ahalf \xi$ are in $C_c^\infty(\Omega,
-\mathbb{R}^2)$, the space of $C^\infty$ vector fields with compact support.
-In the following we let $\norm{\xi}_A = \sup_x \sqrt{\xi^T A \xi}$ and
-$\norm{\eta}_A^* = \sup_x \sqrt{\eta^T A^{-1} \eta}$, and with that we
+\mathbb{R}^2)$, the space of smooth vector fields with compact support.
+In the following we define the norms $\lnorm{\xi}_A = \sup_x (\xi^T A
+\xi)^{\sfrac{1}{2}}$ and
+$\lnorm{\eta}_A^* = \sup_x (\eta^T A^{-1} \eta)^{\sfrac{1}{2}}$, and with that we
present the formal definition of the anisotropic total variation.
\begin{definition}[Anisotropic total variation]
For a function $u \in BV(\Omega)$ and a continuous symmetric
Section~\ref{sec:diffusion_filtering}: $f_\sigma = K_\sigma * \tilde{f}$,
where $\tilde{f}$ is the symmetric extension of the initial image $f$ in
$\mathbb{R}^2$. The smoothing parameter $\sigma$ is called the
-\emph{noise scale}, and with this we aim to avoid being too sensitive to
-noise in $f$.
+\emph{noise scale}, and it controls the scale at which details are
+considered to be noise.
As seen in Figure~\ref{fig:edges}, the edge detector is fine for
detecting edges, but it can not give us information about larger
$\lambda_2 = 0$ with corresponding eigenvectors $v_1$ and $v_2$ parallel
and perpendicular to $\nabla f_\sigma(x)$ respectively.
-We are also interested in identifying and being sensitive to features in
-a neighborhood around the point $x$, such as corners or curved edges and
-\fixme{coherent structures}. This is why we introduce the component-wise
-convolution with $K_\rho$ such that
+To detect features in a neighborhood around the point $x$, such as
+corners, curved edges and coherent structures we introduce the
+component-wise convolution with $K_\rho$ such that
\begin{equation}
S_\rho(x) := K_\rho * \big( \nabla f_\sigma(x) \otimes
\nabla f_\sigma(x) \big)(x).
\end{equation}
The parameter $\rho$, called
the \emph{integration scale}, controls the size of the neighborhood
-which affects the structure tensor. It controls the size of the
+which affects the structure tensor. Thus it defines the size of the
structures we want our anisotropy tensor to be sensitive to.
\begin{figure}
\centering{}
- \includegraphics[width=0.4\textwidth]{fig/edges.png}
+ \begin{subfigure}[b]{0.45\textwidth}
+ \includegraphics[width=\textwidth]{fig/factory/finger/n_q200.png}
+ \end{subfigure}
+ ~
+ \begin{subfigure}[b]{0.45\textwidth}
+ \includegraphics[width=\textwidth]{fig/factory/finger/r_p2_n32_b15000_g150_r20_s3_edge.png}
+ \end{subfigure}
\caption[The edge detector $\abs{\nabla f_\sigma(x)}$]{%
- The largest eigenvalue of the structure tensor is $\abs{\nabla
- f_\sigma(x)}^2$ and can as we see function as an edge detector.
+ A noisy fingerprint on the left, and the largest eigenvalue of
+ the structure tensor is $\abs{\nabla f_\sigma(x)}^2$ on the
+ left, which---as we can see---functions as an edge
+ detector.
}
\label{fig:edges}
\end{figure}
The smoothed tensor $S_\rho(x)$ can easily be verified to be symmetric
-positive semi-definite, just like $S_0(x)$. We order the two real
+positive semi-definite, just like $S_0(x)$. In addition, when $\rho >
+0$, the elements of $S_\rho$ are smooth maps from $\Omega$ to
+$\mathbb{R}$.
+
+We order the two real
eigenvalues $\lambda_1 \geq \lambda_2$ with corresponding eigenvectors
$v_1$ and $v_2$. From the characteristic polynomial of $S_\rho(x) =
(\begin{smallmatrix} s_{11} & s_{12} \\ s_{12} & s_{22}
The vector
$v_1$ will then indicate the direction of most variation in the
-neighborhood. From this we see that an edge would give $\lambda_1 \gg
-\lambda_2 \approx 0$, while smooth areas would give $\lambda_1 \approx
+neighborhood. An edge will give $\lambda_1 \gg
+\lambda_2 \approx 0$, while smooth areas will give $\lambda_1 \approx
\lambda_2 \approx 0$. In corners we have variation in the direction of
$v_1$ but also perpendicular to $v_1$, so we will have $\lambda_1
\approx \lambda_2 \gg 0$. Thus the quantity $(\lambda_1 - \lambda_2)^2$
seen in Figure~\ref{fig:tensor_viz} where the two tensors are shown at
an edge in the image.
+\begin{figure}
+ \input{fig/tensor_viz}
+\end{figure}
+
In smooth areas, $\sigma_1 \approx 1$ and
-$A(x)$ will be close to the identity matrix, while $A(x)$ will reduce
-the effect of $\nabla u$ across edges.
+$A(x)$ will be close to the identity matrix. At or around edges,
+$\sigma_1$, which corresponds to the eigenvector perpendicular to the
+edge, will be small.
Around corners $A(x)$ will be close to the identity matrix, which gives
regularization similar to smooth areas. This is one possible down-side
The parameter $\omega$ controls the amount of anisotropy in the method,
such that if it is very large we are left with the identity matrix. Note
-that changing the parameter $\omega$ implicitly affects how much
-restoration is done. For an image $u$, decreasing $\omega$ will, all
+that changing the parameter $\omega$ implicitly affects the amount of
+regularization applied. For an image $u$, decreasing $\omega$ will, all
else being equal, decrease the lowest eigenvalue of $A(x)$ and in turn
decrease the anisotropic total variation $\TVA(u)$.
$\Sigma(x)$ will be the identity matrix, so any orthogonal matrix will
suffice for $U(x)$.
+Note that the eigenvalues of $S_\rho$ are continuous, and so are the
+eigenvectors (ignoring their sign) except possibly when $\lambda_1 =
+\lambda_2$. Thus $A$ is also continuous except possibly in these points.
+When $\lambda_1 = \lambda_2$ however, the eigenvalues of $A$ will both
+be 1, and $A$ is the identity matrix. Thus we can argue that if
+$S_\rho(x) \to \lambda I$ then $A(x) \to I$ and $A$ is continuous in all
+of $\Omega$.
+
See \cite{weickert1999coherence} for a different tensor construction,
made to enhance flow structures in the image, relevant in for example
fingerprint analysis.
-\fixme{argue that $A$ is continuous?}
-
-\begin{figure}
- \input{fig/tensor_viz}
-\end{figure}
-
\section{Well-posedness}
\fixme{Rating: 7/10}
dx
\end{equation}
can be thought of as---and has the properties of---a semi-norm, and is
-therefore convex. The sum of the two is thus strictly convex, which,
-given the existence of a minimizer, implies uniqueness.
+therefore convex. The sum of the fidelity and regularization terms is
+thus strictly convex, which, given the existence of a minimizer, implies
+uniqueness.
\subsection{Coercivity}
Coercivity relates to how the functional behaves when the norm of the
image $u$ tends to infinity. What we need in order to conclude with existence
-is sequential coercivity. Thus we need all level sets $F^\alpha = \{ u \in
-L^2(\Omega) : F(u) \leq \alpha \}$ to be \emph{sequentially
+is weak sequential coercivity. Thus we need all level sets $F^\alpha = \{ u \in
+L^2(\Omega) : F(u) \leq \alpha \}$ to be \emph{weakly sequentially
pre-compact}, meaning that all sequences in the set contain a
-subsequence converging to an element of the closure of the set.
+subsequence weakly converging to an element of the closure of the set.
It is obvious from the fidelity term that for some fixed $f \in
L^2(\Omega)$, if $\lnorm{u}_{L^2} \to \infty$ then $F(u) \to \infty$.
This implies that all the level sets $F^\alpha$ are bounded.
-Since $L^2(\Omega)$ is a Hilbert space, all bounded sequences contain a
+Since $L^2(\Omega)$ is a Hilbert space all bounded sequences contain a
weakly convergent subsequence. Thus all the level sets $F^\alpha$ are
-weakly sequentially pre-compact, since all sequences in these sets have
-a subsequence weakly converging to a point in the closure of the set.
+weakly sequentially pre-compact.
\begin{figure}
\input{fig/lower_semicont}
\label{eq:positive_int}
\end{equation}
Note that \eqref{eq:positive_int} only holds for non-negative images,
-which complicates the following proof a little.
+which complicates the proof of the anisotropic coarea formula a little.
\begin{figure}
\input{fig/eta_r}
\end{cases}
\end{aligned}
\end{equation}
- visualized in Figure~\ref{fig:eta_r} and \ref{fig:eta_r_diff}. By
+ visualized in Figure~\ref{fig:eta_r_both}. By
composing the function $\eta_r$ with our image $u$ and using Green's
formula, for example from \cite[Corollary
9.32]{grasmair2010anisotropic} we obtain
for all vector fields $\xi \in C_c^\infty(\Omega, \mathbb{R}^2)$.
The measure of $\{ x : u(x) = \lambda \text{ and } \nabla u(x) \neq
0\}$ is zero for all $\lambda$ following from \cite[Corollary \rom{1},
- Section 3.1.2]{evans1991measure}. Thus we can ignore the sets $\{ u
- = s\}$ and $\{ u = s + r \}$.
+ Section 3.1.2]{evans1991measure}, thus we can ignore the sets $\{ u
+ = s\}$ and $\{ u = s + r \}$ in the integrals.
Assuming that $\norm{\xi}_A^* \leq 1$ we obtain from \eqref{eq:mdef}
and \eqref{eq:eta_greens} that
\begin{equation}
As the limit when $r \to 0$ of the left-hand side exists almost
everywhere, suppose it exists at $s \in \mathbb{R}$. The integrand
on the right-hand side $-\eta_r(u) \diver \xi$ approaches $-u^s
- \diver \xi$ pointwise almost everywhere. We apply
- Lebesgue's dominated convergence theorem
- using that $\abs{\eta_r(u) \diver \xi} \leq
- \abs{u^s \diver \xi}$ and that $\xi \in C^\infty_c(\Omega,
- \mathbb{R}^2)$ is bounded by the extreme value theorem. Thus
+ \diver \xi$ pointwise almost everywhere.
+ We also have the bound
+ $\abs{\eta_r(u) \diver \xi} \leq
+ \abs{u^s \diver \xi}$ and know that
+ $\xi \in
+ C^\infty_c(\Omega, \mathbb{R}^2)$ is bounded following from the extreme value
+ theorem. Thus we can apply
+ Lebesgue's dominated convergence theorem, giving that
$\abs{u^s \diver \xi}$ is integrable and
\begin{equation}
\int_\Omega -\eta_r(u)
\end{equation}
As this holds for any
$\norm{\xi}_A^* \leq 1$, we get from the extended total variation
- definition in \eqref{eq:extended_tv} that $m'(s) \geq \TVA(u_s)$
+ definition in \eqref{eq:extended_tv} that $m'(s) \geq \TVA(u^s)$
almost everywhere and conclude using \eqref{eq:tva_geq_mder} that
\begin{equation}
\TVA(u) \geq \int_{-\infty}^\infty m'(t) \, dt \geq
Combining the upper and lower bounds just proved, we have equality.
\end{proof}
This coarea formula is our first step in transforming the anisotropic
-total variation into an easily discretizable expression.
+total variation into an easily discretizable expression. It allows us to
+consider each level $\lambda$ separately when calculating the total
+variation.
The anisotropic total variation of the thresholded images occurring in
the anisotropic coarea formula is very much related to the size of the
\PerA(\{u > \lambda\}; \Omega) \, d\lambda.
\label{eq:functional_per}
\end{equation}
-The transformation is motivated by our upcoming Cauchy--Crofton
-integration formula, and discretization, where an approximation of the
+The transformation is motivated by our upcoming anisotropic Cauchy--Crofton
+integration formula, and the discretization, where an approximation of the
perimeter will be computed using a graph cut machinery.
\section{Cauchy--Crofton formulas}
category often referred to as \emph{Cauchy--Crofton style formulas}, and
give ways to measure geometric objects using the set of all lines in the
plane. The formulas presented here will give a way to measure the length
-of a curve by counting the times it intersects line in the set of all
+of a curve by counting the times it intersects lines in the set of all
lines. The first formula will be for the isotropic case, and we will
-use this to prove the anisotropic formula following it.
+use it to prove the anisotropic formula following it.
We write $\mathcal{L}$ for the set of all lines in the plane,
and parametrize them as shown in Figure~\ref{fig:line_param}. A
\label{thm:euclidean_cauchy_crofton}
\end{theorem}
\begin{proof}
- See \cite[Theorem 3, Section 1-7.]{do1976differential}.
+ See \cite[Theorem 3, Section 1-7]{do1976differential}.
\end{proof}
If our space is equipped with a metric tensor $M(x)$ such that the inner
product of two vectors $a$ and $b$ in a point $x$ is calculated as
perimeter calculation. The set of lines $\mathcal{L}$ is then discretized
in a reasonable way, and the length of the curve $C$ can be approximated
by a sum over all these lines.
-\begin{theorem}[The Riemannian Cauchy--Crofton formula]
+\begin{theorem}[The anisotropic Cauchy--Crofton formula]
Assume that our space $\Omega$ is equipped with a continuous
positive definite metric
tensor $M(x)$, whose eigenvalues are bounded by $0 < k \leq
\end{equation}
\label{thm:riemannian_cauchy_crofton}
\end{theorem}
-\begin{proof}[Proof of the Riemannian Cauchy--Crofton formula]
+\begin{proof}[Proof of the anisotropic Cauchy--Crofton formula]
Assume first that our space is equipped with a constant metric
tensor $M$. The length of a curve in this space can be
calculated by transforming the curve and applying the Euclidean
\nu)^{\sfrac{3}{2}}}.
\end{equation}
We are interested in the Jacobian of the whole transformation
- $J_{\Sigma^2}(U^T \ell_{\theta, \rho})$, so
+ $J_{\Sigma^2}(U^T \ell_{\nu, \rho})$, so
all that is left to do is insert $U^T \ell_{\nu, \rho}$ to obtain
\begin{equation}
\abs{J_M(\ell_{\nu, \rho})} = \abs{J_{\Sigma^2}(U^T \ell_{\nu, \rho})}
\right)^{\sfrac{3}{2}}} \, d\mathcal{L} (\ell_{\nu, \rho}).
\label{eq:riemannian_const_m}
\end{equation}
- We now argue that the similar formula in
+ Further we argue that the similar formula in
\eqref{eq:riemannian_cauchy_crofton} holds for a non-constant but
continuous metric tensor $M(x)$. By partitioning the domain into
disjoint sets $U_i$ such that $\Omega = \cup_i U_i$, we make a
\xi^T M_\pi(x) \xi \leq K \, \norm{\xi}^2, \quad \forall \xi.
\end{equation}
Thus the integrand of \eqref{eq:piecewise_curve_length} is bounded
- by $g(t) = \sqrt{K \cdot \dot{C}(t)^T \dot{C}(t)}$. We know that
+ by $g(t) = (K \cdot \dot{C}(t)^T \dot{C}(t))^{\sfrac{1}{2}}$. We know that
$g(t)$ is integrable as its integral is exactly $\sqrt{K} \abs{C}$
and we have assumed that the curve is of finite length. This means
we can apply Lebesgue's dominated convergence theorem to see that
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}.
+each level set of the image. We will now see how the anisotropic
+Cauchy--Crofton formula can help us calculate the perimeters of the
+level sets.
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
\emph{normal vector} around the boundary, while the anisotropic curve
length in \eqref{eq:riemannian_length} is the integral of the norm of
the \emph{tangent vector} of the curve. Thus a 90\textdegree{} rotation
-is needed.
+separates the two.
If $P$ is a 90\textdegree{} rotation matrix we have
\begin{equation}
\label{eq:per_to_length2}
\end{equation}
-Now we make sure that all the assumptions of the Riemannian
+Now we make sure that all the assumptions of the anisotropic
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}.
Section~\ref{sec:anisotropy_tensor}
\begin{equation}
S_\rho(x) = \left(K_{\rho} * \left( \nabla f_{\sigma}
- \otimes \nabla f_{\sigma} \right)\right)(x),
+ \otimes \nabla f_{\sigma} \right)\right)(x).
\end{equation}
Because of the convolutions with the Gaussian function, this is a smooth
map from $\bar{\Omega}$ to $\mathbb{R}^{2 \times 2}$. As we
the structure tensor $S_\rho(x)$. The extreme value theorem states that
a continuous real-valued function on a nonempty compact
space is bounded above. Thus the eigenvalues $s_1 \geq s_2$ of $S_\rho(x)$
-are bounded above and by the construction in
+are bounded above. Moreover, by the construction in
\eqref{eq:sigma_construction}, there exists uniform bound $k$ such that
-the smallest eigenvalue $\sigma_2$ of the anisotropy tensor $A(x)$ is bounded away
+the smallest eigenvalue $\sigma_1$ of the anisotropy tensor $A(x)$ is bounded away
from zero, as
\begin{equation}
- \sigma_2 = \left(1
+ \sigma_1 = \left(1
+ \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$ is continuous and positive
-definite and thus the curve length calculation in
+definite with bounded eigenvalues and thus the curve length calculation in
\eqref{eq:per_to_length2} fulfills all the assumptions of the
anisotropic Cauchy--Crofton formula in
Theorem~\ref{thm:riemannian_cauchy_crofton}.
Within the restrictions that these theorems put on the tensor $M(x)$, we
have chosen a construction where one eigenvalue is always 1, while the
other varies from 1 in smooth areas towards 0 around edges, with the
-corresponding eigenvector perpendicular to the edge. \fixme{meh}
+corresponding eigenvector perpendicular to the edge.
As mentioned previously, these transformations are motivated by the
discretization that follows in the next chapter.
Assume that our discrete images are given on a uniform grid $\mathcal{G}$,
where each grid point is called a pixel. The image is a function giving
each pixel a value in the set of levels $\mathcal{P} = \{0, \hdots,
-L-1\}$. This is a reasonable assumption for digital grayscale
-images.
-
+L-1\}$, such that $u : \mathcal{G} \to \mathcal{P}$. This is a
+reasonable assumption for digital grayscale images.
Thus, when discretizing the functional in \eqref{eq:final_cont_rep}, we
have to consider that our images now have both discrete domain and
co-domain.
\end{aligned}
\end{equation}
The set of lines $\mathcal{L}$ has been discretized into the set
-$\mathcal{L}_D$. Note that we are approximating the length of the
-\emph{differentiable} curve $C$ in $\Omega$. Being a difference in the
+$\mathcal{L}_D$. Note that $C$ is still a differentiable curve, not yet
+discretized.
+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
be the average of the distance to the two neighboring line families as
shown in Figure~\ref{fig:line_neigh}, and thus also depends on $\phi$.
-In the discrete setting our domain $\Omega$ is discretized as a uniform
-grid $\mathcal{G}$. Our image is then reduced to a function $u :
-\mathcal{G} \to \mathcal{P}$.
-Moreover, the level sets $\{ u > \lambda\}$ will be functions taking the
-value of 0 or 1 on this grid, as shown in Figure \fixme{ref}.
-
The choice of our discrete set of lines $\mathcal{L}_D$ is important, as
it will decide the accuracy of our approximation.
We need some sensible restrictions on
the set $\mathcal{L}_D$ to simplify the further discussion. All lines
intersect at least two grid points, and from the periodicity of our grid
they thus intersect an infinite number of grid points. This puts some
-restrictions on the angles we can choose. For each angle, we include all
+restrictions on the angles we can choose. For each angle included, we include all
possible lines of that family, meaning there are no grid points without
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.
+shown in Figure~\ref{fig:line_neigh}. We write $\mathcal{N}(x)$ for the
+neighborhood of grid point $x$. Thus if there is an edge $e_{ab}$ in our
+discretization $\mathcal{L}_D$, then $a \in \mathcal{N}(b)$. Extending the edges shown in
+Figure~\ref{fig:line_neigh} 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$.
calculations of these points will not fit into our graph cut framework
later, and thus for an edge $e$ we will only consider the question of
``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
+terminals of $e$ lie on different sides of the curve $C$. This
approximation is exact for zero or one intersection points, but will, as
we see in Figure~\ref{fig:curve_edge}, be wrong when we have more.
\label{lem:delta_rho}
\end{lemma}
\begin{proof}
- Consider a line $\ell$ intersecting the point $(p, q) \in
- \mathbb{Z}^2$ in the grid.
+ Consider a line $\ell$ intersecting the point in the grid given by
+ the indices $(p, q) \in
+ \mathbb{Z}^2$.
The distance $\Delta \rho$ from this line $\ell$ to the neighboring
lines can then be calculated as a minimum over the distance to all
other grid points.
\end{equation}
where the sum is over all edges crossing the curve.
-The curve length we initially wanted to calculate was the perimeter
+The curve length we initially wanted to calculate was the level set perimeter
$\PerA(\{u > \lambda\}; \Omega)$ in \eqref{eq:per_approx1}. To find
-edges that crosses the curve, we identify the edges that have one
+edges that crosses this boundary curve, we identify the edges that have one
terminal inside the level set, and the other outside.
Thus we rewrite the sum over $e \cap C$ such that
\begin{equation}
The absolute value $\abs{u^\lambda_a - u^\lambda_b}$ is 1 if one of
$a$ and $b$ lie inside the level set and the other lies outside, and
0 otherwise.
-
In other words the absolute value is one if $e_{ab}$ crosses the
perimeter of $\{ u > \lambda\}$ an odd number of times, and zero
otherwise.
$F_\lambda$. However, it is not guaranteed that the obtained thresholded
images $u^\lambda$ can be combined to make an output image $u$. They
were defined as $u^\lambda = \idfun_{u > \lambda}$, so we need them to
-be monotonically decreasing \fixme{(?)} in increasing level values, i.e.\
+be monotonically decreasing in increasing level values, i.e.\
\begin{equation}
u_x^\lambda \geq u_x^\mu, \quad \forall \lambda \leq \mu, \quad
\forall x \in \mathcal{G}.
images minimizing each level, \emph{while guaranteeing that they meet
this requirement.}
-\fixme{maybe with a $w_{xy}$ definition}
-
\subsubsection{Consistency}
\fixme{Rating: 6/10}
in \eqref{eq:tensor_approx} converges to the tensor in the intersection
point.
-For each $\phi$ parameter, our discretization in the $\rho$ dimension
-can be regarded as a midpoint rule as shown in
+Consider now the discretization in $\rho$. 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.
+part of the discretization is consistent.
\begin{figure}
\input{fig/line_midpoint}
\end{figure}
We include all grid points inside the square in our neighborhood, except
-for multiple points that lie on the same line from the origin. Thus if
+for multiple points that lie on the same line from the origin. If
two or more grid points lie on the same line, we include only the one
closest to the origin. This implies that for each grid point along the
outer edge of this square, we include in our neighborhood a grid point
when one of $\phi_k$ or $\phi_{k-1}$ is horizontal or vertical, shown in
Figure~\ref{fig:square_cons} as the angle $a$.
Thus the
-largest $\Delta \phi_k$ will then be around the horizontal or
-vertical edge.
-These angles can be calculated to be
+largest $\Delta \phi_k$ will be when $\phi_k = \sfrac{m \pi}{2}$ for $m
+\in \mathbb{Z}$, so
+around vertical and horizontal edges.
+The supremum can then be calculated to be
\begin{equation}
- \sup \Delta \phi_k = \frac{2}{2} \arctan \frac{1/n}{n/2} = \arctan
+ \sup \Delta \phi_k = 2 \cdot \sup \frac{\phi_{k+1} - \phi_k}{2} = \arctan \frac{1/n}{n/2} = \arctan
\frac{2}{n^2} \to 0.
\end{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
-$\delta^2 = \Delta \rho \norm{e}$ and the fact that $\norm{e} \geq
-\delta$. Thus for the inter-line distance $\Delta \rho$ we have
+$\delta^2 = \Delta \rho \norm{e}$ and $\norm{e} \geq
+\delta$. Thus for the inter-line distance $\Delta \rho_k$ we have
\begin{equation}
- \sup \Delta \rho = \sup \frac{\delta^2}{\norm{e}} \leq
+ \sup \Delta \rho_k = \sup \frac{\delta^2}{\norm{e}} \leq
\frac{\delta^2}{\delta} = \delta \to 0.
\end{equation}
-And thus the approximation has been shown to be equivalent to
+Hence the approximation has been shown to be equivalent to
well-known, and consistent integral approximations, where the
summand converges to the integrand, and the differences $\Delta \phi$
and $\Delta \rho$ go to zero. Thus the perimeter approximation in
In this section we will look at how these graphs are constructed such
that their minimum cuts correspond to the minimizers of the functional
-$F^\lambda$. The description is taken from my project work
+$F^\lambda$. The description is taken with some small adjustments from my project work
\cite{project}, and is included here for completeness.
\subsection{Graphs}
will let $c(u, v) = 0$ for any pair $(u, v) \not\in E$, and we do not
allow self-loops in our graph. When a directed graph $G$ is equipped
with capacity function $c$, one might call it a capacitated graph or a
-graph, but as all our graphs will be capacitated from this point, we
-will just call them graphs and we write $G = (V, E, c)$.
+flow network, but as all our graphs will be capacitated from this point,
+we will just call them graphs and we write $G = (V, E, c)$.
There are two special vertices in the graph, the source $s$ and the
sink $t$. Contrary to other vertices, which can neither produce nor
We will construct a graph in such a way that if a vertex
$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
+= 1$, and if it ends up on the sink side we set $u^\lambda_x = 0$, as in
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.
\begin{table}[b]
\centering
- \caption[Different combinations for the term
- $F_\lambda^x(u^\lambda_x)$]{Each row represents one of the two possible values of
+ \caption[Graph construction for $F_\lambda^x(u^\lambda_x)$]{Each row
+ represents one of the two possible values of
$u^\lambda_x \in \{0,1\}$. The functional $F_\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 \subref{fig:norm_subgraph_pos} cut cap. & Graph
\subref{fig:norm_subgraph_neg} cut cap. \\
\hline
- $0$ & $0$ & $(\{s, u^\lambda_x\}, \{t\})$ & $0$ & $-F_\lambda^x(1)$ \\
- $1$ & $F_\lambda^x(1)$ & $(\{s\}, \{u^\lambda_x, t\})$ & $F_\lambda^x(1)$ & $0$
+ $0$ & $0$ & $(\{s\}, \{u^\lambda_x, t\})$ & $0$ &
+ $-F_\lambda^x(1)$ \\
+ $1$ & $F_\lambda^x(1)$ & $(\{s, u^\lambda_x\}, \{t\})$ &
+ $F_\lambda^x(1)$ & $0$
\end{tabular}
\label{tab:fid_energy}
\end{table}
\label{eq:fidelity_energy1}
\end{align}
where $F_\lambda^x(1)$ might be positive or negative depending on
-$\lambda$ and the pixel value $v_x$.
+$\lambda$ and the pixel value $f_x$.
Figure~\ref{fig:norm_subgraph} shows how graphs can be constructed to
represent this part of functional. The construction differs
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.
+See Table~\ref{tab:neigh_energy} for an overview of how the values
+of $(u^\lambda_x, u^\lambda_y)$ corresponds to cuts in the graph.
\begin{table}[t]
\centering
- \caption[Different combinations for the term $F^{x,y}(u^\lambda_x,
+ \caption[Graph construction for $F^{x,y}(u^\lambda_x,
u^\lambda_y)$]{An overview of the possible configurations of the variables
in the term $F^{x,y}(u^\lambda_x, u^\lambda_y)$. For each
configuration the corresponding functional value and the cut yielding this
& Alt.\ \subref{fig:neigh_subgraph_alt1} cut cap. & Alt.\
\subref{fig:neigh_subgraph_alt2} cut cap. \\
\hline
- $(0, 0)$ & $0$ & $(\{s, u^\lambda_x, u^\lambda_y\}, \{t\})$ & $w_{xy}$ & $0$ \\
- $(0, 1)$ & $w_{xy}$ & $(\{s, u^\lambda_x\}, \{u^\lambda_y, t\})$ & $2w_{xy}$ & $w_{xy}$ \\
- $(1, 0)$ & $w_{xy}$ & $(\{s, u^\lambda_y\}, \{u^\lambda_x, t\})$ & $2w_{xy}$ & $w_{xy}$ \\
- $(1, 1)$ & $0$ & $(\{s\}, \{u^\lambda_x, u^\lambda_y, t\})$ & $w_{xy}$ & $0$
+ $(0, 0)$ & $0$ & $(\{s\}, \{u^\lambda_x, u^\lambda_y, t\})$
+ & $w_{xy}$ & $0$ \\
+ $(0, 1)$ & $w_{xy}$ & $(\{s, u^\lambda_y\}, \{u^\lambda_x, t\})$ & $2w_{xy}$ & $w_{xy}$ \\
+ $(1, 0)$ & $w_{xy}$ & $(\{s, u^\lambda_x\}, \{u^\lambda_y, t\})$ & $2w_{xy}$ & $w_{xy}$ \\
+ $(1, 1)$ & $0$ & $(\{s, u^\lambda_x, u^\lambda_y\}, \{t\})$ & $w_{xy}$ & $0$
\end{tabular}
\label{tab:neigh_energy}
\end{table}
}
}
- \draw[blue] (0,0) rectangle (1,1);
- \draw[blue] (1,0) rectangle (2,1);
- \draw[blue] (2,0) rectangle (3,1);
- \draw[blue] (3,0) rectangle (4,1);
- \draw[blue] (4,0) rectangle (5,1);
+ \draw[blue, thick] (0,0) rectangle (1,1);
+ \draw[blue, thick] (1,0) rectangle (2,1);
+ \draw[blue, thick] (2,0) rectangle (3,1);
+ \draw[blue, thick] (3,0) rectangle (4,1);
+ \draw[blue, thick] (4,0) rectangle (5,1);
- \draw[blue] (0,1) rectangle (1,2);
- \draw[blue] (1,1) rectangle (2,2);
- \draw[blue] (2,1) rectangle (3,2);
- \draw[blue] (3,1) rectangle (4,2);
- \draw[blue] (4,1) rectangle (5,2);
+ \draw[blue, thick] (0,1) rectangle (1,2);
+ \draw[blue, thick] (1,1) rectangle (2,2);
+ \draw[blue, thick] (2,1) rectangle (3,2);
+ \draw[blue, thick] (3,1) rectangle (4,2);
+ \draw[blue, thick] (4,1) rectangle (5,2);
- \draw[blue] (0,2) rectangle (1,3);
- \draw[blue] (1,2) rectangle (2,3);
- \draw[blue] (2,2) rectangle (3,3);
- \draw[blue] (3,2) rectangle (4,3);
- \draw[blue] (4,2) rectangle (5,3);
+ \draw[blue, thick] (0,2) rectangle (1,3);
+ \draw[blue, thick] (1,2) rectangle (2,3);
+ \draw[blue, thick] (2,2) rectangle (3,3);
+ \draw[blue, thick] (3,2) rectangle (4,3);
+ \draw[blue, thick] (4,2) rectangle (5,3);
\clip (0,0) rectangle (5,3);
- \path[nedge=4cm,red] (0,3) -- (2,2);
- \path[nedge=4cm,red] (0,2) -- (2,1);
- \path[nedge=4cm,red] (0,1) -- (2,0);
+ \path[nedge=4cm,red, thick] (0,3) -- (2,2);
+ \path[nedge=4cm,red, thick] (0,2) -- (2,1);
+ \path[nedge=4cm,red, thick] (0,1) -- (2,0);
- \path[nedge=4cm,red] (-1,1) -- (1,0);
+ \path[nedge=4cm,red, thick] (-1,1) -- (1,0);
- \path[nedge=4cm,red] (1,3) -- (3,2);
- \path[nedge=4cm,red] (1,2) -- (3,1);
- \path[nedge=4cm,red] (1,1) -- (3,0);
+ \path[nedge=4cm,red, thick] (1,3) -- (3,2);
+ \path[nedge=4cm,red, thick] (1,2) -- (3,1);
+ \path[nedge=4cm,red, thick] (1,1) -- (3,0);
- \path[nedge=4cm,red] (2,3) -- (4,2);
- \path[nedge=4cm,red] (3,3) -- (5,2);
- \path[nedge=4cm,red] (4,3) -- (6,2);
+ \path[nedge=4cm,red, thick] (2,3) -- (4,2);
+ \path[nedge=4cm,red, thick] (3,3) -- (5,2);
+ \path[nedge=4cm,red, thick] (4,3) -- (6,2);
- \path (5,1) [edge,-,red] -- ($(6,0)!(5,1)!(4,1)$);
- \path (4,1) [edge,-,red] -- ($(5,0)!(4,1)!(3,1)$);
- \path (3,1) [edge,-,red] -- ($(4,0)!(3,1)!(2,1)$);
- \path (2,1) [edge,-,red] -- ($(3,0)!(2,1)!(1,1)$);
- \path (1,1) [edge,-,red] -- ($(2,0)!(1,1)!(0,1)$);
+ \path (5,1) [edge,-,red, thick] -- ($(6,0)!(5,1)!(4,1)$);
+ \path (4,1) [edge,-,red, thick] -- ($(5,0)!(4,1)!(3,1)$);
+ \path (3,1) [edge,-,red, thick] -- ($(4,0)!(3,1)!(2,1)$);
+ \path (2,1) [edge,-,red, thick] -- ($(3,0)!(2,1)!(1,1)$);
+ \path (1,1) [edge,-,red, thick] -- ($(2,0)!(1,1)!(0,1)$);
- \path (5,2) [edge,-,red] -- ($(6,1)!(5,2)!(4,2)$);
- \path (4,2) [edge,-,red] -- ($(5,1)!(4,2)!(3,2)$);
- \path (3,2) [edge,-,red] -- ($(4,1)!(3,2)!(2,2)$);
- \path (2,2) [edge,-,red] -- ($(3,1)!(2,2)!(1,2)$);
- \path (1,2) [edge,-,red] -- ($(2,1)!(1,2)!(0,2)$);
+ \path (5,2) [edge,-,red, thick] -- ($(6,1)!(5,2)!(4,2)$);
+ \path (4,2) [edge,-,red, thick] -- ($(5,1)!(4,2)!(3,2)$);
+ \path (3,2) [edge,-,red, thick] -- ($(4,1)!(3,2)!(2,2)$);
+ \path (2,2) [edge,-,red, thick] -- ($(3,1)!(2,2)!(1,2)$);
+ \path (1,2) [edge,-,red, thick] -- ($(2,1)!(1,2)!(0,2)$);
- \path (5,3) [edge,-,red] -- ($(6,2)!(5,3)!(4,3)$);
- \path (4,3) [edge,-,red] -- ($(5,2)!(4,3)!(3,3)$);
- \path (3,3) [edge,-,red] -- ($(4,2)!(3,3)!(2,3)$);
- \path (2,3) [edge,-,red] -- ($(3,2)!(2,3)!(1,3)$);
- \path (1,3) [edge,-,red] -- ($(2,2)!(1,3)!(0,3)$);
+ \path (5,3) [edge,-,red, thick] -- ($(6,2)!(5,3)!(4,3)$);
+ \path (4,3) [edge,-,red, thick] -- ($(5,2)!(4,3)!(3,3)$);
+ \path (3,3) [edge,-,red, thick] -- ($(4,2)!(3,3)!(2,3)$);
+ \path (2,3) [edge,-,red, thick] -- ($(3,2)!(2,3)!(1,3)$);
+ \path (1,3) [edge,-,red, thick] -- ($(2,2)!(1,3)!(0,3)$);
\end{tikzpicture}
\caption[A visual argument for $\delta^2 = \Delta \rho \norm{e}$]{
A visual argument showing that $\delta^2 = \Delta \rho \norm{e}$. If
extended to the whole plane, there will be the same amount of blue
- squares as red rectangles, so their areas must be equal.
+ squares as red rectangles as each grid point is the upper left
+ corner of both a blue and a red rectangle. Thus their areas must be equal.
}
\label{fig:area_proof}
\node[vertex] (u) at (0, -1) {$u^\lambda_x$};
\node[vertex] (t) at (0, -2) {t};
- \path[edge] (s) -- node[weight noslope] {$\max\{F_{L-1}^x(1), 0\}$} (u);
- \path[edge] (u) -- node[weight noslope] {$\max\{F_{L-1}^x(1), 0\} - F_\lambda^x(1)$} (t);
+ \path[edge] (s) -- node[weight noslope] {$\max\{F_{L-1}^x(1), 0\} - F_\lambda^x(1)$} (u);
+ \path[edge] (u) -- node[weight noslope] {$\max\{F_{L-1}^x(1), 0\}$} (t);
\end{tikzpicture}
\caption[Modified fidelity subgraph for the Boykov--Kolmogorov
algorithm]{%
% Title: glps_renderer figure
% Creator: GL2PS 1.3.8, (C) 1999-2012 C. Geuzaine
% For: Octave
-% CreationDate: Fri Jan 30 17:04:57 2015
+% CreationDate: Tue Feb 3 19:56:32 2015
\begin{pgfpicture}
\pgfsetlinewidth{0.01pt}
\color[rgb]{1.000000,1.000000,1.000000}
\node[tiny vertex] (\x\y) at (\x, \y) {};
}
}
- \node[tiny vertex,draw] (a) at (1, 2) {};
- \node[tiny vertex,draw] (b) at (1, 3) {};
- \node[tiny vertex,draw] (c) at (2, 2) {};
- \node[tiny vertex,draw] (d) at (2, 3) {};
- \node[tiny vertex,draw] (e) at (2, 4) {};
- \node[tiny vertex,draw] (f) at (3, 3) {};
- \node[tiny vertex,draw] (g) at (3, 1) {};
- \node[tiny vertex,draw] (h) at (2, 1) {};
- \node[tiny vertex,draw] (i) at (2, 0) {};
+ \node[tiny vertex,draw] (a) at (2, 2) {};
+ \node[tiny vertex,draw] (b) at (2, 3) {};
+ \node[tiny vertex,draw] (c) at (3, 2) {};
+ \node[tiny vertex,draw] (d) at (3, 3) {};
+ \node[tiny vertex,draw] (e) at (3, 4) {};
+ \node[tiny vertex,draw] (f) at (4, 3) {};
+ \node[tiny vertex,draw] (g) at (4, 1) {};
+ \node[tiny vertex,draw] (h) at (3, 1) {};
+ \node[tiny vertex,draw] (i) at (3, 0) {};
+
+ \node[tiny vertex,draw] (j) at (2, 1) {};
+ \node[tiny vertex,draw] (k) at (1, 2) {};
+ \node[tiny vertex,draw] (l) at (1, 3) {};
+ \node[tiny vertex,draw] (m) at (1, 4) {};
+ \node[tiny vertex,draw] (n) at (1, 0) {};
+ \node[tiny vertex,draw] (o) at (1, 1) {};
+ \node[tiny vertex,draw] (p) at (0, 3) {};
+ \node[tiny vertex,draw] (q) at (0, 1) {};
\path[edge,-] (a) -- (b);
\path[edge,-] (a) -- (c);
\path[edge,-] (a) -- (h);
\path[edge,-] (a) -- (i);
- \path[dash] (1, 2) -- +(35.785:4cm);
- \path[dash] (1, 2) -- +(54.215:4cm);
+ \path[edge,-] (a) -- (j);
+ \path[edge,-] (a) -- (k);
+ \path[edge,-] (a) -- (l);
+ \path[edge,-] (a) -- (m);
+ \path[edge,-] (a) -- (n);
+ \path[edge,-] (a) -- (o);
+ \path[edge,-] (a) -- (p);
+ \path[edge,-] (a) -- (q);
- \draw (1,2) +(35.785:3.1cm) arc (35.785:54.215:3.1cm);
- \path (1,2) ++(45:3.4cm) node {$\Delta \phi$};
+ \path[dash] (2, 2) -- +(35.785:4cm);
+ \path[dash] (2, 2) -- +(54.215:4cm);
+
+ \draw (2, 2) +(35.785:3.1cm) arc (35.785:54.215:3.1cm);
+ \path (2, 2) ++(45:3.4cm) node {$\Delta \phi$};
%\path (00) [edge, <->, anchor=center] -- node[anchor=south] {h} (10) ;
\begin{scope}
\clip (0, 0) rectangle (3, 2);
- \path[nedge=4cm] (0,5) -- (2,4);
- \path[nedge=4cm] (0,4) -- (2,3);
- \path[nedge=4cm] (0,3) -- (2,2);
- \path[nedge=4cm] (0,2) -- (2,1);
- \path[nedge=4cm] (0,1) -- (2,0);
+ \path[nedge=4cm, thick] (0,5) -- (2,4);
+ \path[nedge=4cm, thick] (0,4) -- (2,3);
+ \path[nedge=4cm, thick] (0,3) -- (2,2);
+ \path[nedge=4cm, thick] (0,2) -- (2,1);
+ \path[nedge=4cm, thick] (0,1) -- (2,0);
- \path[nedge=4cm] (-1,1) -- (1,0);
+ \path[nedge=4cm, thick] (-1,1) -- (1,0);
- \path[nedge=4cm] (1,5) -- (3,4);
- \path[nedge=4cm] (1,4) -- (3,3);
- \path[nedge=4cm] (1,3) -- (3,2);
- \path[nedge=4cm] (1,2) -- (3,1);
- \path[nedge=4cm] (1,1) -- (3,0);
+ \path[nedge=4cm, thick] (1,5) -- (3,4);
+ \path[nedge=4cm, thick] (1,4) -- (3,3);
+ \path[nedge=4cm, thick] (1,3) -- (3,2);
+ \path[nedge=4cm, thick] (1,2) -- (3,1);
+ \path[nedge=4cm, thick] (1,1) -- (3,0);
- \path[nedge=4cm] (2,5) -- (4,4);
- \path[nedge=4cm] (3,5) -- (5,4);
- \path[nedge=4cm] (4,5) -- (6,4);
+ \path[nedge=4cm, thick] (2,5) -- (4,4);
+ \path[nedge=4cm, thick] (3,5) -- (5,4);
+ \path[nedge=4cm, thick] (4,5) -- (6,4);
\end{scope}
- \path[dash] (2.5,1.5) -- (3.5,1);
- \path[dash] ($(3.5,0.5)!(2.5,1.5)!(1.5,1.5)$) -- ($(4.5,0.0)!(3.5,1.0)!(2.5,1.0)$);
+ \path[dash, thick] (2.5,1.5) -- (3.5,1);
+ \path[dash, thick] ($(3.5,0.5)!(2.5,1.5)!(1.5,1.5)$) -- ($(4.5,0.0)!(3.5,1.0)!(2.5,1.0)$);
- \path (3.5,1.0) [edge, <->, anchor=center] -- node[anchor=north
+ \path (3.5,1.0) [edge, thick, <->, anchor=center] -- node[anchor=north
west,yshift=+3.0pt,xshift=+1pt, fill=white, font=\small, inner sep=1pt]
{$\Delta \rho$}
($(4.5,0.0)!(3.5,1.0)!(2.5,1.0)$) ;
\centering
-\begin{subfigure}[t]{0.4\textwidth}
- \centering
- \begin{tikzpicture}
- \draw[->] (0,0) -- (4,0) node[right] {$\lambda$};
- \draw[->] (0,-2) -- (0,2) node[above] {$F_\lambda^x(1)$};
- \draw[line,domain=0:1] plot ({\x},{-1});
- \draw[line,domain=1:2,dashed] plot ({\x},{2*\x-3});
- \draw[line,domain=2:4] plot ({\x},{1});
+\begin{tikzpicture}
+ \draw[->] (0,0) -- (4,0) node[right] {$\lambda$};
+ \draw[->] (0,-2) -- (0,2) node[above] {$F_\lambda^x(1)$};
+ \draw[line,domain=0:4] plot ({\x},{.8*\x - 1.6});
- \draw (1,2pt) -- (1,-2pt) node[
- font=\small,
- anchor=south
- ]
- {$f_x-1$};
+ \draw (2.5,2pt) -- (2.5,-2pt) node[
+ font=\small,
+ anchor=north
+ ]
+ {$f_x$};
- \draw (2,2pt) -- (2,-2pt) node[
- font=\small,
- anchor=north
- ]
- {$f_x$};
+ \draw (1.5,2pt) -- (1.5,-2pt) node[
+ font=\small,
+ anchor=south
+ ]
+ {$f_x-1$};
- \end{tikzpicture}
- \caption{$L^1$ fidelity term.}
- \label{fig:l1_norm_evolution}
-\end{subfigure}
-~
-\begin{subfigure}[t]{0.4\textwidth}
- \centering
- \begin{tikzpicture}
- \draw[->] (0,0) -- (4,0) node[right] {$\lambda$};
- \draw[->] (0,-2) -- (0,2) node[above] {$F_\lambda^x(1)$};
- \draw[line,domain=0:4] plot ({\x},{.8*\x - 1.6});
-
- \draw (2.5,2pt) -- (2.5,-2pt) node[
- font=\small,
- anchor=north
- ]
- {$f_x$};
-
- \draw (1.5,2pt) -- (1.5,-2pt) node[
- font=\small,
- anchor=south
- ]
- {$f_x-1$};
-
- \end{tikzpicture}
- \caption{$L^2$ fidelity term.}
- \label{fig:l2_norm_evolution}
-\end{subfigure}
+\end{tikzpicture}
\caption[Fidelity energy $F_\lambda^x(1)$ as a function of $\lambda$]{%
- Two figures showing how the fidelity energy term $F_\lambda^x(1)$
+ A sketch of how the fidelity energy term $F_\lambda^x(1)$
in \eqref{eq:final_discretization} increases monotonically with $\lambda$.
}
\label{fig:norm_evolution}
\node[vertex] (u) at (0, -1) {$u^\lambda_x$};
\node[vertex] (t) at (0, -2) {t};
- \path[edge] (s) -- node[weight noslope] {$F_\lambda^x(1)$} (u);
+ \path[edge] (u) -- node[weight noslope] {$F_\lambda^x(1)$} (t);
\end{tikzpicture}
\caption{%
The graph when $F_\lambda^x(1) > 0$, with constant equal to 0.
\node[vertex] (u) at (0, -1) {$u^\lambda_x$};
\node[vertex] (t) at (0, -2) {t};
- \path[edge] (u) -- node[weight noslope] {$-F_\lambda^x(1)$} (t);
+ \path[edge] (s) -- node[weight noslope] {$-F_\lambda^x(1)$} (u);
\end{tikzpicture}
\caption{Graph when $F_\lambda^x(1) < 0$, with constant equal to
$-F_\lambda^x(1)$.
%the result section will focus on looking at how these parameters affect
%the restoration performance.
-\chapter{Methods in image restoration}
-
-\fixme{Rating: 7/10}
-
-There are numerous methods in image restoration, and since we do not
-have time nor space to discuss them all, this short overview will focus on the ones
-related to the anisotropic total variation method considered later in
-this thesis. See \cite{jahne} and \cite{aubert2006proc} for more
-background on image processing in general.
-
-In this chapter, and also in the rest of the thesis we will assume that
-we are given an image $f : \Omega \to \mathbb{R}$ where $\Omega$ is a
-rectangular, open domain. Because of limitations in the numerical method
-used, the codomain is $\mathbb{R}$ and we are thus restricted to
-monochrome, or grayscale images. Such images are produced in large
-numbers by for example ultrasound, X-ray and MRI machines.
-
-The space in which the image $f$ resides in will vary, but since we are
-looking at image restoration methods, we assume that it includes some
-kind of noise. Depending on the application and how the image is
-obtained, one might construct different models describing different
-types of noise.
-
-We will assume that the given image $f$ is a combination of an
-underlying, actual image $u^*$, and some noise $\delta$. The simplest
-model is additive noise where the assumption is that $f = u^* + \delta$.
-There is also multiplicative noise where $f = u^* \cdot \delta$. An
-other much seen noise type is salt and pepper noise, which is when black
-and white pixels randomly appear in the image, independently of the
-actual image $u^*$.
-
-These are only models, and in the real world the noise might be more
-complex, and even come from a combination of sources. Depending on the
-application, the goal might not even be to recover $u^*$, but rather to
-obtain an output which fulfills certain smoothness or regularity
-properties. In any case, we will continue using $f$ as the input and $u$
-as the output in the description of our methods.
-
-\section{Diffusion filtering}
-\label{sec:diffusion_filtering}
-
-\fixme{Rating: 8/10}
-
-Diffusion filtering is a broad group of filtering and restoration methods
-based on physical diffusion processes. The basic idea is to take the
-noisy image as the initial value of some diffusion process, and then let it
-evolve for some time. The best known method is probably the Gaussian
-filter or Gaussian blur, in which one convolves the image with the
-Gaussian function
-\begin{equation}
- K_\sigma(x,y) = \frac{1}{2\pi\sigma^2} \exp \left( - \frac{x^2 +
- y^2}{2\sigma^2} \right).
- \label{eq:gaussian_function}
-\end{equation}
-\nomenclature{$K_\sigma(x,y)$}{The Gaussian kernel.}%
-In the discrete setting where the image consists of a grid of pixels,
-the Gaussian blur amounts to calculating each pixel in the output image
-as a weighted average of its neighboring pixels in the input image.
-
-The Gaussian function happens to be the fundamental solution of the heat
-equation $\partial_t u = \Delta u$. Convolving $K_\sigma(x,y)$ with the
-original image $f$ is therefore equivalent to solving the heat equation
-with $f$ as initial value, until some time $T > 0$ depending on
-$\sigma$. Boundary conditions have to be specified of course, and one
-common choice is to symmetrically extend the image in $x$ and $y$
-directions, which corresponds to zero flux boundary conditions.
-
-By basic Fourier analysis it is possible to show that the Gaussian
-filter is a low-pass filter which attenuates high frequencies. Further
-theory can be found in Weickert's book on anisotropic diffusion
-\cite{weickert1998anisotropic}.
-
-The main concern with the Gaussian filter is that it will, in addition
-to smoothing out possible noise, remove details from the image, which
-motivates the next set of methods, where the amount of diffusion can
-vary for different parts of the image.
-
-\subsection{Non-linear diffusion filtering}
-
-\fixme{Rating: 8/10}
-
-In the theory of the heat equation one can introduce a \emph{thermal
-diffusivity} $\alpha$ such that the equation becomes
-\begin{equation}
- \left\{
- \begin{aligned}
- \partial_t u &= \diver \big( \alpha(\nabla u) \nabla u\big), \\
- u |_{t=0} &= f.
- \end{aligned}
- \right.
-\end{equation}
-\nomenclature{$\alpha(\nabla u)$}{Scalar thermal diffusivity.}%
-The thermal diffusivity $\alpha(\nabla u) = \alpha(x, \nabla u)$ is material dependent, and can also
-vary throughout the object. It specifies how well heat travels through
-the specific point in the object. We can make use of this in the image
-restoration context by controlling the diffusivity in different parts of
-the image, where our goal is to reduce noise without loosing image
-detail. Optimally, we would like there to be a lot of diffusion in
-smooth parts of the image, and not so much in areas with a lot of
-details.
-
-One much-studied non-linear diffusion equation is the Perona--Malik
-equation
-\begin{equation}
- \partial_t u = \diver \left( \frac{\nabla u}{1 +
- \frac{\abs{\nabla u}^2}{\lambda^2}} \right).
- \label{eq:perona_malik}
-\end{equation}
-The thermal diffusivity $\alpha(\nabla u) = (1 + \sfrac{\abs{\nabla u}^2}{
-\lambda^2})^{-1}$ varies from $1$ in smooth areas to $0$ as the
-norm of the gradient $\abs{\nabla u}$ grows.
-
-This particular form of the thermal diffusivity has been shown to be
-related to how brightness is perceived by the human visual system. The
-model has some theoretical problems related to well-posedness, for more
-information see \cite{weickert1998anisotropic}.
-
-A different kind of non-linear diffusion model is the total variation
-flow which can be formulated as
-\begin{equation}
- \partial_t u = \diver \frac{\nabla u}{\abs{\nabla u}},
-\end{equation}
-where the diffusivity has a similar effect of reducing the diffusion in
-areas of high variation. As the name suggests this model can be related
-to the total variation formulation presented later. One forward Euler
-time-step in the solution of this partial differential equation
-corresponds to the Euler--Lagrange equation of the variational
-formulation.
-
-Note that we follow Weickert's terminology when it comes to the
-distinction between non-linear and anisotropic diffusion methods. The
-Perona--Malik equation, and other diffusion equations with
-non-homogenous diffusivities, are often by others called anisotropic, as
-the diffusivity depends on the location. We will name these methods
-non-linear and spare the anisotropy term for the ``real'' anisotropic
-methods. These are methods where the diffusivity is a tensor, and thus
-both location and direction dependent.
-
-\subsection{Anisotropic diffusion}
-\label{sec:anisotropic_diffusion}
-
-The diffusivity is made directionally dependent by introducing a
-diffusion \emph{tensor} $A(u)$ such that the initial boundary value
-problem becomes
-\begin{equation}
- \left\{
- \begin{alignedat}{10}
- \partial_t u &{}={} \diver \big(A(u) \nabla u\big)
- &\quad \text{ on }& \Omega \times (0, \infty),\\
- u|_{t=0} &{}={} f
- &\text{ on }& \Omega,\\
- A(u) \nabla u \cdot \nu &{}={} 0
- &\text{ on }& \partial \Omega \times (0, \infty),
- \end{alignedat}
- \right.
- \label{eq:aniso_diff}
-\end{equation}
-\nomenclature{$A(u)$}{Thermal diffusivity tensor, or anisotropy
-tensor.}%
-where $\nu$ is the outer normal of $\Omega$. The tensor $A(u)$
-is constructed such as to diminish the effect of $\nabla
-u$ across what we believe to be edges in the image. This way, there will
-also be less diffusion through these edges. Weickert
-\cite{weickert1998anisotropic} suggests constructing $A(u)$ based on the
-edge estimator $\nabla u_\sigma$ where
-\nomenclature{$\nabla u_\sigma$}{Edge estimator.}%
-\begin{equation}
- u_\sigma := K_\sigma * \tilde{u}
-\end{equation}
-\nomenclature{$u_\sigma$}{Image $u$ smoothed with a Gaussian of
-parameter $\sigma$.}%
-and $\tilde{u}$ is an extension of $u$ from $\Omega$ to $\mathbb{R}^2$
-made by symmetrically extending $u$ across the edges.
-Assuming we are at an edge in the image, the
-direction of $\nabla u_\sigma$ should be perpendicular to the edge,
-while its magnitude will provide information on the steepness of the
-edge.
-
-To extract this information, and also to identify features on a larger
-scale, the \emph{structure tensor} is introduced
-\begin{equation}
- S_\rho(x) := K_\rho * (\nabla u_\sigma \otimes \nabla u_\sigma),
-\end{equation}
-\nomenclature{$S_\rho(x)$}{Structure tensor.}%
-where the convolution with the Gaussian function $K_\rho$ is done
-component-wise. The anisotropy tensor $A(u)$ can then be constructed
-based on the eigenvectors and eigenvalues of $S_\rho(x)$.
-The structure tensor and its properties will be discussed further when
-we introduce our anisotropic total variation functional.
-
-Assuming some
-smoothness, symmetry and uniform positive definiteness on $A(u)$ one can
-prove well-posedness, regularity and an extremum principle of the
-problem \eqref{eq:aniso_diff} as done in \cite{weickert1998anisotropic}.
-
-Even if the diffusivity tensor was introduced to reduce the amount of
-smoothing across edges, the solution of \eqref{eq:aniso_diff} will still be
-infinitely differentiable \cite{weickert1998anisotropic}, i.e.\ $u(T) \in
-C^\infty(\Omega)$ for $T > 0$. Thus we have no real discontinuities, and
-no real edges in our solution.
-
-Further, the anisotropic diffusion may introduce structure based on
-noise, when there really was no structure to begin with. This is a
-problem we hope to avoid in our anisotropic total variation method.
-
-\begin{figure}
- \centering
- \begin{subfigure}[t]{0.4\textwidth}
- \centering
- \includegraphics[width=\textwidth]{fig/fnoise.png}
- \caption{Noisy, $\sigma > 9000$}
- \end{subfigure}
- ~
- \begin{subfigure}[t]{0.4\textwidth}
- \centering
- \includegraphics[width=\textwidth]{fig/fanisodiff.png}
- \caption{Regularized, parameters: bleep bloop}
- \end{subfigure}
- \caption[bleapofsijdosij]{\fixme{UPDATE THIS}}
- \label{fig:ansi_diff}
-\end{figure}
-
-\section{Total variation}
-\label{sec:total_variation}
-
-\fixme{Rating: 8/10, numerical methods: 5/10}
-
-The total variation image restoration method is usually formulated as a
-minimization problem
-\begin{equation}
- \begin{gathered}
- \min_{u \in L^p(\Omega)} F(u), \\
- F(u) = \underbrace{\int_\Omega \abs{u - f}^p \, dx}_{\text{fidelity
- term}} {}+{} \beta \underbrace{\int_\Omega
- \abs{\nabla u} \, dx}_{\mathclap{\text{regularization
- term}}},
- \label{eq:first_min_presentation}
- \end{gathered}
-\end{equation}
-where $p$ is normally taken to be 1 or 2.
-The fidelity term penalizes images $u$ that are far from the original
-image $f$, and thus controls the \emph{fidelity} of our solution. The
-regularization term is the total variation, and minimizing it will
-reduce the variation and thus regularize the image. The $\beta$ parameter
-controls the strength of the regularization.
-
-Since we do not only want to consider differentiable images $u \in
-C^1(\Omega)$ for which the gradient exists, we introduce the total
-variation using the distributional derivative.
-\begin{definition}[Total variation]
- Given a function $u \in L^1(\Omega)$, the total variation of $u$,
- often written $\int_\Omega \abs{Du}\, dx$, where the $D$ is the
- gradient taken in the distributional sense, is
- \begin{equation}
- \TV(u)
- = \int_\Omega \abs{Du} \, dx
- = \sup \left\{ \int_\Omega u \cdot \diver \varphi \, dx :
- \varphi \in C^\infty_c\left(\Omega, \mathbb{R}^2\right),
- \norm{\varphi}_{L^\infty(\Omega)} \leq 1
- \right\}.
- \end{equation}
- The test functions $\varphi$ are taken from
- $C^\infty_c\left(\Omega, \mathbb{R}^2\right)$, the space of smooth
- functions from $\Omega$ to $\mathbb{R}^2$ with compact support in
- $\Omega$.
- \label{def:tv}
-\end{definition}
-\nomenclature{$\TV(u)$}{Total variation of image $u$.}%
-
-Note that since $\Omega$ is open and bounded, the test functions
-$\varphi$ vanish on the boundary of $\Omega$. Thus no variation is
-measured at the boundary. As this restoration method is the one
-which will be extended later in this thesis, we will look a little bit
-more deeply into the background and the numerical methods relating to
-it.
-
-As we are searching for an image with low total variation, it is useful
-to introduce the space of functions of bounded variation.
-\begin{definition}[Functions of bounded variation]
- The space of functions of bounded variation $\BV(\Omega)$ is the
- space of functions $u \in L^1(\Omega)$ for which the total variation
- is finite, i.e.,
- \begin{equation}
- \BV(\Omega) = \left\{ u \in L^1(\Omega) : \TV(u) < \infty
- \right\}.
- \end{equation}
-\end{definition}
-\nomenclature{$\BV(\Omega)$}{Functions of bounded variation in
-$\Omega$.}%
-Our optimization problem has thus become
-\begin{equation}
- \min_{u \in \BV(\Omega)} \int_\Omega \abs{u - v}^p \, dx + \beta \,
- \TV(u).
- \label{eq:second_min_presentation}
-\end{equation}
-
-As with any restoration method, the total variation method has its
-strengths and weaknesses. Its main strength is its ability to recover
-edges in the input image. The total variation of a section only takes
-the absolute change into account, and does not favor gradual changes
-like the diffusion methods.
-
-There is also a theoretical result stating that the set of edges in the
-solution $u$ is contained in the set of edges in the original image $f$,
-thus no new edges are created \cite{caselles2011total}. However, in the
-presence of noise, the method may introduce or rather ``find'' new edges
-that were not in the original image, since flat sections of zero
-variation are encouraged by the functional. This effect is called the
-\emph{stair-casing effect}, and can be seen in Figure~\ref{fig:grad_tv}
-where a noisy gradient has been restored using total variation
-restoration.
-
-\begin{figure}
- \centering
- \begin{subfigure}[t]{0.4\textwidth}
- \centering
- \includegraphics[width=\textwidth]{fig/factory/grad/n_q40.png}
- \caption{Noisy gradient}
- \end{subfigure}
- ~
- \begin{subfigure}[t]{0.4\textwidth}
- \centering
- \includegraphics[width=\textwidth]{fig/factory/grad/r_p2_n16_b5000_g1000000000_r5_s3.png}
- \caption{Total variation restoration}
- \end{subfigure}
- \caption[The stair-casing effect of regular total variation
- restoration]{%
- Although the original gradient was smooth, the total variation
- method manages to find structure in the noise, and create edges
- in the restored image.
- }
- \label{fig:grad_tv}
-\end{figure}
-
-Thin objects and corners may also suffer from contrast loss since
-bringing them closer to their surroundings reduces the
-total variation. An example of this is shown in
-Figure~\ref{fig:tv_example}, where a not particularly noisy fingerprint
-image has been strongly regularized. The original black and white levels
-have been brought closer to yield a lower total variation in the
-regularized image.
-
-\begin{figure}
- \centering
- \begin{subfigure}[t]{0.4\textwidth}
- \centering
- \includegraphics[width=\textwidth]{fig/finger.png}
- \end{subfigure}
- ~
- \begin{subfigure}[t]{0.4\textwidth}
- \centering
- \includegraphics[width=\textwidth]{fig/ftv.png}
- \end{subfigure}
- \caption[Heavily regularized fingerprint showing contrast loss]{%
- A fingerprint heavily regularized using the total variation
- method. The originally white and black ridges have been brought
- closer in value, to reduce the total variation.
- }
- \label{fig:tv_example}
-\end{figure}
-
-\subsection{Numerical methods}
-
-See \cite{chan2011numerical} for an overview of some of the numerical
-methods relating to total variation image restoration. Amongst others it
-describes some dual and primal-dual methods, as well as the graph cut
-approach we take in this thesis.
-
-\subsubsection{Graph cut approach}
-
-\fixme{Could probably be shortened a lot, but I thought it could be
-useful to at least give an idea how it is done.}
-
-Using graph cuts is the approach we will be taking later when
-considering the anisotropic total variation regularization, and it is
-therefore valuable to briefly look into how graph cuts are used in the
-case of regular total variation.
-
-In the theory of graphs, a graph cut is a set of edges that when
-removed will separate the graph into two disconnected parts. A minimum
-cut is a cut such that the sum of the weights is minimal. It
-has been shown that using special graph constructions, the graph cuts
-can be used to minimize a special type of energy functions.
-
-In the discrete setting our image consists of pixels, and is represented
-by a function $u : \mathcal{G} \to \mathcal{P}$ where $\mathcal{G}$ is a
-regular grid over $\Omega$, and $\mathcal{P} = \{0, \hdots, L-1\}$ is
-the discrete set of pixel values, or \emph{levels}.
-
-The idea of the graph cut approach is to decompose the minimization
-problem into one minimization problem for each level of the image, and
-then solve them separately before combining the results.
-
-For an image $u$ and a level $\lambda$ we denote the \emph{level set} by
-$\{ u > \lambda\}$, defined as the set $\{ x \in \Omega : u_x > \lambda
-\}$. The thresholded image $u^\lambda$, an indicator function, is then
-defined as
-\begin{equation}
- u^\lambda = \idfun_{u > \lambda}.
-\end{equation}
-Here, $\idfun_E$ signifies the characteristic function of the set $E$,
-\nomenclature{$\idfun_E$}{Characteristic function of the set $E$.}%
-the function which is equal to one in every point in $E$, and zero
-elsewhere. The graph cut we find later will for each level $\lambda$
-give us the thresholded image $u^\lambda$ which then can be combined to
-form the complete image $u$.
-
-Through careful manipulation of the continuous functional in
-\eqref{eq:second_min_presentation} it is possible to obtain a discrete
-functional decomposed as a sum over all the level values on the form
-\begin{equation}
- F(u) =
- \sum_{\lambda=0}^{L-2} \sum_x F^x_\lambda(u^\lambda_x)
- + \beta \sum_{\lambda = 0}^{L-2} \sum_{(x, y)} F^{x,y}(u^\lambda_x,
- u^\lambda_y)
- =: \sum_{\lambda=0}^{L-2} F_\lambda(u^\lambda)
- \label{eq:old_total_energy}
-\end{equation}
-where the sum over $(x, y)$ is over all pixel pairs $(x, y)$ in a neighbor
-relation, i.e.\ where the pixels are ``close'' to each other.
-
-The actual form of the functional, and the steps to construct it will be
-presented later, but hopefully, it does not seem too infeasible. And for
-those with some knowledge of graph theory, it might be valuable to
-consider the intuition behind the graph cut approach.
-
-When constructing the graph used to find the thresholded image
-$u^\lambda$, we have two special vertices, one representing the set $\{
-u > \lambda\}$, while the other represents the set $\{ u \leq
-\lambda\}$. The pixels are then connected to these vertices with a
-weight representing how strongly they are related to the corresponding
-set. This weight will be based on the value of $F^x_\lambda$.
-
-Then there are connections between pixels in a neighborhood relation,
-representing the energy $F^{x,y}$. Thus when finding a cut, we partition
-the pixels into the sets $\{ u > \lambda \}$ and $\{ u \leq \lambda \}$.
-And if in addition the cut is minimal, we know that the edges cut have
-minimal weight, and can prove that the $u^\lambda$ produced minimizes
-the functional in \eqref{eq:total_energy}.
-
\usepackage{algpseudocode}
\newtheorem{theorem}{Theorem}
-\newtheorem{proposition}{Proposition}
-\newtheorem{lemma}{Lemma}
-\newtheorem{definition}{Definition}
-
+\numberwithin{theorem}{chapter}
+\newtheorem{proposition}[theorem]{Proposition}
+\newtheorem{lemma}[theorem]{Lemma}
+\newtheorem{definition}[theorem]{Definition}
+
% Indeksgenerering er kommentert ut her. Ta bort prosenttegnene
% hvis du vil ha en indeks:
%\usepackage{makeidx}
% Inkluder alle de andre kildefilene
\input{introduction}
+\input{methods}
\input{continuous}
\input{discrete}
\input{maxflow}
this and the next section see how such a minimum cut can be found by
sending flow through the graph and trying to identify the
``bottleneck''. This chapter, except for the description of the
-Boykov--Kolmogorov algorithm is taken from my project work
-\cite{project} and is included here for completeness.
+Boykov--Kolmogorov algorithm is taken with some adjustments from my
+project work \cite{project} and is included here for completeness.
\section{Flow graphs}
%graph flow theory, stating that the capacity of a minimum $s$-$t$-cut
%in a graph, is equal to the maximum flow from the source to the sink.
-But how does this help us? We know that if we know the maximum flow
+But how does this help us? If we know the maximum flow
value, and we have an $s$-$t$-cut with capacity equal to the maximum
flow, we actually have a minimum cut. The question is then, how do we
find a maximum flow, and how do we find a minimum cut?
\section{Augmenting path algorithms}
+
The family of augmenting flow algorithms represent a popular approach to
the maximum flow problem. The idea is simply to look for paths from the
source to the sink that can carry additional flow, so-called augmenting
the source to the sink, and a maximum flow has been reached.
\subsection{Residual graph}
+
When further discussing approaches to solving the maximum flow problem
we will need the notion of a residual graph $G_f = (V_f, E_f, c_f)$,
which is derived from the original graph $G$ and contains the edges
\begin{equation}
e(u) = \sum_{v \in V} f(v, u) - \sum_{v \in V} f(u, v),
\end{equation}
-which represents the amount of flow which \emph{disappears} in vertex
+which represents the amount of flow which disappears in vertex
$u$. Equivalent to the preflow conservation constraint is stating that
$e(u) \geq 0$ for all vertices $u \in V - \{s,t\}$.
complexity of $O(\abs{V}^3)$, which can be very good if you have a dense
graph.
-Using highest level selection rule one always discharges the active
+Using the highest level selection rule one always discharges the active
vertex with the largest distance label. Goldberg and Tarjan state that
this rule also gives a complexity of $O(\abs{V}^3)$ while this bound is
improved to $O(\abs{V}^2 \sqrt{\abs{E}})$ in an article by Cheriyan and
relabeling procedure, as only the labels $d$ are changed. For the
labeling $d$ one has to verify that $d(u) \leq d(v) + 1$ for every edge
$(u,v) \in E_f$ in the residual network. If none, both, or only $v$ is
-relabeled, this is trivial. It is not possible that only $v$ would be
+relabeled, this is trivial. It is not possible that only $u$ would be
relabeled, as this would imply that $d(u) \geq d(v) + 2$ which is not a
valid initial labeling.
satisfied in the new graph.
Through this procedure we have easily created the graph for $\lambda =
-k-1$, and the distance labels remain the same. As these labels always
+k-1$, and the distance labels remain the same. As the distance labels always
increase monotonically, we have a head start compared to the case where
we reset the flow and labels.
\{ u \in V : d^\lambda(u) \geq N \} \subseteq \{ u \in V :
d^{\lambda-1}(u) \geq N \}.
\end{equation}
-For a pixel $x \in S$ we will set $u^\lambda_x = 0$, which together with
+For a pixel $x \in S$ we will set $u^\lambda_x = 1$, which together with
the previous inclusion property implies
\begin{equation}
- u^\lambda_x \geq u^{\lambda-1}_x
+ u^\lambda_x \leq u^{\lambda-1}_x
\end{equation}
for all $x \in \mathcal{G}$. Being equivalent with the inequality in
\eqref{eq:decreasing_sets}, this means our algorithm produces stackable
We then construct our output image $u$ by giving each pixel the value
\begin{equation}
- u_x = \min \{ \lambda \in \{0, \ldots, L-1\} : u^\lambda_x = 1 \}.
+ u_x = \min \{ \lambda \in \{0, \ldots, L-1\} : u^\lambda_x = 0 \}.
\end{equation}
This marks the end of the description of the implemented algorithm, but
-we will further discuss some possible improvements, and also look at the
-results when using the method on different kinds of noisy images.
+we will further discuss some possible improvements.
\subsection{Divide and conquer}
+
The possibility of re-using the graph between separate level
is a very nice property of the push-relabel algorithm, but there are
further room for improvements. Consider one pixel $x$ with value $u_x$,
\subsection{Graph reuse}
-As in the push-relabel algorithm described earlier, it is not necessary
+As with the push-relabel algorithm described earlier, it is not necessary
to completely restart the Boykov--Kolmogorov algorithm for every level
of the image, but we can reuse the graph and the trees $S$ and $T$ in
successive runs.
stay non-negative. The edge from $u_\lambda^x$ to $t$ is non-decreasing
with decreasing $\lambda$ parameter.
+\begin{table}[b]
+ \centering
+ \caption[Graph construction for $F_\lambda^x(u^\lambda_x)$ in the
+ Boykov--Kolmogorov algorithm]{
+ Each row represents one of the two possible values of
+ $u^\lambda_x \in \{0,1\}$. The functional $F_\lambda^x(u^\lambda_x)$
+ and minimum cut obtaining this configuration is shown. The last
+ column shows the capacities of the cut for the graph
+ construction in Figure~\ref{fig:bk_norm_subgraph}. We verify
+ that the cut capacities are equal to the functional value, plus
+ a constant.
+ }
+ \begin{tabular}{cccl}
+ \hline
+ $u^\lambda_x$ & $F_\lambda^x(u^\lambda_x)$ & Min.\ cut $(S, T)$
+ & Graph cut cap. \\
+ \hline
+ $0$ & $0$ & $(\{s\}, \{u^\lambda_x, t\})$ & $\max \{
+ F_{L-1}^x(1), 0\} - F_\lambda^x(1)$ \\
+ $1$ & $F_\lambda^x(1)$ & $(\{s, u^\lambda_x\}, \{t\})$ &
+ $\max \{F_{L-1}^x(1), 0\}$
+ \end{tabular}
+ \label{tab:fid_energy}
+\end{table}
+
This construction allows us to update the capacities of the edges, while
retaining the flow and the trees $S$ and $T$. A nice property that
-follows is that the partition $(V - T, T)$ is after each run a minimum
-cut since there is no more residual capacity between the trees $S$
-and $T$. Further, no vertex already in the tree $T$ will leave $T$ when
+follows is that the partition $(S, V-S)$ is after each run a minimum
+cut since the tree $S$ has been grown as much as possible.
+Further, no vertex already in the tree $S$ will leave $S$ when
the trees are updated and algorithm is run for lower $\lambda$ values.
-This can be seen from the fact that if $u$ is in $T$ and the algorithm
+This can be seen from the fact that if $u$ is in $S$ and the algorithm
has terminated, all paths from $s$ to $t$ going
-through $u$ have a saturated edge somewhere before $u$. For all these
-paths, the capacity will only change for edges after $u$, as we only
-change the capacity of edges $(v, t)$. Thus these saturated edges will
-stay saturated, and $u$ will stay in $T$. \fixme{shaky}
+through $u$ have a saturated edge somewhere after $u$. For all these
+paths, the capacity will only change for edges before $u$, as we only
+change the capacity of edges $(s, v)$. Thus these saturated edges will
+stay saturated, and $u$ will stay in $S$.
\subsection{Performance improvements}
--- /dev/null
+\chapter{Methods in image restoration}
+
+\fixme{Rating: 7/10}
+
+There are numerous methods in image restoration, and since we do not
+have time nor space to discuss them all, this short overview will focus on the ones
+related to the anisotropic total variation method considered later in
+this thesis. See \cite{jahne} and \cite{aubert2006proc} for more
+background on image processing in general.
+
+In this chapter, and also in the rest of the thesis we will assume that
+we are given an image $f : \Omega \to \mathbb{R}$ where $\Omega$ is a
+rectangular, open domain. Because of limitations in the numerical method
+used, the codomain is $\mathbb{R}$ and we are thus restricted to
+monochrome, or grayscale images. Such images are produced in large
+numbers by for example ultrasound, X-ray and MRI machines.
+
+The space in which the image $f$ resides in will vary, but since we are
+looking at image restoration methods, we assume that it includes some
+kind of noise. Depending on the application and how the image is
+obtained, one might construct different models describing different
+types of noise.
+
+We will assume that the given image $f$ is a combination of an
+underlying, actual image $u^*$, and some noise $\delta$. The simplest
+model is additive noise where the assumption is that $f = u^* + \delta$.
+There is also multiplicative noise where $f = u^* \cdot \delta$. An
+other much seen noise type is salt and pepper noise, which is when black
+and white pixels randomly appear in the image.
+
+These are only models, and in the real world the noise might be more
+complex, and even come from a combination of sources. Depending on the
+application, the goal might not even be to recover $u^*$, but rather to
+obtain an output which fulfills certain smoothness or regularity
+properties. In any case, we will continue using $f$ for the noisy input
+image and use $u$ for the output image in the description of the
+restoration methods.
+
+\section{Diffusion filtering}
+\label{sec:diffusion_filtering}
+
+\fixme{Rating: 8/10}
+
+Diffusion filtering is a broad group of filtering and restoration methods
+based on physical diffusion processes. The basic idea is to take the
+noisy image as the initial value of some diffusion process, and then let it
+evolve for some time. The best known method is probably the Gaussian
+filter or Gaussian blur, in which one convolves the image with the
+Gaussian function
+\begin{equation}
+ K_\sigma(x,y) = \frac{1}{2\pi\sigma^2} \exp \left( - \frac{x^2 +
+ y^2}{2\sigma^2} \right).
+ \label{eq:gaussian_function}
+\end{equation}
+\nomenclature{$K_\sigma(x,y)$}{The Gaussian kernel.}%
+In the discrete setting where the image consists of a grid of pixels,
+the Gaussian blur amounts to calculating each pixel in the output image
+as a weighted average of its neighboring pixels in the input image.
+
+The Gaussian function happens to be the fundamental solution of the heat
+equation $\partial_t u = \Delta u$. Convolving $K_\sigma(x,y)$ with the
+original image $f$ is therefore equivalent to solving the heat equation
+with $f$ as initial value, until some time $T > 0$ depending on
+$\sigma$. Boundary conditions have to be specified of course, and one
+common choice is to symmetrically extend the image in $x$ and $y$
+directions, which corresponds to zero flux boundary conditions.
+
+By basic Fourier analysis it is possible to show that the Gaussian
+filter is a low-pass filter which attenuates high frequencies. Further
+theory can be found in Weickert's book on anisotropic diffusion
+\cite{weickert1998anisotropic}.
+
+The main concern with the Gaussian filter is that it will, in addition
+to smoothing out possible noise, remove details from the image, which
+motivates the next set of methods, where the amount of diffusion can
+vary between different parts of the image.
+
+\subsection{Non-linear diffusion filtering}
+
+\fixme{Rating: 8/10}
+
+In the theory of the heat equation one can introduce a \emph{thermal
+diffusivity} $\alpha$ such that the equation becomes
+\begin{equation}
+ \left\{
+ \begin{aligned}
+ \partial_t u &= \diver \big( \alpha(\nabla u) \nabla u\big), \\
+ u |_{t=0} &= f.
+ \end{aligned}
+ \right.
+\end{equation}
+\nomenclature{$\alpha(\nabla u)$}{Scalar thermal diffusivity.}%
+The thermal diffusivity $\alpha(\nabla u) = \alpha(x, \nabla u)$ is material dependent, and can also
+vary throughout the object. It specifies how well heat travels through
+the specific point in the object. We can make use of this in the image
+restoration context specifying different diffusivity in different parts of
+the image, in an effort to reduce noise without loosing image
+detail. Optimally, we would like there to be a lot of diffusion in
+smooth parts of the image, and not so much in areas with a lot of
+details.
+
+One much-studied non-linear diffusion equation is the Perona--Malik
+equation
+\begin{equation}
+ \partial_t u = \diver \left( \frac{\nabla u}{1 +
+ \frac{\abs{\nabla u}^2}{\lambda^2}} \right).
+ \label{eq:perona_malik}
+\end{equation}
+The thermal diffusivity $\alpha(\nabla u) = (1 + \sfrac{\abs{\nabla u}^2}{
+\lambda^2})^{-1}$ varies from $1$ in smooth areas to $0$ as the
+norm of the gradient $\abs{\nabla u}$ grows.
+
+This particular form of the thermal diffusivity has been shown to be
+related to how brightness is perceived by the human visual system. The
+model has some theoretical problems related to well-posedness, for more
+information see \cite{weickert1998anisotropic}.
+
+A different kind of non-linear diffusion model is the total variation
+flow which can be formulated as
+\begin{equation}
+ \partial_t u = \diver \frac{\nabla u}{\abs{\nabla u}},
+\end{equation}
+where the diffusivity has a similar effect of reducing the diffusion in
+areas of high variation. As the name suggests this model can be related
+to the total variation formulation presented later. One forward Euler
+time-step in the solution of this partial differential equation
+corresponds to the Euler--Lagrange equation of the variational
+formulation.
+
+Note that we follow Weickert's terminology when it comes to the
+distinction between non-linear and anisotropic diffusion methods. The
+Perona--Malik equation, and other diffusion equations with
+non-homogenous diffusivities, are often by others called anisotropic, as
+the diffusivity depends on the location. We will name these methods
+non-linear and spare the anisotropy term for the ``real'' anisotropic
+methods. These are methods where the diffusivity is a tensor, and thus
+both location and direction dependent.
+
+\subsection{Anisotropic diffusion}
+\label{sec:anisotropic_diffusion}
+
+The diffusivity is made directionally dependent by introducing a
+diffusion \emph{tensor} $A(u)$ such that the initial boundary value
+problem becomes
+\begin{equation}
+ \left\{
+ \begin{alignedat}{10}
+ \partial_t u &{}={} \diver \big(A(u) \nabla u\big)
+ &\quad \text{ on }& \Omega \times (0, \infty),\\
+ u|_{t=0} &{}={} f
+ &\text{ on }& \Omega,\\
+ A(u) \nabla u \cdot \nu &{}={} 0
+ &\text{ on }& \partial \Omega \times (0, \infty),
+ \end{alignedat}
+ \right.
+ \label{eq:aniso_diff}
+\end{equation}
+\nomenclature{$A(u)$}{Thermal diffusivity tensor, or anisotropy
+tensor.}%
+where $\nu$ is the outer normal of $\Omega$. The tensor $A(u)$
+is constructed such as to diminish the effect of $\nabla
+u$ across what we believe to be edges in the image. This way, there will
+also be less diffusion through these edges. Weickert
+\cite{weickert1998anisotropic} suggests constructing $A(u)$ based on the
+edge estimator $\nabla u_\sigma$ where
+\nomenclature{$\nabla u_\sigma$}{Edge estimator.}%
+\begin{equation}
+ u_\sigma := K_\sigma * \tilde{u}
+\end{equation}
+\nomenclature{$u_\sigma$}{Image $u$ smoothed with a Gaussian of
+parameter $\sigma$.}%
+and $\tilde{u}$ is an extension of $u$ from $\Omega$ to $\mathbb{R}^2$
+made by symmetrically extending $u$ across the edges.
+Assuming we are at an edge in the image, the
+direction of $\nabla u_\sigma$ should be perpendicular to the edge,
+while its magnitude will provide information on the steepness of the
+edge.
+
+To extract this information, and also to identify features on a larger
+scale, the \emph{structure tensor} is introduced
+\begin{equation}
+ S_\rho(x) := K_\rho * (\nabla u_\sigma \otimes \nabla u_\sigma),
+\end{equation}
+\nomenclature{$S_\rho(x)$}{Structure tensor.}%
+where the convolution with the Gaussian function $K_\rho$ is done
+component-wise. The anisotropy tensor $A(u)$ can then be constructed
+based on the eigenvectors and eigenvalues of $S_\rho(x)$.
+The structure tensor and its properties will be discussed further when
+we introduce our anisotropic total variation functional.
+
+Assuming some
+smoothness, symmetry and uniform positive definiteness on $A(u)$ one can
+prove well-posedness, regularity and an extremum principle of the
+problem \eqref{eq:aniso_diff} as done in \cite{weickert1998anisotropic}.
+
+Even if the diffusivity tensor was introduced to reduce the amount of
+smoothing across edges, the solution of \eqref{eq:aniso_diff} will still be
+infinitely differentiable \cite{weickert1998anisotropic}, i.e.\ $u(T) \in
+C^\infty(\Omega)$ for $T > 0$. Thus there are no real discontinuities, and
+no real edges in the solution.
+
+Further, the anisotropic diffusion may introduce structure based on
+noise, when there really was no structure to begin with. This is a
+problem we aim to avoid in our anisotropic total variation method.
+
+%\begin{figure}
+% \centering
+% \begin{subfigure}[t]{0.4\textwidth}
+% \centering
+% \includegraphics[width=\textwidth]{fig/fnoise.png}
+% \caption{Noisy, $\sigma > 9000$}
+% \end{subfigure}
+% ~
+% \begin{subfigure}[t]{0.4\textwidth}
+% \centering
+% \includegraphics[width=\textwidth]{fig/fanisodiff.png}
+% \caption{Regularized, parameters: bleep bloop}
+% \end{subfigure}
+% \caption[bleapofsijdosij]{\fixme{UPDATE THIS}}
+% \label{fig:ansi_diff}
+%\end{figure}
+
+\section{Total variation}
+\label{sec:total_variation}
+
+\fixme{Rating: 8/10, numerical methods: 5/10}
+
+Total variation was initially introduced to the field of image
+restoration by Rudin, Osher and Fatemi in \cite{rudin1992nonlinear} and
+is usually formulated as a minimization problem
+\begin{equation}
+ \begin{gathered}
+ \min_{u \in L^p(\Omega)} F(u), \\
+ F(u) = \underbrace{\int_\Omega \abs{u - f}^p \, dx}_{\text{fidelity
+ term}} {}+{} \beta \underbrace{\int_\Omega
+ \abs{\nabla u} \, dx}_{\mathclap{\text{regularization
+ term}}},
+ \label{eq:first_min_presentation}
+ \end{gathered}
+\end{equation}
+where $p$ is normally taken to be 1 or 2.
+The fidelity term penalizes images $u$ that are far from the original
+image $f$. The
+regularization term is the total variation of the image, and minimizing it will
+reduce the variation and thus regularize the image. The $\beta$ parameter
+controls the strength of the regularization. Note that $u = f$ is a
+minimizer of the fidelity term, while a constant image $u = c$ is a
+minimizer of the regularization term.
+
+As this restoration method is the one
+which will be extended later in this thesis, we will look a little bit
+more deeply into the background and the numerical methods relating to
+it.
+
+Since we do not only want to consider differentiable images $u \in
+C^1(\Omega)$ for which the gradient exists, we introduce the total
+variation using the distributional derivative.
+\begin{definition}[Total variation]
+ Given a function $u \in L^1(\Omega)$, the total variation of $u$,
+ often written $\int_\Omega \abs{Du}\, dx$, where the $D$ is the
+ gradient taken in the distributional sense, is
+ \begin{equation}
+ \TV(u)
+ = \int_\Omega \abs{Du} \, dx
+ = \sup \left\{ \int_\Omega u \cdot \diver \varphi \, dx :
+ \varphi \in C^\infty_c\left(\Omega, \mathbb{R}^2\right),
+ \norm{\varphi}_{L^\infty(\Omega)} \leq 1
+ \right\}.
+ \end{equation}
+ The test functions $\varphi$ are taken from
+ $C^\infty_c\left(\Omega, \mathbb{R}^2\right)$, the space of smooth
+ functions from $\Omega$ to $\mathbb{R}^2$ with compact support in
+ $\Omega$.
+ \label{def:tv}
+\end{definition}
+\nomenclature{$\TV(u)$}{Total variation of image $u$.}%
+
+Note that since $\Omega$ is open and bounded, the test functions
+$\varphi$ vanish on the boundary of $\Omega$. Thus no variation is
+measured at the boundary.
+
+As we are searching for an image with low total variation, it is useful
+to introduce the space of functions of bounded variation.
+\begin{definition}[Functions of bounded variation]
+ The space of functions of bounded variation $\BV(\Omega)$ is the
+ space of functions $u \in L^1(\Omega)$ for which the total variation
+ is finite, i.e.,
+ \begin{equation}
+ \BV(\Omega) = \left\{ u \in L^1(\Omega) : \TV(u) < \infty
+ \right\}.
+ \end{equation}
+\end{definition}
+\nomenclature{$\BV(\Omega)$}{Functions of bounded variation in
+$\Omega$.}%
+Our optimization problem has thus become
+\begin{equation}
+ \min_{u \in \BV(\Omega)} \int_\Omega \abs{u - v}^p \, dx + \beta \,
+ \TV(u).
+ \label{eq:second_min_presentation}
+\end{equation}
+
+As with any restoration method, the total variation method has its
+strengths and weaknesses. Its main strength is its ability to recover
+edges in the input image. The total variation of a section only takes
+the absolute change into account, and does not favor gradual changes
+like the diffusion methods.
+
+There is also a theoretical result stating that the set of edges in the
+solution $u$ is contained in the set of edges in the original image $f$,
+thus no new edges are created \cite{caselles2011total}. However, in the
+presence of noise, the method may introduce or rather ``find'' new edges
+that were not in the original image, since flat sections of zero
+variation are encouraged by the functional. This effect is called the
+\emph{stair-casing effect}, and can be seen in Figure~\ref{fig:grad_tv}
+where a noisy gradient has been restored using total variation
+restoration.
+
+\begin{figure}
+ \centering
+ \begin{subfigure}[t]{0.4\textwidth}
+ \centering
+ \includegraphics[width=\textwidth]{fig/factory/grad/n_q40.png}
+ \caption{Noisy gradient}
+ \end{subfigure}
+ ~
+ \begin{subfigure}[t]{0.4\textwidth}
+ \centering
+ \includegraphics[width=\textwidth]{fig/factory/grad/r_p2_n16_b5000_g1000000000_r5_s3.png}
+ \caption{Total variation restoration}
+ \end{subfigure}
+ \caption[The stair-casing effect of regular total variation
+ restoration]{%
+ Although the original gradient was smooth, the total variation
+ method manages to find structure in the noise, and create edges
+ in the restored image.
+ }
+ \label{fig:grad_tv}
+\end{figure}
+
+Fine details, thin objects and corners may suffer from contrast loss since
+bringing them closer to their surroundings reduces the
+total variation. An example of this is shown in
+Figure~\ref{fig:tv_example}, where a not particularly noisy fingerprint
+image has been strongly regularized. The original black and white levels
+have been brought closer to yield a lower total variation in the
+regularized image.
+
+\begin{figure}
+ \centering
+ \begin{subfigure}[t]{0.4\textwidth}
+ \centering
+ \includegraphics[width=\textwidth]{fig/finger.png}
+ \end{subfigure}
+ ~
+ \begin{subfigure}[t]{0.4\textwidth}
+ \centering
+ \includegraphics[width=\textwidth]{fig/ftv.png}
+ \end{subfigure}
+ \caption[Heavily regularized fingerprint showing contrast loss]{%
+ A fingerprint heavily regularized using the total variation
+ method. The originally white and black ridges have been brought
+ closer in value, to reduce the total variation.
+ }
+ \label{fig:tv_example}
+\end{figure}
+
+\subsection{Numerical methods}
+
+See \cite{chan2011numerical} for an overview of some of the numerical
+methods relating to total variation image restoration. Amongst others it
+describes some dual and primal-dual methods, as well as the graph cut
+approach we take in this thesis.
+
+\subsubsection{Graph cut approach}
+
+\fixme{Could probably be shortened a lot, but I thought it could be
+useful to at least give an idea how it is done.}
+
+Using graph cuts is the approach we will be taking later when
+considering the anisotropic total variation regularization, and it is
+therefore valuable to briefly look into how graph cuts are used in the
+case of regular total variation.
+
+In the theory of graphs, a graph cut is a set of edges that when
+removed will separate the graph into two disconnected parts. A minimum
+cut is a cut such that the sum of the weight of the edges cut is minimal. It
+has been shown that for some energy functionals, it is possible to
+construct graphs for which the minimum cuts correspond to minimizers of
+the functional.
+
+In the discrete setting our image consists of pixels, and is represented
+by a function $u : \mathcal{G} \to \mathcal{P}$ where $\mathcal{G}$ is a
+regular grid over $\Omega$, and $\mathcal{P} = \{0, \hdots, L-1\}$ is
+the discrete set of pixel values, or \emph{levels}. We denote the value
+in pixel $x$ as $u(x) = u_x$.
+
+For an image $u$ and a level $\lambda$ we denote the \emph{level set} by
+$\{ u > \lambda\}$, defined as the set $\{ x \in \Omega : u_x > \lambda
+\}$. The thresholded image $u^\lambda$, an indicator function, is then
+defined as
+\begin{equation}
+ u^\lambda = \idfun_{u > \lambda}.
+\end{equation}
+Here, $\idfun_E$ signifies the characteristic function of the set $E$,
+\nomenclature{$\idfun_E$}{Characteristic function of the set $E$.}%
+the function which is equal to one in every point in $E$, and zero
+elsewhere.
+
+The idea of the graph cut approach is to decompose the minimization
+problem into one minimization problem for each level of the image, and
+then solve them separately before combining the results.
+
+Through careful manipulation of the continuous functional in
+\eqref{eq:second_min_presentation} it is possible to obtain a discrete
+functional decomposed as a sum over all the level values on the form
+\begin{equation}
+ F(u) =
+ \sum_{\lambda=0}^{L-2} \sum_x F^x_\lambda(u^\lambda_x)
+ + \beta \sum_{\lambda = 0}^{L-2} \sum_{(x, y)} F^{x,y}(u^\lambda_x,
+ u^\lambda_y)
+ =: \sum_{\lambda=0}^{L-2} F_\lambda(u^\lambda)
+ \label{eq:old_total_energy}
+\end{equation}
+where the sum over $(x, y)$ is over all pixel pairs $(x, y)$ in a neighbor
+relation, i.e.\ where the pixels are ``close'' to each other.
+The actual form of the functional, and the steps to construct it will be
+presented later.
+
+The graph cut we find will for each level $\lambda$ give us the
+thresholded image $u^\lambda$, and they can then be combined to form the
+complete image $u$.
+
+When constructing the graph used to find the thresholded image
+$u^\lambda$, we have two special vertices, one representing the set $\{
+u > \lambda\}$, while the other represents the set $\{ u \leq
+\lambda\}$. The pixels are then connected to these vertices with a
+weight representing how strongly they are related to the corresponding
+set. This weight will be based on the value of $F^x_\lambda$.
+
+Additionally there are connections between pixels in a neighborhood relation,
+representing the energy $F^{x,y}$. Thus when finding a cut, we partition
+the pixels into the sets $\{ u > \lambda \}$ and $\{ u \leq \lambda \}$.
+And if in addition the cut is minimal, we know that the edges cut have
+minimal weight, and can prove that the $u^\lambda$ found minimizes
+the functional in \eqref{eq:old_total_energy}.
+
Figure~\ref{fig:lena_process} we show the different stages of our image
restoration algorithm. We start out with a section of the much-used Lena
test image in Figure~\ref{fig:lena_process_orig}. The first step in the
-construction of the structure tensor is to blur the image with parameter
+construction of the structure tensor is to apply a Gaussian blur with parameter
$\sigma$, and the result is shown in Figure~\ref{fig:lena_process_blur}.
The blurring is done so that the edge detector $\nabla f_\sigma$,
visualized as $\abs{\nabla f_\sigma}^2$ in Figure~\ref{fig:lena_process_edge}, is not too
have been scaled by the corresponding eigenvalue.
Another way of visualizing the tensor is shown in
-Figure~\ref{fig:lena_process_color}. The color in each point is decided
-by the size of the smallest eigenvalue, and the direction of its
-corresponding eigenvector, using the color wheel in
+Figure~\ref{fig:lena_process_color}. The brightness and color in each point is decided
+by the size of the smallest eigenvalue and the direction of its
+corresponding eigenvector using the color wheel in
Figure~\ref{fig:color_wheel}. The direction decides the color,
-while the brightness, here the radius, is set to $1 / \sigma_1 - 1$,
+while the brightness, which is the radius in the color wheel, is set to $1 / \sigma_1 - 1$,
where $\sigma_1$ is the smallest eigenvalue.
Thus the stronger the anisotropy, the brighter the color, while
we expect uniform areas in the original image to be black in the
\end{overpic}
%\includegraphics[width=\textwidth]{fig/factory/lines/r_p2_n8_b10000_g100_r4_s2_color.png}
\caption{%
- Restored with relatively high noise scale of $\sigma = 2$.
+ High noise scale of $\sigma = 2$ means thinnest lines are
+ regularized.
}
\label{fig:scale_comp_high}
\end{subfigure}
\end{overpic}
%\includegraphics[width=\textwidth]{fig/factory/lines/r_p2_n8_b10000_g100_r15_s0.2_color.png}
\caption{%
- Restored with low noise scale $\sigma = 0.2$, but large integration
- scale $\rho = 15$.
+ Low noise scale $\sigma = 0.2$, but large integration
+ scale $\rho = 15$ means inner structure is ignored.
}
\label{fig:scale_comp_int}
\end{subfigure}
An example constructed to show the effects of the parameters
$\sigma$ and $\rho$ in the anisotropy tensor.
The noise scale $\sigma$ controls the smoothing done before edge
- detection, and $\sigma=2$ means the thinnest lines are
- not considered to be edges by the tensor, and they are thus
- regularized.
+ detection.
The integration scale $\rho$ controls the size of the structures
- considered by the tensor, and does in this case allow us to
- ignore the smaller structures. Parameters: $\abs{\mathcal{N}} =
+ considered by the tensor.
+ Parameters: $\abs{\mathcal{N}} =
8$, $\beta = 10000$, $\omega = 100$.
}
\label{fig:scale_comp}
variation in Figure~\ref{fig:finger_contrast_tv}. By visual inspection,
the last image seems to have somewhat higher contrast. This can be
confirmed by calculating the standard deviation of the images, one of
-many possible contrast measures, which gives 22.9, 37.9 and 26.6
+many possible contrast measures, which gives 22.9, 37.9 and 27.2
respectively.
To further compare the three restored images, a single row has been
Figure~\ref{fig:finger_contrast_tv}, some details in the singularity in
the upper left is actually lost. The metric tensor is approximately the
identity matrix there, and thus an increased restoration parameter
-$\beta$ leads to increased smoothing.
+$\beta$ leads to increased smoothing compared to the regular total
+variation of Figure~\ref{fig:finger_contrast_tv}.
\begin{figure}
\centering
So far we have just inspected the restored image visually to judge the
results of the method. Depending on the application, we might have
different wishes for the results. In many cases we just want to remove the
-noise $\delta$ in \fixme{ref} where the captured image $f$ is assumed to
-consist of an actual image $u$ and an additional noise term $\delta$
-such that $f = u + \delta$. The noise removed by the restoration method,
+noise $\delta$ in $f = u^* + \delta$ where the captured image $f$ is assumed to
+consist of an actual image $u^*$ and an additional noise term $\delta$.
+The noise removed by the restoration method,
$\hat{\delta} = f - u^*$ is called the method noise. Any assumptions of
properties we have on $\delta$, we optimally want $\hat{\delta}$ to
fulfill as well. Thus if we assumed independently distributed Gaussian
can in that case be if the method noise contains a lot of the details of
the image.
+\fixme{fuck}
+
%Figure~\ref{fig:lena_method_noise} shows the results after restoring the
%noisy Lena picture in Figure \fixme{noisylena}, with different
%parameters, and their method noise. We note that the method noise of the
projects / master / commitdiff