Method for reconstructing an object by x-ray tomography, comprising a calculation of regularization terms
The method addresses the challenges of noisy and incomplete data in X-ray tomography by calculating voxel-by-voxel regularization terms from tomographic incompleteness maps, enhancing image quality and simplifying hyperparameter tuning for diverse X-ray systems.
Patent Information
- Application Number
- PCT/EP2024/068894
- Authority / Receiving Office
- WO · WO
- Patent Type
- Applications
- Current Assignee / Owner
- Priority Date
- 2023-07-07
- Filing Date
- 2024-07-04
- Publication Date
- 2025-01-16
AI Technical Summary
Current X-ray tomography reconstruction algorithms face challenges with noisy and incomplete data, leading to inconsistent solutions, and require manual tuning of hyperparameters, which is difficult due to the need for reference images and diverse data, especially when dealing with limited angles and tomographic incompleteness.
A method that calculates voxel-by-voxel regularization terms based on tomographic incompleteness, adapting to each scanner geometry without requiring experimental data or human intervention, by pre-calculating directions, computing incompleteness maps, and using these to adjust regularization terms for improved reconstruction.
This approach allows for precise, locally adapted regularization, improving image quality in scenarios with missing data and limited angles, and simplifies hyperparameter adjustment, applicable to various X-ray tomography systems, including medical imaging and baggage screening.
Smart Images

Figure EP2024068894_16012025_PF_FP_ABST
Abstract
Description
DESCRIPTION Title of the invention: Method for reconstructing an object by X-ray tomography comprising a calculation of regularization terms.
[0001] The invention relates to the field of X-ray imaging, and more particularly to the field of X-ray tomographic reconstruction. It relates to a method for reconstructing an object by X-ray tomography comprising a calculation of regularization terms.
[0002] The invention finds its application in the field of security, in particular for baggage searches taking place in sensitive sites (airports, train stations, etc.), during events or on any occasion of baggage inspection. It also finds its application in the medical field, for example when using X-ray medical imaging equipment. More generally, the invention finds its application when it is desired to reconstruct an object by tomography.
[0003] X-ray computed tomography provides a two- or three-dimensional image from a series of tomographic projections via a reconstruction algorithm. Reconstruction algorithms used in tomography have weaknesses when the measured data are highly noisy and / or missing. The reconstruction algorithm then faces the lack of data and can then lead to very inconsistent solutions.
[0004] In order to reconstruct an object, the reconstruction algorithm then requires introducing, in addition to input data, a priori on the object to be reconstructed. The a priori on the object are introduced via a regularization, thus giving rise to a regularized reconstruction algorithm.
[0005] Regularized reconstruction algorithms thus solve an inverse problem from projection data. This family of algorithms requires the adjustment of regularization terms. These regularization terms aim to provide a compromise between the information in the data and the a priori on the object to be reconstructed. Regularization makes it possible to introduce an a priori on the statistical distribution of the imaged object and / or take into account the acquisition geometry. This regularization depends on parameters called here hyper-parameters. The hyper-parameters are calculated before the object reconstruction process, as detailed later.
[0006] Most of the time, regularization terms are constructed to limit the number of hyperparameters to be tuned. Indeed, the larger the number of hyperparameters, the more difficult it is to tune them manually.
[0007] To automatically tune hyper-parameters, some approaches rely on cross-validation methods [Kohavi, 1995], or more recently on machine learning methods [Shen, 2018], [Ding, 2021]. In 2021, work [Zhang, 2021] showed that a Directional Total Variation type regularization is well suited for reconstruction problems with limited angles. This regularization nevertheless requires tuning more hyper-parameters than a conventional Total Variation [Rudin, 1992].
[0008] Current methods for adjusting regularization terms require reference images and / or a large amount of diverse data. While it is possible to obtain reference images through simulations, these methods remain difficult to obtain under experimental conditions.
[0009] Furthermore, ideally, the choice of regularization terms should incorporate the geometry of the scanner, i.e. the position of the source(s) and the detector(s).
[0010] 3D computed tomography theory has some tools to provide local information about the object to be reconstructed according to the scanner geometry. Assuming untruncated projections, [Tuy, 1983] gives a condition to check whether a continuous X-ray source trajectory is sufficient to reconstruct the object. The condition can be written as follows: "every plane that intersects the imaged object must intersect the source trajectory at least once". In practice, all scanners are limited by a finite number of source positions and therefore Tuy's condition is never met. We speak of tomographic incompleteness when Tuy's condition is not met.
[0011] Figure 1 shows examples of source-detector trajectories T of a scanner used for object reconstruction. Example a) shows a circular source-detector trajectory. Such a trajectory does not allow for exact reconstruction of the object outside the plane of the trajectory. Indeed, a plane intersecting the object and which is parallel, for example, to the circle of the trajectory T does not intersect the said source-detector trajectory. A helical geometry of the trajectory T, as represented in example b), allows for an exact reconstruction to be obtained. Indeed, whatever plane is taken, it necessarily intersects the helix of trajectory T.
[0012] Several metrics have been studied to quantify the impact of tomographic incompleteness. [Metzler, Bowsher, & Jaszczak, 2003], [Lin & Meikle, 2010] and [Liu, et al., 2012] numerically evaluate the percentage of planes that are intersected by the source trajectory by sampling the unit sphere. However, this type of measurement is only possible for a continuous source trajectory. However, in practice, acquisitions are performed at discrete points along the source trajectory. Thus, planes necessarily pass between the sources. Furthermore, the authors of these publications indicate that it is difficult to predict the quality of the reconstruction from this metric.
[0013] [Stopp, Winne, Jank, & Keeve, 2010] proposed to measure locally the average of the minimal angles, that is, the angles defined by the planes passing through the voxel under study and the associated optimal source. [Clackdoyle & Noo, 2020] also quantify tomographic incompleteness to measure how far a voxel is from respecting the Tuy condition in a given direction. These methods can be applied to discrete source trajectories.
[0014] Prior art methods do not exploit the theoretical local reconstruction conditions of tomography. The resulting adjustment is identical for any object to be reconstructed, regardless of local sampling differences. Indeed, in the case of a regularization that would be adapted locally, the number of parameters or terms to be adjusted would be too large. For example, a non-local regularization such as Total Directional Variation has three hyper-parameters to adjust, whereas a local total directional variation applied to a volume of one million voxels would require adjusting three million hyper-parameters. This makes hyper-parameter selection impossible with current methods.
[0015] It should be noted that in this document, voxel means either a pixel or a voxel depending on whether we are considering a two-dimensional or three-dimensional reconstruction.
[0016] The invention aims to propose for the reconstruction of an object by X-ray tomography, a voxel-by-voxel regularization adapting to each system acquisition or scanner to which it is applied while being inexpensive in data.
[0017] The invention proposes for this purpose a method for reconstructing an object by X-ray tomography using X-ray sources and detectors, said object belonging to an imaged space defined by a set of voxels. The method comprises the following steps: - obtaining the positions of the X-ray sources and detectors; - pre-calculating a set of directions; - calculating a tomographic incompleteness vector map comprising for each voxel a vector along a direction ^ ∗ among said set of directions; - calculating a tomographic incompleteness map of the metrics defining for each voxel a loss of information along the direction ^ ∗; - calculating regularization terms from the tomographic incompleteness vector map and the tomographic incompleteness map of the metrics; and - solving, iteratively, a regularized tomographic reconstruction problem using the calculated regularization terms.
[0018] The solution provided by the invention thus makes it possible to exploit the information from the tomographic incompleteness to adjust the terms of a regularization, for example by adjusting the hyper-parameters of a local Total Variation type regularization. The objective is to introduce a voxel-by-voxel specific regularization which will exploit the directional information as well as the regularization strength (i.e. the weight given to the regularization) estimated by the tomographic incompleteness.
[0019] Since the incompleteness information is local and specific to each acquisition system or scanner, this leads to a regularization that adapts to each acquisition system. The advantages obtained by the present invention are: - taking into account the geometry of the scanner; - introduction of a local regularization allowing a more precise reconstruction; - possibility of reconstructing a better quality image, particularly for all acquisition geometries that lead to image reconstruction problems with missing data (limited angle problems, ill-posed inverse problems, ill-conditioned inverse problems); - possibility of removing the lock on hyper-parameter adjustment for local Total Directional Variation regularization by using an analytical criterion exploiting the geometry of the acquisition system without the slightest need to generate experimental data, simulations or reconstructions, nor human intervention or artificial intelligence; - invention applicable to known geometries as well as new geometries (stationary multi-sources, etc.); - invention applicable to all X-ray tomography systems, in the fields of medical imaging, non-destructive testing or even baggage control in airports; - invention opening new perspectives for new reconstruction methods in X-ray tomography (dynamic reconstruction with 3D Total Directional Variation + local time).
[0020] Particularly convenient preferred features of the reconstruction method according to the invention are presented below.
[0021] The calculation of the tomographic incompleteness vector map comprises for each voxel of the imaged space the sub-steps of: - selection from said sources of a source for each direction from said set of directions; and - choice for each selected source of the direction ^ ∗, the loss of information defined by the tomographic incompleteness map of the metrics corresponding to the loss of information according to the chosen direction ^ ∗ by said detector associated with said selected source.
[0022] The step of calculating the tomographic incompleteness vector map comprises a prior sub-step of selecting, from among said sources, sources projecting the voxel of the object onto their associated detectors, the source being selected from the sources projecting the object voxel onto their associated detectors.
[0023] The selection of said source is carried out by the calculation of a directional local tomographic incompleteness ^^^, ^^ ∈ ℝ ^ according to the expression ^^^, ^^ = where ^ ^ ∈ ℝ ! , is the projection of a source on the plane defined by the point ^ of the voxel and the normal direction ^.
[0024] The chosen direction is chosen based on the local directional tomographic incompleteness ^^^, ^^.
[0025] The chosen direction ^ ∗ is chosen from the calculation of a norm of several directions for which the directional local tomographic incompleteness ^^^, ^^ is calculated.
[0026] The norm is the Lp norm for p ∈ [1;+∞]. Preferably, p is equal to +∞.
[0027] The chosen direction ^ ∗ is chosen from the calculation of a weighted average of several directions for which the directional local tomographic incompleteness ^^^, ^^ is calculated.
[0028] The regularization terms are calculated for a regularization function ‖% ‖ &é()*+,-*.next: where ^u8^ 1,2,3 corresponds to the finite differences in different directions, j, k, l corresponding to the indexing for each voxel, and d to the direction of the finite difference.
[0029] The regularization function ‖%‖&é()*+,-*.is of the directional total variation type.
[0030] Other features and advantages of the invention will also appear in the description below with reference to the appended drawings, given as non-limiting examples: - figure 1 represents examples of source-detector trajectories of a scanner used for the reconstruction of an object; - Figure 2 represents a diagram of a method for reconstructing an object by X-ray tomography according to an embodiment of the invention; - Figure 3 represents an example of calculating a direction set; - Figure 4 represents an example of calculating the directional local tomographic incompleteness; and - Figure 5 represents examples of 3D maps of the tomographic incompleteness.
[0031] Figure 2 represents a method for reconstructing an object imaged by X-ray tomography according to one embodiment. The reconstruction algorithm used in the invention for reconstructing the object is a regularized reconstruction algorithm. This method firstly involves a calculation of regularization terms, the steps of which are represented, according to an exemplary embodiment, in Figure 2.
[0032] X-ray tomography is performed using a scanner. The scanner includes X-ray sources and detectors. We denote {" ; , "9, … , "_=} ∈ ℝ !×@ and {A ; , A9, … , A_=} ∈ ℝ !×@ the set of source-detector pair positions in the scanner geometry. There are N source-detector pairs. Each source has an associated detector.
[0033] The calculation of regularization terms and the reconstruction of the object are implemented by computer.
[0034] The reconstruction method comprises a step of obtaining E0 the positions of the X-ray sources and detectors and defining an imaged space Ω containing the object to be reconstructed. The object to be reconstructed is not taken into account in the method. In particular, whatever the shape of the object, this does not influence the calculation of the incompleteness and therefore the regularization. Only the dimensions of the imaged space Ω containing the object to be reconstructed are taken into account in the method.
[0035] The image space Ω is meshed into several voxels. The image space Ω can have any type of shape and dimensions. The image space Ω is for example a cube of dimensions 256DD × 256DD × 256DD.
[0036] The reconstruction method also includes a step of pre-computing a set of directions. The pre-computing step allows the directions to be discretized in space. A finite number of directions in space is thus predefined. The greater the number of pre-computed or predefined directions, the more precise the calculated incompleteness direction. However, the greater the number of pre-computed directions, the greater the number of calculations. The number of directions is chosen so as to describe the entire space well while limiting the calculations to a reasonable number.
[0037] The pre-calculation step consists, for example, of calculating a set of homogeneous unit vectors on the sphere E 9, as illustrated in Figure 3. For example, the Fibonacci lattice method [Gonzalez, 2009] which places the points on a spherical spiral uniformly is used. Iteratively, each new point is placed between the largest gap of the previous points.
[0038] The reconstruction method then comprises steps E2, E3 of calculating tomographic incompleteness maps. In particular, the method comprises a step E2 of calculating a vector map of tomographic incompleteness and a step E3 of calculating a tomographic incompleteness map of the metrics.
[0039] According to an exemplary embodiment, the calculation step E2 of the tomographic incompleteness vector map comprises, for each voxel of the imaged space or point ^ ∈ Ω, a first sub-step E20 of selection of the sources projecting the voxel of the imaged space onto their associated detectors. In other words, the sources such as the direction ^" G, ^^ intersects detector A G are selected. The set of sources projecting the voxel of the imaged space Ω onto their associated detectors is noted
[0040] The number of sources is thus reduced. For each voxel in the imaged space Ω, only the sources that can see the said voxel are selected.
[0041] The calculation step E2 of the tomographic incompleteness vector map includes for each voxel of the imaged space or point ^ ∈ Ω, a second sub-step E21 of selection of the source, called best source, best seeing the voxel.
[0042] The selection of the best source is performed by calculating a directional local tomographic incompleteness ^^^, ^^ ∈ ℝ ^ for a point ^ ∈ Ω and a pre-calculated direction ^ ∈ E 9 . Local directional tomographic incompleteness is calculated according to is the projection of a source " ^ on the plane defined by the point ^ of the voxel and the normal direction ^ : ^ ^ = " ^ ^^ ⋅ ^M^. This expression is defined by [Clackdoyle & Noo, 2020].
[0043] The above calculation of the directional local tomographic incompleteness ^^^, ^^ evaluates the minimum tangent of the angles defined by the plane ^^, ^^ and the X-ray lines, i.e. the lines passing through ^ and the source points along the trajectory. The directional local tomographic incompleteness ^^^, ^^ thus evaluates to what extent a voxel or point ^ is reconstructable in the direction ^.
[0044] The higher the value ^ ^ ^, ^ ^ is large, the further the point ^ is from satisfying the local Tuy condition and lacks tomographic information in the ^ direction. Indeed, the higher the ratio is small, the more it means that the source “approaches” the normal plane. ^ ^^, ^ ^ is 0 when the source is in the normal plane. Conversely, the more the direction ^ is aligned with the source lines seeing the point ^, the higher the ratio ^ ^ tends towards infinity. Thus, the higher the ratio is large, the more difficult it is to satisfy the local Tuy condition and the less well we can see the point ^ in the plane defined by the codirection ^ (or normal plane).
[0045] By calculating the minimum of the ratio ||^^^^^^^||, the direction that is most orthogonal to the source-point vector is sought. The best source is thus obtained.
[0046] The best source is selected here from the sources projecting the voxel of the imaged space Ω onto their associated detectors. The first sub-step of selection E20 of the sources projecting the voxel of the imaged space Ω onto their associated detectors is optional. According to an exemplary embodiment, the set of sources {" ; , "9, … , "_^} is considered. The selection of the best source is carried out among the set of sources {" ; , "9, … , "_^}.
[0047] Figure 4 shows an illustration of the calculation of directional local tomographic incompleteness ^ ^ ^, ^ ^ for the selection of the best source for the reconstruction of a point x. For the sake of simplification, only two sources s1, s2 are represented. We note the angle between the line passing through the source s1 and the point x, and the plane orthogonal to ^ . The local tomographic incompleteness directional ^^^, ^^ is equal to tan(ψ G ^. Calculating the tangent tan(ψ ;^ for source s1 gives 0.28 while the calculation of the tangent tan(ψ9^ for source s2 gives 0.65. The local directional tomographic incompleteness ^^^, ^^ is therefore 0.28 here.
[0048] According to another exemplary embodiment, the source selected in step E21 may be a source different from the best source. For example, the second best source may be selected.
[0049] The calculation step E2 of the tomographic incompleteness vector map includes for each voxel of the imaged space or point ^ ∈ Ω, a third sub-step of choice E22 for the best selected sources of a direction ^ ∗.
[0050] The chosen direction ^ ∗ can be chosen or calculated in different ways. The chosen direction ^ ∗ is chosen or calculated based on the local directional tomographic incompleteness ^^^, ^^.
[0051] The chosen direction ^ ∗ can be chosen from the calculation of a norm. The chosen direction ^ ∗ is for example chosen from the calculation of the norm of several directions (v1, …., vn), said directions (called considered directions) corresponding to the set of pre-calculated directions or to a subset of said set of pre-calculated directions. The norm is for example the Lp norm of the considered directions, with p ∈ [1 ;+∞].
[0052] A preferred example is to choose the worst direction ^ ∗. The worst direction ^ ∗ corresponds to the direction in which each voxel is least well seen by the detector associated with the best selected source.
[0053] The worst direction ^ ∗ ∈ E 9 at point ^ ∈ Ω is calculated by the expression: ^ ∗ ^^^ = argmax S∈TU{^^^, ^^}. In other words, the worst direction ^ ∗ corresponds to the Lp norm of the directions considered for which p tends towards +∞.
[0054] The worst direction ^ ∗ is thus the one which maximizes the local directional tomographic incompleteness ^ ^ ^, ^ ^ This is the direction in which, even when choosing the best source, we obtain the worst possible result, i.e. we see point ^ the least well.
[0055] In a simplified example for explanation in which the direction set includes three directions { ^1 , ^2, ^3} , the tomographic incompleteness local directional ^^^, ^^ is calculated for each of the three directions. We thus obtain three values ^^^, ^1^, ^^^, ^2^, ^^^, ^3^. The worst direction ^ ∗ corresponds to the direction among the three directions { ^1 , ^2, ^3} having the largest value of local directional tomographic incompleteness.
[0056] According to another exemplary embodiment, the chosen direction ^ ∗ can be chosen by calculating a weighted average of the directions considered. The directional local tomographic incompleteness is calculated for each of the directions considered. A weight is associated with each of the directions considered, the weight being a function of the local incompleteness of the direction considered. The weight is all the more important as the value of the local incompleteness of the direction considered is large.
[0057] Thus, in the same simplified example where a set of directions comprising three directions { ^1 , ^2, ^3} is considered, the directional local tomographic incompleteness ^^^, ^^ is calculated for each of the three directions. We thus obtain three values ^^^, ^1^, ^^^, ^2^, ^^^, ^3^. Depending on the values obtained ^ ^ ^, ^1 ^ , ^ ^ ^, ^2 ^ , ^ ^ ^, ^3 ^, weights p1, p2, p3 are respectively associated with each of the three directions { ^1 , ^2, ^3} (weight p1 is associated with ^ ^ ^, ^1 ^ , the weight p2 at ^ ^ ^, ^2 ^ , and the weight p3 ^ ^ ^, ^3 ^ ). The chosen direction is calculated by calculating the weighted average of the three directions { ^1, ^2, ^3}.
[0058] Non-limiting examples of calculating the chosen direction ^ ∗ , some of which are mentioned above, are listed below: - weighted average: ^ ∗^ ^ ^ , with N the number of pre- directions calculated; [ ^ weighted average over the =^ worst directions ∗^ ^ ∑Y\] - : ^ ^ = , avec _ ` the i ème largest incompleteness value (or smallest incompleteness value) and =^ chosen less than or equal to N the number of pre-calculated directions; - average - the best direction
[0059] However you calculate or choose the direction ^ ∗ , at each point ^ of the imaged space Ω, we obtain a vector giving ^ ∗ the chosen direction. We thus obtain the 3D vector map of the tomographic incompleteness of the directions ^ ∗ : ℝ ! → E 9 , i.e. the 3D spatial distribution of tomographic incompleteness in the imaged space Ω.
[0060] The calculation step E3 of the tomographic incompleteness map of the metrics is carried out after this choice of the direction ^ ∗ . The tomographic incompleteness of metrics is a function of the directional local tomographic incompleteness ^^^, ^^. In other words, the 3D map of the tomographic incompleteness of metrics ^0: ℝ ! → ℝ ^ is defined by ^0^^^ = ℎ^^^^, ^^^, where h is a function.
[0061] The 3D map of tomographic incompleteness of metrics ^0: ℝ ! → ℝ ^ is for example defined by: ^0 ^ ^ ^ = ^^^, ^ ∗^ ^ ^ ^ . In other words, the tomographic incompleteness of the IM metrics for a given voxel corresponds to the incompleteness function applied to the chosen direction ^ ∗ .
[0062] The tomographic incompleteness map of the metrics defines for each voxel a loss of information according to the chosen direction ^ ∗ by the detector associated with the best selected source. The tomographic incompleteness of the metrics indicates to what extent the point ^ of the imaged space Ω is not clearly seen in the chosen direction ^ ∗ The larger the value of ^0^^^, the more incomplete the point ^ is in the sense of the tomographic reconstruction conditions.
[0063] The 3D map of tomographic incompleteness of ^0 metrics can be defined by another function. For example, the 3D map of tomographic incompleteness of ^0 metrics: ℝ ! → ℝ ^ can be defined by: ^0^^^ = ^²^^, ^ ∗ ^^^^, or by
[0064] Once calculated as explained above, the 3D map of the tomographic incompleteness of the ^0 metrics can be further reduced to a segment [a,b] with a and b ∈ ℝ ^ . We thus obtain: ^0^^^ = f + h^i I i^kZl^^^ . ^0^^^.
[0065] We thus obtain two 3D maps of tomographic incompleteness directly linked to the intrinsic geometry of the scanner: - the vector map of tomographic incompleteness indicates at each point ^ ∈ Ω, a direction ^ ∗ ^^^ that there is a lack of information (by example the worst direction ^ ∗^^^ in which we have the least information on this point in terms of projection); and - the tomographic incompleteness map of the metrics quantifies this lack
[0066] Figure 5 shows examples of 3D maps of tomographic incompleteness where the imaged space Ω is a cube. For each point ^ ∈ Ω, the maps display a cone with: the chosen direction ^ ∗ ^^^ and a value corresponding to the lack of information ^0 ^ ^ ^ . For the sake of readability of the representation, the mesh of the imaged space Ω is here 3^3 (i.e. 27 cones). A finer mesh, for example 9^9 is preferred.
[0067] The variation of lack of information ^0 ^ ^ ^ is represented by means of point clouds, said points being arranged on the cones. The higher the value of ^0 ^ ^ ^ The higher the density of points, the greater the density of points.
[0068] In a) of Figure 5, a case is represented where the source-detector trajectory T is circular. The figure shows that only the points on the plane of the circle of the trajectory respect the Tuy condition. The value of the tomographic incompleteness increases as one moves away from the plane of the circle. In other words, the voxel is less and less reconstructable tomographically as one moves away from the plane of the circle.
[0069] In b) of Figure 5, a case is shown where the source-detector trajectory T is helical. The figure shows that the tomographic incompleteness is more homogeneous and very low at all points. All directions (even the worst direction) remain satisfactory for tomographic reconstruction.
[0070] The reconstruction method comprises a step E4 of calculating the regularization terms from the tomographic incompleteness vector map and the tomographic incompleteness map of the metrics. The information ^ ∗ ^^^ and ^0 ^ ^ ^ are used to adjust the regularization terms of the regularization function.
[0071] Regularization makes it possible to give a statistical a priori on the imaged object and / or take into account the acquisition geometry of the scanner.
[0072] The regularization is achieved here by an iterative approach. In this approach, a high-dimensional minimization problem is solved, iteratively, for tomographic reconstruction. The problem here is to find an approximation % m of f such that: where A is the tomographic projection matrix that allows the transition from the object domain to the projection domain. In other words, matrix A is the transformation that allows the transition from image space to projection data space. b are the acquired data (also called sinogram). ‖ p% − q ‖ is called the data attachment term. The data attachment term ensures that the solution remains close to the measurements made. ‖ f ‖ &é()*+,-*.is called a regularization term where r is a parameter allowing to adjust the importance of the regularization with respect to the data attachment term. The regularization used in the example presented is of the local total variation type. The regularization term can be noted in this case
[0073] The above equation imposes a regularity condition to control the divergences of the reconstruction algorithm due to inconsistencies in the input data. The regularity condition constitutes a priori information.
[0074] The data attachment term can be for example ‖ p% − q ‖ . This is a data attachment making the a priori that the acquisition noise is Gaussian. The data attachment can of course be different.
[0075] The regularization function ‖%‖&é()*+,-*.can be in particular the following: where ^u^ w,x,y corresponds to the finite differences in different directions (e.g. the three directions (x, y, z)) and j, k, l correspond to the indexing for each voxel. g is a function which can for example be 4^^ ∗ ^ = z. The differences are for example calculated in the following way (subtraction of two neighboring voxels): Of course, finite differences exploiting a different number of neighboring voxels are conceivable.
[0076] Other non-limiting examples of the g function different from the aforementioned standard are conceivable and listed below: ^ ^ ^5 ∗ , 1,2^;,3 + ^ ^ ^5 ∗ , 1^;,2,3 with {^ ; , ^9, …} weights chosen, for example according to a Gaussian distribution.
[0077] The terms ∗ and ^ 5,1,2,3 allow weights to be introduced into the regularization. These weights are adapted in the directions of space (index d) as well as for each voxel (indices j, k, l). ^5 ∗ , 1,2,3 gives in which direction we want to regulate more or less. gives the desired regularization strength.
[0078] Let's take the example of a voxel for which there is a difference between the terms ^ ^ ∗ , ^ ^ ∗ , ^ ^ ∗ along each of the three directions x, y, z. If the local directional tomographic incompleteness ^^^, ^^ of the voxel is low, although there is a directional effect (the terms ^ ^ ∗ , ^ ^ ∗ , ^ ^ ∗ being different), the overall regularization strength ^0 1,2,3 is set to be small since the incompleteness ^ ^ ^, ^ ^ of the voxel is weak. Conversely, for a very poorly reconstructable voxel (i.e. for which the local directional tomographic incompleteness ^^^, ^^ is strong), the overall regularization strength is fixed so as to be strong and the ∗ e term ^ 5,1,2,3allows you to put a weighting in the different directions.
[0079] The regularization carried out here is therefore not only directional but also local, that is to say specific to each voxel.
[0080] The regularization function ‖%‖&é()*+,-*.defined above is given as a non-limiting example. The regularization function may be defined differently. For example, the regularization function may be an adaptation of [Bayram, 2012] to the local regularization of the invention. The regularization function ‖%‖&é()*+,-*.may be in particular the following: with ^ = 4^^5 ∗ , 1,2,3 ^ being the direction of the regularization (for example 4^^5 ∗ , 1,2,3 ^ = and r being the stretch factor at
[0081] The reconstruction method according to the invention further comprises a step E5 of solving the minimization problem using calculated regularization terms.
[0082] The reconstruction method according to the invention is implemented by computer. The computer may comprise various computing, storage and communications units configured to interact with each other. The computer may comprise a processor, one or more storage devices, an input / output interface and a human-machine interface. The computer here refers to any computer device or system capable of implementing the reconstruction method according to the invention.
[0083] The invention proposes a criterion for automatic or semi-automatic adjustment of local regularization terms. This criterion is based on tomographic incompleteness calculated from the geometry of the scanner. The regularization according to the invention is particularly suitable in the case of a directional total variation type regularization. However, any reconstruction algorithm with a regularization term (or terms) that can adapt locally is conceivable.
[0084] The main benefit is better image reconstruction quality both visually and quantitatively.
[0085] The invention further allows both to take into account the geometry of the scanner and to build a regularization analytically. Moreover, this regularization works for all scanner geometries. Finally, the a priori on the object used in the tomographic reconstruction algorithm is adjusted locally, which makes it more accurate.
[0086] The invention allows the introduction of a voxel-by-voxel specific regularization that will exploit the directional information as well as the regularization strength. The regularization strength allows the definition of the weight of the regularization in relation to the attachment to the data that it is data-driven. For high incompleteness, a high weight is given to regularization and a lower weight to data. Conversely, for low incompleteness (meaning that the object is highly reconstructable only with the acquired data), a high weight is given to data attachment while a lower weight is given to regularization.
[0087] Bibliographie [Bayram, 2012] : Bayram and Kamasak 2012, « Directional Total Variation » , IEEE Signal Processing Letters, Vol.19, No.12 [Clackdoyle & Noo, 2020] : Clackdoyle, R., & Noo, F. (2020, January). Quantification of Tomographic Incompleteness in Cone-Beam Reconstruction.4, 63– 80. doi:10.1109 / trpms.2019.2918222. [Ding, 2021] : Ding, Q. N. (2021). Deep learning with adaptive hyper- parameters for low-dose CT image reconstruction. IEEE Transactions on Computational Imaging, 7, 648-660. [Gonzalez, 2009] : Gonzalez, A. (2009, November). Measurement of Areas on a Sphere Using Fibonacci and Latitude–Longitude Lattices. Mathematical Geosciences, 42, 49–64. doi:10.1007 / s11004-009-9257-x. [Kohavi, 1995] : Kohavi, R. (1995). A study of cross-validation and bootstrap for accuracy estimation and model selection.14, 1137-1145. [Lin & Meikle, 2010] : Lin, J., & Meikle, S. R. (2010, October). Truncated pinhole SPECT: Sufficient sampling criteria and applications.IEEE Nuclear Science Symposuim: Medical Imaging Conference. IEEE. doi:10.1109 / nssmic.2010.5874140. [Liu, et al., 2012] : Liu, B., Bennett, J., Wang, G., Man, B. D., Zeng, K., Yin, Z., . . . Yu, H. (2012, April). Completeness map evaluation demonstrated with candidate next-generation cardiac CT architectures. Medical Physics, 39, 2405–2416. doi:10.1118 / 1.3700172. [Metzler, Bowsher, & Jaszczak, 2003] : Metzler, S. D., Bowsher, J. E., & Jaszczak, R. J. (2003, October). Geometrical similarities of the Orlov and Tuy sampling criteria and a numerical algorithm for assessing sampling completeness. IEEE Transactions on Nuclear Science, 50, 1550–1555. doi:10.1109 / tns.2003.817385. [Rudin, 1992] : Rudin, L. I. (1992). Nonlinear total variation based noise removal algorithms. Physica D: nonlinear phenomena, 60.1(4), 259-268. [Shen, 2018] : Shen, C. G. (2018). Intelligent parameter tuning in optimization- based iterative CT reconstruction via deep reinforcement learning. IEEE Transactions on Medical Imaging, 37.6, 1430-1439. [Stopp, Winne, Jank, & Keeve, 2010] : Stopp, F., Winne, C., Jank, E., & Keeve, E. (2010, February). Quality evaluation of image recording strategies for limited angle tomography. Tsinghua Science and Technology, 15, 25–29. doi:10.1016 / s1007- 0214(10)70004-3. [Tuy, 1983] : Tuy, H. K. (1983, June). An Inversion Formula for Cone-Beam Reconstruction. SIAM Journal on Applied Mathematics, 43, 546–552. doi:10.1137 / 0143035. [Zhang, 2021] : Zhang, Z. C. (2021). Directional-TV algorithm for image reconstruction from limited-angular-range data. Medical image analysis, 70.
Claims
CLAIMS 1. Method for reconstructing an object by X-ray tomography using sources ^" ; , "9, … , " @ ^ and detectors ^A ; , A9, … , A @^ of X-rays, said object belonging to an imaged space (Ω) defined by a set of voxels, said method comprising the following steps: - obtaining (E0) the positions of the sources and the X-ray detectors; - pre-calculation (E1) of a set of directions; - calculation (E2) of a tomographic incompleteness vector map comprising for each voxel a vector in a direction ^ ∗ among said set of directions; - calculation (E3) of a tomographic incompleteness map of the metrics ^0^^^ defining for each voxel a loss of information in the direction ^ ∗; - calculation (E4) of regularization terms from the tomographic incompleteness vector map and the tomographic incompleteness map of the metrics; and - resolution (E5), iteratively, of a regularized tomographic reconstruction problem using the calculated regularization terms. 2.Reconstruction method according to claim 1, in which the calculation (E2) of the tomographic incompleteness vector map comprises for each voxel of the imaged space ^Ω) the sub-steps of: - selection (E21) from among said sources ^". ; , "9, … , " @ ^ of a source for each direction among said set of directions; and - choice (E22) for each selected source of the direction ^ ∗, the loss of information defined by the tomographic incompleteness map of the metrics ^0^^^ corresponding to the loss of information according to the direction chosen ^ ∗ by said detector associated with said selected source.
3. Reconstruction method according to claim 2, wherein the step of calculating (E2) the tomographic incompleteness vector map comprises a preliminary selection sub-step (E20), among said sources ^" ; , "9, … , " @^, sources projecting the object voxel onto their associated detectors, the source being selected from the sources projecting the object voxel onto their detectors ^A ; , A9, … , A @ ^ associated.
4. Reconstruction method according to claim 2 or claim 3, in which the selection (E21) of said source is carried out by calculating a directional local tomographic incompleteness ^^^, ^^ ∈ ℝ ^ according to the expression ∈ ℝ ! , is the projection of a source " ^ onto the plane defined by the point ^ of the voxel and the normal direction ^.
5. Reconstruction method according to claim 4, in which said chosen direction is chosen as a function of the directional local tomographic incompleteness ^ ^ ^, ^ ^.
6. Reconstruction method according to claim 4 or claim 5, in which said chosen direction ^ ∗ is chosen from the calculation of a norm of several directions (v1, …., vn) for which the directional local tomographic incompleteness ^ ^ ^, ^ ^ is calculated.
7. Reconstruction method according to claim 6, wherein said norm is the Lp norm for p ∈ [1;+∞], p preferably being equal to +∞.
8. Reconstruction method according to claim 4 or claim 5, wherein said chosen direction ^ ∗ is chosen from the calculation of a weighted average of several directions (v1, …., vn) for which the directional local tomographic incompleteness ^ ^ ^, ^ ^ is calculated.
9. Reconstruction method according to one of claims 1 to 8, in which the regularization terms are calculated for a regularization function ‖%‖ &é()*+,-*. next: where ^u8^ 1,2,3corresponds to the finite differences in different directions, j, k, l corresponding to the indexing for each voxel, and d to the direction of the finite difference.
10. A reconstruction method according to claim 9, wherein the regularization function ‖%‖ &é()*+,-*. is of the directional total variation type.
