EP4264553A1 - Methode de reconstruction d'une image a partir d'une acquisition comprimee de mesures - Google Patents

Methode de reconstruction d'une image a partir d'une acquisition comprimee de mesures

Info

Publication number
EP4264553A1
EP4264553A1 EP21823852.5A EP21823852A EP4264553A1 EP 4264553 A1 EP4264553 A1 EP 4264553A1 EP 21823852 A EP21823852 A EP 21823852A EP 4264553 A1 EP4264553 A1 EP 4264553A1
Authority
EP
European Patent Office
Prior art keywords
matrix
image
vector
rows
reconstructing
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
EP21823852.5A
Other languages
German (de)
English (en)
Inventor
Antoine Souloumiac
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Commissariat a lEnergie Atomique et aux Energies Alternatives CEA
Original Assignee
Commissariat a lEnergie Atomique CEA
Commissariat a lEnergie Atomique et aux Energies Alternatives CEA
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Commissariat a lEnergie Atomique CEA, Commissariat a lEnergie Atomique et aux Energies Alternatives CEA filed Critical Commissariat a lEnergie Atomique CEA
Publication of EP4264553A1 publication Critical patent/EP4264553A1/fr
Pending legal-status Critical Current

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/00—2D [Two Dimensional] image generation
    • G06T11/003—Reconstruction from projections, e.g. tomography
    • G06T11/006—Inverse problem, transformation from projection-space into object-space, e.g. transform methods, back-projection, algebraic methods

Definitions

  • the invention relates to the field of imaging, for example imaging by X-rays or by magnetic resonance or any other imaging method making it possible to reconstruct a 2D or 3D image of an object.
  • the invention applies in particular to the field of medical imaging or astronomical imaging.
  • the invention also relates to the field of compressed acquisition or "compressed sensing" in English which aims to reconstruct a signal from a limited number of acquisitions of measurements of this signal.
  • the source and the detectors are positioned on either side of the object to be imaged. and perform a measurement of the X-ray absorptions of a voxel of the object along a direction defined between the source and the detectors.
  • the source and the detectors rotate around the object in order to modify, for each acquisition, the direction, that is to say the angle of view.
  • each acquisition has a cost in terms of acquisition time by the imaging device, energy consumed by the device or even X-ray dose absorbed by the patient (which must be limited to limit the side effects harmful effects of ionizing radiation).
  • the invention proposes a new compressed acquisition method based on the search for linear combinations of maximum kurtosis to determine a parsimonious solution to the aforementioned underdetermined linear system.
  • the compressed acquisition method according to the invention makes it possible to reduce the number of acquisitions of an imaging device in order to reconstruct an image without loss of quality.
  • the subject of the invention is a method for reconstructing an image I of an object, the image I having a number d of pixels, the method comprising the steps of:
  • said linear combination is obtained by searching, for at least one row of the second matrix , a linear combination of a pair of components formed by the vector x M c and this line, having a maximum kurtosis.
  • said linear combination is obtained by searching for a Givens rotation matrix, defined by a rotation angle ⁇ , for which the product of said pair of components and of this Givens rotation matrix exhibits a maximum kurtosis for one of its two components or exhibits a maximum difference between the kurtosis of the two components.
  • the search for a linear combination of maximum kurtosis is iterated for several rows of the second matrix
  • the dictionary matrix D consists of a base of elementary signals, for example a base of wavelet functions, or is generated by automatic learning from learning images.
  • the vector m is acquired by means of an X-ray detector and being made up of the concatenation of several measurement vectors, each measurement vector comprising several measurements carried out by the detector for a different angle of view of the object.
  • the projection matrix P contains projection vectors respectively defining directions associated with an angle of view.
  • the image I is reconstructed in two dimensions or in three dimensions.
  • the invention also relates to a system for reconstructing an image I of an object, the image I having a number d of pixels, the system comprising an imaging device for acquiring, in a vector m , a number p, strictly less than d, of object imaging measurements and a processing unit configured for:
  • the processing unit is configured to execute the steps of the image reconstruction method according to the invention.
  • the imaging device is an X-ray detector or a magnetic resonance imaging unit.
  • the system comprises a display device for displaying the reconstructed image I.
  • Figure 1 shows a diagram illustrating the principle of an X-ray imager
  • FIG. 2 represents a flowchart detailing the steps for implementing an image reconstruction method according to the invention
  • FIG. 3 represents a flowchart detailing the steps for implementing the search for a sparse solution of an underdetermined linear system according to an embodiment of the invention.
  • Figure 1 schematically illustrates the operating principle of an X-ray imager or X-ray scanner comprising an SRX X-ray source, for example an X-ray tube and a DRX detection device consisting of several elementary detectors for example arranged linearly or in the form of a matrix.
  • the object to be imaged is placed between the SRX source and the DRX detection device and is irradiated by RX X-rays emitted by the source towards the detectors in a direction predetermined by an angle ⁇ .
  • Each pixel of the object is characterized by an X-ray absorption value and it is sought to reconstruct the image, denoted l(x,y), of these values.
  • the measurement made by a detector corresponds to the combination of these absorption values along the ray connecting the SRX source to the detector.
  • the elementary detectors of the DRX device therefore measure the projection of an image l(x, y) on these detectors in the direction defined by the angle Q.
  • the assembly consisting of the SRX source and the DRX detectors pivots around the object in order to carry out several successive acquisitions at different angles Q.
  • the object to be imaged corresponds to a thin anatomical slice or section.
  • the SRX source and the DRX detectors are constrained to be positioned in the plane of the anatomical slice and rotate around the slice.
  • the two-dimensional image of the anatomical slice is reconstructed from the different projections measured by the detectors for different angles Q.
  • the same principle of acquisition can be achieved in 3D tomography by generalizing this principle to a body volume.
  • the SRX source and detectors XRDs pivot in this case around the volume to be imaged according to optimized trajectories ensuring the most complete coverage possible of the different angles of view.
  • the directions of projection are, in the case of 3D imaging, defined by two angles which determine the point of view of the imaging device with respect to the volume to be imaged.
  • One objective of the invention is to provide a method making it possible to reconstruct the image I (in 2D or in 3D) from a reduced number of acquisitions.
  • Figure 2 details the steps for implementing the image reconstruction method according to one embodiment of the invention.
  • the first step 201 consists in acquiring several imaging measurements by means of an imaging device.
  • each measurement consists of a vector of N d values corresponding to the number of elementary detectors of the ray detection device X XRD.
  • N a measurements are carried out for various different viewing angles ⁇ as described above.
  • One objective of the invention is to reconstruct an image I of the object to be imaged, the image having a number d of pixels strictly greater than the dimension p of the vector m.
  • any other imaging device capable of producing a vector of imaging measurements m of an object, according to different projections can be used as a replacement for the X-ray scanner.
  • a magnetic resonance imaging (MRI) unit can be used to perform the measurement step 201.
  • the projection angle variation described above for the X-ray scanner is replaced by a variation of a magnetic field created by currents passing through coils called field gradient coils. This variable field is added to the permanent field, fixed in time and useful space, created by a superconducting magnet.
  • the MRI signal is collected by radio frequency antennas placed close to the patient inside the superconducting magnet.
  • the vector m corresponds to different measurements taken by the antennas of the MRI unit, these measurements being associated with a given position in "k-space" space.
  • the vector m corresponds to the result of the projection of the desired image I on the projection base defined by the matrix P.
  • the image I can be expressed as a sparse linear combination of reference images taken from a dictionary D defining an image decomposition base:
  • l D.x with D a matrix of dimensions d by N (N>d) and x a sparse vector of size N, that is to say comprising a large proportion of zero values.
  • the matrix D contains the components of a signal decomposition base, for example a decomposition base into wavelets or into ridgelet functions or into curvelet functions or any other signal base which makes it possible to decompose the image I .
  • a signal decomposition base for example a decomposition base into wavelets or into ridgelet functions or into curvelet functions or any other signal base which makes it possible to decompose the image I .
  • the matrix A T PD (A T ap rows and N columns) is determined from the projection matrix P and the dictionary D, these two matrices P and D being determined as a function of the intended application.
  • FIG. 3 schematizes the detail of step 203 making it possible to determine the solution x.
  • the search for the most parsimonious possible solution x is equivalent to the search for a linear combination of x MC and the rows of A which is the most parsimonious possible.
  • step 203 of the method according to the invention is carried out by means of the sub-steps described in FIG. 3:
  • Step 302 determination of the matrix A ⁇ whose rows form an orthonormal basis of the complementary subspace orthogonal to the rows of the matrix A T .
  • Step 303 search for a linear combination of the vector x MC and the rows of the matrix which exhibits maximum kurtosis.
  • the step 302 does not depend on the measurement vector m and can be performed initially once and for all and not for each image to be reconstructed.
  • the vector x MC is normalized so as to obtain a matrix X with orthonormal rows.
  • Step 303 is performed by determining different linear combinations of the first row of the matrix X with the other rows.
  • step 303 is performed by determining a square rotation matrix U, of dimensions (N-p+1)x(N-p+1) and by multiplying X on the left by this matrix U, or U x X.
  • the matrix X has dimension N-p+1 by N.
  • the orthogonal matrix U such that the first row of U x X has maximum kurtosis cannot be calculated in a single step; it is constructed iteratively as a product (the product matrices of orthogonal matrices are orthogonal) of particular elementary orthogonal matrices called Givens rotations (which correspond to plane rotations) denoted G(i,j,0) and j ⁇ N — p + 1 and 0 ⁇ ⁇ ⁇ ⁇ .
  • X(i,:) denotes the line of index i of the matrix X.
  • the term “:” is used to designate all the values of the index.
  • the step 303 of the method thus consists in optimizing the angle of rotation 0 for several pairs (l,j) of rows of the matrix X in order to maximize the kurtosis of the first row.
  • Givens rotations thus decomposes the global optimization problem into a series of 2-dimensional optimization problems that are easier to perform.
  • Two variant embodiments can be envisaged for determining the angle of rotation 0 which makes it possible to obtain a maximum empirical kurtosis for the first row of the matrix.
  • a first embodiment variant consists in finding the angle of rotation
  • the angle of rotation 0 which maximizes the difference of the empirical kurtosis K u -K v is then sought.
  • u corresponds to the line of highest kurtosis
  • v corresponds to the line of lowest kurtosis.
  • a second embodiment variant consists in directly seeking the angle of rotation ⁇ which makes it possible to maximize the empirical kurtosis K u of the line u which is equal to:
  • the expression (2) corresponds to a linear and quadratic criterion under quadratic constraint.
  • the angle 0 which maximizes this expression is obtained by numerical resolution.
  • One of the two variants described above is therefore applied to each pair of rows (1,j) of the matrix X according to a predetermined sequence.
  • Different line scanning structures can be envisaged for this purpose.
  • an iteration loop is performed by varying j from 2 to p+1 until convergence.
  • the iterative procedure is organized as follows.
  • the invention can be implemented by means of an imaging device of the type described in FIG. 1 coupled to a processing unit capable of executing the steps of the image reconstruction method, according to the invention, from the measurements m made by the imaging device.
  • the invention can be implemented as a computer program comprising instructions for its execution.
  • the computer program may be recorded on a processor-readable recording medium.
  • references to a computer program which, when executed, performs any of the functions previously described, is not limited to an application program running on a single host computer. Rather, the terms computer program and software are used herein in a general sense to refer to any type of computer code (e.g., application software, firmware, microcode, or other form of computer instruction) that can be used to program one or more processors to implement aspects of the techniques described herein.
  • the IT means or resources can in particular be distributed (“Cloud computing”), possibly using peer-to-peer technologies.
  • the software code may be executed on any suitable processor (e.g., a microprocessor) or processor core or set of processors, whether provided in a single computing device or distributed among multiple computing devices (e.g.
  • the executable code of each program allowing the programmable device to implement the processes according to the invention can be stored, for example, in the hard disk or in ROM.
  • the program or programs can be loaded into one of the storage means of the device before being executed.
  • the central unit can control and direct the execution of the instructions or portions of software code of the program or programs according to the invention, instructions which are stored in the hard disk or in the ROM or else in the other aforementioned storage elements.

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Algebra (AREA)
  • Image Processing (AREA)
  • Image Analysis (AREA)
  • Compression Of Band Width Or Redundancy In Fax (AREA)
  • Apparatus For Radiation Diagnosis (AREA)
  • Length Measuring Devices By Optical Means (AREA)

Abstract

Méthode de reconstruction d'une image I d'un objet comprenant les étapes de : - Acquérir (201), dans un vecteur m, un nombre p de mesures d'imagerie de l'objet, - Déterminer (202) une première matrice AT égale au produit d'une matrice de projection P définissant un ensemble de vecteurs de projection et d'une matrice dictionnaire D constituée des composantes d'une base de décomposition de l'image, - Rechercher (203) une solution parcimonieuse x, de dimension N strictement supérieure à d, au système linéaire sous-déterminé ATx=m - Reconstruire (204) l'image I par le produit de la matrice dictionnaire D et du vecteur x, solution dudit système.

Description

DESCRIPTION
Titre de l’invention : Méthode de reconstruction d’une image à partir d’une acquisition comprimée de mesures
[0001] L’invention concerne le domaine de l’imagerie, par exemple l’imagerie par rayons X ou par résonance magnétique ou toute autre méthode d’imagerie permettant de reconstruire une image en 2D ou 3D d’un objet. L’invention s’applique notamment au domaine de l’imagerie médicale ou de l’imagerie astronomique.
[0002] L’invention concerne également le domaine de l’acquisition comprimée ou « compressed sensing >> en anglais qui vise à reconstruire un signal à partir d’un nombre limité d’acquisitions de mesures de ce signal.
[0003] Pour pouvoir reconstruire une image d’un objet donné avec une résolution précise, cela nécessite de réaliser un grand nombre d’acquisitions au moyen d’un dispositif d’imagerie, selon différents angles de prise de vue.
[0004] Par exemple, pour le cas d’imageurs par rayons X constitués d’une source de rayons X et d’un ensemble de détecteurs, la source et les détecteurs sont positionnés de part et d’autre de l’objet à imager et réalisent une mesure des absorptions de rayons X d’un voxel de l’objet selon une direction définie entre la source et les détecteurs. Pour reconstruire une image de l’objet, la source et les détecteurs tournent autour de l’objet afin de modifier, pour chaque acquisition, la direction, c’est-à-dire l’angle de prise de vue.
[0005] Plus le nombre d’acquisitions selon différents angles de prise de vue est important, meilleure sera la résolution de l’image reconstruite. Cependant, chaque acquisition présente un coût en termes de temps d’acquisition par le dispositif d’imagerie, d’énergie consommée par le dispositif ou encore de dose de rayons X absorbée par le patient (qu’il faut limiter pour limiter les effets secondaires nocifs des rayonnements ionisants).
[0006] Il existe donc un besoin pour réduire le nombre d’acquisitions tout en conservant une qualité d’image satisfaisante. [0007] Les techniques connues d’acquisition comprimée (« compressed sensing >>) sont basées sur la recherche d’une solution parcimonieuse à un système linéaire sous-déterminé.
[0008] Une difficulté algorithmique réside dans le calcul de cette solution parcimonieuse. Différentes solutions de l’état de l’art ont été proposées pour résoudre ce problème : des algorithmes de type « matching pursuit >> ou « orthogonal matching pursuit », des techniques de seuillage ou encore faisant intervenir la norme L1 Ces solutions sont par exemple décrites dans les références [1 ], [2], [3] et [4],
[0009] Toutes ces méthodes présentent des inconvénients liés à des problèmes de convergence, d’exactitude de la solution ou de complexité de mise en oeuvre.
[0010] L’invention propose une nouvelle méthode d’acquisition comprimée basée sur la recherche de combinaisons linéaires de kurtosis maximum pour déterminer une solution parcimonieuse au système linéaire sous-déterminé précité.
[001 1 ] La méthode d’acquisition comprimée selon l’invention permet de diminuer le nombre d’acquisitions d’un dispositif d’imagerie pour reconstruire une image sans perte de qualité.
[0012] L’invention a pour objet une méthode de reconstruction d’une image I d’un objet, l’image I ayant un nombre d de pixels, la méthode comprenant les étapes de :
- Acquérir, dans un vecteur m, un nombre p, strictement inférieur à d, de mesures d’imagerie de l’objet, au moyen d’un dispositif d’imagerie,
- Déterminer une première matrice AT égale au produit d’une matrice de projection P définissant un ensemble de vecteurs de projection et d’une matrice dictionnaire D constituée des composantes d’une base de décomposition de l’image,
- Rechercher une solution parcimonieuse x, de dimension N strictement supérieure à d, au système linéaire sous-déterminé ATx=m en : i. Déterminant une solution xMc de norme quadratique minimale dudit système linéaire et normaliser le vecteur xMc ii. Déterminant une seconde matrice dont les lignes forment une base orthonormée du sous-espace complémentaire orthogonal aux lignes de la première matrice AT, iii. Recherchant une combinaison linéaire x du vecteur xMc et des lignes de la seconde matrice ayant un kurtosis maximal.
- Reconstruire l’image I par le produit de la matrice dictionnaire D et du vecteur x, solution dudit système.
[0013] Selon un aspect particulier de l’invention, ladite combinaison linéaire est obtenue en recherchant, pour au moins une ligne de la seconde matrice , une combinaison linéaire d’un couple de composantes formé du vecteur xMc et de cette ligne, ayant un kurtosis maximal.
[0014] Selon un aspect particulier de l’invention, ladite combinaison linéaire est obtenue en recherchant une matrice de rotation de Givens, définie par un angle de rotation θ , pour laquelle le produit dudit couple de composantes et de cette matrice de rotation de Givens présente un kurtosis maximal pour l’une de ses deux composantes ou présente une différence maximale entre les kurtosis des deux composantes.
[0015] Selon un aspect particulier de l’invention, la recherche d’une combinaison linéaire de kurtosis maximal est itérée pour plusieurs lignes de la seconde matrice
[0016] Selon un aspect particulier de l’invention, la matrice dictionnaire D est constituée d’une base de signaux élémentaires, par exemple une base de fonctions ondelettes, ou est générée par apprentissage automatique à partir d’images d’apprentissage.
[0017] Selon un aspect particulier de l’invention, le vecteur m est acquis au moyen d’un détecteur à rayons X et étant constitué de la concaténation de plusieurs vecteurs de mesures, chaque vecteur de mesure comprenant plusieurs mesures réalisées par le détecteur pour un angle de prise de vue de l’objet différent.
[0018] Selon un aspect particulier de l’invention, la matrice de projection P contient des vecteurs de projection définissant respectivement des directions associées à un angle de prise de vue.
[0019] Selon un aspect particulier de l’invention, l’image I est reconstruite en deux dimensions ou en trois dimensions. [0020] L’invention a aussi pour objet un système de reconstruction d’une image I d’un objet, l’image I ayant un nombre d de pixels, le système comprenant un dispositif d’imagerie pour acquérir, dans un vecteur m, un nombre p, strictement inférieur à d, de mesures d’imagerie de l’objet et une unité de traitement configurée pour :
- Déterminer une première matrice AT égale au produit d’une matrice de projection P définissant un ensemble de vecteurs de projection et d’une matrice dictionnaire D constituée des composantes d’une base de décomposition de l’image,
- Rechercher une solution parcimonieuse x, de dimension N strictement supérieure à d, au système linéaire sous-déterminé ATx=m en : i. Déterminant une solution xMc de norme quadratique minimale dudit système linéaire et normaliser le vecteur xMc ii. Déterminant une seconde matrice dont les lignes forment une base orthonormée du sous-espace complémentaire orthogonal aux lignes de la première matrice AT, iii. Recherchant une combinaison linéaire x du vecteur xMc et des lignes de la seconde matrice A^ ayant un kurtosis maximal.
- Reconstruire l’image I par le produit de la matrice dictionnaire D et du vecteur x, solution dudit système.
[0021] Selon un aspect particulier de l’invention, l’unité de traitement est configurée pour exécuter les étapes de la méthode de reconstruction d’image selon l’invention.
[0022] Selon un aspect particulier de l’invention, le dispositif d’imagerie est un détecteur à rayons X ou une unité d’imagerie par résonance magnétique.
[0023] Selon un aspect particulier de l’invention, le système comprend un dispositif d’affichage pour afficher l’image I reconstruite.
[0024] D’autres caractéristiques et avantages de la présente invention apparaîtront mieux à la lecture de la description qui suit en relation aux dessins annexés suivants. [0025] La figure 1 représente un schéma illustrant le principe d’un imageur à rayons X,
[0026] la figure 2 représente un organigramme détaillant les étapes de mise en oeuvre d’une méthode de reconstruction d’image selon l’invention,
[0027] la figure 3 représente un organigramme détaillant les étapes de mise en oeuvre de la recherche d’une solution parcimonieuse d’un système linéaire sous- déterminé selon un mode de réalisation de l’invention.
[0028] La figure 1 illustre schématiquement le principe de fonctionnement d’un imageur à rayons X ou scanner à rayons X comportant une source de rayons X SRX, par exemple un tube à rayons X et un dispositif de détection DRX constitué de plusieurs détecteurs élémentaires par exemple agencés linéairement ou sous forme d’une matrice. L’objet à imager est placé entre la source SRX et le dispositif de détection DRX et est irradié par des rayons X RX émis par la source vers les détecteurs selon une direction prédéterminée par un angle 9. Chaque pixel de l’objet est caractérisé par une valeur d’absorption de rayons X et on cherche à reconstruire l’image, notée l(x,y), de ces valeurs. La mesure réalisée par un détecteur correspond à la combinaison de ces valeurs d’absorption le long du rayon reliant la source SRX au détecteur.
[0029] Les détecteurs élémentaires du dispositif DRX mesurent donc la projection d’une image l(x, y) sur ces détecteurs selon la direction définie par l’angle Q.
[0030] Pour imager complètement l’objet, l’ensemble constitué de la source SRX et des détecteurs DRX pivote autour de l’objet afin de réaliser plusieurs acquisitions successives selon différents angle Q.
[0031] Dans le domaine de l’imagerie médicale, l’objet à imager correspond à une tranche anatomique de faible épaisseur ou coupe. La source SRX et les détecteurs DRX sont contraints à être positionnés dans le plan de la tranche anatomique et tournent autour de la tranche.
[0032] L’image en deux dimensions de la tranche anatomique est reconstruite à partir des différentes projections mesurées par les détecteurs pour différents angles Q.
[0033] Le même principe d’acquisition peut être réalisé en tomographie 3D en généralisant ce principe à un volume corporel. La source SRX et les détecteurs DRX pivotent dans ce cas autour du volume à imager selon des trajectoires optimisées assurant une couverture la plus complète possible des différents angles de vue. Les directions de projection sont, dans le cas d’une imagerie 3D, définies par deux angles qui déterminent le point de vue du dispositif d’imagerie par rapport au volume à imager.
[0034] Un objectif de l’invention est de fournir une méthode permettant de reconstruire l’image I (en 2D ou en 3D) à partir d’un nombre réduit d’acquisitions.
[0035] La figure 2 détaille les étapes de mise en oeuvre de la méthode de reconstruction d’image selon un mode de réalisation de l’invention.
[0036] La première étape 201 consiste à acquérir plusieurs mesures d’imagerie au moyen d’un dispositif d’imagerie.
[0037] Dans le cas où le dispositif d’imagerie utilisé est un scanner à rayons X du type décrit à la figure 1 , chaque mesure est constituée par un vecteur de Nd valeurs correspondant au nombre de détecteurs élémentaires du dispositif de détection de rayons X DRX. On réalise Na mesures pour différents angles 9 de prise de vue différents comme décrit ci-dessus.
[0038] Les mesures acquises sont numérisées et concaténées dans un vecteur m de taille p=Nd x Na. Un objectif de l’invention est de reconstruire une image I de l’objet à imager, l’image ayant un nombre d de pixels strictement supérieur à la dimension p du vecteur m.
[0039] Sans sortir du cadre de l’invention, tout autre dispositif d’imagerie apte à produire un vecteur de mesures m d’imagerie d’un objet, selon différentes projections, peut être utilisé en remplacement du scanner à rayons X. En particulier, une unité d’imagerie par résonance magnétique (IRM) peut être utilisée pour réaliser l’étape de mesures 201. Dans ce cas, la variation de l’angle de projection décrit ci-dessus pour le scanner à rayons X est remplacée par une variation d’un champ magnétique créé par des courants passant dans des bobines appelées bobines de gradient de champ. Ce champ variable s’ajoute au champ permanent, fixe dans le temps et l’espace utile, créé par un aimant supraconducteur.
[0040] Ces variations de champs sont modélisées comme une trajectoire dans un espace abstrait appelé « k-space >> en anglais. Une séquence IRM consiste en une certaine trajectoire pour échantillonner cet espace, plus on a d’échantillons meilleure est l’image. La durée nécessaire pour acquérir une image est très variable : de l’ordre d’une seconde pour l’IRM fonctionnelle à plusieurs dizaines de minutes pour une image en T2 mais dans tous les cas il est utile de raccourcir cette durée sans dégrader la qualité image.
[0041] Pour chaque position donnée dans l’espace « k-space », le signal IRM est recueilli par des antennes radiofréquences disposées à proximité du patient à l’intérieur de l’aimant supraconducteur.
[0042] Ainsi, pour une application d’imageries par résonance magnétique, le vecteur m correspond à différentes mesures réalisées par les antennes de l’unité IRM, ces mesures étant associées à une position donnée dans l’espace « k-space >>.
[0043] En revenant à l’exemple initial d’un scanner à rayons X, si on note P la matrice de projection qui contient les vecteurs associés aux directions données par les différents angles 9 de prise de vue, on peut écrire :
[0044] m= P.l,
[0045] En effet, le vecteur m correspond au résultat de la projection de l’image recherchée I sur la base de projection définie par la matrice P.
[0046] Par ailleurs, l’image I peut s’exprimer comme une combinaison linéaire parcimonieuse d’images de référence prises dans un dictionnaire D définissant une base de décomposition d’images :
[0047] l= D.x avec D une matrice de dimensions d par N (N>d) et x un vecteur de taille N parcimonieux, c’est-à-dire comportant une grande proportion de valeurs nulles.
[0048] La matrice D contient les composantes d’une base de décomposition de signaux, par exemple une base de décomposition en ondelettes ou en fonctions de ridgelet ou en fonctions de curvelet ou toute autre base de signaux qui permet de décomposer l’image I.
[0049] Dans le cas où l’image I à acquérir appartient à un domaine d’applications particulier, la matrice dictionnaire D peut être obtenue à l’issue d’un apprentissage automatique de dictionnaire à partir d’images d’apprentissage, par exemple via la méthode décrite dans la référence [5]. [0050] A partir des relations précédentes on peut donc écrire m= P.D.x= AT.x.
[0051] Dans une étape 202, on détermine la matrice AT = P.D (AT a p lignes et N colonnes), à partir de la matrice de projection P et du dictionnaire D, ces deux matrices P et D étant déterminées en fonction de l’application visée.
[0052] Dans une étape 203, on recherche ensuite une solution x la plus parcimonieuse possible au système linéaire sous-déterminée AT.x = m.
[0053] Enfin, on reconstruit (étape 204) l’image l=D.x à partir de la solution x déterminée à l’étape 203.
[0054] La figure 3 schématise le détail de l’étape 203 permettant de déterminer la solution x.
[0055] Le calcul d’une solution x aussi parcimonieuse que possible du système linéaire sous-déterminé précédent peut être ramené à la maximisation d’un kurtosis de la manière suivante. On peut toujours considérer que la solution parcimonieuse x est égale à la somme de :
- une solution de norme quadratique minimale xMC = A ATA')~1m (cette solution est une combinaison linéaire des lignes de la matrice AT et n’est, en général, pas parcimonieuse) et
- d’un vecteur orthogonal aux lignes de la matrice AT que l’on note x±. Si on note A une matrice à N-p lignes, et N colonnes et dont les lignes forment une base orthonormée du sous-espace complémentaire orthogonal aux p lignes de AT, alors est une combinaison linéaire des colonnes de 4±.
[0056] On peut vérifier que x = xMC + x± est bien une solution du système ATx = m quel que soit x± puisque
[0057] Ainsi, la recherche d’une solution x la plus parcimonieuse possible est équivalente à la recherche une combinaison linéaire de xMC et des lignes de A qui soit la plus parcimonieuse possible.
[0058] Cela revient à chercher la combinaison linéaire des lignes de kurtosis maximum de la matrice X égale à [0059] En effet, un kurtosis très élevé caractérise une variable aléatoire très proche de 0 avec une grande probabilité et non nulle avec une très faible probabilité, autrement dit une variable aléatoire parcimonieuse.
[0060] Autrement dit, l’étape 203 de la méthode selon l’invention est réalisée au moyen des sous-étapes décrites à la figure 3 :
- Etape 301 : détermination de la solution de norme quadratique minimale XMC = A(AT A')-1m
- Etape 302, détermination de la matrice A^ dont les lignes forment une base orthonormée du sous-espace complémentaire orthogonal aux lignes de la matrice AT.
- Etape 303 : recherche d’une combinaison linéaire du vecteur xMC et des lignes de la matrice qui présente un kurtosis maximum.
[0061] L’étape 302 ne dépend pas du vecteur de mesures m et peut être réalisée initialement une fois pour toutes et non pour chaque image à reconstruire.
[0062] Le vecteur xMC est normalisé de sorte à obtenir une matrice X aux lignes orthonormées.
[0063] On définit le kurtosis empirique d’une ligne (x1,x2, ...,xn, ...xN) de la matrice X par la formule Cette formule se simplifie du fait des lignes orthonormées car on a
[0064] La formule simplifiée du kurtosis empirique devient
[0065] La recherche d’une valeur maximale du kurtosis empirique Kx est équivalente à la recherche d’une valeur maximale de la quantité simplifiée
[0066] L’étape 303 est réalisée en déterminant différentes combinaisons linéaires de la première ligne de la matrice X avec les autres lignes.
[0067] Autrement dit, l’étape 303 est réalisée en déterminant une matrice carrée U de rotation, de dimensions (N-p+1)x(N-p+1) et en multipliant X à gauche par cette matrice U, soit U x X. La matrice X est de dimension N-p+1 par N. On remarque qu’il est avantageux de ne multiplier à gauche que par des matrices orthogonales (matrices U qui vérifient UTU = UUT = IN-P+1 = la matrice identité de dimension (/V - p + 1) x (/V - p + 1)) car elles conservent la normalisation des lignes c’est- à-dire la validité de la formule simplifiée du kurtosis empirique.
[0068] La matrice orthogonale U telle que la première ligne de U x X soit de kurtosis maximal ne peut être calculée en une seule étape ; elle est construite itérativement comme produit (les matrices produits de matrices orthogonales sont orthogonales) de matrices orthogonales élémentaires particulières appelées rotations de Givens (qui correspondent à des rotations planes) notées G(i,j, 0) et j ≤ N — p + 1 et 0 ≤ θ ≤ π.
[0069] On remarque que le produit G(i,j, 0) x X = X’ remplace les lignes X(i, : ) et X(j, : ) de X par les lignes et X'(j, : ) égales aux combinaisons linéaires suivantes :
[0070] X'(i, : ) = cos θ X(i, : ) — sin θ X(j, : ) et X'(j, : ) = sin θ X(i, : ) + cos θ X(j, : ) .
[0071 ]X(i, : ) désigne la ligne d’indice i de la matrice X. Le terme « : >> est utilisé pour désigner toutes les valeurs de l’indice.
[0072] L’étape 303 de la méthode consiste ainsi à optimiser l’angle de rotation 0 pour plusieurs couples (l,j) de lignes de la matrice X afin de maximiser le kurtosis de la première ligne.
[0073] L’utilisation des rotations de Givens décompose ainsi le problème global d’optimisation en une série de problèmes d’optimisation de dimension 2 plus faciles à réaliser.
[0074] Deux variantes de réalisation sont envisageables pour déterminer l’angle de rotation 0 qui permet d’obtenir un kurtosis empirique maximum pour la première ligne de la matrice.
[0075] Une première variante de réalisation consiste à rechercher l’angle de rotation
0 qui maximise l’écart des kurtosis entre les deux lignes du couple. En effet, en maximisant cet écart, cela permet d’obtenir un kurtosis approximativement maximal pour la ligne de kurtosis le plus élevé.
[0076] Notons X(1, : ) = (x1, ...xn, ... xN) et X(j, : ) = (y1, ...yn, ...yN') deux lignes de la matrice X dont on cherche une combinaison linéaire de kurtosis maximum.
[0077] On définit les deux combinaisons linéaires u et v par :
[0079] En supposant les lignes X(1, : ) et X(j, : ) orthonormées, on obtient des vecteurs u et v orthonormés, leurs kurtosis empiriques respectifs sont définis par
[0080] On recherche alors l’angle de rotation 0 qui maximise la différence des kurtosis empiriques Ku - Kv. Dans la suite on considère par convention que « u >> correspond à la ligne de kurtosis le plus élevé et « v >> correspond à la ligne de kurtosis le plus faible. Ainsi, on cherche à maximiser la différence des kurtosis empiriques Ku - Kv.
[0082] A partir de l’expression (1 ), on détermine les valeurs de cos 2θ et sin 2θ qui maximisent cette expression puis on en déduit les valeurs de cos θ et sin θ et enfin la valeur de l’angle de rotation θ qui permet d’obtenir un kurtosis empirique maximal pour la ligne u résultant de la rotation des lignes X(1, : ) et X(j, : ).
[0085] Un avantage de cette première variante est sa faible complexité de résolution numérique. [0086] Une seconde variante de réalisation consiste à rechercher directement l’angle de rotation θ qui permet de maximiser le kurtosis empirique Ku de la ligne u qui vaut :
[0087]
[0088] L’expression (2) correspond à un critère linéaire et quadratique sous contrainte quadratique. L’angle 0 qui maximise cette expression est obtenu par résolution numérique.
[0089] L’une des deux variantes décrites ci-dessus est donc appliquée à chaque couple de lignes (1,j) de la matrice X selon une séquence prédéterminée. Différentes structures de balayage des lignes peuvent être envisagées à cet effet.
[0090] Par exemple, on réalise une boucle d’itérations en faisant varier j de 2 à p+1 jusqu’à convergence.
[0091] A chaque itération (numérotée par l’indice « k »), les produits des matrices de rotation de Givens obtenues sont accumulés dans une matrice orthogonale Uk.
[0092] La procédure itérative est organisée de la manière suivante.
[0093] On initialise U0 = IN-P+1 et X0 = X , avec IN-P+1 la matrice identité de dimension (N - p + 1) x (N - p + 1).
[0094] Puis à chaque itération, pour chaque matrice de rotation déterminée, on calcule
[0095] L’orthogonalité des rotations de Givens (c’est-à-dire la propriété G (i,j, θ)T x G(i,j, θ) = IN-P+1) garantit la stabilité du produit Uk x Xk au cours des itérations, en effet
[0096] On a donc Xk = Uk x X tout au long de la procédure d’optimisation. Une fois la convergence atteinte après un certain nombre K d’itérations, la première ligne de est aussi parcimonieuse que possible, le vecteur x, parcimonieux et solution du système
[0097] L’invention peut être implémentée au moyen d’un dispositif d’imagerie du type de celui décrit à la figure 1 couplé à une unité de traitement apte à exécuter les étapes de la méthode de reconstruction d’image, selon l’invention, à partir des mesures m réalisées par le dispositif d’imagerie.
[0098] Ainsi, l’invention peut être mise en oeuvre en tant que programme d’ordinateur comportant des instructions pour son exécution. Le programme d’ordinateur peut être enregistré sur un support d’enregistrement lisible par un processeur.
[0099] La référence à un programme d'ordinateur qui, lorsqu'il est exécuté, effectue l'une quelconque des fonctions décrites précédemment, ne se limite pas à un programme d'application s'exécutant sur un ordinateur hôte unique. Au contraire, les termes programme d'ordinateur et logiciel sont utilisés ici dans un sens général pour faire référence à tout type de code informatique (par exemple, un logiciel d'application, un micro logiciel, un microcode, ou toute autre forme d'instruction d'ordinateur) qui peut être utilisé pour programmer un ou plusieurs processeurs pour mettre en oeuvre des aspects des techniques décrites ici. Les moyens ou ressources informatiques peuvent notamment être distribués ("Cloud computing"), éventuellement selon des technologies de pair-à-pair. Le code logiciel peut être exécuté sur n'importe quel processeur approprié (par exemple, un microprocesseur) ou cœur de processeur ou un ensemble de processeurs, qu'ils soient prévus dans un dispositif de calcul unique ou répartis entre plusieurs dispositifs de calcul (par exemple tels qu’éventuellement accessibles dans l’environnement du dispositif). Le code exécutable de chaque programme permettant au dispositif programmable de mettre en œuvre les processus selon l'invention, peut être stocké, par exemple, dans le disque dur ou en mémoire morte. De manière générale, le ou les programmes pourront être chargés dans un des moyens de stockage du dispositif avant d'être exécutés. L'unité centrale peut commander et diriger l'exécution des instructions ou portions de code logiciel du ou des programmes selon l'invention, instructions qui sont stockées dans le disque dur ou dans la mémoire morte ou bien dans les autres éléments de stockage précités. [01001 Références
[0101 ] [1 ] D. L. Donoho, "Compressed sensing," in IEEE Transactions on Information Theory, vol. 52, no. 4, pp. 1289-1306, April 2006.
[0102] [2] S. G. Mallat, S. G.; Zhang, Z. (1993). "Matching Pursuits with Time- Frequency Dictionaries". IEEE Transactions on Signal Processing. 1993 (12): 3397-3415.
[0103] [3] A. Chambolle, R. A. DeVore, N. Y. Lee, and B. J. Lucier, Nonlinear wavelet image processing: Variational problems, compression, and noise removal through wavelet shrinkage, IEEE Trans. Image Process., 7 (1998), pp. 319-335.
[0104] [4] LASSO (norme L1 ) : Tibshirani, Robert (1996). "Regression Shrinkage and Selection via the lasso". Journal of the Royal Statistical Society. Series B (methodological). Wiley. 58 (1 ): 267-88.
[0105] [5] M. Aharon, M. Elad and A. Bruckstein, "K-SVD: An algorithm for designing overcomplete dictionaries for sparse representation," in IEEE Transactions on Signal Processing, vol. 54, no. 1 1 , pp. 431 1 -4322, Nov. 2006, doi: 10.1 109/TSP.2006.881 199.

Claims

REVENDICATIONS Méthode de reconstruction d’une image I d’un objet à partir d’une acquisition comprimée de mesures, l’image I ayant un nombre d de pixels, la méthode comprenant les étapes de :
- Acquérir (201 ), dans un vecteur m, un nombre p, strictement inférieur à d, de mesures d’imagerie de l’objet, au moyen d’un dispositif d’imagerie,
- Déterminer (202) une première matrice AT égale au produit d’une matrice de projection P définissant un ensemble de vecteurs de projection et d’une matrice dictionnaire D constituée des composantes d’une base de décomposition de l’image,
- Rechercher (203) une solution parcimonieuse x, de dimension N strictement supérieure à d, au système linéaire sous-déterminé ATx=m en : i. Déterminant (301 ) une solution xMC de norme quadratique minimale dudit système linéaire et normaliser le vecteur xMc ii. Déterminant (302) une seconde matrice dont les lignes forment une base orthonormée du sous-espace complémentaire orthogonal aux lignes de la première matrice AT, iii. Recherchant (303) une combinaison linéaire x du vecteur xMc et des lignes de la seconde matrice A ayant un kurtosis maximal.
- Reconstruire (204) l’image I par le produit de la matrice dictionnaire D et du vecteur x, solution dudit système. Méthode de reconstruction d’image selon la revendication 1 dans laquelle ladite combinaison linéaire est obtenue en recherchant, pour au moins une ligne de la seconde matrice A , une combinaison linéaire d’un couple de composantes formé du vecteur xMc st de cette ligne, ayant un kurtosis maximal. Méthode de reconstruction d’image selon la revendication 2 dans laquelle ladite combinaison linéaire est obtenue en recherchant une matrice de rotation de Givens, définie par un angle de rotation 9, pour laquelle le produit dudit couple de composantes et de cette matrice de rotation de Givens présente un kurtosis maximal pour l’une de ses deux composantes ou présente une différence maximale entre les kurtosis des deux composantes. Méthode de reconstruction d’image selon la revendication 3 dans laquelle la recherche d’une combinaison linéaire de kurtosis maximal est itérée pour plusieurs lignes de la seconde matrice AT. Méthode de reconstruction d’image selon l’une quelconque des revendications précédentes dans laquelle la matrice dictionnaire D est constituée d’une base de signaux élémentaires, par exemple une base de fonctions ondelettes, ou est générée par apprentissage automatique à partir d’images d’apprentissage. Méthode de reconstruction d’image selon l’une quelconque des revendications précédentes dans laquelle le vecteur m est acquis au moyen d’un détecteur à rayons X et étant constitué de la concaténation de plusieurs vecteurs de mesures, chaque vecteur de mesure comprenant plusieurs mesures réalisées par le détecteur pour un angle de prise de vue de l’objet différent. Méthode de reconstruction d’image selon la revendication 6 dans laquelle la matrice de projection P contient des vecteurs de projection définissant respectivement des directions associées à un angle de prise de vue. Méthode de reconstruction d’image selon l’une quelconque des revendications précédentes dans laquelle l’image I est reconstruite en deux dimensions ou en trois dimensions. Système de reconstruction d’une image I d’un objet à partir d’une acquisition comprimée de mesures, l’image I ayant un nombre d de pixels, le système comprenant un dispositif d’imagerie pour acquérir, dans un vecteur m, un nombre p, strictement inférieur à d, de mesures d’imagerie de l’objet et une unité de traitement configurée pour :
Déterminer une première matrice AT égale au produit d’une matrice de projection P définissant un ensemble de vecteurs de projection et d’une matrice dictionnaire D constituée des composantes d’une base de décomposition de l’image,
- Rechercher une solution parcimonieuse x, de dimension N strictement supérieure à d, au système linéaire sous-déterminé ATx=m en : i. Déterminant une solution xMc de norme quadratique minimale dudit système linéaire et normaliser le vecteur xMc ii. Déterminant une seconde matrice ATj_ dont les lignes forment une base orthonormée du sous-espace complémentaire orthogonal aux lignes de la première matrice AT, iii. Recherchant une combinaison linéaire x du vecteur xMC et des lignes de la seconde matrice AT ayant un kurtosis maximal.
- Reconstruire l’image I par le produit de la matrice dictionnaire D et du vecteur x, solution dudit système. Système de reconstruction d’une image selon la revendication 9 dans lequel l’unité de traitement est configurée pour exécuter les étapes de la méthode de reconstruction d’image selon l’une quelconque des revendications précédentes 1 à 8. Système de reconstruction d’une image selon la revendication 9 dans lequel le dispositif d’imagerie est un détecteur à rayons X ou une unité d’imagerie par résonance magnétique. Système de reconstruction d’une image selon l’une des revendications 9 à 11 comprenant un dispositif d’affichage pour afficher l’image I reconstruite.
EP21823852.5A 2020-12-17 2021-12-01 Methode de reconstruction d'une image a partir d'une acquisition comprimee de mesures Pending EP4264553A1 (fr)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
FR2013539A FR3118251B1 (fr) 2020-12-17 2020-12-17 Méthode de reconstruction d’une image à partir d’une acquisition comprimée de mesures
PCT/EP2021/083710 WO2022128465A1 (fr) 2020-12-17 2021-12-01 Methode de reconstruction d'une image a partir d'une acquisition comprimee de mesures

Publications (1)

Publication Number Publication Date
EP4264553A1 true EP4264553A1 (fr) 2023-10-25

Family

ID=75438904

Family Applications (1)

Application Number Title Priority Date Filing Date
EP21823852.5A Pending EP4264553A1 (fr) 2020-12-17 2021-12-01 Methode de reconstruction d'une image a partir d'une acquisition comprimee de mesures

Country Status (3)

Country Link
EP (1) EP4264553A1 (fr)
FR (1) FR3118251B1 (fr)
WO (1) WO2022128465A1 (fr)

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8922210B2 (en) * 2011-03-31 2014-12-30 General Electric Company Method and apparatus for performing diffusion spectrum imaging
US9208587B2 (en) * 2014-04-25 2015-12-08 General Electric Company Systems and methods for compressed sensing for multi-shell magnetic resonance imaging
WO2016187148A1 (fr) * 2015-05-15 2016-11-24 New York University Système, procédé et support accessible par ordinateur pour l'estimation de bruit, la suppression de bruit et la suppression de phénomène de gibbs
CN109993105A (zh) * 2019-03-29 2019-07-09 北京化工大学 一种改进的自适应稀疏采样故障分类方法

Also Published As

Publication number Publication date
FR3118251B1 (fr) 2022-12-09
FR3118251A1 (fr) 2022-06-24
WO2022128465A1 (fr) 2022-06-23

Similar Documents

Publication Publication Date Title
Ongie et al. Deep learning techniques for inverse problems in imaging
Darestani et al. Accelerated MRI with un-trained neural networks
WO2019128660A1 (fr) Procédé et dispositif d'entraînement de réseau neuronal, procédé de traitement d'image, et dispositif et support de stockage
US9542761B2 (en) Generalized approximate message passing algorithms for sparse magnetic resonance imaging reconstruction
US9801591B2 (en) Fast iterative algorithm for superresolving computed tomography with missing data
Zhang et al. Energy preserved sampling for compressed sensing MRI
JP6573538B2 (ja) 信号を再構成する方法および装置
US9147267B2 (en) Reconstruction of image data
US8897529B2 (en) Apparatus, system, and method for non-convex prior image constrained compressed sensing
US11941787B2 (en) Denoising depth data of low-signal pixels
WO2019115381A1 (fr) Procédé, dispositif et programme de traitement d'images de diffraction d'un matériau cristallin
WO2022179386A1 (fr) Procédé, appareil et système d'imagerie quantique distribuée, et support de stockage lisible par ordinateur
Rai et al. An unsupervised deep learning framework for medical image denoising
Barbieri et al. Circumventing the curse of dimensionality in magnetic resonance fingerprinting through a deep learning approach
Blocker et al. Low-rank plus sparse tensor models for light-field reconstruction from focal stack data
EP4264553A1 (fr) Methode de reconstruction d'une image a partir d'une acquisition comprimee de mesures
JP2023525468A (ja) 複屈折測定値のバックグラウンド補正
JP2023525464A (ja) 複屈折データの効率的な読み取り
US20180365810A1 (en) Object image recovery from digital holograms
Faisan et al. Joint filtering estimation of Stokes vector images based on a nonlocal means approach
Burger et al. Reconstruction methods in thz single-pixel imaging
Thiébaut Image reconstruction with optical interferometers
Gao et al. Positron emission tomography image reconstruction using feature extraction
Dinh et al. A minimalist approach to 3D photoemission orbital tomography: algorithms and data requirements
Paleo Iterative Methods in regularized tomographic reconstruction

Legal Events

Date Code Title Description
STAA Information on the status of an ep patent application or granted ep patent

Free format text: STATUS: UNKNOWN

STAA Information on the status of an ep patent application or granted ep patent

Free format text: STATUS: THE INTERNATIONAL PUBLICATION HAS BEEN MADE

PUAI Public reference made under article 153(3) epc to a published international application that has entered the european phase

Free format text: ORIGINAL CODE: 0009012

STAA Information on the status of an ep patent application or granted ep patent

Free format text: STATUS: REQUEST FOR EXAMINATION WAS MADE

17P Request for examination filed

Effective date: 20230623

AK Designated contracting states

Kind code of ref document: A1

Designated state(s): AL AT BE BG CH CY CZ DE DK EE ES FI FR GB GR HR HU IE IS IT LI LT LU LV MC MK MT NL NO PL PT RO RS SE SI SK SM TR

DAV Request for validation of the european patent (deleted)
DAX Request for extension of the european patent (deleted)