·Table of Contents ·Computer Processing and Simulation | Dynamic Image Reconstruction: General estimation Principles for Dynamic tomographyV. M. Artemiev, A. O. NaumovInstitute of Applied Physics, Akademicheskaya str.16, 220072, Minsk, Belarus. G.-R.Tillack Federal Institute for Materials Research and Testing (BAM), Unter den Eichen 87,12205 Berlin, Germany. Contact |
A dynamic Gaussian random field with the above properties can be represented by the following S-dimensional linear stochastic difference equation (SDE)
x(k+1)=Ax(k)+u+w(k) | (2.1) |
with the S ´ S dimensional constant matrix A and the S dimensional constant vector u . The S dimensional random vector w(k) represents Gaussian distributed white noise with zero mean and constant covariance matrix
given by the probability density distribution (PDD) function
(2.2) |
The SDE image model (2.1) define the image field as a multidimensional Markovian random sequence. Assuming Gaussian initial conditions the process x(k) remains Gaussian for all time. Additionally it is supposed that the process
x(k) tends to a stationary regime, i.e. the field PDD function tends to a stationary Gaussian distribution for k(r)¥
The Markovian random sequence (2.1) is a priori completely defined by two functions [3,4]:
(2.3) |
The S dimensional vector m_{o} gives the image field prior mean value and the S ´ S matrix p_{o} describes the prior covariance matrix in the space domain. The corresponding TP function is given by
(2.4) |
In the next step the matrices A,Q as well as the vector u of the SDE image model (2.1) have to be determined for the given prior m_{o} and P_{o}.
From eq. (2.1) the following expression for the image mean value m_{o}(k)=áx(k)ñ is found
(2.5) |
(2.6) |
Hence the stationary field mean value m_{o} is determined by the choice of the constant vector u if the matrix A is known.
The SDE for the centralized image vector follows from eq.(2.1) and eq.(2.5)
(2.7) |
The space-time S´S covariance matrix P_{0}(k+1)|k)can be introduced by the expression
(2.8) |
with the time delay parameter l =0,1,2,.... For l=0 the covariance matrix (2.8) transfers to the space covariance matrix P_{o}(k)=P_{o}( k | k ).
As for eq. (2.7) the equation for the image vector are given by
(2.9) |
By substituting the equations for the image vectors into eq. (2.9) the following recursive equation is obtained
(2.10) |
In the steady state and .The S´S matrix F_{o}(l) is the stationary space-time covariance matrix that satisfies
(2.11) |
For l=0 the matrix F_{o}(0) transforms to P_{o} with
(2.12) |
After multiplying matrix A^{1} by eq.(2.12) and substituting the result into eq. (2.11) the final equation for the stationary space-time covariance matrix is obtained
(2.13) |
Eqs.(2.6), (2.12) and (2.13) define the mean value m_{o}, space covariance matrix P_{o} and stationary space-time covariance matrix F_{o}(l) for given A,Q and u .
The following inverse equations are valid to obtain the model parameters A , Q and u for given P_{o},F_{o}(l) and m_{o} :
(2.14) |
It is important to note that the solution of eq. (2.14) is unique for arbitrary mean image vector m_{o} and space covariance matrix P_{o} but only for one chosen time delay parameter l .
The stochastic boundary conditions for the field model (2.1) can be specified by choosing the mean image vector m_{o}and the space covariance matrix P_{o} on the boundary of the hyper parallelepiped under consideration.
(3.1) |
with the 1´S row matrix H_{1}(k) having S-S_{k} elements equal to zero. For the following time k + 1 the observation is given by another set of S_{k+1} < S components of the vector x(k+1) with z_{1}(k +1)=H_{1}(k +1)x(k + 1). The observations z_{1}(k) for k=1,2,... are scalars and the corresponding observation matrices have the properties
(3.2) |
Here the coefficient c_{k} give the square sum of the matrix elements of matrix H_{1}(k). If these elements hold only the values zero or one c_{k} gives the number of non-vanishing elements of H_{1}(k). The coefficients c_{kl} describe the degree of the matrices' intersection. If the coefficient c_{kl} is equal to zero the intersection of the matrices is the empty set.
The m-dimensional vector z_{m}(k) holds the simultaneous observation of m ray sums(3.1) defined by
z_{m}(k)=H_{m}(k)x(k) | (3.3) |
with the m ´ S observation matrix H_{m}(k)which is built from the row matrices H_{li}(k) ,i=1,2,...,m , i.e.
H_{m}(k)=[H_{11}(k),...,H_{1i}(k),...,H_{1m}(k)]^{T} | (3.4) |
(3.5) |
(3.6) |
The dynamic image reconstruction process is fully observable if Rank(H_{1}^{k})=S.
It is supposed that the number of projections is equal toN and that each projection consists of M ray sums. Then the rank of the observation matrix (3.6) can not be greater than NM. In practical applications it usually turns out that Rank(H_{1}^{N})<S and the dynamic image reconstruction is an ill-posed problem [3,4]. This is the main difficulty one has to face compared to standard estimation problems.
Usually the observation is corrupted by additive noise v(k)and linear image blurring given by y(k)=D(k)x(k) with the blurring matrix D(k). As a result the observation takes the form
(3.7) |
(3.8) |
(3.9) |
(3.10) |
(4.1) |
(4.2) |
(4.3) |
(4.4) |
(4.5) |
(4.6) |
(4.7) |
(4.8) |
P(k+1|k)=AP(k)A^{T}+Q | (4.9) |
(4.10) |
(4.11) |
The static random image is a partial case of the above dynamic image(2.1) with A =I,Q = 0 ,u=0 given by the SDE
x(k +1)=x(k) The Bayesian approach for the reconstruction of random static images was studied for a long period of time after the publication[2] where the one-step optimal linear reconstruction algorithm was found which is similar to proposed KF algorithm(4.11).
Fig 1: Geometrical setup. |
The image vector is supposed to haveS=250 x 250 image elements, the number of detector elements is chosen to M = 400, and the total number of source positions is N_{o}=50. The image is assumed to be a Markovian random sequence with zero-mean value, axial symmetric covariance matrix with covariance radius of 30 image elements. The relaxation time depends on the motion velocity and is chosen to 200 time intervals.
The image reconstruction is carried out sequentially applying the KF algorithm (4.5)-(4.9) for the case m = M. Two reconstruction variants are investigated. First it is supposed that the image reconstruction procedure is applied to the quasi-static regime using 50 projections acquired during the time interval between the reconstruction steps kand k+1. The reconstruction performance, given by the criteria(4.10), is shown in Fig.2 with representing the number of projections actually used for the reconstruction. Each reconstruction step is performed without taking into account any dynamic properties of the image, i.e. the time covariance matrix (2.11) is neglected. This procedure yields the reconstruction error of at the end of every reconstruction step.
Fig 2: Reconstruction error for quasi-static regime | Fig 3: Reconstruction error for dynamic regime |
In the second variant the dynamic reconstruction technology is employed. Here only 10 projections are included in each reconstruction step within the same time interval as used in the quasi-static regime. To account for the dynamic properties of the image in this case the time covariance matrix (2.11) is employed. Fig.3 shows the performance of the dynamic reconstruction. At the end of intervals the same accuracy is obtained as in the previous case. Thus the dynamic reconstruction requires a smaller number of projections per reconstruction step to reach the same accuracy because of considering the time covariance matrix as a kind of prior knowledge about the dynamic properties of the image to complete the previously introduced spatial prior.
The reconstruction results are shown in fig.4. Here fig. 4a gives a part of the original image represented over a total time of 300 intervals:k=1,...,300. . Fig4b represents the result of the dynamic reconstruction for the same period of time considering a total number of 3000 projections.
Fig 4: Reconstruction result for dynamic regime (b) in comparison to original image (a) |
© AIPnD , created by NDT.net | |Home| |Top| |