TAKING INTO ACCOUNT A PRIORI INFORMATION IN THE ITERATIVE RECONSTRUCTION OF IMAGES OF FOUNDRY PRODUCTS

. Methods of restoring images and properties of non-destructive testing objects based on solving inverse problems (problems of restoring distribution functions of unknown characteristics of an object based on the results of indirect measurements) are considered. Management methods are based on solving inverse problems and allow you to get the most complete information about the distributed properties of an object. The need to attract additional information imposes serious restrictions on the development of universal applied algorithms for solving incorrectly set tasks. As a rule, individual additional information is available for each specific non-destructive testing task. An effective numerical algorithm for solving an incorrectly posed problem should be focused on taking this information into account at each stage of the solution search. When solving an applied problem, it is also necessary that the algorithm corresponds to both the measuring capabilities and the capabilities of available computing tools. The problem of low-projection X-ray tomography is always associated with a lack of initial data and can only be solved using a priori information. To introduce the necessary additional information into the numerical algorithm, the methods of iterative reconstruction of tomographic images are identified as the most suitable. One of the approaches to the presentation of this kind of information is described. A practical solution to this problem will expand the scope of the X-ray tomography method.


Introduction.
Recently, a number of new directions have emerged in the field of iterative image reconstruction, for example. in [1], the method of iterative reconstruction is considered, which uses graphics processors to speed up calculations.At the bottom, a method based on model iterative reconstruction (Model-Based Iterative Reconstruction) is also used, including all influencing factors: the geometric scheme of data collection, taking into account different types of a priori information about the object, regularization of the reconstruction process using total variation.All the improvements proposed by the authors are intended for breast image reconstruction, but can also be applied to industrial tomography.The aim of this paper [2] was to propose a new TV-based optimization framework for the reconstruction of DBT images.The framework, described in a rigorous numerical setting, includes both constrained and unconstrained models, thus it is a flexible tool easily enabling the use of different data-fitting or regularization terms as well as the addition of further box constraints, to reconstruct reliable volumes from subsampled noisy data.They are also interested in finding an automatic strategy to set the regularization parameter, which strongly affects the quality of the reconstructions, in order to avoid its manually tuning (which is infeasible in a clinical setting).Recently [3], an extended family of power-divergence measures with two parameters was proposed together with an iterative reconstruction algorithm based on minimization of the divergence measure as an objective function of the reconstructed images for computed tomography.Numerical experiments on the reconstruction algorithm illustrated that it has advantages over conventional iterative methods from noisy measured projections by setting appropriate values of the parameters.In this paper, they present a novel neural network architecture for determining the most appropriate parameters depending on the noise level of the projections and the shape of the target image.Through experiments, they show that the algorithm of the architecture, which has an optimization sub-network with multiplicative connections rather than additive ones, works well.
Let's take into account that, among the methods of non-destructive testing, a special place is occupied by the method of X-ray tomography.The advantage of this method is that its information content about each elementary volume of the object under study is many times higher than in other known methods.High efficiency of X-ray tomography method was first demonstrated on examples of its use in medicine and biology.
Technical implementation and scope of the X-ray tomography method are rigidly connected with the mathematical apparatus, which is used to reconstruct the spatial distribution of the coefficient of linear attenuation of X-ray radiation from X-ray projections, and the capabilities of computer technology.The first X-ray tomographs were based on a two-dimensional inverse Radon transformation, which showed that the distribution of the linear attenuation coefficient μ(x, y) in an infinitely thin layer of an object is determined by the totality of all linear integrals [4]: where p(l, θ) is the integral of the function μ(x, y) along the line, which is at a distance l from the origin, and its perpendicular makes an angle θ with the axis x.
Based on transformation (1), the principle of modern computed tomography is as follows.The object under study is irradiated at various angles, the recorded results of measurements of the characteristics of radiation passing through a thin layer of the object are converted into numerical codes.These codes are entered into a computer, where the spatial distribution of the linear attenuation coefficient is determined by numerical inverse Radon transformation, which is then visualized using appropriate devices [5, p. 12].It should be noted that, on the one hand, the use of reverse Radon conversion is most effective when using narrowly collimated radiation sources (in this case, the nonlinear contribution of scattered photons is small) and electronic radiation detectors (due to the large volume of measured data), on the other hand, it requires a strictly defined viewing scheme of the object, requiring a significant number of angles (usually several hundred) for a satisfactory approximation of the function p(l, θ) when the detector-source system moves relative to the object at an angle of 180 degrees or 360 degrees.
Further improvement of X-ray tomographs due to the development and use of the mathematical apparatus of three-dimensional generalizations of the Radon transformation [6,7] allowed the use of a conical radiation source and a two-dimensional matrix detector, which significantly increased the speed of information collection for three-dimensional reconstruction.But practical applications have been found by tomographs that use some approximate three-dimensional transformations.
Iterative Algebraic Methods for Reconstructing Images of Objects.Small amount of projection data, which is typical for industrial low-projection tomography systems, makes it necessary to develop reconstruction algorithms that allow the use of additional a priori information about the reconstructed object and are easily adaptable to various non-traditional scanning schemes.The choice of a specific reconstruction algorithm is determined, first of all, by the features of the tomography system, more precisely, by the features of the projection data recorded in the system.
There are situations when, due to the features of the registration system, either the number of projections is very small, or they are not fully known, or obtained in a limited range of angles.Restoration problems with such a specification of projection data become strongly ill-posed and the uniqueness and stability theorems are violated for them.Therefore, the direct application of analytical algorithms to solve them does not give acceptable results, and it is necessary to develop other, more suitable algorithms.In these cases, the most preferable methods are expansions into finite series of orthogonal functions.
Restoring images from their projections using methods of decomposition into finite series differs fundamentally from the methods of integral transformations.The fact is that the methods of integral transformations use operators acting on a set of functions that are given on the entire continuum of real numbers.And only to implement the solution on computers and workstations, continuous operators are replaced by discrete operators that act on functions defined on finite sets, but this is done at the very end of the recovery algorithm execution procedure.
Decomposition methods into finite series imply image discretization prior to the start of the restoration algorithm, which reduces the main restoration problem to solving systems of linear (or non-linear) algebraic equations, i.e. to problems of computational linear algebra.It is known that the system of linear algebraic equations (SLAE) can be solved using both direct and iterative methods.For medium-sized systems, direct methods are often the most preferable.Iterative methods are used mainly for problems of very large dimensions, which, as a rule, are problems of computed tomography.
To apply finite series expansion algorithms to solving problems of reconstructive computed tomography, it is first necessary to construct a complete discrete model of the restoration problem.To do this, we choose a system of basis functions b j (x), j = 1,2,.., J and consider the discretized form f  of the func- tion f with respect to this system: where f j = const.Let R i , i = 1,2,.., I be a set of linear continuous functionals, each of which associates the function f(x) with a real number R i f, which we denote by p i .Then, applying the operators R ito equality (2), taking into account their continuity and linearity, and denoting we get the desired SLAE Af = p, where is the image vector, p= (p 1 ,..., p I ) -projection vector.Thus, the solution of the main problem of restoring an image from a given set of projections p i reduces to solving the SLAE of the form (4), while the vector p is known to be given with some error.There are several different sources of noise for the vector p: 1) instead of the function f(x), we consider its approximate expression as a discretized function f  (x); 2) the projection vector p is obtained in practice from experimental or practical physical measurements, which are inevitably noisy; 3) when calculating elements a ij of matrix A, we assume the idealization of natural rays, as will be discussed below.
Assume that the function f(x) is given in the spatial domain V. Let's divide the region V into some finite number of sub-regions, which we will call voxels.(Voxel is a 3D image element.)Let's renumber all voxels in some convenient way from 1 to N. We will assume that the restored function f(x) takes a constant value f j inside the j-th voxel, i.e. the function f(x) is replaced by its discretized expression: -th be a system of basis functions defined as follows: ( ) ( ) -th Let us assume that as linear functionals R choose the direct Radon transform along a set of some lines L s : geometrically represents the length of the intersection of the i-th ray with the j-th voxel.
From the foregoing, the flexibility of the proposed algebraic approach to solving restoration problems is obvious.First, it should be noted the freedom of choice of basis functions , which is associated with the choice of the partition grid of the domain D (in addition to the Cartesian square grid, one can choose a polar grid in the two-dimensional case and, in general, a grid of arbitrary geometric configuration).
Moreover, the formulation of the problem does not depend on the geometry of the rays.In this case, only the elements a ij change, which, as a rule, are calculated in advance.It should also be noted that all of the above can be transferred to the three-dimensional case.
Thus, the algebraic method for a fixed choice of continuous linear operators R i and basis functions b j reduces the recovery problem to the solution of SLAE (4), i.e., it would seem, to a standard problem of computational linear algebra.However, as applied to the problem of image reconstruction, this problem has a number of characteristic features: the dimension of the system is extremely large: as a rule, the number of equations and unknowns is on the order of 10 7 -10 10 ; the projection matrix A = (a ij ) is very sparse, since each ray intersects a very small number of pixels, so the overwhelming number of its elements is zero (approximately only 1 % of the matrix A coefficients are non-zero).In this case, the non-zero elements of the matrix A, as a rule, do not form any specific ordered structure in order to be able to apply any of the known methods for solving systems of linear algebraic equations with sparse matrices; matrix A is a rectangular matrix m × n, as a rule, m < n (here m is the number of equations, n is the number of unknowns).In the latter case, the system is underdetermined; the system of equations ( 4) is unstable with respect to the initial data, i.e. small changes in the vector р = (p 1 ,..., p m ) can respond to arbitrarily large changes in the vector of unknowns f = These characteristic features turn a seemingly simple problem of computational linear algebra into a computational procedure that is very difficult to implement.In this case, for each specific problem, it is necessary to select the appropriate reconstruction algorithm.
Classical Algorithm for Algebraic Reconstruction.Here we consider the solution of the problem of reconstructive computational tomography using classical iterative algebraic restoration methods based on expansion into finite orthogonal series [8].For the case of specifying systems of basis functions in the form of formulas (6), this leads to the need to solve the SLAE of the form Ax = p, (8) where А = (a ij ) - projection matrix, х = . .
-image vector, р = (p 1 ,..., p m ) -projection vector.For the first time, the application of an iterative algebraic algorithm for computed tomography was described in 1970 in [9].In this algorithm, an arbitrary value x (0) is chosen as the initial approximation of the image vector R n , if the system is joint, the (k+1)-st iteration is obtained from the k-th iteration in some additive way.In this case, only one ray, for example, the i-th ray, is considered sequentially in the enumeration order, and only those components of the vector x (k) that correspond to the pixels intersected by this ray are subject to change.The value of the discrepancy between the measured value p i and the approximate value of the projection ( ) k ij j j a x ∑ obtained by substituting the k-th iteration х (k) is redistributed between the pixels located along the i-th beam in proportion to their weights a ij in the beam.Thus, in one cycle of the k-th iteration, the values of only those pixels that are along the given i-th ray are changed, and the remaining values of the function remain unchanged.This algorithm is known in the literature as ART and can be defined as follows: 1) Initial approximation x (0) of R n is chosen arbitrarily; 2) (k+1)-st iteration is calculated by the formula where a i = (a ij ) 1 n j= is the i-th row of the matrix A, the relaxation parameters λ (k) are a sequence of real numbers, and i = i k = k(mod m ) + 1, i.e. rays move cyclically.This algorithm is characterized by efficient use of memory, since the n -dimensional vector x (k+1) can be stored in the same block of RAM as the vector x (k) , the need to store which disappears after the k-th iteration.Researchers have proven that it belongs to the class of iterative algorithms based on the method of projecting onto convex sets.
Image Reconstruction Algorithms with Using a Priori Information.As noted above, in practice it is often necessary to solve the so-called restoration problems with an incomplete data set.This data set is usually not enough to solve the recovery problem.However, for many of these problems, in particular for the problems of studying physical processes, it happens that some additional information about the original object being restored is known.This various a priori information can be effectively taken into account in algebraic iterative algorithms.The use of this a priori additional information often provides a significant improvement in the rate of convergence of iterative processes.In this subsection, various combinations of iterative algebraic reconstruction algorithms considered in the previous sections are considered, using the most important and available a priori information about the object being restored, and the convergence properties of such algorithms are established.
Let us assume that the iterative algebraic algorithm for solving SLAE ( 4) is given in the form of the following recurrent formula where Q : R n → R n is some operator.Additional information about the object in the general case can also be specified using some operator S : R n → R n as then the general iterative process will be represented as a superposition of these two operators: ( . It should be noted that the convergence of sequence (10) to the solution of system (4) does not always imply the convergence of sequence ( 12) and vice versa.The operators Q and S must satisfy certain properties in order to guarantee this convergence.Suppose it is known a priori that the desired solution x belongs to some subspace L ⊂ R n .Then as the operator S we can take the operator of orthogonal pro- jection Р = Р L , onto the subspace L. For example, if the boundaries of the values of the desired solution vector x are known a priori, i.e. it is known that the i-th component of x i is in the interval [α i , β i ], then as the set L the direct product of these intervals can be taken: ] In this case, the orthogonal projection operator As a priori knowledge about an object, the boundaries of its extent can also be known in advance, i.e. known set S ⊂ R n , outside of which the restored function is equal to 0. In this case, the projection operator onto the subset S defined as We will apply the principle (13) described above of introducing a priori information into an iterative process for the reconstruction of images of aluminum cases.The first set of 32 X-ray projections (problem for low angle tomography) was provided by M. Simon [10].All these projections were obtained in the angular range of 180 degrees, that is, with one-sided access to the reconstruction object.The number of projection views was several dozen, and all of them were located in an incomplete angular range of less than 360 degrees.
Construction of the Initial Approximation for Calculating the Image of an Aluminum Case.To reduce the number of voxels for which it is necessary to estimate the approximate value of the linear attenuation coefficient, it is advisable to reduce their number, based on the principle of observability of the reconstructed object in all projections.
Define the three-dimensional area of the estimated position of the object where η(x) is the threshold function: This means that the voxel j belongs to Volume if the values of its corresponding pixels Pixel jn , into which it is projected, are greater than some minimum value on all projections.
According to (15) voxel j belongs to Volume if the values of the corresponding pixels exceed the measurement error ε jn on all N projections.
Perspective images of the reconstructed object for a different number of iterations, when the entire reconstruction area is taken as an initial approximation, as well as the volume of the area constructed in accordance with rule (15), are given in [11].
Volume region is taken as the initial approximation.It is clear that in this case we take into account a priori information about the spatial position of the reconstruction object.The role of subset S in relations (14) plays a subset of the reconstruction zone Volume.It can be seen that when using the second method of choosing the initial approximation, we more accurately and efficiently restore the shape and internal structure of the object (Figure 1).

Image Reconstruction with Using a Priori Information about the Boundaries of the Solution
Values.Since the attenuation coefficient of the aluminum body is uniform, we can use condition (13) to improve the quality of the image reconstruction.By multiplying all projections by some constant coefficient, one can calibrate them so that the density of the reconstructed object becomes equal to one.Yu.B. Denkevich in [12] to improve the fulfillment of conditions (13) for voxel values whose boundaries belong to the segment 0 = < x j < = 1, j = 1…N suggested simultaneously minimizing the functional B(x) = (4.0*x*x]-6.0*x + 2.0)).A numerical algorithm for solving our problem can be constructed using the variational problem of minimizing the functional: where x  is an approximate solution of the equation ˆ, x p = A ||Ax -p|| 2 is the discrepancy, B(x) is the stabilizing functional, α is the regularization parameter.
Figure 2 shows an image of an aluminum block, which is obtained by solving the variational problem (16).
The second set of 2,000 X-ray projections was provided by I. Georgiev, Research Fellow at the Institute of Information and Communication Technologies (IICT) of the Bulgarian Academy of Sciences.They were obtained on an industrial tomograph XT H 225 from Nikon Metrology.The angular step between projections was 0.18 degrees.The pixel size is 0.127 mm.The resolution of the detector matrix Figure 3 shows X-ray projections of an aluminum part for angles of 0 degrees, 45 degrees and 90 degrees.First, a virtual space was built according to relations (15), shown in Figure 4.
Figure 5 shows a 3D perspective view of an aluminum part.The image reconstruction was carried out by the iterative SART method using a GPU [5], which made it possible to carry out several high-dimensional reconstructions with areas of both 512 × 512 × 512 voxels and 1024 × 1024 × 1024 voxels.