Shape reconstruction in X-ray tomography from a small number of projections using deformable models
Ali Mohammad--Djafari$^1$ \quad and \quad Ken Sauer$^2$
[12pt] $^1$ Laboratoire des Signaux et Syst\`emes (CNRS--SUPELEC--UPS)
\'Ecole Sup\'erieure d'\'Electricit\'e
Plateau de Moulon
91192 Gif--sur--Yvette Cedex
France
[12pt] $^2$ Department of Electrical Engineering
Notre Dame
IN 46556
USA
} \def\eabs{
latexdef.tex
%LATEXDEF.LATEX
%LaTeX document
%Author: Ali Mohammad-djafari
%Date  : 22 Avril 1991.
%
\def\bdoc{\begin{document}}
\def\edoc{\end{document}}

\def\babs{\begin{abstract}}
\def\eabs{\end{abstract}}

\def\bbib{\begin{thebibliography}{999}}
\def\ebib{\end{thebibliography}}

\def\beq{\begin{equation}}
\def\eeq{\end{equation}}

\def\beqn{\begin{eqnarray}}
\def\eeqn{\end{eqnarray}}

\def\beqnn{\begin{eqnarray*}}
\def\eeqnn{\end{eqnarray*}}

\def\barr{\begin{array}}
\def\earr{\end{array}}

\def\bqu{\begin{quote}}
\def\equ{\end{quote}}

\def\bqun{\begin{quotation}}
\def\equn{\end{quotation}}

\def\bit{\begin{itemize}}
\def\eit{\end{itemize}}

\def\ben{\begin{enumerate}}
\def\een{\end{enumerate}}

\def\bpic{\begin{picture}}
\def\epic{\end{picture}}

\def\bfig{\protect \begin{figure}}
\def\efig{\protect \end{figure}}

\def\bcc{\begin{center}}
\def\ecc{\end{center}}

\def\brr{\begin{flushright}}
\def\err{\end{flushright}}

\def\bll{\begin{flushleft}}
\def\ell{\end{flushleft}}

\def\btab{\begin{tabular}}
\def\etab{\end{tabular}}

\def\ul{\underline}
\def\ol{\overline}


month-bi.tex
% TEX 7(ascii) bits
%--------------------- month-fr.tex ---------------------------
% Fichier de d\'efinition pour les mois, en bilingue
%
\def\jan{janvier\xspace}
\def\feb{f\'evrier\xspace}
\def\mar{mars\xspace}
\def\apr{avril\xspace}
\def\may{mai\xspace}
\def\jun{juin\xspace}
\def\jul{juillet\xspace}
\def\aug{ao\^ut\xspace}
\def\sep{septembre\xspace}
\def\oct{octobre\xspace}
\def\nov{novembre\xspace}
\def\dec{d\'ecembre\xspace}
%
\def\Jan{January\xspace}	
\def\Feb{February\xspace}
\def\Mar{March\xspace}
\def\Apr{April\xspace}
\def\May{May\xspace}
\def\Jun{June\xspace}
\def\Jul{July\xspace}
\def\Aug{August\xspace}
\def\Sep{September\xspace}
\def\Oct{October\xspace}
\def\Nov{November\xspace}
\def\Dec{December\xspace}
% plus
\def\sub{soumis \`a\xspace}
paper1.tex
%paper1.tex
\documentclass{article}
\usepackage{kenmax,xspace,graphicx}
\graphicspath{{Eps/}}

\input{latexdef.tex}
\input{trdef.tex}
\input{spgraph.tex}

\title{Shape reconstruction in X-ray tomography from a small number of
projections using deformable models}

\author{Ali Mohammad--Djafari$^1$ \quad and \quad Ken Sauer$^2$\\[12pt]
$^1$ Laboratoire des Signaux et Syst\`emes (CNRS--SUPELEC--UPS) \\
\'Ecole Sup\'erieure d'\'Electricit\'e \\
Plateau de Moulon, 91192 Gif--sur--Yvette Cedex, France \\[12pt]
$^2$ Department of Electrical Engineering \\
University of Notre Dame \\
Notre Dame, IN 46556, USA
}
\date{}

\begin{document}
\maketitle

\index{Tomography}
\index{Computed tomography}
\index{Deformable templates}
\index{Shape reconstruction}

\begin{abstract}
X-ray tomographic image reconstruction consists of determining an
object function from its projections. In many applications such as
non-destructive testing, we look for a fault region (air) in a
homogeneous, known background (metal). The image reconstruction problem
then becomes the determination of the shape of the default region.
Two approaches can be used: modeling the image as a binary Markov
random field and estimating the pixels of the image, or modeling
the shape of the fault and estimating it directly from the
projections. In this work we model the fault shape by a deformable
polygonal disc or a deformable polyhedral volume
and propose a new method for directly estimating the coordinates
of its vertices from a very limited number of its projections.
The basic idea is not new, but in other competing methods, in general,
the fault shape is modeled by a small number of parameters (polygonal
shapes with very small number of vertices, snakes and deformable
templates) and these parameters are estimated either by least squares
or by maximum likelihood methods. We propose modeling the shape
of the fault region by a polygon with a large number of vertices,
allowing modeling of nearly any shape and estimation of its vertices'
coordinates directly from the projections by defining the solution as
the minimizer of an appropriate regularized criterion. This formulation
can also be interpreted as a maximum a posteriori (MAP) estimate in a
Bayesian estimation framework. To optimize this criterion we use either a
simulated annealing or a special purpose deterministic algorithm based
on iterated conditional modes (ICM).
The simulated results are very encouraging, especially when the number
and the angles of projections are very limited.
\\�~\\
{\bf key words:}~Computed tomography, Shape reconstruction, non destructive
testing, Bayesian MAP estimation
\end{abstract}

\section{Introduction}
Tomographic image reconstruction in non destructive testing (NDT) is recent and
consists in determining an object $f(x,y)$ from its projects $p(r,\phi)$:

\beq \label{eq.1}
p(r,\phi)=\intd f(x,y) \delta(r-x\cosf-y\sinf) \d{x}\d{y}
\eeq
In many image reconstruction applications, especially in NDT, we
know that $f(x,y)$ has a constant value $c_1$ inside a region $P$ (fault
region) and another constant value $c_2$ outside that region (safe or
background region), e.g. metal \& air.

\beq \label{eq.2}
f(x,y)=\left\{\barr{ll}
c_{1} & \hbox{if}\quad (x,y)\in P, \\
c_{2} & \hbox{elsewhere} \earr \right.
\eeq
The image reconstruction problem becomes then the determination of the shape
of the fault region $P$.
In this work, without loss of generality, we assume that $c_1=1$ and $c_2=0$
and model the shape of the object by its contour.

There has been many works in image reconstruction dealing with this problem.
To emphasis the originality and the place of this work, we give here a
summary of the different approaches for this problem:

\medskip\noindent$\bullet$
The first approach consists in discretizing the equation (\ref{eq.1})
to obtain:

\beq \label{eq.3}
\pb=\Hb\fb+\nb
\eeq
where,
$\fb$ is the discretized values of the object $f(x,y)$ (the pixel values of
the image),
$\pb$ is values of the projection data $p(r,\phi)$,
$\nb$ is a vector to represent the modeling and measurement errors (noise) and
$\Hb$ the discretized Radon operator.
Then the solution is defined as the minimizer of a compound criterion

\beq \label{eq.4}
J(\fb)=Q(\pb - \Hb \fb) + \lambda \Omega(\fb),
\eeq
where $\lambda$ is the regularization parameter and $Q$ and $\Omega$
has to be chosen appropriately to reflect our prior knowledge on the
noise and on the image. This is the classical regularization approach
of general image reconstruction problem.
One can also interpret $J(\fb)$ as the maximum a posteriori (MAP)
criterion in the Bayesian estimation framework where
$\expf{-Q(\fb)}$ represents the likelihood term and $\expf{-\Omega(\fb)}$
the prior probability law.

This approach has been used with success in many applications
(e.g. \cite{Herman87,Kak87,Geman87,Demoment89})
but the cost of its calculation is huge due to the great dimension of $\fb$.
Many works have been done on choosing appropriate regularization functionals or
equivalently appropriate prior probability laws for $\fb$ to enforce some
special properties of the image such as smoothness, positivity or piecewise
smoothness
\cite{Brette94a,Djafari96a,Bouman93,Bedini94,Nikolova95d}.
Among these, one can mention mainly two types of functions for $\Omega(\fb)$:

Entropic laws:
\beq \label{eq.5}
\Omega(\fb)=\sum_{j=1}^N \phi(f_j)
\quad\hbox{with}\quad
\phi(x)=\left\{|x|^p, \, -x\log x, \, \log x, \, \cdots\right\}
\eeq

Markovian laws:
\beq \label{es.6}
\Omega(\fb)=\sum_{j=1}^N \sum_{i\in{\cal N}_j} \phi(f_j,f_i)
\hspace*{5mm}
\eeq
with convex potential functions:
\beq \label{eq.7}
\phi(x,y)=\left\{|x-y|^p,\,  -|x-y|\log\frac{x}{y}, \, \log\cosh|x-y|, \cdots\right\}
\eeq
or non convex potential functions:
\beq \label{eq.8}
\phi(x,y)=\left\{\min\{|x-y|^2,1\}, \, \frac{-1}{1+|x-y|^2}, \cdots\right\}
\eeq
See for example \cite{Djafari96a} for the entropic laws,
\cite{Bouman93,Brette94a} for scale invariant markovian laws with
convex potential functions and
\cite{Bedini94,Nikolova94,Nikolova95d,Nikolova96a} for markovian laws with
non convex potential functions and other specific choices.

\medskip\noindent$\bullet$
The second approach consists in modeling directly the closed contour of the
object as the zero-crossing of a smooth function $u(x,y)$:

\beq
\partial D=\left\{(x,y): u(x,y)=0\right\}
\quad\hbox{and}\quad
f(x,y)=\left\{\barr{ll}
c_{1} & \hbox{if~} u(x,y) > 0, \\
c_{2} & \hbox{if~} u(x,y) < 0
\earr\right.,
\eeq
and in defining a time evolution for $u$
(consequently the corresponding contour $\partial D$ and function $f$)
such that

\beq
\partial D(t)=\left\{(x,y): u(x,y,t)=0\right\}
\eeq
and such that as time
$t\mapsto \infty$ we arrive at a function $u(x,y)$ such that the
associated $f(x,y)$ is a solution to the inverse problem
in a least square (LS) sense. This means that the evolution of $u$ and
consequently the corresponding contour $\partial D$ and object $f$ is such that
the LS criterion

\beq \label{LS_crit}
J(f)=\|p-R(f)\|^{2}
\eeq
decreases during the evolution. In this approach the function $u$
is assimilated to a surface (of heat front wave) and the evolution of
the contours $\partial D$ is constraint to be perpendicular to the
surface. This means that the variation of the interior region
$(\delta x, \delta y)$ is such that
\beq
(\delta x, \delta y)=\alpha(x,y,t) \frac{\nabla u}{|\nabla u|}
\eeq
This is the Level-Set approach originally developed by Osher and Sethian
\cite{Osher88} for
problems involving the motion of curves and surfaces and then adapted,
used and referred as {\em snakes} or {\em active contour models} by many
authors in computer vision \cite{Catte92,Malladi95} and recently for
inverse problems \cite{Santosa96}. 
See also \cite{Hanson97a,Hanson97b} for the application of the deformable 
templates in tomographic image reconstruction. 

This approach also needs pixel or voxel representation of the image and
its calculation cost is huge, specially in 3D imaging systems.
We are presently working on this approach trying to extend it for
minimizing a regularized criterion in place of the the LS criterion
(\ref{LS_crit}) and implementing it in 2D and 3D tomographic image
reconstruction.

\medskip\noindent$\bullet$
The third approach starts by giving a parametric model for the object
and then tries to estimate these parameters using least squares (LS)
or maximum likelihood (ML) methods.
In general, in this approach one chooses a parametric model such as
superposition of circular or elliptical homogeneous regions to be able to relate
analytically the projections to these parameters.
For example, for a superposition of elliptical homogeneous regions we have:

\beq \label{eq.9}
f(x,y)=\sum_{k=1}^K d_k f_k(x-\alpha_k,y-\beta_k)
\eeq

\beq \label{eq.10}
\hbox{with~}
f_k(x,y)=\left\{\barr{ll}
1 & \hbox{if}\quad
(x/a_{k})^2+(y/b_{k})^2<1, \\
0 & \hbox{elsewhere} \earr \right.
\eeq
where $\thetab=\{d_k, \alpha_k, \beta_k, a_k, b_k, k=1,\cdots,K\}$
is a vector of parameters defining the parametric model of the image
(density values, coordinates of the centers and the two diameters of
the ellipses). It is then easy to calculate analytically the
projections and the relation between the data and the unknown
parameters becomes:

\beq \label{eq.11}
p(r,\phi)=h(r,\phi;\thetab)+n(r,\phi)
\eeq
where $h(r,\phi;\thetab)$ has an analytic expression in $\thetab$.
The LS or the ML estimate when the noise is assumed to be zero mean,
white and Gaussian, is then given by:

\beq \label{eq.12}
\wh{\thetab}=\argmins{\thetab}{\|p(r,\phi)-h(r,\phi;\thetab)\|^2}
\eeq
This approach has also been used with success
in image reconstruction \cite{Shepp82,Rossi84,Prince90,Prince91}.
But, the range of applicability of these methods is
limited to the cases where the parametric models are actually appropriate.

\medskip\noindent$\bullet$
The fourth approach which is more appropriate to our problem of shape
reconstruction, consists in modeling directly the contour of the object
by a function, say $g(\theta)$ in a cylindrical coordinates
$(\rho, \theta)$ such as:

\beq \label{eq.13}
D=\left\{(x,y): x^2+y^2 = g^2(\theta) \right\}.
\eeq
The next step is then to relate the projections $p(r,\phi)$ to $g(\theta)$
which, in this case is:

\beq \label{eq.14}
p(r,\phi) = \izdpi\int_0^{g(\theta)} \delta(r-\rho\cos(\phi-\theta))
\rho\d{\rho}\d{\theta}.
\eeq
and finally to discretize this relation to obtain:

\beq \label{eq.15}
\pb=\hb(\gb)+\nb
\eeq
where $\gb$ represents the discretized values of $g(\theta)$ defining the
contour of the object and $\hb(\gb)$ represents the discretized version of
the nonlinear operator (\ref{eq.14}) relating projection data $\pb$ and $\gb$.
Then, one defines the solution as the argument which minimizes

\beq \label{eq.16}
J(\gb)=||\pb - \hb(\gb)||^2 + \lambda \Omega(\gb),
\eeq
where $\Omega(\gb)$ has to be chosen appropriately to reflect some regularity
property of the object's contour.

In this case also one can consider $J(\gb)$ as the MAP criterion with
$Q(\gb)=||\pb - \hb(\gb)||^2$ as the likelihood term and $\Omega(\gb)$ as
the prior one. \\
This approach has been used in image restoration \cite{Friedland92},
but it seems to be new in image reconstruction applications and
the proposed method in this work is in this category.
The originality of our work is to model the contour of the object by a
piecewise linear function which means that the object is modeled as a
polygonal region whose vertices are estimated directly from the
projection data.

\section{Proposed method}
In this paper we propose to model the contour of the object (fault region)
as a periodic piecewise linear function or equivalently to model the shape
of the object as a polygonal region with a great number $N$ of vertices to be
able to approximate any shape.
Then we propose to estimate directly the coordinates
$\{(x_j,y_j), j=1,\cdots,N\}$
of the vertices of this polygonal region from the projection data.

\bfig[htb]
\setlength{\unitlength}{.4mm}
\hfil\input{tomox5.pic}\hfil
\caption{Proposed shape reconstruction modeling.}
\efig


The idea of modeling the shape of the object as a polygonal region is not
new and some works have been done in image reconstruction applications.
See for example Milanfar, Karl \& Willsky
\cite{Milanfar94,Milanfar94b,Milanfar95},
but in general, in these works, either a hypothesis of convexity of the
polygonal region has been used which is very restrictive in real applications
or the number of vertices of the polygon is restricted to a very small
number.
In our work we do not make neither of these hypothesis; the polygon
can be convex or not and we choose $N$ sufficiently
great to to be able to approximate appropriately any shape.

The solution is then defined as the minimizer the following criterion

\beq \label{eq.17}
J(\zb)=||\pb - \hb(\zb)||^2 + \lambda \Omega(\zb),
\eeq
where $\zb=\xb+i\yb$ is a complex vector whose real and imaginary
parts represent the $x$ and the $y$ coordinates of the polygon vertices,
$\hb(\zb)$ represents the direct operator which calculates the projections
for any given $\zb$ and $\Omega(\zb)$ is chosen to be a function which
reflects the regularity of the object contour.
One can for example choose the following function:

\beq \label{eq.18}
\Omega(\zb)=\left(\sum_{j=1}^N |z_{j+1}-z_j|-2\pi R_0\right)^2
\eeq
which favors a shape whose contour length is near a prior known value
$2\pi R_0$.
In this work we used the following function:

\beq \label{eq.19}
\Omega(\zb)=\sum_{j=1}^N |z_{j-1}-2 z_{j}+z_{j+1}|^2
=4 \sum_{j=1}^N \left|z_{j}-(z_{j-1}-z_{j+1})/2\right|^2,
\eeq
which favors a shape whose local curvature is limited.
Note that
$\left|z_{j}-(z_{j-1}-z_{j+1})/2\right|$
is just the Euclidian distance between
the point $z_{j}$ and the midpoint of the line segment passing through
$z_{j-1}$ and $z_{j+1}$ and so this choice favors a shape whose local
curvature is limited.
We can also give a probabilistic interpretation to this choice. In
fact we can consider $z_{j}$ as random variables with the following
Markovian law:

\beq  \label{eq.20}
p(z_{j}|\zb)=p(z_{j}|z_{j-1},z_{j+1})
\propto \expf{-\frac{1}{2\sigma^{2}}
\left|z_{j}-(z_{j-1}-z_{j+1})/2\right|^2}.
\eeq
Other functions are possible and are studied in this work.

In both cases, the criterion $J(\zb)$ is multi-modal essentially due to the
fact that $\hb(\zb)$ is a nonlinear function of $\zb$.
Calculating the optimal solution corresponding to the global minimum
of (\ref{eq.17}) needs then carefully designed algorithms.
For this we propose the two following strategies:

\medskip\noindent$\bullet$
The first is to use a global optimization technique such as simulated
annealing. This technique has given satisfactory results as it can be seen
from the simulations in the next section.
However, this algorithm needs a great number of iterations and some skills
for choosing the first temperature and cooling schedule, but the overall
calculations is not very important due to the fact that, in this algorithm,
at each iteration only one of the vertices $z_j$ is changed. So, at each step
we need to calculate the variation of the criterion due to this change which
can be done with a reasonable cost.

\medskip\noindent$\bullet$
The second is to find an initial solution in the attractive region of
the global optimum and to use a local descent type algorithm such as
the ICM (Iterated conditional modes) of Besag
\cite{Besag86,Besag89,Besag93} to find the solution.

The main problem here is how to find this initial solution.
For this, we used a moment based method proposed by
Milanfar, Karl \& Willsky \cite{Milanfar94,Milanfar94b,Milanfar95}
which is accurate enough to obtain an initial solution which is not very
far from the optimum.
The basic idea of this method is to relate the moments of the projections
to the moments of a class of polygonal regions obtained by an affine
transformation of a regular polygonal region, and so to estimate a polygonal
region whose corners are on an ellipse and whose moments
up to the second order matches those of the projections.

However, there is no theoretical proof that this initial solution will be
in the attractive region of the global optimum.
The next section will show some results comparing the performances of these
two methods as well as a comparison with some other classical methods.

\section{Simulation results}
To measure the performances of the proposed method and keeping the objective
of using this method for NDT applications where the number of projections
are very limited, we simulated two cases where the objects have polygonal
shapes with $N=40$ corners (hand-made) and calculated their projections
for only 5 directions

\[
\phi=\{-45,-22.5,0,22.5,45\} \quad \hbox{~degrees}.
\]
Then, we added some noise (white, Gaussian and centered) on them to simulate
the measurement errors. The signal to noise ratio (SNR) was chosen 20dB.
We define SNR as follows:
\[
SNR=10 \log \frac{\sum_i (p_i -\bar{p})^2}{\sum_i (n_i -\bar{n})^2}
\]
where $p_i$ and $n_i$ are respectively data and noise samples and $\bar{p}$
and $\bar{n}$ their respective mean values.

Finally, from these data we estimated the solutions by the
proposed method using either the simulated annealing (SA) or the
iterated conditional modes (ICM) algorithms.

\smallskip
Figure~2 shows these two objects and
their relative simulated projections data.

\smallskip
In Figures 3 and 4, we give the reconstruction results obtained by the
proposed method using either the SA algorithm (Figure~3)
or the ICM algorithm (Figure~4).
In these figures we show the original objects,
the initial solutions, the intermediate solutions during the iterations
and the final reconstructed objects obtained after 200 iterations.

Note that, the SA algorithm is theoretically independent of initialization
while the ICM is not.
However, in these figures, for the purpose of the comparison,
we show the results obtained by the two algorithms with the
same initialization.

\smallskip
To show that the method is not very sensible on the prior knowledge of the
vertices number, we give in Figure~5, the reconstruction results of the
object 2 in 4 cases: $N=10$, $N=20$ and $N=30$ and $N=40$.
As we can remark all the reconstructed results seem satisfactory.

\smallskip
In Figure~6 we show a comparison between the results obtained by the
proposed method and those obtained either by a classical backprojection
method or by some other methods in the first approach using (\ref{eq.3})
and (\ref{eq.4}) with different regularization functionals
$\Omega(\fb)$, more specifically:
\bit
\item {\em Gaussian Markov models} (\ref{eq.5}) with
the potential function $\phi(x,y)=|x-y|^2$
which can also be considered as a quadratic regularization method; and
\item {\em Compound Markov models} with non convex potential functions
$\phi(x,y)=\min\left\{|x-y|^2,1\right\}$ which is a truncated quadratic
potential function.
\eit
In the first case the criterion to optimize is convex and we used a conjugate
gradient (CG) algorithm to find the optimized solution.
In the second case the criterion is not convex and we used a
Graduated non convexity (GNC) based
optimization algorithm developed in \cite{Nikolova94,Nikolova95d,Nikolova96a}
to find the solution.

Note that, these results are given here to show the relative performances
of these methods in a very difficult situation where we have only five
projections. In fact, in more comfortable situations (more projections
uniformly distributed around the object and high SNR) all these methods, even
the simplest one such as the classical backprojection will give similar
and satisfactory results.
Here, we compare the results obtained from the same restricted set of
data.

%\clearpage\newpage
\bfig[htb]
\bcc
\includegraphics[width=9cm,height=9cm]{proj1} \\
\includegraphics[width=9cm,height=9cm]{proj2}
\ecc
\caption{Original images and simulated projections data.}
\efig

%\clearpage\newpage
\bfig[htb]
\bcc
\btab{@{}c@{}c@{}}
\includegraphics[width=6.5cm,height=6.5cm]{fig21}
\includegraphics[width=6.5cm,height=6.5cm]{fig22}
\etab
\ecc
\caption[Reconstruction results using simulated annealing.]{
Reconstruction results using simulated annealing.\\
o) Original objects, \quad
+) Initializations, \quad
.) Evolution of the solutions during the iterations and \\
 $\star$) Final reconstructed objects.
}
\efig

%\clearpage\newpage
\bfig[hbt]
\bcc
\btab{@{}c@{}c@{}}
  \includegraphics[width=6.5cm,height=6.5cm]{fig71} &
  \includegraphics[width=6.5cm,height=6.5cm]{fig72}
\etab
\ecc
\caption[Reconstruction using a local descent (ICM) minimizer.]{
Reconstruction results using a moment-based initialization and a local
descent (ICM) minimizer.\\
o) Original objects, \quad
+) Initializations, \quad
.) Evolution of the solutions during the iterations and \\
 $\star$) Final reconstructed objects.
}
\efig

%\clearpage\newpage
\bfig[htb]
\bcc
\btab{@{}c@{}c@{}}
  a & b \\
  \includegraphics[width=6.5cm,height=6.5cm]{fig2a} &
  \includegraphics[width=6.5cm,height=6.5cm]{fig2b} \\
  c & d \\
  \includegraphics[width=6.5cm,height=6.5cm]{fig2c} &
  \includegraphics[width=6.5cm,height=6.5cm]{fig2d} \\
\etab
\ecc
\caption[Reconstruction results with different number of vertices.]{
Reconstruction results with different number of vertices $N$ of the same
object. The number of vertices of the original object is $N=40$. \\
\centerline{
a) $N=10$, \qquad b) $N=20$ \qquad c) $N=30$ \qquad and \qquad d) $N=40$.}
}
\efig

%\clearpage\newpage
\bfig[htb]
\bcc
\btab{@{}c@{}c@{}c@{}}
a & b & c \\ ~ \\
\includegraphics[width=4.2cm,height=4.2cm]{original1}&
\includegraphics[width=4.2cm,height=4.2cm]{recuit11} &
\includegraphics[width=4.2cm,height=4.2cm]{icm11}
\\ ~\\ ~\\
d & e & f \\ ~ \\
\includegraphics[width=4.2cm,height=4.2cm]{back1} &
\includegraphics[width=4.2cm,height=4.2cm]{regq11}&
\includegraphics[width=4.2cm,height=4.2cm]{gnc11}
\etab
\ecc
\caption[A comparison with some classical methods.]{
A comparison with backprojection and some other classical methods:\\
a) Original objects,  \\
b) Results obtained by the proposed method using the SA
   optimization algorithm,  \\
c) Results obtained by the proposed method using the ICM algorithm,\\
d) Backprojection, \\
e) Gaussian Markov modeling MAP estimation and  \\
f) Compound Markov modeling and GNC optimization algorithm.
}
\efig

\bfig[htb]
\bcc
\btab{@{}c@{}c@{}c@{}}
a & b & c \\
\includegraphics[width=4.2cm,height=4.2cm]{original2}&
\includegraphics[width=4.2cm,height=4.2cm]{recuit12} &
\includegraphics[width=4.2cm,height=4.2cm]{icm12}
\\
d & e & f \\
\includegraphics[width=4.2cm,height=4.2cm]{back2} &
\includegraphics[width=4.2cm,height=4.2cm]{regq12}&
\includegraphics[width=4.2cm,height=4.2cm]{gnc12}
\etab
\ecc
\caption[A comparison with some classical methods.]{
A comparison with backprojection and some other classical methods: \\
a) Original objects,  \\
b) Results obtained by the proposed method using the SA
   optimization algorithm,  \\
c) Results obtained by the proposed method using the ICM algorithm,\\
d) Backprojection, \\
e) Gaussian Markov modeling MAP estimation and  \\
f) Compound Markov modeling and GNC optimization algorithm.
}
\efig

\clearpage\newpage
\section{Conclusions}
A new method for tomographic image reconstruction of a compact binary
object from a small number of its projections is proposed.
The basic idea of the proposed method is to model the compact binary
object as a polygonal region whose vertices coordinates are estimated
directly from the projections using the Bayesian MAP estimation
framework or equivalently by optimizing a regularized criterion.

Unfortunately, this criterion is not unimodal.
To find the optimized solution two algorithms have been proposed: \\�
-- a global optimization algorithm based on simulated annealing (SA) and \\
-- a local descent-based method based on the Iterated Conditional Modes
(ICM) algorithm proposed originally by Besag,
with a good initialization obtained by using a moment based method.

The first algorithm seems to give entire satisfaction.
The second can also give satisfaction, but it may also be plugged in a
local minimum.
In both algorithms the main cost calculation is due to the calculus of
the variation of the criterion when one of the vertices coordinates
is changed. We have written an efficient program to do this
\cite{Djafari96c,Djafari96e}.

An extension of this work to 3D image reconstruction with small number of
conic projections is in preparation \cite{Djafari97a,Djafari97f}.
The final objective of the proposed method is for non destructive
testing (NDT) image reconstruction applications where we can use not only
X-rays but also ultrasound or Eddy currents or a combination of them
\cite{Gautier95b,Premel95,Nikolova96a} to localize and to characterize more
accurately any anomalies (air bulbs) in metallic structures.

\input{revuedef.tex}
\input{month-bi.tex}
\def\UP#1{\uppercase{#1}}
\def\CGVIP{CGVIP}
\bibliographystyle{ieeetr}

%\bibliography{\bibdir hanson,\bibdir gpibase,\bibdir gpipubli,\bibdir 
%amd,\bibdir amd96,\bibdir amd97}

\begin{thebibliography}{10}

\bibitem{Herman87}
G.~Herman, H.~Tuy, H.~Langenberg, and P.~Sabatier, {\em Basic Methods of
  Tomography and Inverse Problems}.
\newblock Adams Hilger, 1987.

\bibitem{Kak87}
A.~Kak and M.~Slaney, {\em Principles of Computerized Tomographic Imaging}.
\newblock New York, NY: IEEE Press, 1987.

\bibitem{Geman87}
S.~Geman and D.~McClure, ``Statistical methods for tomographic image
  reconstruction,'' in {\em Proc. of the 46-th Session of the ISI, Bulletin of
  the ISI}, vol.~52, pp.~22--26, 1987.

\bibitem{Demoment89}
G.~Demoment, ``Image reconstruction and restoration : Overview of common
  estimation structure and problems,'' {\em \ieeeASSP}, vol.~ASSP-37,
  pp.~2024--2036, \Dec 1989.

\bibitem{Brette94a}
S.~Brette, J.~Idier, and A.~Mohammad-Djafari, ``Scale invariant \up{M}arkov
  models for linear inverse problems,'' in {\em Proc. of the Section on
  \UP{B}ayesian Statistical Sciences}, (Alicante, Spain), pp.~266--270,
  American Statistical Association, 1994.

\bibitem{Djafari96a}
A.~Mohammad-Djafari and J.~Idier, {\em A scale invariant \UP{B}ayesian method
  to solve linear inverse problems}, pp.~121--134.
\newblock MaxEnt93, Santa Barbara, U.S.A.: \KAP, \up{G}. \up{H}eidbreder~ed.,
  1996.

\bibitem{Bouman93}
C.~Bouman and K.~Sauer, ``A generalized \up{G}aussian image model for
  edge-preserving {\sc map} estimation,'' {\em \ieeeIP}, vol.~IP-2,
  pp.~296--310, \Jul 1993.

\bibitem{Bedini94}
L.~Bedini, I.~Gerace, and A.~Tonazzini, ``A deterministic algorithm for
  reconstructing images with interacting discontinuities,'' {\em \CVGIP},
  vol.~56, pp.~109--123, \Mar 1994.
\newblock AMD.

\bibitem{Nikolova95d}
M.~Nikolova, J.~Idier, and A.~Mohammad-Djafari, ``Inversion of large-support
  ill-posed linear operators using a piecewise \up{G}aussian {\sc mrf},'' tech.
  rep., {\sc gpi--lss}, to appear in \ieeeIP, Gif-sur-Yvette, France, 1997.

\bibitem{Nikolova94}
M.~Nikolova, A.~Mohammad-Djafari, and J.~Idier, ``Inversion of large-support
  ill-conditioned linear operators using a \up{M}arkov model with a line
  process,'' in {\em \ICASSP}, vol.~V, (Adelaide, Australia), pp.~357--360,
  1994.

\bibitem{Nikolova96a}
M.~Nikolova and A.~Mohammad-Djafari, ``Eddy current tomography using a binary
  \up{M}arkov model,'' {\em Signal Processing}, vol.~49, pp.~119--132, May
  1996.

\bibitem{Osher88}
S.~Osher and J.~Sethian, ``Fronts propagating with curvature-dependent speed:
  Algorithms based on hamilton-jacobi formulations,'' {\em Journal of
  Computational Physics}, vol.~79, pp.~12--49, 1988.

\bibitem{Catte92}
F.~Catt\'e, P.~Lions, J.~Morel, and T.~Coll, ``Image selective smoothing and
  edge detection by nonlinear diffusion,'' {\em SIAM J. Num. Anal.}, vol.~29,
  pp.~182--193, 1992.

\bibitem{Malladi95}
R.~Malladi, J.~Sethian, and B.~Vemuri, ``Shape modelling with front
  propagation: A level-set approch,'' {\em \ieeePAMI}, vol.~17, 1995.

\bibitem{Santosa96}
F.~Santosa, ``A level-set approch for inverse problems involving obstacles,''
  {\em Control, Optimisation and Calculus of Variations}, vol.~1, pp.~17--33,
  1996.

\bibitem{Hanson97a}
K.~M. Hanson, G.~S. Cunningham, and R.~J. McKee, ``Uncertainties in tomographic
  reconstructions based on deformable models,'' {\em {\rm in} Image
  Processing{\rm , K. M. Hanson, ed.,} Proc. SPIE}, vol.~3034, pp.~276--286,
  1997.

\bibitem{Hanson97b}
K.~M. Hanson, G.~S. Cunningham, and R.~J. McKee, ``Uncertainty assessment for
  reconstructions based on deformable models,'' {\em Int. J. Imaging Systems
  and Technology}, vol.~8, pp.~506--512, 1997.

\bibitem{Shepp82}
L.~A. Shepp and Y.~Vardi, ``Maximum likelihood reconstruction for emission
  tomography,'' {\em IEEE Trans. on Medical Imaging}, vol.~MI-1, pp.~113--122,
  1982.

\bibitem{Rossi84}
D.~Rossi and A.~Wilsky, ``Reconstruction from projections based on detection
  and estimation of objects,'' {\em ieeeASSP}, vol.~ASSP-32, no.~4,
  pp.~886--906, 1984.

\bibitem{Prince90}
J.~Prince and A.~Wilsky, ``Reconstructing convex sets from support line
  measurements,'' {\em ieeePAMI}, vol.~12, no.~3, pp.~377--389, 1990.

\bibitem{Prince91}
J.~Prince and A.~Wilsky, ``Convex set reconstruction using prior shape
  information,'' {\em CVGIP}, vol.~53, no.~5, pp.~413--427, 1991.

\bibitem{Friedland92}
N.~Friedland and A.~Rosenfeld, ``Compact object recognition using
  energy-function-based optimization,'' {\em \ieeePAMI}, vol.~14, no.~7,
  pp.~770--777, 1992.

\bibitem{Milanfar94}
P.~Milanfar, W.~Karl, and A.~Wilsky, ``A moment-based variational approach to
  tomographic reconstruction,'' {\em \ieeeIP}, vol.~25, no.~9, pp.~772--781,
  1994.

\bibitem{Milanfar94b}
P.~Milanfar, W.~Karl, and A.~Wilsky, ``Reconstructing binary polygonal objects
  from projections: A statistical view,'' {\em \CGVIP}, vol.~56, no.~5,
  pp.~371--391, 1994.

\bibitem{Milanfar95}
P.~Milanfar, G.~Vergese, W.~Karl, and A.~Wilsky, ``Reconstructing polygons from
  moments with connections to array processing,'' {\em \ieeeSP}, vol.~43,
  no.~2, pp.~432--443, 1995.

\bibitem{Besag86}
J.~E. Besag, ``On the statistical analysis of dirty pictures (with
  discussion),'' {\em \JRSSB}, vol.~48, no.~3, pp.~259--302, 1986.

\bibitem{Besag89}
J.~E. Besag, ``Digital image processing : Towards \up{B}ayesian image
  analysis,'' {\em \JAS}, vol.~16, no.~3, pp.~395--407, 1989.

\bibitem{Besag93}
J.~E. Besag and P.~Green, ``Spatial statistics and \up{B}ayesian computation,''
  {\em J. R. Statist. Soc. B}, vol.~55, pp.~25--37, 1993.

\bibitem{Djafari96c}
A.~Mohammad-Djafari, ``Image reconstruction of a compact object from a few
  number of projections,'' in {\em \UP{IASTED}, Int. Conf. on Signal and Image
  Processing (SIP'96)}, (Orlando, Florida, USA), pp.~325--329, 1996.

\bibitem{Djafari96e}
A.~Mohammad-Djafari, ``A \up{B}ayesian approach to shape reconstruction of a
  compact object from a few number of projections,'' in {\em \UP{M}aximum
  \UP{E}ntropy and \UP{B}ayesian \UP{M}ethods}, (MaxEnt96, Aug. 96,
  Berg-en-Dal, South Africa), \KAP, 1996.

\bibitem{Djafari97a}
A.~Mohammad-Djafari and Y.~Khayi, ``Array processing techniques and shape
  reconstruction in tomography,'' in {\em \ICASSP}, vol.~II, (Munich, Germany),
  pp.~3785--3789, 1997.

\bibitem{Djafari97f}
A.~Mohammad-Djafari, K.~Sauer, Y.~Khayi, and E. Cano, ``Reconstruction of the
  shape of a compact object from a few number of projections,'' in {\em IEEE
  Int. Conf. on Image Processing (ICIP)}, vol.~1, (IEEE ICIP 97, Santa Barbara,
  CA, USA), pp.~165--169, 1997.

\bibitem{Gautier95b}
S.~Gautier, G.~Le~Besnerais, A.~Mohammad-Djafari, and B.~Lavayssi\`ere, {\em
  Data fusion in the field of non destructive testing}.
\newblock Santa Fe, U.S.A.: \KAP, \up{K}.~\up{H}anson~ed., 1995.

\bibitem{Premel95}
D.~Pr\'emel and A.~Mohammad-Djafari, ``Eddy current tomography in cylindrical
  geometry,'' {\em \ieeeM}, vol.~M-31, pp.~2000--2003, \May 1995.

\end{thebibliography}

\end{document}


revuedef.tex
% TEX 7(ascii) bits
% fichier : revuedef.tex	Auteur: AMD, GLB	09/92
% Modifie et completee par	L.Lemitre le 05 10 94
%				et al...
%
% Ce fichier contient des alias sur les noms de revues
% pour faciliter la saisie de fiches biblio.
%
% En commentaire, le facteur d'impact (nombre moyen de citation
% des articles, normalise a la moyenne mondiale sur l'ensemble
% des revues indexees par ISI) pour l'annee 1994 (J. I.). 

\def\AsAs{Astrononmy and Astrophysics}					% 2.328
\def\AAP{Advances in Applied Probability}				% 0.549
\def\ABE{Annals of Biomedical Engineering}				% 0.709
\def\AISM{Annals of Institute of Statistical Mathematics}		% 0.175
\def\AMS{Annals of Mathematical Statistics}			% ANNALS OF MATHEMATICS STUDIES ? 0.8            
\def\AO{Applied Optics}							% 1.033
\def\AP{The Annals of Probability}					% 0.68
\def\ARAA{Annual Review of Astronomy and Astrophysics}			% 9.865
\def\AST{The Annals of Statistics}					% 0.78
\def\AT{Annales des T\'el\'ecommunications}				% 0.155
\def\BMC{Biometrics}							% 1.207
\def\BMK{Biometrika}							% 0.832
\def\CPAM{Communications on Pure and Applied Mathematics}		% 1.282
\def\EMK{Econometrica}							% 2.54
\def\CRAS{Compte-rendus de l'acad\'emie des sciences}			% ?
\def\CVGIP{Computer Vision and Graphics and Image Processing}		% ?
\def\GJRAS{Geophysical Journal of the Royal Astrononomical Society}	% 1.371
\def\GSC{Geoscience}						% GEOSCIENCE CANADA ? 0.273
\def\GPH{Geophysics}							% 0.824
\def\GRETSI#1{Actes du #1$^{\mbox{e}}$ Colloque GRETSI} 		% ?
\def\CGIP{Computer Graphics and Image Processing}			% ?
\def\ICASSP{Proceedings of IEEE ICASSP}					% ?
\def\ICEMBS{Proceedings of IEEE EMBS}					% ?
\def\ICIP{Proceedings of the International Conference on Image Processing}% ?
\def\ieeP{Proceedings of the IEE}					% ?
\def\ieeeAC{IEEE Transactions on Automatic and Control}			% 0.867
\def\ieeeAES{IEEE Transactions on Aerospace and Electronic Systems}	% 0.459
\def\ieeeAP{IEEE Transactions on Antennas and Propagation}		% 0.806
\def\ieeeAPm{IEEE Antennas and Propagation Magazine}			% ?
\def\ieeeASSP{IEEE Transactions on Acoustics Speech and Signal Processing}% ?
\def\ieeeBME{IEEE Transactions on Biomedical Engineering}		% 1.061
\def\ieeeCS{IEEE Transactions on Circuits and Systems}			% 0.54
\def\ieeeCT{IEEE Transactions on Circuit Theory}			% 0.732
\def\ieeeC{IEEE Transactions on Communications}				% 0.969
\def\ieeeGE{IEEE Transactions on Geoscience and Remote Sensing}		% 1.356
\def\ieeeGEE{IEEE Transactions on Geosciences Electronics}		% ?
\def\ieeeIP{IEEE Transactions on Image Processing}			% ?
\def\ieeeIT{IEEE Transactions on Information Theory}			% 1.971
\def\ieeeMI{IEEE Transactions on Medical Imaging}			% 1.372
\def\ieeeMTT{IEEE Transactions on Microwave Theory and Technology}	% 1.004
\def\ieeeM{IEEE Transactions on Magnetics}				% 0.758
\def\ieeeNS{IEEE Transactions on Nuclear Sciences}			% 1.183
\def\ieeePAMI{IEEE Transactions on Pattern Analysis and Machine Intelligence}% 2.006
\def\ieeeP{Proceedings of the IEEE}					% ?
\def\ieeeRS{IEEE Transactions on Radio Science}				% ?
\def\ieeeSMC{IEEE Transactions on Systems, Man and Cybernetics}		% 0.649
\def\ieeeSP{IEEE Transactions on Signal Processing}			% 1.234
\def\ieeeSSC{IEEE Transactions on Systems Science and Cybernetics}	% ?
\def\ieeeSU{IEEE Transactions on Sonics and Ultrasonics}		% ?
\def\ieeeUFFC{IEEE Transactions on Ultrasonics Ferroelectrics and Frequency Control}% 0.927
\def\IJC{International Journal of Control}				% 0.617
\def\IJCV{International Journal of Computer Vision}			% 1.153
\def\IJIST{International Journal of Imaging Systems and Technology}	% ?
\def\IP{Inverse Problems}						% 0.98
\def\ISR{International Statistical Review}				% 0.75
\def\IUSS{Proceedings of International Ultrasonics Symposium}		% ?
\def\JAPH{Journal of Applied Physics}					% 1.658
\def\JAP{Journal of Applied Probability}				% 0.432
\def\JAS{Journal of Applied Statistics}					% ?
\def\JASA{Journal of Acoustical Society America}			% 1.273
\def\JASAS{Journal of American Statistical Association}			% 1.244
\def\JBME{Journal of Biomedical Engineering}				% 0.7
\def\JCAM{Journal of Computational and Applied Mathematics}		% 0.349
\def\JEWA{Journal of Electromagnetic Waves and Applications}		% 0.514
\def\JMO{Journal of Modern Optics}					% 1.005
\def\JNDE{Journal of Nondestructive Evaluation}				% ?
\def\JMP{Journal of Mathematical Physics}				% 0.969
\def\JOSA{Journal of the Optical Society of America}			% 1.425 (A-OPTICS IMAGE SCIENCE AND VISION)
\def\JP{Journal de Physique}						% 1.773 JOURNAL DE PHYSIQUE I
									% 1.923 JOURNAL DE PHYSIQUE II
									% 0.475 JOURNAL DE PHYSIQUE III
									% 0.115 JOURNAL DE PHYSIQUE IV
\def\JRSSA{Journal of the Royal Statistical Society A}			% 0.976
\def\JRSSB{Journal of the Royal Statistical Society B}			% 2.538
\def\JRSSC{Journal of the Royal Statistical Society C}			% ?
\def\JSPI{Journal of Statistical Planning and Inference}  		% 0.145
\def\JTSA{Journal of Time Series Analysis}                   		% ?
\def\JVCIR{Journal of Visual Communication and Image Representation} 	% ?
	\def\MMAS{???} % Trouve dans gpi base 				% ?
\def\KAP{Kluwer \uppercase{A}cademic \uppercase{P}ublishers}							% 
\def\MNAS{Mathematical Methods in Applied Science}			% 0.247
\def\MNRAS{Monthly Notices of the Royal Astronomical Society}		% 3.089
\def\MP{Mathematical Programming}					% 0.763
	\def\NSIP{NSIP}  % Trouve dans gpi base				% ?
\def\OC{Optics Communication}						% 1.205
\def\PRA{Physical Review A}						% 2.292
\def\PRB{Physical Review B}						% 3.187
\def\PRC{Physical Review C}						% 1.842
\def\PRD{Physical Review D}						% 3.233
\def\PRL{Physical Review Letters}					% 6.626
\def\RGSP{Review of Geophysics and Space Physics}			% 4.314
\def\RPA{Revue de Physique Appliqu\'e}							% 0.753
\def\RS{Radio Science}							% 0.753
\def\SP{Signal Processing}						% 0.44
\def\siamAM{SIAM Journal of Applied Mathematics}			% 0.743
\def\siamCO{SIAM Journal of Control}					% 0.968
\def\siamJO{SIAM Journal of Optimization}				% ?
\def\siamMA{SIAM Journal of Mathematical Analysis}			% 0.765
\def\siamNA{SIAM Journal of Numerical Analysis}				% 1.021
\def\siamR{SIAM Review}							% 1.216
\def\SSR{Stochastics and Stochastics Reports}       			% ?
\def\TPA{Theory of Probability and its Applications}			% 0.045
\def\TMK{Technometrics}							% 1.585
\def\TS{Traitement du Signal}						% ?
\def\UCMMP{U.S.S.R. Computational Mathematics and Mathematical Physics}	% ?
\def\UMB{Ultrasound in Medecine and Biology}				% 1.252
\def\US{Ultrasonics}							% 1.046
\def\USI{Ultrasonic Imaging}						% 1.409

spgraph.tex
% SpGpraph Version 0.0  < April 6, 1988 >
%        Written by Conrad Kwok
%
% The macros are in public domain.
% You may distribute or modify it in any ways you like.
%
% You need an dvi-to-printer program that supports the \special
% commands used by the macros. The commands are subset of the commands
% used by tpic program. For more information, read the READ.ME file
%
% Usage:
%    For LaTeX
%  Include the file after the \documentstyle command by saying
%  \include spgraph
%
%    For TeX
% Include the file after including the file latexpicobjs.tex
%  \include latexpicobj
%  \include spgraph
% in the beginning of the TeX file.
%
% The macros are for use with LaTeX picture environment
% By including the macro file, you can draw
%   1) \line in any slope
%   2) \circle with any radius
%   3) \ellipse with any x-axis and y-axis length
%   4) Circular \arc by giving starting and ending angle (in radius)
%
% Furthermore you may draw lines in any thickness.
%
%
\makeatletter

% Create a new counter to store the line thickness in thousandth of an inch
\newcount\@gphlinewidth
\@gphlinewidth\@wholewidth \divide\@gphlinewidth 4736

% Redefine \thinlines, \thicklines and define \Thickline
% See also latex.tex
\def\thinlines{\let\@linefnt\tenln \let\@circlefnt\tencirc
    \@wholewidth\fontdimen8\tenln \@halfwidth .5\@wholewidth
    \@gphlinewidth\@wholewidth \divide\@gphlinewidth 4736\relax}
\def\thicklines{\let\@linefnt\tenlnw \let\@circlefnt\tencircw
    \@wholewidth\fontdimen8\tenlnw \@halfwidth .5\@wholewidth
    \@gphlinewidth\@wholewidth \divide\@gphlinewidth 4736
    \advance\@gphlinewidth\@ne   % Make the output looks better
    \relax}
\def\Thicklines{\let\@linefnt\tenlnw \let\@circlefnt\tencircw
    \@wholewidth\fontdimen8\tenlnw \@wholewidth 1.5\@wholewidth
    \@halfwidth .5\@wholewidth
    \@gphlinewidth\@wholewidth \divide\@gphlinewidth 4736\relax}
%
% Generate the \special command for drawing arc
\def\@circlespecial#1#2#3#4{%
       \special{pn \the\@gphlinewidth}%
       \special{ar 0 0 #1 #2 #3 #4}
}
%
% Command for drawing an arc. Elliptical arc command can be generated
% but all iptex program I saw so far does not support that.
\def\@arc#1#2#3#4{%
% convert TeX dimension to length in terms thousandth of an inch
 \@tempdima #1\unitlength
 \@tempdimb #2\unitlength
        \@tempcnta\@tempdima \advance\@tempcnta 4736 \divide\@tempcnta 9473
 \@tempcntb\@tempdimb \advance\@tempcntb 4736 \divide\@tempcntb 9473
 \setbox\@tempboxa\hbox{%
     \@circlespecial{\the\@tempcnta}{\the\@tempcntb}{#3}{#4}}%
        \wd\@tempboxa\z@ \box\@tempboxa
}
%
% Command for drawing Circle and Ellipse in terms of \@arc
% replace original \@circle
\def\@circle#1{\@arc{#1}{#1}{0}{6.2832}}
\def\ellipse#1#2{\@arc{#1}{#2}{0}{6.2832}}
\def\arc#1#2#3{\@arc{#1}{#1}{#2}{#3}}
%
% Generate \special commands for drawing line
\def\@linespecial#1#2{%
       \special{pn \the\@gphlinewidth}%
       \special{pa 0 0}%
       \special{pa #1 #2}%
       \special{fp}
}
%
% Replace original \@sline
% For description, see latex.tex
\def\@sline{%
 \ifnum\@xarg< 0
   \@negargtrue \@xarg -\@xarg \@tempdima -\@linelen
 \else
   \@negargfalse \@tempdima\@linelen
 \fi
% truncation is used in arithmetic
 \@tempcnta\@linelen \divide\@tempcnta 4736
        \@yyarg -\@yarg \multiply\@yyarg \@tempcnta \divide\@yyarg\@xarg
  \if@negarg
     \@tempcnta -\@tempcnta
 \fi
 \setbox\@linechar\hbox{\@linespecial{\the\@tempcnta}{\the\@yyarg}}%
 \wd\@linechar\@tempdima
 \@clnht\@linelen
        \multiply\@clnht\@yarg
        \divide\@clnht\@xarg
 \ifnum\@yarg< 0
   \@clnht -\@clnht
   \ht\@linechar\z@ \dp\@linechar\@clnht
 \else
   \ht\@linechar\@clnht \dp\@linechar\z@
 \fi
 \box\@linechar
% Codes below (till end of the command) is only necessary 
% when used with \@svector
 \if@negarg
   \@yyarg -\@yarg
 \else
   \@yyarg \@yarg
 \fi
 \setbox\@linechar\hbox{\@linefnt\@getlinechar(\@xarg,\@yyarg)}%
 \ifnum\@yarg> 0
   \let\@upordown\raise
   \advance\@clnht -\ht\@linechar
 \else
   \let\@upordown\lower
 \fi
 \if@negarg \kern\wd\@linechar \fi
}
%
\makeatother


trdef.tex
%TRDEF.LATEX
%LaTeX document
%Author: Ali Mohammad-djafari
%Date  : 22 Augoust 1991.

%BOLD MATH

\def\bm#1{\mbox{\boldmath $#1$}}
%\def\bm#1{\mbox{$\boldmath{#1}$}}
%\def\bm#1{\mathbf{#1}}

%Numbers
\def\zerob{{\bm 0}}
\def\oneb{{\bm 1}}

%Small letters
\def\ab{{\bm a}}
\def\bb{{\bm b}}
\def\cb{{\bm c}}
\def\db{{\bm d}}
\def\eb{{\bm e}}
\def\fb{{\bm f}}
\def\gb{{\bm g}}
\def\hb{{\bm h}}
\def\ib{{\bm i}}
\def\jb{{\bm j}}
\def\kb{{\bm k}}
\def\lb{{\bm l}}
\def\mb{{\bm m}}
\def\nb{{\bm n}}
\def\ob{{\bm o}}
\def\pb{{\bm p}}
\def\qb{{\bm q}}
\def\rb{{\bm r}}
\def\sb{{\bm s}}
\def\tb{{\bm t}}
\def\ub{{\bm u}}
\def\vb{{\bm v}}
\def\wb{{\bm w}}
\def\xb{{\bm x}}
\def\yb{{\bm y}}
\def\zb{{\bm z}}

%Capital letters
\def\Ab{{\bm A}}
\def\Bb{{\bm B}}
\def\Cb{{\bm C}}
\def\Db{{\bm D}}
\def\Eb{{\bm E}}
\def\Fb{{\bm F}}
\def\Gb{{\bm G}}
\def\Hb{{\bm H}}
\def\Ib{{\bm I}}
\def\Jb{{\bm J}}
\def\Kb{{\bm K}}
\def\Lb{{\bm L}}
\def\Mb{{\bm M}}
\def\Nb{{\bm N}}
\def\Ob{{\bm O}}
\def\Pb{{\bm P}}
\def\Qb{{\bm Q}}
\def\Rb{{\bm R}}
\def\Sb{{\bm S}}
\def\Tb{{\bm T}}
\def\Ub{{\bm U}}
\def\Vb{{\bm V}}
\def\Wb{{\bm W}}
\def\Xb{{\bm X}}
\def\Yb{{\bm Y}}
\def\Zb{{\bm Z}}

%Small Greek letters
\def\alphab{\bm{\alpha}}
\def\betab{\bm{\beta}}
\def\deltab{\bm{\delta}}
\def\epsilonb{\bm{\epsilon}}
\def\gammab{\bm{\gamma}}
\def\omegab{\bm{\omega}}
\def\thetab{\bm{\theta}}
\def\xib{\bm{\xi}}
\def\lambdab{\bm{\lambda}}
\def\taub{\bm{\tau}}
\def\phib{\bm{\phi}}
\def\mub{\bm{\mu}}
\def\psib{\bm{\psi}}
\def\chib{\bm{\chi}}
\def\sigmab{\bm{\sigma}}

%Capital Greek letters
\def\Deltab{\bm{\Delta}}
\def\Lambdab{\bm{\Lambda}}
\def\Phib{\bm{\Phi}}
\def\Psib{\bm{\Psi}}
\def\Sigmab{\bm{\Sigma}}

%Caligraphic capital letters
\def\Ac{{\cal A}}
\def\Bc{{\cal B}}
\def\Cc{{\cal C}}
\def\Dc{{\cal D}}
\def\Ec{{\cal E}}
\def\Fc{{\cal F}}
\def\Gc{{\cal G}}
\def\Hc{{\cal H}}
\def\Ic{{\cal I}}
\def\Jc{{\cal J}}
\def\Kc{{\cal K}}
\def\Lc{{\cal L}}
\def\Mc{{\cal M}}
\def\Nc{{\cal N}}
\def\Oc{{\cal O}}
\def\Pc{{\cal P}}
\def\Qc{{\cal Q}}
\def\Rc{{\cal R}}
\def\Sc{{\cal S}}
\def\Tc{{\cal T}}
\def\Uc{{\cal U}}
\def\Vc{{\cal V}}
\def\Wc{{\cal W}}
\def\Xc{{\cal X}}
\def\Yc{{\cal Y}}
\def\Zc{{\cal Z}}

%widetilde and widehat
\def\wt#1{\widetilde{#1}}
\def\wh#1{\widehat{#1}}
%

%widehat
\def\xh{\widehat{x}}
\def\thetah{\widehat{\theta}}
\def\betah{\widehat{\beta}}

\def\xbh{\widehat{\xb}}
\def\thetabh{\widehat{\thetab}}
\def\betabh{\widehat{\betab}}

\def\xbhk{\widehat{\xb}^{k}}
\def\thetahk{\widehat{\theta}^{k}}
\def\betahk{\widehat{\beta}^{k}}

\def\xbhkp{\widehat{\xb}^{k+1}}
\def\thetahkp{\widehat{\theta}^{k+1}}
\def\betahkp{\widehat{\beta}^{k+1}}

\def\thetabhk{\widehat{\thetab}^{k}}
\def\betabhk{\widehat{\betab}^{k}}

\def\thetabhkp{\widehat{\thetab}^{k+1}}
\def\betabhkp{\widehat{\betab}^{k+1}}

\def\thetamin{\theta_{\mbox{\tiny min}}}
\def\thetamax{\theta_{\mbox{\tiny max}}}
\def\betamin{\beta_{\mbox{\tiny min}}}
\def\betamax{\beta_{\mbox{\tiny max}}}

%Arrows
\def\ra{\rightarrow}
\def\la{\leftarrow}
\def\da{\downarrow}
\def\ua{\uparrow}

\def\Ra{\Rightarrow}
\def\La{\Leftarrow}
\def\Da{\Downarrow}
\def\Ua{\Uparrow}

\def\lra{\longrightarrow}
\def\lla{\longleftarrow}
\def\Lra{\Longrightarrow}
\def\Lla{\longleftarrow}

\def\lrarr{\leftrightarrow}
\def\Lrarr{\Leftrightarrow}
\def\udarr{\updownarrow}
\def\Uparr{\Updownarrow}

%
\def\d#1{\,\mbox{d}#1}
\def\dxdy{\d{x}\d{y}}
\def\dwxdwy{\d{\omega_x}\d{\omega_y}}
\def\dxdydz{\d{x}\d{y}\d{z}}

%DEFINITIONS
\def\disp#1{{\displaystyle #1}}
\def\diag#1{\mbox{diag}\left\{#1\right\}}

\def\Prob#1{\mbox{Pr}\left\{#1\right\}}
\def\var#1{\mbox{Var}\left\{#1\right\}}
\def\cov#1{\mbox{Cov}\left\{#1\right\}}
\def\corr#1{\mbox{Corr}\left\{#1\right\}}
\def\trace#1{\mbox{Tr}\left\{#1\right\}}
\def\rang#1{\mbox{rang}\left\{#1\right\}}
\def\det#1{\mbox{d\'et}\left\{#1\right\}}

\def\cosf{\cos \phi}
\def\sinf{\sin \phi}
\def\cost{\cos \theta}
\def\sint{\sin \theta}

\def\sgn{\mbox{sgn}}
\def\sinc{\mbox{sinc}}
\def\rect{\mbox{rect}}
\def\sincf#1{\mbox{sinc}\left(#1\right)}
\def\rectf#1{\mbox{rect}\left(#1\right)}
\def\trif#1{\mbox{tri}\left(#1\right)}
%
%VECTEURS
\def\xvec#1#2#3{\left\{#1_#2,\ldots,#1_#3\right\}}

\def\vx{\left[x_1,\ldots, x_n\right]^t}
\def\vz{\left[z_1,\ldots, z_n\right]^t}
\def\vw{\left[\omega_1,\ldots, \omega_n\right]^t}
\def\vxi{\left[\xi_1,\ldots, \xi_n\right]^t}
%
%INTEGRALES
\def\iii{\int_{-\infty}^{+\infty}}
\def\izi{\int_{0}^{\infty}}
\def\izpi{\int_{0}^{\pi}}
\def\izdpi{\int_{0}^{2\pi}}
\def\intd{\int\kern-.8em\int}
\def\intt{\int\kern-.8em\int\kern-.8em\int}
\def\intg{\int\kern-1.1em\int}
%
\def\sumd{\mathop{\sum\sum}}

\def\sumi{\sum_{i=1}^{M}}
\def\sumj{\sum_{i=1}^{N}}
\def\sumk{\sum_{k=1}^{K}}
\def\sumn{\sum_{n=1}^{N}}
\def\summ{\sum_{m=1}^{M}}

%
%TRANSFORMS
\def\TA#1{{\cal A}\left\{ {#1} \right\}}
\def\TH#1{{\cal H}\left\{ {#1} \right\}}
\def\TP#1{{\cal P}\left\{ {#1} \right\}}
\def\TR#1{{\cal R}\left\{ {#1} \right\}}
\def\TRa#1{{\cal R}^{\dag}\left\{ {#1} \right\}}
\def\BR#1{{\cal B}\left\{ {#1} \right\}}
\def\TF#1{{\cal F}\left\{ {#1} \right\}}
\def\TFI#1{{\cal F}^{-1}\left\{ {#1} \right\}}
\def\TFn#1#2{{\cal F}_{#1}\left\{ {#2} \right\}}
\def\TFnI#1#2{{\cal F}_{#1}^{-1}\left\{ {#2} \right\}}
%
\def\Im#1{{\cal I}\mbox{m}\left(#1\right)}
\def\Ker#1{{\cal K}\mbox{er}\left(#1\right)}
%
\def\Imag#1{\mbox{Im}\left(#1\right)}
\def\Re#1{\mbox{Re}\left(#1\right)}
%
\def\expf#1{\exp\left[ {#1} \right]}

%derivee
\def\dfdx#1#2{{\mbox{d} {#1}\over{\mbox{d} {#2}}}}
\def\dfdxd#1#2{{\mbox{d}^2 {#1}\over{\mbox{d} {#2}^2}}}
\def\dfdxt#1#2{{\mbox{d}^3 {#1}\over{\mbox{d} {#2}^3}}}
\def\dfdxn#1#2{{\mbox{d}^n {#1}\over{\mbox{d} {#2}^n}}}
\def\dfdxk#1#2{{\mbox{d}^k {#1}\over{\mbox{d} {#2}^k}}}

%derivee partielle
\def\dpdx#1#2{{{\partial {#1}\over \partial {#2}}}}
\def\dpdxd#1#2{{{\partial^2 {#1}}\over{\partial {#2}^2}}}
\def\dpdxdy#1#2#3{{{\partial ^2 {#1}}\over{\partial {#2} \partial {#3}}}}

\def\arg{\mbox{arg}}
\def\argmins#1#2{\mbox{arg}\min_{#1}\left\{{#2}\right\}}
\def\argmaxs#1#2{\mbox{arg}\max_{#1}\left\{{#2}\right\}}
\def\argmin#1#2{\mathop{\mbox{arg}\min}_{#1}\left\{{#2}\right\}}
\def\argmax#1#2{\mathop{\mbox{arg}\max}_{#1}\left\{{#2}\right\}}

\def\esp#1{\mbox{E}\left\{ #1 \right\}}
\def\espx#1#2{\mbox{E}_{#1}\left\{ #2 \right\}}

\def\wth#1{\widehat{\widetilde{\phantom{#1}}}\!\!\!\! #1}

%OTHERS
\def\lrf{L_{r,\phi}}
\def\fw{\widehat{f}(\omegab)}
\def\fthwxi{\wth{f}(\Omega,\xib)}
\def\fthwfi{\wth{f}(\Omega,\phi)}
\def\ftrfi{\widetilde{f}(r,\phi)}

\def\fwxwy{\widehat{f}(\omega_x, \omega_y)}
\def\wxpwy{(\omega_x \, x + \omega_y \, y)}

\def\wtx{\omegab^t \cdot \xb}
\def\ejwtx{\exp\left[j \omegab^t \cdot \xb\right]}
\def\xitx{\xib^t \cdot \xb}

%
\def\ftrxi{\widetilde{f}(r,\xib)}
\def\ent{-\int p(x) \, \ln p(x) \d{x}}

%YOURS

\def\mean#1{\left< #1 \right>}
\def\slnhn{\sum_{n=1}^N \lambda_n h_n(\rb)}
\def\slngn{\sum_{n=1}^N \lambda_n g_n(\rb)}
\def\smngn{\sum_{n=1}^N \mu_n g_n(\rb)}
\def\slmhm{\sum_{m=1}^N \lambda_m h_m(\rb)}
\def\vlambda{\bm{\lambda} = [\lambda_1,\ldots,\lambda_n]}

\def\apriori{{\em a priori} }
\def\aposteriori{{\em a posteriori} }

\def\titre#1{\bcc{\Large\bf #1}\ecc}

% My definitions
\def\AMD{Ali Mohammad--Djafari}
\def\LSSa{Laboratoire des Signaux et Syst\`emes 
(CNRS--ESE--UPS) \\ 
\'Ecole Sup\'erieure d'\'Electricit\'e \\ 
Plateau de Moulon, 91192 Gif sur Yvette Cedex, France.}

\def\ME{maximum entropy}
\def\pdf{probability distribution function}
\def\lm{Lagrange multipliers}
%%
\def\fix#1{\phi _#1(x)}
\def\fin{\fix n}
\def\fik{\fix k}
\def\fiz{\fix 0}
\def\sfinz{\sum_{n=0}^N \lambda_n \, \fin}
\def\sfinu{\sum_{n=1}^N \lambda_n \, \fin}
\def\bl{\bm{\lambda}}
\def\bd{\bm{\delta}}
\def\blz{\bl ^0}
\def\gnl{G _n(\bl)}
\def\gnlz{G _n(\blz)}
\def\un{n=1,\dots, N}
\def\nn{n=0,\dots, N}

\def\finn{\fin , \nn}
\def\esfinz{\exp\,\left[ -\sfinz \right] }
\def\esfinu{\exp\,\left[ -\sfinu \right] }
\def\esxm{\exp\,\left[ -\sum_{m=0}^N \lambda_m \, x^m \right] }
\def\efin{\esp \fin }
\def\zl{Z(\bl)}
\def\finxi{\phi _n(x_i)}
\def\snfinxi{\sum_{n=1}^N \lambda_n \finxi}
\def\esnfinxi{\exp \left[ - \snfinxi \right]}
\def\smfinxi{\sum_{i=1}^M \finxi}

\def\ejnw{\exp \left( -j n \omega_0 x \right) }
\def\eejnw{\mbox{E} \left\lbrace \ejnw \right\rbrace}

\def\signed#1{{\unskip\nobreak\hfil\penalty50\hskip2em\mbox{}
\nobreak\hfil\tt#1\parfillskip=0pt \finalhyphendemerits=0 \par}}

\def\uncatcodespecials{\def\do##1{\catcode`##1=12 }\dospecials}
\def\listing#1{\par\begingroup\setupverbatim\input#1 \endgroup}
\newcount\lineno
\def\setupverbatim{\tt \lineno=0
 \obeylines \uncatcodespecials \obeyspaces
 \everypar{\advance\lineno by1 \llap{\sevenrm\the\lineno\ \ }}}
{\obeyspaces\global\let =\ }

\def\defined{\stackrel{\mbox{def}}{=}}
\def\str{\stackrel}

\def\ER{\mbox{I\kern-.25em R}}
\def\EC{\mbox{C\kern-.8em C}}
\def\EZ{\mbox{Z\kern-.55em Z}}
\def\EN{\mbox{N\kern-.8em N}}

\def\singles{
 \abovedisplayskip 12pt plus 3pt minus 9pt
 \belowdisplayskip 12pt plus 3pt minus 9pt
 \abovedisplayshortskip 0pt plus 3pt
 \belowdisplayshortskip 7pt plus 3pt minus 4pt
 \baselineskip 14.4pt
 \lineskip 1pt
 \lineskiplimit 0pt}
\def\oneandhalf{
 \abovedisplayskip 18pt plus 3pt minus 9pt
 \belowdisplayskip 18pt plus 3pt minus 9pt
 \abovedisplayshortskip 0pt plus 3pt
 \belowdisplayshortskip 9.333pt plus 3pt
 \baselineskip 20pt
 \lineskip 2pt
 \lineskiplimit 1pt}

\def\double{
 \abovedisplayskip 24pt plus 3pt minus 9pt
 \belowdisplayskip 24pt plus 3pt minus 9pt
 \abovedisplayshortskip 0pt plus 3pt
 \belowdisplayshortskip 12pt plus 3pt
 \baselineskip 27pt
 \lineskip 3pt
 \lineskiplimit 2pt}

\def\dadb{\d{\alpha}\d{\beta}}

\def\ffbox#1{\fbox{\mbox{\vbox{#1}}}}

\def\rot{\mbox{rot}}
\def\case#1#2#3#4{
    \left\{
           \begin{array}{ll}
            {\displaystyle #1} & {\displaystyle #2} \cr 
            {\displaystyle #3} & {\displaystyle #4}
           \end{array}
    \right. }

\def\beqnarr#1&#2&#3\\#4&#5&#6\eeqnarr{
    \left\{
           \begin{array}{lcl}
            {\displaystyle #1} & #2 & {\displaystyle #3} \\ 
            {\displaystyle #4} & #5 & {\displaystyle #6} 
           \end{array}
    \right. }

\def\pyx{p(\yb|\xb)}
\def\pxy{p(\xb|\yb)}

%\newtheorem{definition}{D\'efinition}
%\newtheorem{theorem}{Th\'eor\`eme}

\def\ie{{\em i.e.}}
\def\unsdpi{\left(\frac{1}{2\pi}\right)}
\def\unspi{\left(\frac{1}{\pi}\right)}
\def\up{\uppercase}
\def\zjm{z_{j-1}}
\def\zjp{z_{j+1}}
\def\fxyp{f(x,y)=\left\{
\barr{ll} 1 & (x,y)\in P\\ 0 & (x,y)\not\in P\earr
\right.}


paper1.bbl
\begin{thebibliography}{}

\end{thebibliography}
