discretization and graph cut framework used in my project work
\cite{project}. The anisotropy tensor was introduced in hopes of
reducing the amount of regularization applied across edges in the image,
-and by that prevent contrast loss. It was
+and by that preventing contrast loss. It was
constructed based on the structure tensor which contains local information
about edge direction and steepness.
implemented.
A goal was to give readers a good understanding of the inner workings
-of the method, as well as everything necessary for a working
-implementation. Thus we put effort into describing all steps, from the
-initial continuous formulation, to the discretization, and the maximum
-flow algorithms in the end.
+of the method, as well as all details necessary for a working
+implementation. Thus we put effort into describing all steps going from the
+initial continuous formulation, to the discretization, and finally the maximum
+flow algorithms.
%Effort was put into giving a complete overview of the method, describing
%each part needed to go from the initial continuous problem, to the
Further work is possible in the study of the continuous problem, its
well-foundedness, and also the anisotropic coarea and Cauchy--Crofton
formulas, which can be studied on a measure-theoretic foundation. Also,
-the construction of the tensor offers choices that can be explored
+the construction of the anisotropy tensor offers choices that can be explored
further.
Regarding the discretization, the choice of neighborhood stencil also
algorithm performance. One particular possibility would be the
application of non-uniform stencils with varying stencil size depending
on the local level of detail. This could give better performance without
-sacrificing too much in solution accuracy.
+sacrificing too much solution accuracy.
%Because we have spread our attention across all parts of the restoration
%method, it is also possible to go deeper into the theory behind all of
\TVA(u) = \int_\Omega \sqrt{\nabla u(x)^T A(x) \nabla u(x)} \, dx
\label{eq:aniso_tv_sqrt}
\end{equation}
+\nomenclature{$A(x)$}{Anisotropy tensor}%
for all $u \in C^1(\Omega)$. We assume here that $A(x)$ is continuous
and positive definite,
and we will later need the eigenvalues of $A(x)$ to be uniformly bounded
below and above. If $A(x)$ is the identity matrix we get the
regular total variation found in \eqref{eq:first_min_presentation}.
-When reducing the regular total variation, we will also try to reduce
+When minimizing the regular total variation, we will also try to reduce
the variation over known edges in the image. This can lead to unwanted
contrast loss, especially in fine details.
By controlling $A(x)$ such that the contribution of $\nabla u(x)$ is
reduced across known edges, we hope to retain the regularization
-properties of the original method while reducing some of the negative
-effects. If the variation across an edge is ``ignored'' by the
+properties of the original method while reducing this contrast loss.
+If the variation across an edge is ``ignored'' by the
functional, there is no gain in reducing the height of the edge as
before.
in his book on Anisotropic Diffusion \cite{weickert1998anisotropic}, and
briefly mentioned in Section~\ref{sec:anisotropic_diffusion}, a good
starting point is the \emph{edge detector} $\nabla f_\sigma$.
+\nomenclature{$\nabla f_\sigma$}{Edge detector}%
The image is smoothed by a Gaussian filter as described in
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}$.
We order the two real
-eigenvalues $\lambda_1 \geq \lambda_2$ with corresponding eigenvectors
+eigenvalues such that $\lambda_1 \geq \lambda_2$ and denote the 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}
\end{smallmatrix})$ we obtain a closed form
\end{aligned}
\label{eq:sigma_construction}
\end{equation}
-\nomenclature{$\omega$}{Anisotropy parameter}
+\nomenclature{$\omega$}{Anisotropy parameter}%
Thus the eigenvectors of $A(x)$ and $S_\rho(x)$ are equal, while the
eigenvalues are different. A visualization of the two tensors can be
seen in Figure~\ref{fig:tensor_viz} where the two tensors are shown at
Around corners $A(x)$ will be close to the identity matrix, which gives
regularization similar to smooth areas. This is one possible down-side
-of this tensor choice.
+of this tensor choice, as rounded corners may occur.
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
+such that if it is very large we are left with the identity matrix and
+our method becomes the regular total variation method. Note also
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)$.
For the case where $\lambda_1 = \lambda_2$, the $U(x)$ in our
-decomposition is not well-defined. This is not a problem though, since
+decomposition is not well-defined. This is not a problem however, since
$\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
+When $\lambda_1 = \lambda_2$ however, the eigenvalues $\sigma_1$ and
+$\sigma_2$ 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$.
to have a \emph{well-posed} problem are lower semi-continuity
and coercivity for existence, and convexity for uniqueness. We restrict
ourself to $L^2(\Omega)$ which makes sense with our fidelity term,
-assuming that $f \in L^2(\Omega)$ initially.
+assuming that $f \in L^2(\Omega)$.
We consider the weak topology, as it will allow us to arrive at an
existence result relatively easily.
\end{proof}
This coarea formula is our first step in transforming the anisotropic
total variation into an easily discretizable expression. It allows us to
-consider each level $\lambda$ separately when calculating the total
-variation.
+consider each level $\lambda$ separately when calculating the
+anisotropic 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
here, but can for the isotropic case be found in for example
\cite[Section 5.4 and 5.5]{ziemer1989}.
-The anisotropic coarea formula allows us to transform the anisotropic
-total variation such that we are left with minimizing the following
+Using the anisotropic coarea formula and inserting the anisotropic
+perimeter definition we transform the anisotropic
+total variation and are left with the problem of minimizing the following
functional
\begin{equation}
F(u) = \int_\Omega (u - f)^2 \, dx + \beta \, \int_{-\infty}^\infty
product of two vectors $a$ and $b$ in a point $x$ is calculated as
$\langle a, b\rangle_M = \langle a, M(x) b \rangle$, then the length of
a curve $\gamma$ parametrized by some parameter $t$ becomes
-\nomenclature{$M(x)$}{A metric tensor}%
+\nomenclature{$M(x)$}{The metric tensor}%
\nomenclature{$\abs{C}$}{The length of the curve $C$}%
\nomenclature{$\abs{C}_M$}{The length of the curve $C$ calculated using
the metric tensor $M$}%
\end{pmatrix}
\end{align}
As $U$ and $U^T$ correspond to rotations and our measure
- $\mathcal{L}$ is invariant under rigid motions, $U$ and $U^T$ do not
+ $\mathcal{L}$ is invariant under rotations, $U$ and $U^T$ do not
have direct contributions to the Jacobian. They do however affect the
input angle of the operator $\Lambda$ such that $J_M(\ell_{\phi,
\rho}) = J_{\Sigma^2}(U^T \ell_{\phi, \rho})$. Thus we will now
tensor $A$ equal to the identity matrix $I$, the perimeter coincides
nicely with the length of the boundary curve, assuming some regularity
for the boundary.
-
In the general case we need
to be more careful. As can be seen in \eqref{eq:perimeter_calc}, the
anisotropic perimeter is calculated by integrating the norm of the
elements in
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)$
+space is bounded above. Thus the eigenvalues $\lambda_1 \geq \lambda_2$ of $S_\rho(x)$
are bounded above. Moreover, by the construction in
\eqref{eq:sigma_construction}, there exists uniform bound $k$ such that
the smallest eigenvalue $\sigma_1$ of the anisotropy tensor $A(x)$ is bounded away
from zero, as
\begin{equation}
\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.
+ + \frac{(\lambda_1 - \lambda_2)^2}{\omega^2}\right)^{-1} \geq \left(1 +
+ \frac{\lambda_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 with bounded eigenvalues and thus the curve length calculation in
+definite with bounded eigenvalues $k \leq \sigma_1 \leq \sigma_2 \leq K
+= 1$ 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}.
-Thus we can apply the formula to calculate the
+Hence we can apply the formula to calculate the
perimeter in \eqref{eq:per_to_length2} as
\begin{equation}
\PerA(U; \Omega) = \int_\mathcal{L} \sum_{x \in
other varies from 1 in smooth areas towards 0 around edges, with the
corresponding eigenvector perpendicular to the edge.
-As mentioned previously, these transformations are motivated by the
-discretization that follows in the next chapter.
-
The set of lines can then be represented by the neighborhood of a pixel as
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
+neighborhood of grid point $x$.
+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$.
\forall x \in \mathcal{G}.
\label{eq:decreasing_sets}
\end{equation}
-Later we will present a graph cut algorithm that find the thresholded
+Later we will present two graph cut algorithms that find the thresholded
images minimizing each level, \emph{while guaranteeing that they meet
this requirement.}
This construction can be seen in Figure~\ref{fig:square_cons} for $n =
5$.
The maximal angle between two lines $\phi_k - \phi_{k-1}$ will always be
-when one of $\phi_k$ or $\phi_{k-1}$ is horizontal or vertical, shown in
+when one of $\phi_k$ and $\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 be when $\phi_k = \sfrac{m \pi}{2}$ for $m
Therefore
numerous different approaches to denoising exist, each having their own
strengths and weaknesses. Some of these are introduced in
-Chapter~\ref{chap:methods}, all having the weakness that edges in the
-original are smoothed over.
+Chapter~\ref{chap:methods}, and one of the main challenges they all face
+is the
+recovery of edges.
A method well known for recovering edges is the total variation
method, as the total variation does not favor smooth gradients over
-edges. In my project work~\cite{project}, I described how a graph cut
-framework can be used to obtain a numerical
-solution in a total variation based
-image restoration method.
-The method consisted of trying to reduce the total variation of the
+edges.
+I gave an overview of this method in my project work~\cite{project},
+where I used a graph cut
+framework to obtain a numerical
+solution.
+The method consists of trying to reduce the total variation of the
image, while still staying ``close'' to the original.
A problem with the total varation method is that
tensor into the total variation, thus making it directionally dependent.
This means we can control the regularization
applied to the image based on position and direction. The main idea is then to reduce
-the regularization across edges in the image. Information about the
-edges are extracted using the structure tensor.
+the regularization applied across edges in the image, while we still
+regularize along them.
-This continuous problem is a convex minimization problem, and many
-optimization approaches exist. However, we choose to discretize in such
+The variational problem we obtain is a convex minimization problem, and many
+optimization approaches exist. We choose to discretize in such
a way that we can apply the same graph cut framework used in my project
work~\cite{project}. Through the coarea formula, the functional is
decomposed into a sequence of minimization problems, one for each level of the image.
functionals. The push-relabel algorithm is considered to be the fastest
and most versatile for general graphs, while the Boykov--Kolmogorov
algorithm is specially tailored for the type of graphs we find in these
-kinds of imaging applications.
+kinds of imaging applications. We describe every part of the method in
+detail such that it can be easily implemented by the reader. In
+addition, a \cpp{} implementation is attached.
In the end we present numerical results that show how the
different parameters affect the restoration, and we look into and explain some
artifacts caused by approximations in the discretization. Further we look at how the
-introduction of the anisotropy in some cases amend some of the
-weaknesses of the total variation method.
+introduction of the anisotropy in certain cases amend some of the
+weaknesses of the total variation method. We particularly look at how
+contrast loss is reduced in images containing thin structures such as
+fingerprints.
%artifacts caused by approximations in the discretization, and how they
%affect the restoration.
%Options: Sonny, Lenny, Glenn, Conny, Rejne, Bjarne, Bjornstrup
\usepackage[Lenny]{fncychap}
+%\mathtoolsset{showonlyrefs=true}
+
\newminted{c++}{fontsize=\tiny}
\usepackage{polyglossia}
thesis we consider a particular kind of edge-enhancing image restoration method based on
total variation. We want to address the fact that the total variation
method in
-some cases leads to contrast loss in thin structures, and so
+some cases leads to contrast loss in thin structures. To reduce the
+contrast loss
a directional dependence is introduced
through an
-anisotropy tensor. The tensor controls the regularization applied based on
-position and direction, and is
-constructed using edge information extracted from the noisy image. In
-order to minimize the resulting functional, it is discretized in such a
-way that we can construct graphs where minimum cuts correspond to
-minimizers of the discrete functional. This transformation is made
-possible by a coarea and Cauchy--Crofton formula.
-We end with some experiments on the different parameters to see how the
-method can be controlled and we show that the introduced anisotropy in
-some cases improve the results.
+anisotropy tensor.
+The tensor controls the regularization applied based on
+the position in the image and the direction of the gradient. It is
+constructed using edge information extracted from the noisy image.
+We optimize the resulting functional using a graph cut framework;
+%In
+%order to minimize the resulting functional, it is discretized in such a
+%way that we can construct graphs where minimum cuts correspond to
+%minimizers of the discrete functional.
+%This transformation is made
+a discretization which is made
+possible by a coarea and a Cauchy--Crofton formula.
+%We end with some experiments on the different parameters to see how the
+%method can be controlled and we show that the introduced anisotropy in
+%some cases improve the results.
+In the end we perform numerical studies, experiment with the parameters and
+discuss the results.
+
%There are numerous methods for image denoising and restoration. In this
%thesis we consider a total variation based method modified with an
\section*{Sammendrag}
-Lololol.
+I denne master-avhandlingen ser vi på en spesifikk kant-bevarende
+støyfjerningsalgoritme basert på «total variation». Vi vil ta for oss
+at «total variation» i noen tilfeller fører til tap av kontrast i
+detaljer og tynne strukturer. For å redusere kontrast-tapet introduserer
+vi en retningsavhengig anisotrop tensor. Denne tensoren kontrollerer
+støyfjerningen basert på posisjonen i bildet, og retningen til
+gradienten. Den blir konstruert basert på kant-informasjon fra det
+opprinnelige støyete bildet. Vi optimerer funksjonalen ved hjelp av et
+graf-kutt-rammeverk, som er gjort mulig ved hjelp av en coarea- og en
+Cauchy--Crofton-likning. Vi avslutter med en numerisk studie og
+exprimentering på parametrene og diskusjon av resultatene.
\vspace*{\fill}
\section*{Preface}
+\begingroup
+\setlength{\parskip}{0.4cm}
+\setlength{\parindent}{0cm}
This master thesis concludes my study at the Applied Physics and
Mathematics Master's degree programme with specialization in Industrial
Mathematics at the Norwegian University of Science and Technology
work with my project and this thesis.
Finally I would like to thank my family for their support, and Mats,
-Lars, Kine, Hager, Henrik and Edvard for life coaching and productive
+Lars, Kine, Hager, Henrik and Edvard for productive
discussions around the coffee pot.
Bjørn Rustad, \today.
+\endgroup
\vspace*{\fill}
\cleardoublepage
% Bibliografi/referanseliste skal komme før appendiks
+\phantomsection
+\addcontentsline{toc}{chapter}{\bibname}
\bibliography{bib}
\bibliographystyle{unsrt}
\cleardoublepage
+\phantomsection
+\addcontentsline{toc}{chapter}{\listfigurename}
\listoffigures
\cleardoublepage
+\phantomsection
+\addcontentsline{toc}{chapter}{\listtablename}
\listoftables
\cleardoublepage
+\phantomsection
+\addcontentsline{toc}{chapter}{\nomname}
\markboth{\MakeUppercase\nomname}{\MakeUppercase\nomname}
\printnomenclature[2.5cm]
\chapter{Methods in image restoration}
\label{chap:methods}
-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
+There are numerous methods in image restoration, but we do not
+have time nor space to discuss them all. In short overview, which is an
+extension of the one given in my project \cite{project}, we will focus
+on the methods
related to the anisotropic total variation method considered later in
-this thesis. See \cite{jahne} and \cite{aubert2006proc} for more
+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
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
+properties. In any case, we will continue denoting the noisy input
+image $f$ and use $u$ for the output image in the description of the
restoration methods.
\section{Diffusion filtering}
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
+restoration context by 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
\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
+areas of high variation. As the name suggests, this model can be related
+to the variational 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.
\right.
\label{eq:aniso_diff}
\end{equation}
-\nomenclature{$A(u)$}{Thermal diffusivity tensor, or anisotropy
-tensor}%
+\nomenclature{$A(u)$}{Thermal diffusivity 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}
and $\tilde{u}$ is an extension of $u$ from $\Omega$ to $\mathbb{R}^2$
\nomenclature{$\tilde{u}$}{Symmetric extension of $u$ from $\Omega$ to
$\mathbb{R}^2$}%
-made by symmetrically extending $u$ across the edges.
+made by symmetrically extending $u$ across the boundary of $\Omega$.
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
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
+However, 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
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.
+where a noisy gradient has been restored using the total variation
+method.
\begin{figure}
\centering
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
+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
+cut is a cut such that the sum of the weight of the edges in the cut is minimal. It
+has been shown that for some discrete functionals, it is possible to
construct graphs for which the minimum cuts correspond to minimizers of
the functional.
\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.
+relation, i.e.\ pixels that are ``close'' to each other.
The actual form of the functional, and the steps to construct it will be
presented later.
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
+u > \lambda\}$, and one which 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$.