+++ /dev/null
-% Dette eksempelet er laget for article-dokumentklassen. Hvis
-% skriver i 'book'-dokumentklassen vil du kanskje bytte ut
-% \section med \chapter, \subsection med \section, osv...
-
-\section{Total variation}
-The method was initially introduced by L. Rudin, S. Osher, and E. Fatemi
-in \cite{rudin1992nonlinear} in the context of image restoration, but is
-now used in numerous other applications. One of its strong theoretical
-points is its ability to recover edges, compared to other restoration
-algorithms which might smooth over edges in the original image. The
-method can be formulated as minimizing the energy function
-\begin{equation}
- E_v(u) =
- \int_\Omega \abs{u - v}^p
- + \beta \TV(u)
- \label{eq:energy_function}
-\end{equation}
-where $v$ is the original image, and $\Omega$ is the image domain, in
-our case a subspace of $\mathbb{R}^2$. Normally, $p$ is chosen to be 1
-or 2. \fixme{Why? Theoretical niceness? Models the real world? There are
-also other norms describing other noise models?}
-
-The term $\TV(u)$ is the total variation of the image, and is defined as
-follows:
-\begin{definition}[Total variation]
- Given a function $u \in L^1(\Omega)$, the total variation of $u$,
- often written $\int_\Omega \abs{Du}$, where the $D$ is the gradient
- taken in the distributional sense, is
- \begin{equation}
- \TV(u)
- = \int_\Omega \abs{Du}
- = \sup \left\{ \int_\Omega u \,\, \mathrm{div} \, \varphi :
- \varphi \in C^\infty_c\left(\Omega, \mathbb{R}^N\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}^N\right)$, the space of smooth
- functions from $\Omega$ to $\mathbb{R}^N$ with compact support.
-\end{definition}
-
-Minimizing the total variationn will smooth out differences in the
-image. The term $\int_\Omega \abs{u - v}^p$ of
-\eqref{eq:energy_function} will constrain the solution image $u$ to be
-close to the input image $v$. Combining the two will give us an image
-close to the original, but more regular. The parameter $\beta$ controls
-how strongly we want to regularize the image.
-
-Following the notation used in \cite{darbon2006image}, we let $u_s$
-denote the value of the image $u$ at position $s \in \Omega$ and
-introduce the level sets $u^\lambda$ of the image
-\begin{equation}
- u^\lambda = \idfun_{u \leq \lambda}
- \label{eq:level_set}
-\end{equation}
-which will be important when dividing the problem into smaller
-subproblems.
-
-When discretizing the energy function \eqref{eq:energy_function} later
-on, we want to decompose it as a sum over the different levels, or pixel
-values, of the image. This is why we present the coarea formula here,
-and use it to write the total variation $\TV(u)$ as an integral over the
-range level values $\lambda$. The coarea formula and further references
-can be found in \cite{caselles2011total}.
-
-But first we need to introduce the space of functions with bounded
-variation and the perimeter of a set.
-\begin{definition}[Functions of bounded variation]
- The space of functions with bounded variation $\BV(\Omega)$ is the
- space of functions $u \in L^1(\Omega)$ such that 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}
-\begin{definition}[Set perimeter]
- The perimeter of a set $E$ is defined as
- \begin{equation}
- P(E,\Omega) = \int_\Omega \abs{D \idfun_E}.
- \end{equation}
- A measurable set $E \subset \Omega$ is of finite perimeter in
- $\Omega$ if $\idfun_E \in \BV(\Omega)$.
-\end{definition}
-These two definitions allows us to present the coarea formula
-\begin{theorem}[The coarea formula]
- Let $u \in BV(\Omega)$. Then for almost every $\lambda$ the set
- $u^\lambda$ is of finite perimeter, and one has the coarea formula
- \begin{equation}
- \TV(u)
- = \int_\Omega \abs{D u}
- = \int_\mathbb{R} P(u^\lambda) \, d\lambda
- \label{eq:coarea_formula}
- \end{equation}
-\end{theorem}
-This allows us to rewrite \eqref{eq:energy_function} to
-\begin{equation}
- E_v(u) =
- \int_\Omega \abs{u - v}^p
- + \beta \int_\mathbb{R} P(u^\lambda) \, d\lambda,
- \label{eq:energy_perimeter}
-\end{equation}
-which is a good starting point for the discretization. More analysis of
-the continuous problem is given in \cite{caselles2011total}, together
-with different numerical methods for solving it. Results on how the set
-of jumps in the resulting image $u$ is contained in the set of jumps in
-the original image $v$ are also presented.
-
-\fixme{\cite{caselles2011total} has a nice introduction to total
- variation. It describes the foundations, the spaces, and some
- methods, but mostly continuous. Could include some more theory here,
- like the coarea stuff, and the number of discontinuities?}
-
-\fixme{Write something about the different methods that exist?}
-
-\subsection{Discrete problem}
-Since digital images are given on a discrete grid, with values taken
-from a discrete and finite set of levels, we want to discretize
-\eqref{eq:energy_perimeter}, and solve the resulting discrete problem.
-We assume that the image is given on a discrete grid $S$,
-where the value of each pixel is taken from the set $\mathcal{L} = [0,
-\hdots, L-1]$. A reasonable assumption for grayscale images.
-
-\begin{figure}
- \centering
- \begin{tikzpicture}[
- baseline=(M.center),
- bg1/.style={fill=black!10},
- bg2/.style={fill=black!25},
- bg3/.style={fill=black!40},
- bg4/.style={fill=black!55},
- bg5/.style={fill=black!70},
- ]
- \matrix (M) [
- matrix of nodes,
- %column sep={6pt,between origins},
- %row sep={6pt,between origins}
- ]
-{
- \node (A) {}; & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 \\
- & 1 & 1 & 1 & 1 & 1 & 1 & 1 & & 1 & 1 & 1 & 1 \\
- & 1 & 1 & & 1 & 1 & 1 & & & & 1 & 1 & 1 \\
- & 1 & & & & 1 & 1 & & & & & 1 & 1 \\
- \node (L) {$\mathcal{L}$}; & 1 & & & & & & & & & & & 1 \\
- & |[bg5]| 5 & |[bg3]| 3 & |[bg2]| 2 & |[bg3]| 3 & |[bg4]| 4 & |[bg4]| 4 & |[bg2]| 2 & |[bg1]| 1 & |[bg2]| 2 & |[bg3]| 3 & |[bg4]| 4 & |[bg5]| 5 \\
-};
-
-\draw[->] (L) -- (A);
-\end{tikzpicture}
-\caption{The level sets of a one-dimensional image.}
-\label{fig:level_sets}
-\end{figure}
-After discretizing the image, we obtain a finite number of level sets
-$u^\lambda$. Figure \ref{fig:level_sets} shows an example of a
-one-dimensional image divided into its level sets.
-
-\begin{figure}[b]
- \centering
- \begin{subfigure}[b]{0.4\textwidth}
- \centering
- \begin{tikzpicture}[->,>=stealth',shorten >=1pt,auto,node distance=1.5cm,
- thick,main node/.style={circle,draw}]
-
- \node[main node] (1) {u};
- \node[main node] (2) [below of=1] {};
- \node[main node] (3) [right of=1] {v};
- \node[main node] (4) [left of=1] {};
- \node[main node] (5) [above of=1] {};
-
- \path[every node/.style={font=\sffamily\small}]
- (1) edge node [right] {} (2)
- (1) edge node [right] {$w_{uv}$} (3)
- (1) edge node [left] {} (4)
- (1) edge node [left] {} (5);
- \end{tikzpicture}
- \caption{Size four neighborhood.}
- \label{fig:c4_grid}
- \end{subfigure}
- ~
- \begin{subfigure}[b]{0.4\textwidth}
- \centering
- \begin{tikzpicture}[->,>=stealth',shorten >=1pt,auto,node distance=1.5cm,
- thick,main node/.style={circle,draw}]
-
- \node[main node] (1) {u};
- \node[main node] (2) [below of=1] {};
- \node[main node] (3) [right of=1] {};
- \node[main node] (4) [left of=1] {};
- \node[main node] (5) [above of=1] {};
- \node[main node] (6) [below right=1.5cm of 1] {};
- \node[main node] (7) [above right=1.5cm of 1] {};
- \node[main node] (8) [above left=1.5cm of 1] {};
- \node[main node] (9) [below left=1.5cm of 1] {};
-
- \path[every node/.style={font=\sffamily\small}]
- (1) edge node [right] {} (2)
- (1) edge node [right] {} (3)
- (1) edge node [left] {} (4)
- (1) edge node [left] {} (5)
- (1) edge node [right] {} (6)
- (1) edge node [right] {} (7)
- (1) edge node [left] {} (8)
- (1) edge node [left] {} (9);
- \end{tikzpicture}
- \caption{Size eight neighborhood.}
- \label{fig:c8_grid}
- \end{subfigure}
- \caption{Two common neighborhood stencils.}
- \label{fig:common_neighborhoods}
-\end{figure}
-We discretize the total variation in \eqref{eq:coarea_formula} as
-follows
-\begin{equation}
- \TV(u)
- = \sum_{\lambda = 0}^{L-2} P(u^\lambda)
- = \sum_{\lambda = 0}^{L-2} \sum_{(u,v)} w_{uv}
- \abs{u_u^\lambda - u_v^\lambda},
- \label{eq:tv_discrete_int}
-\end{equation}
-where the second sum is over nodes in a neighborhood relation and the
-$w_{uv}$ is the weight of this relation. The sum only goes up to the
-level $L-2$ since $u^\lambda$ is constant equal to 1 in every pixel of
-the image. Figure \ref{fig:common_neighborhoods} shows the two most used
-neighbordhood stencils, which will be discussed more thorougly in
-\fixme{Section ??}. Boykov and Kolmogorov argue in
-\cite{boykov2003computing} that if the weight is chosen as
-\begin{equation}
- w_{uv} = \frac{h^2 \cdot \Delta \phi_{uv}}{2 \cdot \abs{d_{uv}}},
-\end{equation}
-the discrete perimeter in \eqref{eq:tv_discrete_int} converges to the
-continuous perimeter in \eqref{eq:coarea_formula}. Here, $h$ is the grid
-size and $d_{uv}$ is the euclidean distance of the edge. The $\Delta
-\phi_{uv}$ is the difference between the angle of this edge and the
-next edge, if the edges are sorted by increasing angles. These
-parameters are also shown in Figure \ref{fig:common_neighborhoods}.
-Boykov and Kolmogorov prove we have convergence when all of $h$, $\Delta
-\phi_{uv}$, and $d_{uv}$ go to zero.
-\fixme{Maybe move (and expand) this somewhere else? Where the grids are
-discussed for example? Also, the section is a bit messy.}
-
-\fixme{BAD TRANSITION}
-
-For two binary variables $a$ and $b$ we can easily verify that $\abs{a -
-b} = a + b - 2 a b$. Using this we rewrite \eqref{eq:tv_discrete_int} to
-\begin{equation}
- \TV(u) = \sum_{\lambda = 0}^{L-2} \sum_{(u,v)} w_{uv}
- \left(
- \left(1 - 2 u_t^\lambda \right) u_s^\lambda + u_t^\lambda
- \right).
- \label{eq:tv_discrete}
-\end{equation}
-
-We define the following function for some pixel value $x$ and some pixel
-position $u$ in the original image $v$
-\begin{equation}
- N_u(x) = \abs{x - v_u}^p,
-\end{equation}
-such that the $\int_\Omega \abs{u - v}^p$ of \eqref{eq:energy_perimeter}
-can be discretized in the following way \fixme{Give this a name?}
-\begin{equation}
- \sum_u \abs{u_u - v_u}^p = \sum_u N_u(u_u)
- \label{eq:norm_discrete_int}.
-\end{equation}
-We want to write the energy as a sum over the different levels of the
-image, so we use the following decomposition formula which holds for any
-function $N(k)$ taking values $k \in \mathcal{L}$
-\begin{equation}
-\begin{aligned}
- N(k) &= \sum_{\lambda=0}^{k-1} \big(
- N(\lambda + 1) - N(\lambda)
- \big) + N(0) \\
- &= \sum_{\lambda=0}^{L-2} \big(
- N(\lambda + 1) - N(\lambda)
- \big) \mathbbm{1}_{\lambda < k} + N(0) \\
- &= \sum_{\lambda=0}^{L-2} \big(
- N(\lambda + 1) - N(\lambda)
- \big) (1 - u^\lambda) + N(0).
-\end{aligned}
-\end{equation}
-This allows us to rewrite \eqref{eq:norm_discrete_int} and obtain
-\begin{equation}
- \sum_{\lambda=0}^{L-2} \sum_u
- \big(
- N_u(\lambda + 1) -
- N_u(\lambda)
- \big)
- (1 - u^\lambda) + N_u(0)
-\end{equation}
-We have now discretized the energy function and decomposed it into a sum
-over all the levels $\lambda$ of the image. If we minimize each level
-separately it is obvious that we also minimize the total energy. The
-question is if the obtained level sets $u^\lambda$ can be combined to
-make an output image. The level sets were defined as $u^\lambda =
-\mathbbm{1}_{u \leq \lambda}$, so we need them to be monotonely
-increasing in increasing level values, i.e.\
-\begin{equation}
- u^\lambda_u \leq u^\mu_u \quad \forall \lambda \leq \mu,
- \quad \forall u \in S.
-\end{equation}
-
-\fixme{Want to use $u$ and $v$ as pixels/nodes, what to use as image?
-Since $s$ and $t$ normally is the source and sink.}
-
-\fixme{Some like to call one term the fidelity term, and the other the
-regularization term?}
-
-\fixme{\cite{boykov2003computing} motivates the discrete estimation of
- the perimeter of a set, and it makes sense. The size of the
- perimeter should be proportional to the number of edges it crosses,
- which is why we sum up $\abs{u^\lambda_s - u^\lambda_t}$ (which is
- only different from zero at the cut).
-}
-
-\fixme{How do we get to the discrete formulation. Why is it OK to do
- what we do with the gradient. How do we decompose the energy
- function into different levels and why is it OK to optimize each
- level separately. How do we choose the different neighbourhoods, and
- what could this mean for the result. Darbon and Sigelle.
- \cite{chan2011numerical} describes different methods for using total
- variation, and especially how one can use it in the discrete
- setting. It also goes into the min-cut method. Yay. Good book.
-}
-
-\subsection{Graph cut formulation}
-Short introduction to graphs. Notation. What is a cut, what is a minimal
-s-t-cut. How do we get from the discrete energy to the graph and why
-does finding the graph cut result in minimizing the energy of the
-current label. Boykov and Kolmogorov.
-
-\subsubsection{Min-cut max-flow theorem}
-Want to find the minimal s-t-cut, so we introduce flow here maybe? Or
-earlier? Finding the max flow results in also finding the minimal cut,
-and there are «good» algorithms for finding the maximal flow.
-
-\subsection{Maximum flow algorithms}
-Two different ideas. One is to always maintain a valid flow in the
-network, Ford Foulkerson, always looking for new paths from s to t along
-which we can send additional flow. This leads to Edmonds-Karp, and
-Dinic/Dinitz. Second idea is push-relabel, allowing positive excess in
-the nodes. We lift the source, and let as much flow as possible flow
-towards the sink.
-
-\subsubsection{Dinic's algorithm}
-The «best» of the augmenting-path algorithms (citation needed). Makes a
-level graph with BFS, and finds a blocking flow through this graph which
-only goes from nodes with one label to nodes with greater labels. When
-this blocking flow is found we restart at the BFS.
-
-\subsubsection{The push-relabel algorithm}
-The nodes are allowed to have a positive excess but we still follow the
-capacity constraints. The nodes also have a labeling, which has to be
-valid. Push flow from active nodes, and relabel, preferrably in a
-specific order, until it is not possible anymore. The minimal cut is
-actually found before the flow becomes valid (has to be explained).
-
-Because of the way the graph evolves, we can reuse the graph between the
-different labels.
-
-Initially proposed in \cite{goldberg1988new} where one can also find
-proof that the algorithm maintains a valid labeling and that it
-terminates at at maximal flow. Also has the lemma that if you have a
-preflow, then there is no path from $s$ to $t$ in the residual graph.
-But what is a preflow? Also, the distance label never decreases. Also
-propose the global relabeling heuristic.
-
-Shiloach and Vishkin originally proposed the FIFO way of doing things
-(according to Golberg and Tarjan)?
-
-Cherkassky and Goldberg present some good implementation tips in
-\cite{cherkassky1997implementing}, for example the global and gap
-relabeling heuristics. Also has OK presentations of the different ways
-of choosing the order of discharge operations. In the results I should
-include some results with and without the heuristics.
-
-\cite{ahuja1997computational} has some good performance analysis.
-
-\cite{boykov2006graph} is a good source on how graph cuts are used in
-image analysis. It considers both one level graphs, and multilevel
-graphs. Which is nice. It has just a small section on the different
-algorithms that exist and that the special Boykov Kolmogorov algorithm
-is the best.
-
-\cite{derigs1989implementing} is also a good source for the different
-heuristics. And also on the neighbor list that should not be reset after
-each discharge (should try to implement this again, and reset at every
-relabel!). Also has some semi-usable notes on how the minimal cut is
-found after phase 1 of the algorithm is finished. Note that in this
-article the gap relabeling heuristic is also called RELABEL-GLOBAL.
-
-\begin{algorithm}
-\begin{algorithmic}
- \Function{Push}{$u$}
- \ForAll{$v$ neighbour of $u$}
- \If{$res(u, v) > 0$}
- \State $f \gets min(res(u, v), excess[u])$
- \State $flow(u, v) \mathrel{+}= f$
- \State $flow(v, u) \mathrel{-}= f$
- \State $excess[u] \mathrel{-}= f$
- \State $excess[v] \mathrel{+}= f$
- \EndIf
- \EndFor
- \EndFunction
-\end{algorithmic}
-\caption{\sf The push procedure of the Push-Relabel algorithm}
-\end{algorithm}
-
-\begin{algorithm}
-\begin{algorithmic}
- \Function{Relabel}{$u$}
- \If{$u$ is only node at its height}
- \Call{Gap}{$u$}
- \Else
- \State $height[u] \gets min(height[v] \, \forall v \in
- neighbors[u] : res(u, v) > 0) + 1$
- \EndIf
- \EndFunction
-\end{algorithmic}
-\caption{\sf The relabel procedure of the Push-Relabel algorithm}
-\end{algorithm}
-
-\begin{algorithm}
-\begin{algorithmic}
- \Function{Gap}{$u$}
- \State $k \gets height[u]$
- \ForAll{$v$ with height $\geq k$}
- \State $height[v] \gets N$
- \EndFor
- \EndFunction
-\end{algorithmic}
-\caption{\sf The gap procedure of the Push-Relabel algorithm}
-\end{algorithm}
-
-\subsubsection{Other algorithms}
-Special made augmenting flow algorithm for this type of graphs. Used in
-software1, software2, etc. Is fast.
-