the anisotropic total variation
\begin{equation}
\TVA(u) = \int_\Omega \sqrt{\nabla u(x)^T A(x) \nabla u(x)} \, dx
+ \label{eq:aniso_tv_sqrt}
\end{equation}
for all $u \in C^1(\Omega)$. If $A(x)$ is the identity
matrix we get the regular total variation found in
\leq 1 \right\}.
\label{eq:extended_tv}
\end{equation}
+ \label{def:extended_tv}
\end{definition}
+With this extended definition, we have arrived at a minimization problem
+where we seek to find a minimizer of the functional
+\begin{equation}
+ F(u) = \int_\Omega (u - f)^2 \, dx + \beta \, \TVA(u)
+\end{equation}
The question is then how to construct this anisotropy tensor $A(x)$ to
get the improvements we hope for, and how this affects our numerical
solution methods.
-\fixme{Include the whole functional somewhere.}
\fixme{Describe how this is similar to a Euler time step in a related
equation \cite{grasmair2010anisotropic}. This is only without the square
root.}
\subsection{Anisotropy tensor}
+\label{sec:anisotropy_tensor}
\fixme{
This tensor can be visualized, and can be used for edge
implementation-focused chapter.
}
-There are many possible choices for the anisotropy tensor $A(x)$, as
-long as one fulfills the assumptions of \fixme{ref}. Additionally we
-will later need the eigenvalues to be bounded above and bounded away
-from zero.
+There are many possible choices for the anisotropy tensor $A(x)$. Our
+constraints are that we have assumed it to be symmetric positive
+definite, and we have some wishes for its properties. We would first and
+foremost like it to down-weight the effect of $\nabla u$ in
+\eqref{eq:aniso_tv_sqrt} across true edges, while maintaining normal
+regularization properties in smooth sections.
-Our wish for the metric tensor is that it down-weighs the variation
-across true edge in the image, meaning that true edges do not contribute
-as much in the total variation as they normally would. This in turn
-means that these edges will be less regularized. By true edges we mean
-that we do not want the tensor to be sensitive to noise in the image, so
-we somehow want to be sure about the edges we find. \fixme{Maybe refer
-to some figure}
+By true edges we mean that that we do not want the tensor to be
+sensitive to noise in the image, and thus find edges where there are
+none, so we somehow want to be sure about the edges we find.
Edges can be found in many different ways, but as suggested by Weickert
-in his book on Anisotropic Diffusion \cite{weickert1998anisotropic}, a
-good starting point is the structure tensor
+in his book on Anisotropic Diffusion \cite{weickert1998anisotropic}, and
+briefly mentioned in Section \ref{sec:anisotropic_diffusion}, a good
+starting point is the structure tensor
\begin{equation}
- S(x) := K_\rho * \left( \nabla \tilde{f}_\sigma(x) \otimes
+ S_\rho(x) := K_\rho * \left( \nabla \tilde{f}_\sigma(x) \otimes
\nabla \tilde{f}_\sigma(x) \right),
\label{eq:s_def}
\end{equation}
-where $\tilde{f}$ is the \fixme{symmetric} extension of the initial
-image from $\Omega$ to $\mathbb{R}^2$ and the $\sigma$ signifies a
-Gaussian smoothing $\tilde{f}_\sigma = K_\sigma * \tilde{f}$. The
-initial smoothing parameter $\sigma$, called the \emph{noise scale} by
-Weickert, will control how sensitive we are to noise and small details.
-
-The convolution with $K_\rho$ is done componentwise, and $\rho$, the
-\emph{integration scale} thus controls the size of the neighborhood from
-which the structure tensor is computed. It controls the size of the
-structures we want our anisotropy tensor to consider.
-
-But we can not use the structure tensor $S(x)$ directly as our
-anisotropy tensor $A(x)$, because $S(x)$ will be an operator with its
-\fixme{principal direction} in the direction of $\nabla
-\tilde{f}_\sigma$ which is the opposite of what we want. \fixme{more
-inspiration from weickert here, he has some good stuff.} The structure
-tensor $S(x)$ is real and symmetric, so it admits an eigendecomposition
+where $\tilde{f}$ is the symmetric extension of the initial image $f$ in
+$\mathbb{R}^2$. The $\sigma$ signifies a Gaussian smoothing
+$\tilde{f}_\sigma = K_\sigma * \tilde{f}$. This initial smoothing
+parameter $\sigma$, called the \emph{noise scale} by Weickert, will
+control how sensitive the edge detector $\nabla \tilde{f}_\sigma$ is to
+noise and small details.
+
+First consider the tensor $S_0(x) = \nabla \tilde{f}_\sigma(x) \otimes
+\nabla \tilde{f}_\sigma(x)$. It obviously contains the same information
+as the edge detector itself. Its eigenvalues will be 0 and $\labs{\nabla
+\tilde{f}_\sigma(x)}^2$ with corresponding eigenvectors $v_1$ and $v_2$
+perpendicular and parallel to $\nabla \tilde{f}_\sigma(x)$ repectively.
+
+We are also interested in identifying and being sensitive to features in
+a neighborhood around the point, such as corners or curved edges and
+\fixme{coherent structures}. This is why we introduce the convolution
+with $K_\rho$, which is done componentwise. The parameter $\rho$, called
+the \emph{integration scale} thus controls the size of the neighborhood
+which affects the structure tensor. It controls the size of the
+structures we want our anisotropy tensor to be sensitive to.
+
+Being a smoothed combination of symmetric positive semi-definite
+matrices, $S_\rho(x)$ can easily be verified to be symmetric positive
+semi-definite itself. We order its two real eigenvalues $\lambda_1 \geq
+\lambda_2$ with corresponding eigenvectors $v_1$ and $v_2$. 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
+\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$.
+
+These eigenvalues indicate that the product $\sqrt{\nabla u^T S_\rho(x)
+\nabla u}$ would amplify the effect of $\nabla u$ across edges, the
+opposite of what we want.
+
+We begin our anisotropy tensor construction by eigendecomposing the
+structure tensor as
\begin{equation}
- S(x) = U(x) \Lambda(x) U(x)^T,
+ S_\rho(x) = U(x) \Lambda(x) U(x)^T,
\end{equation}
-where $\Lambda(x)$ is a matrix with the eigenvalues of $S(x)$ on the
+where $\Lambda(x)$ is a matrix with the eigenvalues of $S_\rho(x)$ on the
diagonal, while $U(x)$ is a rotation matrix and has the eigenvectors of
-$S(x)$ as its columns.
+$S_\rho(x)$ as its columns.
-If we consider the structure tensor of \eqref{eq:s_def} before
-convolving with the Gaussian function, it is clear that there will be
-maximum one non-zero eigenvalue, whose eigenvector is the vector
-$\nabla u_\sigma(x)$ itself. If $\nabla u_\sigma(x) = (u_x, u_y)^T$ the
-eigenvalue equation becomes
-\begin{equation}
- \begin{pmatrix}
- u_x^2 & u_x u_y \\
- u_x u_y & u_y^2
- \end{pmatrix}
- \begin{pmatrix}
- u_x \\
- u_y
- \end{pmatrix} =
- (u_x^2 + u_y^2)
- \begin{pmatrix}
- u_x \\
- u_y
- \end{pmatrix}
-\end{equation}
-and unless the smoothing parameter $\rho$ is very large, the eigenvalues
-of $S(x)$ will always be dominated by the direction of $\nabla
-u_\sigma(x)$. This means that $S(x)$ is not at all suitable for being
-our metric tensor, since it actually weights the direction of $\nabla
-u_\sigma(x)$ more. We only need to change around the eigenvalues a bit.
-If the eigen-decomposition gave us
+The eigen-decomposition gave us
\begin{equation}
\Lambda(x) = \begin{pmatrix}
\lambda_1 & 0 \\
U(x) \Sigma(x) U(x)^T$ where
\begin{equation}
\Sigma(x) = \begin{pmatrix}
- \sigma_1(x) & 0 \\
- 0 & \sigma_2(x)
+ \sigma_1 & 0 \\
+ 0 & \sigma_2
\end{pmatrix}.
\end{equation}
-To reverse the effect of $S(x)$ such that $M(x)$ actually weighs down
+To reverse the effect of $S_\rho(x)$ such that $A(x)$ actually weighs down
the norm across edges, we need $\sigma_1 \leq \sigma_2$, so we construct
\begin{align}
\sigma_1 &= \left(1 + \frac{(\lambda_1 -
\TVA(u) = \sup_{\norm{\xi}_A^* \leq 1} \int_\Omega u \diver \xi \,
dx
\end{equation}
-can be thought of as, and has the properties of a norm, and is therefore
+can be thought of as - and has the properties of - a norm, and is therefore
convex. The sum of the two is thus strictly convex, which, given the
existence of a minimizer, implies uniqueness.
\subsection{Coercivity}
We need coercivity to show that we cannot go further and further away to
-obtain a better and better solution. This means that $\norm{u}_{L^2} \to
+obtain a better and better solution. This means that $\lnorm{u}_{L^2} \to
\infty$ should imply that $F(u) \to \infty$, which is obvious from the
fidelity term for some fixed $v \in L^2(\Omega)$.
This is again is a
supremum of weakly continuous functionals. Thus the regularization term
is also sequentially weakly lower semi-continuous.
-\fixme{Could consider writing the supremums as a supremum over a set
-with $\{ \}$ instead of all the subscripts!}
-\fixme{It would probably be most correct to say: weakly sequentially
-lower semi-continuous?}
The usual ways of going from coercivity and lower semicontinuity to
existence do not work in infinite dimensions. But with our coercivity
\end{theorem}
\begin{proof}
Assume that $u \in C^1(\Omega) \cap \BV(\Omega)$. The extension to
- all functions $u \in \BV(\Omega)$ complicates things and will not be
- considered here.
+ all functions $u \in \BV(\Omega)$ can be made by approximation
+ arguments but will not be considered here.
\paragraph{First we prove that $\TVA(u) \leq \int_{-\infty}^\infty
\TVA(u^s) \, ds$.}
\PerA(U;\Omega) = \TVA(\idfun_U).
\end{equation}
A measurable set $U \subset \Omega$ is of finite anisotropic
- perimeter in $\Omega$ if $\idfun_U \in \BV(\Omega)$.
+ perimeter in $\Omega$ if $\idfun_U \in \BV(\Omega)$. \fixme{nope?}
\end{definition}
The anisotropic set perimeter is not like the regular set perimeter and
does not measure the length of the boundary of the set, but it can be
unit normal vector. The length of the boundary would normally be
calculated by integrating the norm of the tangent vector. The connection
will be made after the Cauchy--Crofton formulas have been introduced.
+
+The anisotropic coarea formula allows us to transform the anisotropic
+total variation such that we are left with minimizing the following
+functional
+\begin{equation}
+ F(u) = \int_\Omega (u - f)^2 \, dx + \beta \, \int_{-\infty}^\infty
+ \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
+perimeter will be computed using a graph cut machinery.
+
\fixme{well, if we assumed differentiability, we wouldn't need the sup
definition of the TV...}
-\fixme{and then some conclusion on where we are at.}
-
\section{Cauchy--Crofton formulas}
\begin{figure}
\input{fig/line_param}
\end{figure}
-\fixme{Something something curious geometric measure theory interesting
-integral formulas etc, introduction one or two sentences}
+In the fields of integral theory and geometric measure theory there are
+a number of interesting integral formulas. Several of them fall in a
+category often refered 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 a curve by
+counting the times it intersects line in the set of all lines.
+Intuitively, a long curve will intersect more lines.
We write $\mathcal{L}$ for the set of all straight lines in the plane,
and parametrize them as shown in Figure \ref{fig:line_param}. Thus a
-line is parametrized by the angle $\phi$ of the normal going to the
-origin, and the distance $\rho$ from origin to the line. Sometimes it is
+line is parametrized by the angle $\phi \in [0, 2\pi)$ of the normal going to the
+origin, and the distance $\rho \in [0, \infty)$ from origin to the line. Sometimes it is
more convenient to consider a unit vector $\nu$ giving the direction of
-the line instead of the angle parameter $\phi$. We will write a line
+the line instead of the angle parameter $\phi$. We will write a line
$\ell_{\phi, \rho} = \ell_{\nu, \rho}$ where $\nu$ is a unit vector
along the line, i.e.\ $\nu = (-\sin \phi, \cos \phi)^T$. By defining the
measure on this set $d\mathcal{L} = \dpdr$ we are ready to introduce the
\label{thm:euclidean_cauchy_crofton}
\end{theorem}
This elegant formula is very useful when we later will discretize our
-energy function. The set of lines $\mathcal{L}$ is then discretized
-depending on our choice of neighborhood, and the length of the curve $C$
-can be approximated by a sum.
+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.
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
M\big(\gamma(t)\big) \, \dot{\gamma} \rangle} \, dt.
\label{eq:riemannian_length}
\end{equation}
-We will now prove a Cauchy--Crofton formula in this case where we have a
-metric tensor in each point.
+We will now present and prove a Cauchy--Crofton formula in this case
+where our domain is equipped with a metric tensor in each point.
\begin{theorem}[The Riemannian Cauchy--Crofton formula]
Assume that our space $\Omega$ is equipped with a continous metric
tensor $M(x)$, whose eigenvalues are bounded $0 < k \leq
\label{eq:piecewise_curve_length}
\end{equation}
We know that $M_\pi(x)$ converges pointwise to $M(x)$, and thus
- $\abs{\dot{C}(t)}_{M_\pi}$ converges pointwise to
- $\abs{\dot{C}(t)}_M$. We have assumed bounds on the eigenvalues
+ $\labs{\dot{C}(t)}_{M_\pi}$ converges pointwise to
+ $\labs{\dot{C}(t)}_M$. We have assumed bounds on the eigenvalues
of $M(x)$ such that according to the Rayleigh principle
\begin{equation}
K \geq \lambda_1 = \max_\xi \frac{\xi^T M_\pi(x) \xi}{ \xi^T \xi
$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
- $\abs{C}_{M_\pi} \to \abs{C}_M$.
-
- \fixme{ugly line spacing}
+ $\labs{C}_{M_\pi} \to \labs{C}_M$.
We apply the same theorem to show that the right hand side of
\eqref{eq:riemannian_segments} converges. Recall the definition of
\end{equation}
\end{proof}
-\section{Yeeep}
-
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. In the Euclidean case, which here would amount to
-setting the anisotropy tensor $A$ equal to the identity matrix $I$, the
-perimeter coincides nicely with the length of the boundary curve,
-assuming some regularity for the boundary. In the general case we need
+set of the image as introduced in \eqref{eq:functional_per}. In the
+Euclidean case, which here would amount to setting the anisotropy tensor
+$A$ equal to the identity matrix $I$, the perimeter coincides nicely
+with the length of the boundary curve, assuming some regularity for the
+boundary.
+
+In the general case we need
to be more careful. As can be seen in \eqref{eq:perimeter_calc}, the
anisotropic perimeter is calculated by integrating the norm of the
\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. If $P$ is a 90\textdegree{}
-rotation matrix and we assume an arclength parametrisation of our curve
-we have
+the \emph{tangent vector} of the curve. The normal and tangent vector
+are always perpendicular, so all we need is a 90\textdegree{} rotation.
+If $P$ is a 90\textdegree{} rotation matrix we have
\begin{equation}
\begin{aligned}
\PerA(U; \Omega) &= \int_{\partial U} \sqrt{ \langle
\end{aligned}
\label{eq:per_to_length1}
\end{equation}
-Before applying the Cauchy--Crofton formula we simplify the equation by
-defining the metric tensor $M(x) = P A(x) P^T$ and $\gamma = \partial U$
-and observing that a 90\textdegree{} rotation of the normal $\nu_{\partial
-U}$ gives us the tangent $\pm \dot{\gamma}$. This yields
+We simplify the equation by defining the metric tensor $M(x) = P A(x)
+P^T$ and let $\gamma = \partial U$ be an arclength parametrization of
+the boundary of $U$. Observe that a 90\textdegree{} rotation of the
+normal $\nu_{\partial U}$ gives us the tangent $\pm \dot{\gamma}$.
+Inserting this into \eqref{eq:per_to_length1} we get
\begin{equation}
\PerA(U; \Omega) = \int_\gamma \sqrt{ \langle \dot{\gamma}, M(x) \,
\dot{\gamma} \rangle }.
\label{eq:per_to_length2}
\end{equation}
-\fixme{what about the boundary of $\Omega$?}
+Note that since $\Omega$ is open, we assume that $\gamma = \partial U
+\subset \Omega$ and thus $\gamma$ can not include parts of the boundary
+of $\Omega$.
+
Now we make sure that all the assumptions of the Riemannian
Cauchy--Crofton formula in Theorem \ref{thm:riemannian_cauchy_crofton}
-are fulfilled so that it can be applied.
+are fulfilled so that it can be applied to the curve length integral we
+have constructed in \eqref{eq:per_to_length2}.
-The structure tensor is constructed as described in \fixme{ref}
+The structure tensor is constructed as described in Section
+\ref{sec:anisotropy_tensor}
\begin{equation}
- S(x) = \left(K_{\rho} * \left( \nabla \tilde{f}_{\sigma} \otimes \nabla
- \tilde{f}_{\sigma} \right)\right)(x),
+ S_\rho(x) = \left(K_{\rho} * \left( \nabla \tilde{f}_{\sigma}
+ \otimes \nabla \tilde{f}_{\sigma} \right)\right)(x),
\end{equation}
-where $\tilde{f}$ is the symmetric extension of our input image $f$ to
-$\mathbb{R}^2$.
Because of the convolutions with the Gaussian function, this is a smooth
-continuous map from $\bar{\Omega}$ to $\mathbb{R}^{2\times 2}$. As we
+continuous map from $\bar{\Omega}$ to $\mathbb{R}^{2 \times 2}$. As we
can see in \fixme{ref} the eigenvalues are the roots of a monic
polynomial and thus depend continuously on the coefficients of the
polynomial, which in turn are continuous functions of the elements in
-the structure tensor $S(x)$. The extreme value theorem \fixme{ref}
+the structure tensor $S_\rho(x)$. The extreme value theorem \fixme{ref}
states that a continuous real-valued function on a nonempty compact
-space is bounded above. Thus the eigenvalues $s_1$ and $s_2$ of $S(x)$
+space is bounded above. Thus the eigenvalues $s_1$ and $s_2$ of $S_\rho(x)$
are bounded from above and by construction the smallest eigenvalue of
our anisotropy tensor $A(x)$ is bounded away from zero as
\begin{equation}
projects / master / commitdiff